離散フーリエ変換(2) - 類似度と相互相関関数

ある信号Aとある信号Bがどれくらい似ているかを表す指標のひとつに相関関数があります。

相互相関関数(そうごそうかんかんすう、英: Cross correlation function)は、ふたつの信号、配列(ベクトル)の類似性を確認するために使われる。関数の配列の結果がすべて1であれば相関があり、すべてゼロであれば無相関であり、すべて-1であれば負の相関がある。(Wikipedia)

"関数の配列の結果"が何を意味するのか分からないのでなんともいえませんが、2つ目の文はすべて1とかのあたりがなんか間違ってるような気がしますね。

連続な定義域無限の関数だと以下です。
 (f*g)(\tau) = \int^{\infty}_{-\infty}f^{*}(t)g(t+\tau)dt

\tauはfとgをどれだけずらすかをあらわすパラメタです。

無限/連続が入ってくるとややこしいので有限長/離散にします。

長さnの二つの複素数 f_0, \dots, f_{n-1}と、 g_0, \dots, g_{n-1}に対しては、

 C_{\tau} = \sum^{n - 1}_{t = 0} f^{*}_t g_{t+\tau} (0 \le \tau < n - 1 ... (*)

となります。

”おい!てめぇ今さっき長さnっつったばっかなのにgの添え字が2n-1になり得るじゃねーか!”

と、僕も思いました。

前回同様fもgも長さnの周期で無限に続いているという仮定があるようです。

英語版Wikipediaだと

 This is also known as a sliding dot product or sliding inner-product.

と書かれているように、片方の添え字を\tauだけずらしながら内積を計算しているのと同じことになります。

同じ長さの2つの信号があったとき、それらが似ていればその内積(の絶対値)が大きく、似ていなければ小さくなると仮定しましょう。

この仮定が正しいとすれば、長い信号A(長さL)のなかのどこに短い信号B(長さS)かを知るためには、各 0 \le \tau \le L - S - 1 に対して、

 C_{\tau} = \sum^{S-1}_{s = 0} A^{*}_{s + \tau}B_s

を計算して、C_{\tau}がもっとも大きくなった位置\tauが、Bがある可能性が一番高い場所だということになります。

上式は信号Bの後ろに長さL-Sの0をくっつけて(*)式を計算していることと同じことです。

これを素直に計算して見ましょう。

#include <iostream>
#include <vector>

using namespace std;

vector<double> CrossCorr(const vector<double> &a, const vector<double> &b) {
  vector<double> result;
  for(int i = 0; i < a.size() - b.size(); ++i) {
    double C_i = 0.0;
    for(int j = 0; j < b.size(); ++j) {
      C_i += a.at(i+j) * b.at(j);
    }
    result.push_back(C_i);
  }
  return result;
}

int main() {
  vector<double> orig;
  for(int i = 0; i < 30; ++i) {
    if(i >= 10 && i < 20) 
      orig.push_back(static_cast<double>(i % 10));
    else
      orig.push_back(2.0);
  }

  vector<double> tmpl;
  for(int i = 0; i < 10; ++i) {
    tmpl.push_back(static_cast<double>(i));
  }

  vector<double> result = CrossCorr(orig, tmpl);

  for(size_t i = 0; i < result.size(); ++i) {
    cout << i << " : " << result.at(i) << endl;
  }

  return 0;
}

結果

0 : 90
1 : 72
2 : 65
3 : 68
4 : 80
5 : 100
6 : 127
7 : 160
8 : 198
9 : 240
10 : 285
11 : 258
12 : 230
13 : 202
14 : 175
15 : 150
16 : 128
17 : 110
18 : 97
19 : 90

aは長さ30で10~19番目の要素が0, 1, ..., 9, そのほかの要素が2の数列、bは長さ10で0~9番目が0, 1, ..., 9の数列です。

aの10から19がbの0から9と一致しているとき、aの10から19を表す10で結果が最大値をとっていることから、なんとなく類似度としては役に立っていそうですが、相互相関関数には重大な欠点があります。

aの20番目から29番目の要素を2ではなく1000にすると結果は以下のようになります。

0 : 90
1 : 72
2 : 65
3 : 68
4 : 80
5 : 100
6 : 127
7 : 160
8 : 198
9 : 240
10 : 285
11 : 9240
12 : 17196
13 : 24154
14 : 30115
15 : 35080
16 : 39050
17 : 42026
18 : 44009
19 : 45000

全く似ていなくても、一方の信号の絶対値が大きいとそれに引きずられてしまうのです。

じゃあだめじゃねえか、とも確かにいえますが、

- 実用上、データの一部分のみが極端に大きかったり小さかったりしないから引きずられることがあんまりない、って言う分野が存在する
- 次回以降に扱う、畳み込み定理とFFTを使用すると計算が簡単かつ高速に終わる

という性質から、一概に切り捨てて良いという話にはならないのです。

次回はこの相互相関関数とDFTの関係を取り上げます。