【カルマンフィルタ#2】【Python】処理の高速化(この漸化式解けますか?)

今回は前回のPythonでのカルマンフィルタの続きです。というか久々に数学を解いた気がしたので、紹介したいと思って記事にしました。

時間間隔が変化する場合のカルマンフィルタ

以前のカルマンフィルタは時間間隔が一定の場合の処理方法でした。でももしかしたら時間間隔が変化する計測データ(例えば、あまり重要でない時間帯は時間間隔を広く、重要な時間帯は時間間隔を短くしたい場合の計測データ)があるかもしれません。

この場合、この資料の最後の方に書かれている不当間隔時間モデルというのがあるみたいなのですが、正直難しくてわかりませんでした。その他の方法として、一番短い時間間隔に合わせる方法が考えられます。つまり、一番短い時間間隔でフィルタし続け、データがない部分はカルマン予測を行う方法です。

ただし、この場合時間間隔が短い場合に予測に時間がかかってしまう恐れがあります(下図)。

この予測のループ計算を高速化できればパラメータ算出やフィルタリングで時間がかからずに済みます(下図)。

コードでいうと以下の部分、この関数をforループで繰り替えす必要がある。

1期先予測の高速化

1次のトレンドモデルの場合

1次のトレンドモデルの場合は簡単です。

コードで書くと、nn-1になっているのは測定点の間隔分だけ予測するようにしているからです。例えば1分間隔で測定していて、最短の1秒間隔で予測する場合、測定してから59回予測してから次の測定点になるからです。

2次のトレンドモデルの場合

2次のトレンドモデルの場合は難しいです。

そこで以下のようにEnを仮定して漸化式を解きます。これがやりたかった!というか他に解き方があったら教えてほしいです。

こういうの高校で数学勉強して以来だったので解いててワクワクしてしまいました。

以下僕なりの解法です。

Enをnの式で表すことができました!あとは総和を求めます。

G行列のn乗についてはどうしても計算しなければいけませんが、numpyのnp.linalg.matrix_power(G,n)で高速で計算できるのでそれほど時間はかからないと思います。

コードの作成

以上の式を使ってカルマンフィルタのコードを記述します。なお、カルマン平滑化については今回考えていないので、削除しました(ごめんなさい)。

1. データの作成

今回は予め作っておいたcsvを読み込みます。中身は前回の記事で作成したものと同じですが、indexとして時間を1分間隔で作成し、400step以降は1秒間隔にしました(計1000step)。

出力図は以下の通り、見ずらいですが、青点線で時間間隔が変わっていて、赤間隔の間が観測データがないところです。

2. カルマンフィルタ(1期先予測・フィルタ)の作成

1期先予測とフィルタのコードです。

3. パラメータ推定(最尤推定法)

ほぼ変わってないですが、前時刻と今時刻を比較してその差が2秒以上あればその分予測をするようにしています。1秒以下の場合はこの行程をスキップして観測値がなければ予測、あればフィルタの処理を行っています。

4. カルマンフィルタの実行

インデックスを比較して秒刻みで予測とフィルタを行うようにしています。

5. 図化+出力

図化と出力です。平滑化を省きました。

図化したグラフがこれ

ん~微妙かも、予測が多すぎるのも良くないんですかね?まあ計算できるってことでよし

時間の比較

通常の方法と比較して計算してみた結果、尤度推定にかかる時間は通常の場合で0.1426秒、上述の方法では0.0688秒でした。半分くらいの時間ですね。多分もっと長い時間であれば効果がより確認できると思います。

今回の短いステップ数でもこれだけ差がでるので十分だと思います。

まとめ

今回はカルマンフィルタの1期先予測で高速に計算する方法を考えてみました。

今後もカルマンフィルタについて勉強していきたいと思いますので、よろしくお願いします。

参考サイト

コメント

タイトルとURLをコピーしました