長期記憶過程の自己相関関数とパワースペクトルのスケーリング関係：解析的導出

　今回は、長期記憶過程の自己相関関数とパワースペクトルのスケーリング指数の関係について説明します。学生から、導出できないという相談がときどきあるので、やや詳しく書いておきます。連続近似して積分したほうが計算は楽ですが、できるだけ離散の形でで攻める方法も紹介しておきます。

［目標］自己相関がべき乗で減衰するときのパワースペクトルを求める
- Wiener–Khinchineの定理
パワースペクトルの評価
連続近似により積分で計算する方法

　まず、基本事項を整理します．

自己相関関数の定義　離散時間 $\displaystyle{ i = 0,1,2,\dots }$ 上の弱定常確率過程 $\displaystyle{ {x [i]} }$ を考え，平均 $\displaystyle{ \langle x [i] \rangle = 0 }$ ，分散 $\displaystyle{ \langle x [i]^ 2\rangle = \sigma^ 2 }$ を仮定します。ここで、 $\displaystyle{ \langle \, \cdot \, \rangle }$ は、期待値の計算を表します。

　ラグ $\displaystyle{ n }$ のときの自己相関関数を

$\displaystyle{ C(n) \;=\; \langle x [i]\, x [i + n ]\rangle }$

と定義します。定常性の仮定により，自己相関は $\displaystyle{ i }$ に依存せず $\displaystyle{ n }$ だけに依存します。

長期記憶過程 (long memory process)　大きな $\displaystyle{ n }$ に対して

$\displaystyle{ C(n)\;\sim\;n^{-\gamma}, \quad (0 < \gamma < 1) }$

となるとき，この過程を長期記憶過程と呼びます (ここでは、 $\sim$ を比例するの意味で使います)。

このとき， $\displaystyle{ \gamma \in (0,1) }$ の範囲では， $\displaystyle{ \sum _ {n=1}^ {\infty} n^ {-\gamma} }$ が発散するので、

$\displaystyle{ \lim_{N \to \infty} \sum_{n=1}^N C(n) \;=\; \infty }$

となり，過去の影響の総和が無限になる長時間相関があること示しています．

主要項の評価が重要　比例を表す $\sim$ の使い方について補足しておきます。べき指数を推定するときには、両対数プロットの傾きを調べますが、その場合、値が最も大きい項 (主要項)のみに注目すれば十分です。ですので、 $n$ が十分大きいとき ( $n \gg 1$ )、

$\displaystyle{ \begin{aligned} F(n) &= 5 n^{2} +4 n + 3 \\ &\sim n^2 \end{aligned} }$

と書き、べき指数が最も大きい主要項のみを評価します。また、 $f$ が0に近いときは ( $0 \lt f \ll 1/2$ )、

$\displaystyle{ \begin{aligned} G(f) &= 3 f^{-2} + 2 f^{-1} + 1 \\ &\sim f^{-2} \end{aligned} }$

と書き、べき指数が最も小さい主要項のみを評価します。

［目標］自己相関がべき乗で減衰するときのパワースペクトルを求める

　ここでの目標は、 $0 \lt \gamma \lt 1$ のとき、

$\displaystyle{ C(n)\;\sim\;n^{-\gamma} \quad \Longrightarrow \quad S(f)\;\sim\; f^{\gamma - 1} \;=\; f^{-\,\beta} }$

が成り立つことを導き、スケーリング指数に

$\displaystyle{ \beta = 1 - \gamma }$

の関係が成り立つことを示すことです。

Wiener–Khinchineの定理

　Wiener–Khinchine の定理により、自己相関関数とパワースペクトルの関係が与えられているので、計算の流れは明確です。

　つまり、自己相関関数 $\displaystyle{ C(n) }$ とパワースペクトル密度 $\displaystyle{ S(f) }$ は，Wiener–Khinchine の定理によってフーリエ変換ペアの関係にあります．離散時間の場合， $\displaystyle{ -\infty \lt n \lt \infty }$ について

$\displaystyle{ \begin{aligned} S(f) &= \sum_{n=-\infty}^{\infty} C(n)\, e^{-\,i\,2\pi f\,n}\\ &= C(0) + 2\sum_{n=1}^{\infty} C(n)\,\cos(2\pi f\,n) \end{aligned} }$

が成り立ちます。

自己相関関数は、実数の値をとり、軸対称 $\displaystyle{ C(n)=C(-n) }$ なので、フーリエ・コサイン変換を使って、自己相関関数からパワースペクトルを求めることができます。

　では、 $n$ が十分大きい ( $n \gg 1$ )とき、 $\displaystyle{ C(n) \sim n^ {-\gamma} }$ に漸近する場合に，小さなの周波数 $\displaystyle{ f }$ ( $0 \lt f \ll 1/2$ )の領域で、 $\displaystyle{ S(f) }$ はどのようなべき乗則に漸近するでしょうか？この答えを導くことが、今回のゴールです。

パワースペクトルの評価

　Wiener–Khinchine の定理を使えば、自己相関関数がべき的に減衰するときのパワースペクトルは、

$\displaystyle{ \begin{align} S(f) &= C(0) + 2\sum_{n=1}^{\infty} C(n)\,\cos(2\pi f\,n) \\ &= C(0) + 2\sum_{n=1}^{\infty} n^{-\gamma}\,\cos(2\pi f\,n) \\ &\sim \sum_{n=1}^{\infty} n^{-\gamma}\,\cos(2\pi f\,n) \end{align} }$

となりますので、小さな周波数 $\displaystyle{ f }$ ( $\displaystyle{ 0 \lt f \ll 1/2 }$ ) について、

$\displaystyle{ \sum_{n=1}^{\infty} n^{-\gamma}\,\cos(2\pi f\,n)}$

が、どのように振る舞うかを評価すればよいことになります．ここでは，以下のテクニックを使います。

(1) 収束を確保するために，まずは $\displaystyle{ |r|\lt1 }$ を満たすパラメータ $\displaystyle{ r }$ を導入し，

$\displaystyle{ \sum_{n=1}^{\infty} n^{-\gamma}\,\bigl(r e^{i2\pi f}\bigr)^{n} }$

を考えます。この級数は $\displaystyle{ |r|\lt1 }$ なら、

$\displaystyle{ \left| \sum_{n=1}^{\infty} n^{-\gamma}\,\bigl(r e^{i2\pi f}\bigr)^{n} \right| < \left| \sum_{n=1}^{\infty} \bigl(r e^{i2\pi f}\bigr)^{n} \right| < \sum_{n=1}^{\infty} r^{n} }$

なので、絶対収束する形になっています。
(2) そして、 $\displaystyle{ r \to 1 - 0 }$ の極限をとることによって、

$\displaystyle{ \lim_{r \to 1 - 0} \sum_{n=1}^{\infty} n^{-\gamma}\,\bigl(r e^{i2\pi f}\bigr)^{n} = \sum_{n=1}^\infty n^{-\gamma} e^{\,i\,2\pi f\,n} }$

と変形します。さらに、

$\displaystyle{ \begin{align} \sum_{n=1}^\infty n^{-\gamma} e^{\,i\,2\pi f\,n} &= \sum_{n=1}^\infty n^{-\gamma} \cos(2\pi f\,n) + i \sum_{n=1}^\infty n^{-\gamma} \sin(2\pi f\,n) \end{align} }$

ですので、この実部をとれば $\displaystyle{ \cos(2\pi f\,n) }$ の和を得ることができます。

(3) $\displaystyle{ \sum _ {n=1}^ {\infty} n^ {-\gamma} x^ n }$ 型の和は，

$\displaystyle{ \sum_{n=1}^{\infty} n^{\delta} x^n \;\approx\; \Gamma(\delta+1)\,(1-x)^{-\delta-1} \quad (|x|<1, \; {\rm Re} (\delta) > -1) }$

という漸近公式を使って評価出来ます。ここで、 $\displaystyle{ \Gamma(z) }$ はガンマ関数です。後で説明するように、これは、連続近似してガンマ関数の定義式を使う方法や、Taylor 展開やStirling の近似公式を使う方法などで、導くことができます。

　今回の問題では， $\displaystyle{ \delta = -\gamma }$ と置くことで

$\displaystyle{ n^{\delta} = n^{-\gamma} }$

を対応させます．

漸近公式の導出：連続近似 (その1)

　漸近公式

$\displaystyle{ \sum_{n=1}^{\infty} n^{\delta} x^n \;\approx\; \Gamma(\delta+1)\,(1-x)^{-\delta-1} \quad (|x|<1, \; {\rm Re}(\delta) > -1) }$

を積分近似で導出します。

　計算過程を詳しく説明しなければ流れは簡単で、以下の通りです。
(1) 連続近似して、和を積分で近似します。積分にはガンマ関数の定義式を使うだけです。

$\displaystyle{ \sum_{n=1}^\infty n^{\delta} x^n \approx \int_0^\infty t^\delta x^t \, dt = (-\log x)^{-\delta-1} \Gamma(\delta+1) }$

(2) さらに、 $\displaystyle{- \log x}$ を、 $\displaystyle{ x = 1 }$ のまわりで、Taylor 展開します。

$\displaystyle{ -\log x = \sum_{n=1}^{\infty} \frac{1}{n}(1 - x)^n \quad (|1 - x| < 1) }$

$\displaystyle{ x }$ が1に非常に近いとき、

$\displaystyle{ -\log x \approx (1 - x) }$

と近似し、

$\displaystyle{ \begin{aligned} \sum_{n=1}^{\infty} n^{\delta} x^n &\approx (-\log x)^{-\delta-1} \Gamma(\delta+1)\\ &\approx (1-x)^{-\delta-1}\, \Gamma(\delta+1) \end{aligned} }$

が導けます。

計算過程の説明 自分で計算できるとなんだかうれしいし、理解できた気になると思いますので、計算に悩む学生向けの説明を加えておきます。
(1) 離散和を連続の積分で置き換える漸近近似

$\displaystyle{ \sum_{n=1}^\infty n^\delta x^n \approx \int_0^\infty t^\delta x^t\,dt\ }$

がわからなければ、区分求積法について勉強してみてください。

　余談ですが、もっと厳密に計算したければ、Mellin-Barnes積分を使って、

$\displaystyle{ \sum_{n=1}^{\infty} n^\delta x^n = \frac{1}{2\pi i} \int_{c - i\infty}^{c + i\infty} \Gamma(s)\, \zeta(-\delta + s)\, (-\log x)^{-s} \, ds }$

を評価してください (今回の目的には、ここまでは必要ありません)。

(2) 被積分関数の変形は、まず、 $\displaystyle{ x^ t = e^ {t \log x} }$ と変形します (この式の両辺の対数をとれば成り立つことが確認できます)。ここでは、 $\log x \lt 0$ の場合を想定します。

　さらに、 $\displaystyle{ \log x = -\alpha }$ のように置換すると見やすいと思います。

　以上により、

$\displaystyle{ \int_0^\infty t^\delta x^t \,dt = \int_0^\infty t^\delta e^{-\alpha t} \,dt. }$

と変形できます。

(3) ガンマ関数と比較しやすい形になったので、ガンマ関数の定義式に合うように変形します。

　ガンマ関数の定義

$\displaystyle{ \Gamma(z)=\int_0^{\infty} t^{z-1} e^{-t} \mathrm{~d} t }$

を踏まえて、 $\displaystyle{ \alpha t = u }$ と置換します。そうすると、

$\displaystyle{ \int_0^\infty t^\delta e^{-\alpha t}\,dt \;=\; \alpha^{-\delta-1}\, \int_0^\infty u^\delta e^{-u}\,du \;=\; \alpha^{-\delta-1}\,\Gamma(\delta+1). }$

がえられます。

　 $\displaystyle{ \alpha = -\log x }$ を戻せば、

$\displaystyle{ \int_0^\infty t^\delta x^t\,dt = \bigl(-\log x\bigr)^{-\delta-1}\, \Gamma(\delta+1) }$

が導けます。あとは、 $\displaystyle{- \log x}$ を、 $\displaystyle{ x = 1 }$ のまわりで、Taylor 展開するだけです。

漸近公式の導出： Taylor展開＆Stirlingの公式 (その2)

　上の方法とは違う導出を説明します。もう満足した人は、ここは飛ばして、「漸近公式を使ってスケーリング関係を導く」まで進んでください。

　以下の導出方法は、過去に私が読んだ論文を参考にしています。引用は、

Integrated approach to the assessment of long range correlation in time series data | Phys. Rev. E

です。とはいえ、以下の計算は、その論文とちょっと違う方法になってます。また、漸近公式は、この論文に登場する形

$\displaystyle{ \sum_{n=1}^{\infty} n^{\delta} x^n \;\approx\; \left\{(1-x)^{-\delta-1} -1 \right\} \, \Gamma(\delta+1) }$

の方が精度が少し高いです。この式は、上で紹介した厳密な関係式

$\displaystyle{ \sum_{n=1}^{\infty} n^\delta x^n = \frac{1}{2\pi i} \int_{c - i\infty}^{c + i\infty} \Gamma(s)\, \zeta(-\delta + s)\, (-\log x)^{-s} \, ds }$

の主要項を評価すれば導くことができます (上の論文ではそのように説明していませんが)。しかし、漸近的な振る舞いを議論する場合に、定数の補正は本質的でないし、結局、値もほとんど同じなどで、ここでは、この形を採用しません。

　ということで、以下の説明は論文を参考にしつつも、私流の違った方法を使っています。

Taylor 展開 により、

$\displaystyle{ (1-x)^{-\delta-1} \;=\; \sum_{n=0}^{\infty} \frac{(\delta+1)(\delta+2)\cdots(\delta+n)}{n!}\, x^n }$

となります。今回の導出のキーとなる変形は，

$\displaystyle{ \frac{(\delta+1)(\delta+2)\cdots(\delta+n)}{n!} \;\approx\; \frac{n^{\delta}}{\Gamma(\delta+1)} \quad (n\to\infty) }$

と近似することです。

　この近似式は、ガンマ関数の定義の一つとして知られている

$\displaystyle{ \Gamma(z)=\lim _{n \rightarrow \infty} \frac{n^z n!}{\prod_{k=0}^n(z+k)} }$

を、 $n$ が十分大きいときの近似式として使えば導けます。

　あるいは、公式

$\displaystyle{ \Gamma(z+1)=z \Gamma(z) }$

を使えば、

$\displaystyle{ \begin{align} \Gamma(\delta + n +1) &= (\delta + n) \Gamma(\delta + n) \\ &= (\delta + n)(\delta + n - 1) \Gamma(\delta + n - 1) \\ &= (\delta + n)\cdots(\delta + 1)\Gamma(\delta + 1) \end{align} }$

となるので、

$\displaystyle{ (\delta+1)(\delta+2) \cdots(\delta+n)=\frac{\Gamma(n+\delta+1)}{\Gamma(\delta+1)} }$

です。さらに、公式

$\displaystyle{ n!=\Gamma(n+1) }$

を使えば、

$\displaystyle{ \begin{align} \frac{(\delta+1)(\delta+2) \cdots(\delta+n)}{n!} &=\frac{\Gamma(n+\delta+1)}{\Gamma(\delta+1)} \cdot \frac{1}{\Gamma(n+1)}\\ & =\frac{\Gamma(n+\delta+1)}{\Gamma(n+1)} \cdot \frac{1}{\Gamma(\delta+1)} \end{align} }$

と変形できます。この

$\displaystyle{ \frac{\Gamma(n+\delta+1)}{\Gamma(n+1)} }$

の部分に対して、Stirling の近似公式

$\displaystyle{ n! = \Gamma(n + 1) \approx \sqrt{2 \pi} n^{n+\frac{1}{2}} e^{-n} }$

を使うことで、

$\displaystyle{ \begin{align} \frac{\Gamma(n+\delta+1)}{\Gamma(n+1)} &\approx \frac{\sqrt{2 \pi}(n+\delta)^{n+\delta+\frac{1}{2}} e^{-(n+\delta)}}{\sqrt{2 \pi} n^{n+\frac{1}{2}} e^{-n}} \\ &=\left(\frac{n+\delta}{n}\right)^{n+\frac{1}{2}}(n+\delta)^\delta e^{-\delta} \\ &= \left(1+\frac{\delta}{n}\right)^n\left(1+\frac{\delta}{n}\right)^{\frac{1}{2}}(n+\delta)^\delta e^{-\delta} \\ &\approx e^{\delta} \cdot 1 \cdot (n+\delta)^\delta e^{-\delta} \quad (n \gg 1) \\ &= e^{\delta} e^{-\delta} (n+\delta)^\delta \\ &\approx n^\delta \end{align} }$

となります。

　以上から、

$\displaystyle{ \begin{aligned} \frac{(\delta+1)(\delta+2) \cdots(\delta+n)}{n!} & =\frac{\Gamma(n+\delta+1)}{\Gamma(n+1)} \cdot \frac{1}{\Gamma(\delta+1)}\\ &\approx \frac{n^\delta}{\Gamma(\delta+1)} \end{aligned} }$

と近似できるので、

$\displaystyle{ (1-x)^{-\delta-1} \;\approx\; \sum_{n=0}^{\infty} \frac{n^\delta}{\Gamma(\delta+1)}\,x^n, \quad (|x|<1,\;n\to\infty). }$

となります。

したがって，

$\displaystyle{ \sum_{n=0}^{\infty} n^\delta\,x^n \;\approx\; \Gamma(\delta+1)\,(1-x)^{-\delta-1}. }$

と見積れるので、 $\displaystyle{ \delta = -\gamma }$ を代入すれば，

$\displaystyle{ \sum_{n=1}^\infty n^{-\gamma}\,x^n \;\approx\; \Gamma(1-\gamma)\,(1-x)^{\gamma -1}. }$

となります。 $\displaystyle{ n=0 }$ のとき、値は0なので、 $\displaystyle{ n=1 }$ から和の計算をはじめるようにしました。

　以上で準備は完了です。

漸近公式を使ってスケーリング関係を導く

上で導いた漸近公式、

$\displaystyle{ \sum_{n=1}^\infty n^{-\gamma}\,x^n \;\approx\; \Gamma(1-\gamma)\,(1-x)^{\gamma -1}. }$

に、 $\displaystyle{ x = r e^ {\,i\,2\pi f} }$ を代入して、

$\displaystyle{ \sum_{n=1}^{\infty} n^{-\gamma}\,\bigl(r e^{\,i\,2\pi f}\bigr)^n \;\approx\; \Gamma(1-\gamma)\,\bigl(1 - r e^{\,i\,2\pi f}\bigr)^{\gamma-1} }$

と書けます。

　 $0 \lt f \ll 1/2$ のとき (0に近いとき)，

$\displaystyle{ 1 - r e^{\,i\,2\pi f} \;\approx\; 1 - \Bigl(r\bigl(1 + i\,2\pi f\bigr)\Bigr) \;=\; 1 - r - i\,2\pi r\,f. }$

さらに $\displaystyle{ r\to 1-0 }$ とすると，

$\displaystyle{ 1 - r \;\to\; 0, \quad -\,i\,2\pi r\,f \;\to\; -\,i\,2\pi f. }$

なので、

$\displaystyle{ 1 - r e^{\,i\,2\pi f} \;\approx\; -\,i\,2\pi f }$

と近似できます。

　よって、

$\displaystyle{ \bigl(1 - r e^{\,i\,2\pi f}\bigr)^{\gamma-1} \;\approx\; \bigl(-\,i\,2\pi f\bigr)^{\gamma-1} \;=\; (2\pi f)^{\gamma-1} \,e^{-\,i\,\frac{\pi}{2}(\gamma-1)}. }$

ここから実部を取りだすと， $\displaystyle{ \cos(2\pi f\,n) }$ に対応する和が得られます。

以上をまとめると， $0 \lt f \ll 1/2$ で

$\displaystyle{ \begin{aligned} \sum_{n=1}^{\infty} n^{-\gamma}\,\cos(2\pi f\,n) &\approx {\rm Re} \Bigl[ \Gamma(1-\gamma)\,(2\pi f)^{\gamma-1} \,e^{-\,i\,\frac{\pi}{2}(\gamma-1)} \Bigr] \\ &= \Gamma(1-\gamma)\,(2\pi f)^{\gamma-1} \cos \left(\frac{\pi}{2}(\gamma-1) \right)\\ &= (2\pi)^{\gamma-1} \, \Gamma(1-\gamma) \sin \left(\frac{\pi \gamma}{2} \right)\, f^{\gamma-1} \\ &\sim f^{\gamma-1}\\ &= \,f^{-(1-\gamma)} \end{aligned} }$

となります。

　したがって、正の周波数 $\displaystyle{ f }$ が0に近い領域で、パワースペクトルが $\displaystyle{ f^ {-\,\beta} }$ 型になるとすれば、

$\displaystyle{ \beta = 1-\gamma }$

が成り立ちます。

連続近似により積分で計算する方法

　上記のような、離散時系列であることにこだわった方法以外にも，

$\displaystyle{ S(f) \sim \sum_{n=1}^{\infty} n^{-\gamma} \cos(2\pi f n) \approx \int n^{-\gamma} \cos(2\pi f n)\,dn }$

と近似して計算する方法もあります。

　いちいち公式を導ていると話が長くなるので、ここでは「岩波数学公式I 微分積分・平面曲線」 (岩波書店、1987)の、255ページにのっている公式、

$\displaystyle{ \int_0^{\infty} x^{a-1} \cos b x\, d x=\frac{\Gamma(a)}{b^a}\, \cos \frac{a \pi}{2} \quad [0 < a < 1, b > 0] }$

を使います。これは、ガンマ関数の定義式

$\displaystyle{ \Gamma(z) = \int_0^{\infty} t^{z-1} e^{-t} \, dt \quad (\Re(z) > 0) }$

を使って、変数変換すれば、

$\displaystyle{ \begin{aligned} \int_0^{\infty} x^{a-1} e^{i b x}\, d x &= \frac{\Gamma(a)}{(-i b)^a}\\ &= \frac{\Gamma(a)}{b^a}\left[\cos \left(\frac{a \pi}{2}\right)+i \sin \left(\frac{a \pi}{2}\right)\right] \end{aligned} }$

となるので、この実部

$\displaystyle{ \frac{\Gamma(a)}{b^a} \cos \left(\frac{a \pi}{2}\right) }$

を残せば導けます。説明を省くつもりが、結局、すべて導いてしまいました。

　この公式を使えば、

$\displaystyle{ \begin{aligned} S(f) &\sim \sum_{n=1}^{\infty} n^{-\gamma} \cos(2\pi f n)\\ &\approx \int_0^{\infty} n^{-\gamma} \cos(2\pi f n)\,dn \\ &= \frac{\Gamma(1-\gamma)}{(2 \pi f)^{1-\gamma}} \left\{\cos \frac{\pi }{2}(1-\gamma)\right\} \\ &= (2 \pi)^{\gamma-1} \Gamma(1-\gamma) \left\{\sin \frac{\pi \gamma}{2}\right\} f^{-(1-\gamma)} \\ &\sim f^{-(1-\gamma)} \quad \left(0 < f \ll \frac{1}{2} \right) \end{aligned} }$