おい、FDAってなんなんだ？

関数データ解析（Functional Data Analysis; FDA）って知ってますか？
時系列データを用いての回帰分析や分類について頭を悩ませていた今日このごろ、いつもお世話になっている先生にポロッと相談したのがきっかけでした。自分の中では、課題と解決策が同時に提示されたような感覚でワクワクが止まりませんでした。
さらにPythonのscikit-fdaなるパッケージに丁寧にチュートリアルが用意されている・・・。試さずにいられようか。いや、いられない。ということで、scikit-fdaのチュートリアルをやったので、備忘録です。

Functional data analysis（関数データ解析）について

計時測定データ（時間の経過や位置の変化に応じて繰り返し計測されたデータ）に対して、得られたデータと着目する事象を統計モデルによって明らかにするために、多変量データとして扱い回帰分析を適用する場合以下のような問題が生じます。

1. 観測データは一般的にはノイズが混入されており、経時測定データのばらつきが大きい場合に本質構造が捉えにくい。
2. 観測時点数の増加はそのままデータの次元の増加につながるため、推定量が不安定になる。
3. 欠損などの理由で個体ごとに観測時点数が不均一であったり、観測時点が不揃いである状況もあり、欠損値補完では補いきれず、一般的には古典的な多変量解析手法を直接適用することは困難である。

これらの問題を解決するため、個々のデータを関数化し、得られた関数をデータとして扱う方法、関数データ解析（Functional Data Analysis; FDA）が確立されています（Ramsay and Silverman, 2005）。
つまり、FDAを用いることで下記の問題を解消することができる！！ということです。・・・そんなに簡単じゃないか。

1. 観測ノイズの除去
2. 高次元化の抑制
3. 観測時点や観測時点数のずれの補完

参考：
1. 関数データに基づく統計的モデリング
2. 情報教育課題合格ログデータによる受講生の類型化

scikit-fda

ドキュメントが充実しているパッケージっていいですよね。こちらもIntroductionから始まり、Getting the data、Basis representation、Scikit-fda and scikit-learnと段階的にパッケージの使用方法を理解できます。手元の環境で動かした結果を貼り付けてくので、興味を持った方はぜひやってみてください。本当は実データ結果をお見せしたかったのですが・・・。

今回は、チュートリアルで用いられている古典的なデータセットの一つである『Berkeley Growthデータセット』（93人の少年少女の身長を誕生~18歳の誕生日までをいくつかの時点で測定したもの）を使っています。出てくるコードはPythonです。参考までに載せています。

## データセットのプロット
まずは概要を確認しましょう。
​

?? skfda.datasets.fetch_growth

データセットの概要について

Berkeley Growth Studyのデータセットを読み込みます。このデータはRパッケージの'fda'から取得しています。Berkeley Growth Studyです。
Berkeley Growth Study (Tuddenham and Snyder, 1954)は、1歳から18歳までの女の子54人男の子39人の身長を記録したものです。身長は各子供の31の年齢で測定され、その標準誤は約3mmで、幼児期に大きく、晩年に低くなる傾向がありました。
参考文献：Tuddenham, R. D., and Snyder, M. M. (1954) "Physical growth of California boys　and girls from birth to age 18", University of California Publications in Child Development, 1,　183-364.

あと関数について

引数（Args）：
return_X_y：データとターゲットのみをタプルで返します。
as_frame： データをPandas DataframeまたはSeriesで返します。
ソースコード（Source）：
def fetch_growth(・・・関数の中身が書いていて長いので割愛)

ふむふむ・・・。では、一旦グラフにしてみましょう。

# デモ用のデータセットを読み込む
X, y = skfda.datasets.fetch_growth(return_X_y=True)
# サンプル数を減らす
X = X[:10]
print(X.plot(marker="^", markersize=0.7, ls=":", lw=0.1).set_size_inches(18, 8))
plt.savefig('./mics/fetch_growth_plot.png', dpi=300)

こういった離散的に得られた経時測定データを解析するのは結構鬱陶しいです。そこで関数化処理します。​

今回は、基底関数はB-スプライン、個数は10として（適当です。）、基底関数展開によって平滑化しました。
パッケージに入っている基底関数については下記のコードで参照ください。

?? skfda.representation.basis

では、ここから何するのか？ということです。平滑化という工程は離散データを複数の基底関数に重みをつけて表すということです。つまり、基底関数を変数とする多変量関数なのです！
だからなにって？なんといっても、素晴らしいのが関数データを変数として統計手法にあてはめることができることです。つまり、以下に示すようなクラスタリングはもちろん、主成分分析とかやりたい放題なのです。

k-nearest neighbor algorithm（k-nn）

パッと思いつくのはやはり分類ですね。ということでやってみました。
まずは、男女別（青色が男子、オレンジが女子）でプロットしてみます。

みた感じで綺麗に分かれているとやる気が出ます。このデータセットを使って性別を目的変数にして、分類器に当てはめてみましょう。今回はチュートリアル通りにk-nnでやってみます。コード簡単だよアピールで載せておきます。

# データセットを分けて
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.25, stratify=y, random_state=0)

# 分類器に多変量データの代わりに関数データFDataGrid（skfdaでのデータの持ち方です。）で渡すだけ。
knn = KNeighborsClassifier()
knn.fit(X_train, y_train)
score = knn.score(X_test, y_test)
print(score)

Score=0.958でaccuracyが約96%の身長の計時計測データから男子と女子を分類するモデルができました。FDAの記事なのでここでストップします。
こんな感じで、関数化した後の拡張解析が充実しています。

主成分分析（PCA）

次は主成分分析してみます。
これはちょっと理解が曖昧なので、結果のみの共有です。申し訳ないです。

まずはPCAなので適当に5つくらいの成分の寄与度を出して、成分の採用個数を決めたいと思います。

二つくらいでほとんどのデータを説明できそうです。（はい、性別です。）
こちらもコードが簡単なので、載せておきます（関数はすでにimportしています。）。

fpca = FPCA(n_components=2)
fpca.fit(X_basis)
print(fpca.components_.plot())

このグラフについては正確な解釈できてないので、師匠に教えてもらおうと思います。こんなんわからんのか！って方で優しい人は是非コメントで教えてください。

微分

最後に・・・。微分はFDAならではの解析方法かと思います。結構データ解釈の時にあーだこーだ言えるのでパッケージに用意してもらえて嬉しいです。

先ほどから使っているデータセットから男子のデータのみを取得し、微分します。コードは少しややこしいので、割愛します。

いわゆる加速度が分かります。自分が使うデータでもやってみたいです。あーだこーだ言えると思います。

アライメント

ついでにアライメントします。このアライメントってのは、関数データの変動傾向から曲線間でピーク位置を揃えてくれます。
おそらく、ユースケースは例1）計測タイミングを揃えたはずだけど、何フレームかずれる、例2）計測条件が一致しているのに何かしらの周期性がある場合とかそんな感じの時に使えると思います。

今回の場合は個人差があるはずの成長タイミングを強引にと揃えた感じになっています。重要な特徴量を一つ消すことになりますので、アライメントの意味はありません。機能の確認としては波がそろっているので、わかりやいのではないでしょうか？

まとめ

以上がFDAの概要とパッケージ使ってみましたってまとめです。実際に使う時には、いろんな問題が出てきそうなのでその都度深掘りしながらチャレンジしたいと思います。続報（成果物）が出せるといいな！！

この記事を読んでまさに！って方やもう知っているよって方、ぜひハートマークください！

※ noteでコード書くのは場違いかなと思いつつも、記事を一箇所に集めておきたかったので無理をしました。分かりづらくて申し訳ないです。