任何時間序列都可以表示為在基波（諧波）頻率上振蕩的余弦和正弦波之和=?j / n，其中j?= 1,2，...，n?/ 2。周期圖給出了關于各種頻率的相對強度的信息，用于解釋時間序列的變化。

周期圖是稱為譜密度的群體函數(shù)的樣本估計，其是群體平穩(wěn)時間序列的頻域表征。譜密度是與自協(xié)方差時域表示直接相關的時間序列的頻域表示。本質上，譜密度和自協(xié)方差函數(shù)包含相同的信息，但以不同的方式表達。

[?回顧注釋：自協(xié)方差是自相關的分子。自相關是自協(xié)方差除以方差。]

假設γ（h）是靜止過程的自協(xié)方差函數(shù)，f（ω）是同一過程的譜密度。在前一句的符號中，h?=時滯和ω=頻率。

?

在高級微積分的語言中，自協(xié)方差和譜密度是傅立葉變換對。我們不會擔心這種情況的微積分。我們將重點關注譜密度的估計 - 一系列的頻域特征。這里僅給出傅里葉變換方程以確定在時域表示和序列的頻域表示之間存在直接鏈接。

在數(shù)學上，頻譜密度是針對負頻率和正頻率定義的。然而，由于函數(shù)的對稱性及其對于-1/2到+1/2范圍之外的頻率的重復模式，我們只需關注0和+1/2之間的頻率。

“總”積分譜密度等于系列的方差。因此，特定頻率區(qū)間內的譜密度可以被視為由那些頻率解釋的方差的量。

估計光譜密度的方法

原始周期圖是人口譜密度的粗略樣本估計。估計是“粗略的”，部分是因為我們只使用離散的基波諧波頻率用于周期圖，而頻譜密度是在連續(xù)的頻率上定義的。

對譜密度的周期圖估計的一種可能的改進是使用居中移動平均值來平滑它。可以使用逐漸減小的方法創(chuàng)建另外的“平滑”，該方法對系列的末端（時間）加權小于數(shù)據(jù)的中心。我們不會在本課中介紹逐漸減少的內容。有興趣的人可以參閱本書第4.4節(jié)和各種互聯(lián)網(wǎng)資源。

平滑周期圖的另一種方法是基于以下事實的參數(shù)估計方法：任何靜止時間序列可以通過某種順序的AR模型來近似（盡管它可能是高階）。在該方法中，找到合適的AR模型，然后將譜密度估計為該估計的AR模型的譜密度。

平滑方法（光譜密度的非參數(shù)估計）

平滑周期圖的常用方法具有如此精美的名稱，聽起來很難。實際上，它只是一個集中的移動平均程序，只有一些可能的修改。對于時間序列，具有參數(shù)m的Daniell核是居中移動平均值，其通過平均時間t-m和t + m（包括）之間的所有值在時間t創(chuàng)建平滑值。例如，m?= 2?的Daniell內核的平滑公式為

^ x t= x t - 2 + x t - 1 + x t + x t + 1 + x t + 2 5x^t=xt?2+xt?1+xt+xt+1+xt+25

在R中，可以使用命令kernel（“daniell”，2）生成m?= 2?的Daniell內核的加權系數(shù)。? 結果是

coef [-2] = 0.2? coef [-1] = 0.2? coef [0] = 0.2? coef [1] = 0.2? coef [2] = 0.2

coef []的下標是指與時間t的平均值中心的時差。因此，在這種情況下的平滑公式是

^ x t=0.2xt-2+0.2xt-1+0.2xt+0.2xt+1+0.2xt+2，x^t=0.2xt?2+0.2xt?1+0.2xt+0.2xt+1+0.2xt+2,

這與上面給出的公式相同。

修改后的Daniell內核使得平均值中的兩個端點接收內部點的重量的一半。對于m?= 2?的修改后的Daniell內核，平滑是

^ x t= x t - 2 +2 x t - 1 +2 x t +2 x t + 1 + x t + 2 8= 0.125 x t - 2 + 0.25 x t - 1 + 0.25 x t + 0.25 x t + 1 + 0.125 x t + 2x^t=xt?2+2xt?1+2xt+2xt+1+xt+28=0.125xt?2+0.25xt?1+0.25xt+0.25xt+1+0.125xt+2

在R中，命令內核（“modified.daniell”，2）將列出剛剛使用的加權系數(shù)。

可以對Daniell內核或修改后的Daniell內核進行卷積（重復），以便將平滑再次應用于平滑值。通過在更寬的時間間隔內取平均值，可以產(chǎn)生更廣泛的平滑效果。例如，為了重復丹尼爾內核米上起因于一個丹尼爾內核與平滑值= 2?米?= 2時，公式將是

^ ^ x t= ^ x t - 2 + ^ x t - 1 + ^ x t + ^ x t + 1 + ^ x t + 2 5x^^t=x^t?2+x^t?1+x^t+x^t+1+x^t+25

這是在任一方向上的兩個時間段t內的平滑值的平均值。

在R中，命令內核（“daniell”，c（2,2））將提供系數(shù)，這些系數(shù)將作為權重應用于對兩個平滑中m?= 2?的回旋Daniell內核的原始數(shù)據(jù)值求平均值。結果是

?coef [-4] = 0.04? coef [-3] = 0.08? coef [-2] = 0.12? coef [-1] = 0.16? coef [0] = 0.20? coef [ 1] = 0.16? coef [2] = 0.12? coef [3] = 0.08? coef [4] = 0.04

這會生成平滑公式

^ x t=0.04xt-4+0.08xt-3+0.12xt-2+0.16xt-1+0.20xt+0.16xt+1+0.12xt+2+0.08xt+3+0.04xt+4。x^t=0.04xt?4+0.08xt?3+0.12xt?2+0.16xt?1+0.20xt+0.16xt+1+0.12xt+2+0.08xt+3+0.04xt+4.

其中端點具有較小重量的修改方法的卷積也是可能的。命令內核（“modified.daniell”，c（2,2））?給出了這些系數(shù)：

coef [-4] = 0.01563? coef [-3] = 0.06250? coef [-2] = 0.12500? coef [-1] = 0.18750? coef [0] = 0.21875? coef [1] = 0.18750? coef [2] = 0.12500? coef [3] = 0.06250? coef [4] = 0.01563

因此，中心值的加權比未修改的Daniell內核略重。

當我們平滑周期圖時，我們在頻率間隔而不是時間間隔上進行平滑。請記住，周期圖是在基本頻率ω確定??=??/ N為??= 1，2，...，??/ 2。讓我（ω???）表示在頻率ω的周期圖值??=??/ N。當我們使用帶參數(shù)m的Daniell內核來平滑周期圖時，平滑值\（\ hat {I}（\ omega_j）\）是范圍（jm）/ n到（j）中頻率的周期圖值的加權平均值+ m）/ n。^ 我（ω?）I^(ωj)

帶寬

?

帶寬應該足以平滑我們的估計，但是如果我們使用太大的帶寬，我們將過多地平滑周期圖并且錯過看到重要的峰值。在實踐中，通常需要一些實驗來找到提供合適平滑的帶寬。

帶寬主要由平滑中平均值的數(shù)量控制。換句話說，Daniell內核的m參數(shù)以及內核是否被卷積（重復）會影響帶寬。

注意：?? 帶有圖表的帶寬R報告與使用上述公式計算的值不匹配。請參閱第12頁的腳注。190您的文字作為解釋。

R代碼

使用Daniell內核對周期圖進行平均/平滑可以使用兩個命令的序列在R中完成。第一個定義了Daniell內核，第二個定義了平滑周期圖。

例如，假設觀察到的序列被命名為x，我們希望使用具有m?= 4的Daniell內核來平滑周期圖。命令是

?mvspec（x，k，log =“no”）

第一個命令創(chuàng)建平滑所需的加權系數(shù)，并將它們存儲在名為k的向量中。（將它稱為k是任意的。它可以被稱為任何東西。）第二個命令要求基于系列x的周期圖的譜密度估計，使用存儲在k中的加權系數(shù)，該圖將是普通的比例，不是對數(shù)刻度。

如果需要卷積，則可以將內核命令修改為類似k = kernel（“daniell”，c（4,4））的內容。

有兩種方法可以實現(xiàn)修改后的Daniell內核。您可以更改內核命令以引用“modified.daniell”而不是“daniell”，也可以跳過使用內核命令并在mvspec命令中使用spans參數(shù)。

spans參數(shù)給出了所需修改的Daniell內核的長度（= 2?m?+1）。例如，m?= 4?的修改后的Daniell內核長度L?= 2?m?+1 = 9，因此我們可以使用該命令

mvspec（x，spans = 9，log =“no”）可以使用每次通過m?= 4?的修改的Daniell內核的兩次通過mvspec（x，spans = c（9,9），log =“no”）

示例：此示例將使用文本中多個位置使用的魚類招募系列，包括第4章中的幾個位置。該系列包含n = 453個月度值，用于衡量南半球位置的魚群數(shù)量。數(shù)據(jù)位于文件recruit.dat中。

可以使用命令創(chuàng)建原始周期圖（或者可以使用第6課中給出的方法創(chuàng)建原始周期圖）。

mvspec（x，log =“no”）

請注意，在剛剛給出的命令中，我們省略了為平滑提供權重的參數(shù)。

原始周期圖如下：

下一個圖是使用具有m?= 4?的Daniell核的平滑周期圖。注意，平滑的一個效果是未平滑版本中的主峰現(xiàn)在是第二高峰。之所以發(fā)生這種情況是因為峰值在未平滑的版本中如此清晰地定義，當我們用幾個周圍值對其進行平均時，高度會降低。

?

下一個圖是一個平滑的周期圖，使用Daniell內核的兩次傳遞，每次傳遞m?= 4。請注意它是如何比以前更平滑。

?

要了解兩個主要峰的位置，請為mvspec輸出指定一個名稱，然后列出它。例如，

specvalues = mvspec（x，k，log =“no”） specvalues $ spec

您可以篩選輸出以查找峰值出現(xiàn)的頻率。頻率和頻譜密度估計單獨列出，但順序相同。確定最大光譜密度，然后找到相應的頻率。

這里，第一個峰值的頻率≈0208。與此周期相關的周期（月數(shù)）= 1 / .0208 = 48個月。第二個峰值出現(xiàn)在頻率≈0.083333。相關時期= 1 / .08333 = 12個月。第一個峰值與厄爾尼諾天氣效應有關。第二個是通常12個月的季節(jié)性影響。

這兩個命令將垂直虛線放在峰值密度的近似位置處的（估計的）譜密度圖上。

abline（v = 1/44，lty =“dotted”） abline（v = 1/12，lty =“dotted”）

這是結果圖：

?

我們已經(jīng)足夠平滑，但出于演示目的，下一個情節(jié)是結果

mvspec（x，spans = c（13,13），log =“no”）

這使用兩次修改的Daniell內核，每次長度L?= 13（所以m?= 6）。情節(jié)有點平滑，但不是很多。順便說一句，峰值與上面的圖中完全相同。

?

絕對有可能過度平滑。假設我們使用總長度= 73（m?= 36）的修改過的Daniell內核。命令是

mvspec（x，spans = 73，log =“no”）

結果如下。高峰消失了！

?

光譜密度的參數(shù)估計

譜密度估計的平滑方法稱為非參數(shù)方法，因為它不使用任何參數(shù)模型用于基礎時間序列過程。另一種方法是參數(shù)方法，該方法需要找到該系列的最佳擬合AR模型，然后繪制該模型的譜密度。該方法由一個定理支持，該定理表明任何時間序列過程的譜密度可以通過AR模型的譜密度（某種順序，可能是高順序）來近似。

在R中，使用命令/功能spec.ar可以輕松完成譜密度的參數(shù)估計。像spec.ar（x，log =“no”）這樣的命令將導致R完成所有工作。同樣，為了識別峰值，我們可以通過執(zhí)行specvalues = spec.ar（x，log =“no”）之類的操作為spec.ar結果指定名稱。

對于魚類補充實例，結果如下。注意，繪制的密度是AR（13）模型的密度。我們當然可以為這些數(shù)據(jù)找到更簡約的ARIMA模型。我們只是使用該模型的譜密度來近似觀察到的序列的譜密度。

?

估計的譜密度的出現(xiàn)與以前相同。估計的厄爾尼諾峰位于略微不同的位置 - 頻率約為0.024，周期約為1 / 0.024 =約42個月。

德趨勢

在進行光譜分析之前，應該對一系列進行去趨勢分析。趨勢將導致在低頻率下的這種主要譜密度，從而不會看到其他峰值。默認情況下，R命令mvspec使用線性趨勢模型執(zhí)行去趨勢。也就是說，使用來自回歸的殘差來估計譜密度，其中y變量=觀測數(shù)據(jù)和x變量=?t。如果存在不同類型的趨勢，例如二次方，則可以使用多項式回歸來在探索估計的譜密度之前使數(shù)據(jù)去趨勢化。需要注意的是R命令spec.ar，但默認情況下不進行解趨勢。

平滑器在原始數(shù)據(jù)中的應用

請注意，此處描述的平滑器也可以應用于原始數(shù)據(jù)。Daniell內核及其修改只是移動平均（或加權移動平均）平滑器。

?

?

還有問題嗎？聯(lián)系我們！

?

大數(shù)據(jù)部落 -中國專業(yè)的第三方數(shù)據(jù)服務提供商，提供定制化的一站式數(shù)據(jù)挖掘和統(tǒng)計分析咨詢服務

統(tǒng)計分析和數(shù)據(jù)挖掘咨詢服務：y0.cn/teradat（咨詢服務請聯(lián)系官網(wǎng)客服）

?QQ：3025393450

?

【服務場景】??

科研項目; 公司項目外包;線上線下一對一培訓;數(shù)據(jù)采集;學術研究;報告撰寫;市場調查。

【大數(shù)據(jù)部落】提供定制化的一站式數(shù)據(jù)挖掘和統(tǒng)計分析咨詢服務

?

轉載于:https://www.cnblogs.com/tecdat/p/10524347.html

總結

以上是生活随笔為你收集整理的用R语言软件估计光谱密度的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內容還不錯，歡迎將生活随笔推薦給好友。