自己相関係数のプログラム

自己相関係数を求めるためのプログラムです。

#include <iostream>
#include <fstream>
#include <cstdlib>
#include <cmath>
using namespace std;
void self_corr (double* x, int imax,double* r){
    // 平均値の算定
    double x_ave = 0.0;
    for(int i = 0; i < imax; i++){
	x_ave+=x[i];
    }
    x_ave/=(double)imax;
    // 平均値を差し引く
    double *y = new double[imax];
    for(int i = 0; i < imax; i++){
	y[i] = x[i]-x_ave;
    }
    // 自己相関係数のうち、分母denomの算定
    double denom=0.0;		// 分母denom
    for(int i = 0; i< imax; i++){
	denom=denom+y[i]*y[i];
    }
    // 自己相関係数のうち、分子numerおよび自己相関係数rの算定
    for(int h = 0; h <imax; h++){
	double numer = 0.0;	// 分子
	for(int i = 0; i < imax-h; i++){
	    numer=numer+y[i]*y[i+h];
	}
	r[h] = numer/denom;	// 自己相関係数
    }
}
//
int main(void) {
    // 入力ファイル
    fstream ifs("in.txt");
    // 時系列の個数
    int imax = 7306;
    // 時系列データxと自己相関係数rの配列確保
    double* x = new double[imax];
    double* r = new double[imax];
    // 入力ファイルからのデータ入力
    for(int i = 0; i < imax; i++) {
	ifs >> x[i];
    }
    // 自己相関係数の算定
    self_corr ( x, imax, r);
    // 自己相関係数の出力
    for(int i = 0; i < imax; i++) {
	cout << r[i] << endl;
    }
}

気象庁のホームページから、日平均気温を引っ張ってきて、
自己相関係数を求めてみました。


気象庁|過去の気象データ検索

日平均気温であれば、365日周期で、相関係数が1に近づきます。
(当たり前ですね。一年前との気温は概ね同じはず。)
ただし、ラグが大きくなればなるほど、周期はそのままに、相関係数は小さくなります。
このように当たり前のことであれば、自己相関を取るまでもなく、その特徴がわかるかと思います。

さらに踏み込むと、二変数に対して相関係数をとる方法に相互相関係数があります。
たとえば、雨と河川流量とかです。なんとなく、雨が降ると流量は増えそうだなと、直感でわかりますが、その時間差は?と聞かれると、不明です。

次回はぜひそれについて。

f:id:nishiru3:20141105184403j:plain