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

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

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

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

?

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

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

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

估計(jì)光譜密度的方法

原始周期圖是人口譜密度的粗略樣本估計(jì)。估計(jì)是“粗略的”，部分是因?yàn)槲覀冎皇褂秒x散的基波諧波頻率用于周期圖，而頻譜密度是在連續(xù)的頻率上定義的。

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

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

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

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

^ 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內(nèi)核的加權(quán)系數(shù)。? 結(jié)果是

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

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

^ 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內(nèi)核使得平均值中的兩個(gè)端點(diǎn)接收內(nèi)部點(diǎn)的重量的一半。對(duì)于m?= 2?的修改后的Daniell內(nèi)核，平滑是

^ 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中，命令內(nèi)核（“modified.daniell”，2）將列出剛剛使用的加權(quán)系數(shù)。

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

^ ^ 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

這是在任一方向上的兩個(gè)時(shí)間段t內(nèi)的平滑值的平均值。

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

?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

這會(huì)生成平滑公式

^ 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.

其中端點(diǎn)具有較小重量的修改方法的卷積也是可能的。命令內(nèi)核（“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

因此，中心值的加權(quán)比未修改的Daniell內(nèi)核略重。

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

帶寬

?

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

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

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

R代碼

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

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

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

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

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

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

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

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

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

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

mvspec（x，log =“no”）

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

原始周期圖如下：

下一個(gè)圖是使用具有m?= 4?的Daniell核的平滑周期圖。注意，平滑的一個(gè)效果是未平滑版本中的主峰現(xiàn)在是第二高峰。之所以發(fā)生這種情況是因?yàn)榉逯翟谖雌交陌姹局腥绱饲逦囟x，當(dāng)我們用幾個(gè)周圍值對(duì)其進(jìn)行平均時(shí)，高度會(huì)降低。

?

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

?

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

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

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

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

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

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

這是結(jié)果圖：

?

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

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

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

?

絕對(duì)有可能過度平滑。假設(shè)我們使用總長(zhǎng)度= 73（m?= 36）的修改過的Daniell內(nèi)核。命令是

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

結(jié)果如下。高峰消失了！

?

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

譜密度估計(jì)的平滑方法稱為非參數(shù)方法，因?yàn)樗皇褂萌魏螀?shù)模型用于基礎(chǔ)時(shí)間序列過程。另一種方法是參數(shù)方法，該方法需要找到該系列的最佳擬合AR模型，然后繪制該模型的譜密度。該方法由一個(gè)定理支持，該定理表明任何時(shí)間序列過程的譜密度可以通過AR模型的譜密度（某種順序，可能是高順序）來近似。

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

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

?

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

德趨勢(shì)

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

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

請(qǐng)注意，此處描述的平滑器也可以應(yīng)用于原始數(shù)據(jù)。Daniell內(nèi)核及其修改只是移動(dòng)平均（或加權(quán)移動(dòng)平均）平滑器。

?

?

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

?

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

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

?QQ：3025393450

?

【服務(wù)場(chǎng)景】??

科研項(xiàng)目; 公司項(xiàng)目外包;線上線下一對(duì)一培訓(xùn);數(shù)據(jù)采集;學(xué)術(shù)研究;報(bào)告撰寫;市場(chǎng)調(diào)查。

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

?

轉(zhuǎn)載于:https://www.cnblogs.com/tecdat/p/10524347.html

總結(jié)

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

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