dtwパッケージで, 時系列間の距離を出す

内容

Dynamic Time Warpingを使って, 時系列間の距離を出してみます. 2つの時系列データのデータ数や間隔が違っても使えるので, 楽しそうです.

DTWのアルゴリズムは, パラメータを1つ含んでいて, これを調整する必要があります. パラメータの値と結果の関係に関しては, Xi, 2007

DTWを使う

Rのパッケージで, dtwを実行できるものがあります. 使いましょう. Giorgino T, 2009

library(data.table)
library(dplyr)

## 
## Attaching package: 'dplyr'
## 
##  以下のオブジェクトはマスクされています (from 'package:data.table') : 
## 
##      last 
## 
##  以下のオブジェクトはマスクされています (from 'package:stats') : 
## 
##      filter, lag 
## 
##  以下のオブジェクトはマスクされています (from 'package:base') : 
## 
##      intersect, setdiff, setequal, union

library(dtw)

## Loading required package: proxy
## 
## Attaching package: 'proxy'
## 
##  以下のオブジェクトはマスクされています (from 'package:stats') : 
## 
##      as.dist, dist 
## 
## Loaded dtw v1.17-1. See ?dtw for help, citation("dtw") for use in publication.

#citation("dtw")
data("aami3a")
plot(aami3a[1:3000], type='l')

## 切り出してみる. 
ref = aami3a %>% window(start = 0, end = 2)
test= aami3a %>% window(start = 2.7, end = 5)
## 形は似ていますが, 
plot(ref)

plot(test)

## 長さが異なります
ref %>% length

## [1] 1441

test %>% length

## [1] 1657

## 距離を計算できます
alignment = dtw(ref, test)
## DTWで計算した距離
alignment$distance

## [1] 15.82

## 対応関係を示します
## 時系列で示すと, こんな感じです. 
## 線で結ばれた点が, 対応関係を表します. 
dtwPlotTwoWay(alignment, xts=ref, yts=test+1)

## 時系列点の対応関係を, 行列の点の位置で表しています.
plot(alignment$index1, alignment$index2, main="Warping function", type='l')

window.sizeの制限

対応する点を探す範囲を制限できます. あまり遠い点が対応しても, 困りますし. window.size で, 探索範囲の上限幅を決めることができます. 当然, 時系列の長さの差より大きくないとダメです.

## 距離を計算できます
alignment = dtw(ref, test, window.type = "sakoechiba", window.size = abs(length(ref) - length(test)))
## DTWで計算した距離
alignment$distance

## [1] 16.51

plot(alignment$index1, alignment$index2, main="Warping function", type='l')

あまり変わりませんね. 最初のところ(50くらい)で, カットオフされています.

open.end と open.begin

最初や最後を捨てることができます. 位相的な性質が違う部分は, 捨てたほうが良さそうですし.

今回だと, 赤い線の最初のほうを捨てると良さそうです. これを実現するのが, open beginです. 最初 -> 最後だと, open endです.

alignment_openEnd = dtw(ref, test, step=asymmetric, open.end=TRUE, open.begin=TRUE)

## 全部対応させたときの距離
alignment$distance

## [1] 16.51

## open endにしたとき
alignment_openEnd$distance

## [1] 5.956

## 対応関係を図示
dtwPlotTwoWay(alignment_openEnd, xts=ref, yts=test+1)

## 時系列点の対応関係を, 行列の点の位置で表しています.
plot(alignment_openEnd$index1, alignment_openEnd$index2, main="Warping function", type='l')

## y軸が200から始まっています
## 最初の200個が捨てられた, という意味でしょう.

いい感じです.

...でもこれって非対称ですよね. \[\mathrm{dtw}(x,y) \neq \mathrm{dtw}(y,x)\] でしょうね. 距離にならないじゃないですか.
確認してみます. 順番を入れ替えます.

alignment_openEnd = dtw(test, ref,step=asymmetric, open.end=TRUE, open.begin=TRUE)
## 対応関係を図示
dtwPlotTwoWay(alignment_openEnd, xts=ref, yts=test+1)

## 時系列点の対応関係を, 行列の点の位置で表しています.
plot(alignment_openEnd$index1, alignment_openEnd$index2, main="Warping function", type='l')

alignment_openEnd$distance

## [1] 10.5

なるほど.

応用先

たとえば, "あいうえお"の5つの波形が正解としてあったとして, 与えられた音波がどの音か, みたいなことができるんでしょうね. 最近は, 隠れマルコフモデルを使うらしいですけど.

正解データがあって, ある程度の長さがあれば, open begin/endでtestデータを埋め込めんでみるのがいいんですかね.

次にやること

1-nearest neiborでクラス分類します.

正解データがない状態なので, 非対称な距離計算だとだめですよね. open begin/endはダメそうです. そのまま突っ込めばいいんですかね. どうしよう.