連続量の多変量解析モデルを分類問題に適用すると起きること

連続量の多変量解析モデルを 分類問題に適用すると 起きること 東京医科大学医療データサイエンス分野 助教 2022/8/7 統計サマーセミナー2022 原田 和治 1

PLS-DAとは？ • クラスラベルを0-1の行列にコーディング • 例）5サンプルのクラスラベルが 1,1,2,3,3 なら 1 1 𝐘= 0 0 0 0 0 1 0 0 0 0 0 1 1 • 𝐘 をターゲット，特徴量（ここでは代謝物の測定量）の行列を 𝐗 として PLS (NIPALS)を実行 • PLSの「得点」に寄与の大きい代謝物がバイオマーカー候補 2022/8/7 統計サマーセミナー2022 2

結果の可視化など分析の枠組みが発達 ※ R package “ropls” をirisに適用．一撃でこれが出力される． 2022/8/7 統計サマーセミナー2022 3

PLSとは？ • 定義の仕方がいろいろあるが， ターゲットが１変量の時は以下が有名（𝑚番目の方向） max Corr 2 𝒚, 𝐗𝛽 Var 𝐗𝛽 𝛽 subject to 𝛽𝑇 𝐒𝑋 𝛽 = 1, 𝛽𝑇 𝐒𝑋 𝛽𝑗 = 0 𝑗 = 1, … , 𝑚 − 1 • 元はWold (1975, 1984) のNIPALSアルゴリズム • Elements of Statistical Learning (Hastie+, 2009) • これを変形すると以下を得る： max Cov 2 𝒚, 𝐗𝛽 Var 𝒚 −1 subject to (上に同じ) 𝛽 • 要するに，PLSは 𝒚 と 𝐗𝛽 の共分散の最大化である 2022/8/7 統計サマーセミナー2022 4

多変量PLS • 以下，𝐘 ∈ ℝ𝑁×𝐾 ， 𝐗 ∈ ℝ𝑁×𝑝 とし，各々中心化済みとする． • 𝐘 が0-1コーディングの場合は後で扱う • ターゲットが多変量の時は，次のように定義する（第一成分） max 𝜃 𝑇 𝐒𝑌𝑋 𝛽 2 𝜃,𝛽 subject to 𝜃 𝑇 𝜃 = 1, 𝛽𝑇 𝐒𝑋 𝛽 = 1 • 明らかに 𝜃 = 𝐒𝑌𝑋 𝛽/ 𝐒𝑌𝑋 𝛽 であるから，次の問題に帰着 max 𝛽𝑇 𝐒𝑋𝑌 𝐒𝑌𝑋 𝛽 𝛽 subject to 𝛽𝑇 𝐒𝑋 𝛽 = 1 • さらに一般化固有値問題 𝐒𝑋𝑌 𝐒𝑌𝑋 𝛽 = 𝜂𝐒𝑋 𝛽 に帰着 ※ 素直に拡張するなら制約は 𝜃 𝑇 𝑺𝑌 𝜃 = 1 では…？という気もするが，NIPALSに基づくものはコレ． 2022/8/7 統計サマーセミナー2022 5

他の多変量解析手法との比較 • Principal Component Analysis max 𝛽𝑇 𝐒𝑋 𝛽 subject to 𝛽𝑇 𝛽 = 1 𝛽 • 固有値問題 𝐒𝑋 𝛽 = 𝜂𝛽 に帰着 • Canonical Correlation Analysis (CCA) max 𝜃 𝑇 𝐒𝑌𝑋 𝛽 subject to 𝜃 𝑇 𝐒𝑌 𝜃 = 𝛽𝑇 𝐒𝑋 𝛽 = 1 𝜃,𝛽 • 2022/8/7 一般化固有値問題に帰着 𝐎 𝐒𝑋𝑌 𝜃 𝐒𝑋 =𝜂 𝐒𝑌𝑋 𝐎 𝛽 𝐎 𝐎 𝜃 𝐒𝑌 𝛽 統計サマーセミナー2022 6

分類問題にPLSを適用するとどうなる？ • 少し計算すると，次がわかる 𝐾 1 𝐒𝑋𝑌 𝐒𝑌𝑋 = 2 ෍ 𝑁𝑘2 𝑋ത𝑘 − 𝑋ത 𝑋ത𝑘 − 𝑋ത 𝑇 𝑁 群間共分散に 似ている 𝑘=1 • さらに，各クラスのサイズが均等の場合には以下を得る 𝐾 1 1 𝑇 ത ത ത ത 𝐒𝑋𝑌 𝐒𝑌𝑋 = ෍ 𝑁𝑘 𝑋𝑘 − 𝑋 𝑋𝑘 − 𝑋 = 𝐒𝑏𝑒𝑡 𝐾𝑁 𝐾 𝑘=1 • 多変量PLSは以下のように書けた max 𝛽𝑇 𝐒𝑋𝑌 𝐒𝑌𝑋 𝛽 𝛽 subject to 𝛽𝑇 𝐒𝑋 𝛽 = 1 • つまり，多変量PLSを分類に使うと，𝐗 の群によらない相関構造を考慮しながら 群間分散（の類似物）は最適化できる． 2022/8/7 統計サマーセミナー2022 7

でも…… • subject to 𝛽𝑇 𝐒𝑋 𝛽 = 1 ← 惜しい • subject to 𝛽𝑇 𝐒𝑤𝑖𝑡ℎ 𝛽 = 1 ← 嬉しい PLS-DA LDA • LDAなら群内分散も考慮できるのに…… • 高次元でもLDA使いたい…… • スパース推定もしたい…… 2022/8/7 統計サマーセミナー2022 8

LDAを正則化しましょう LDAには様々な構成の仕方がある 1. 混合ガウス分布による方法 2. 群間分散と群内分散の比を最大化する方法（Fisher’s LDA） 3. Optimal Scoring による方法 ← これを扱う 2022/8/7 統計サマーセミナー2022 9

Optimal Scoring (Scaling) • Fisherが1936年頃に提案した方法 • 対応分析（複数のカテゴリカル変数の主成分分析と喩えられる）と関係 • OSの問題は次のように定義される 1 min 𝐘𝜃 − 𝐗𝛽 2 𝜃,𝛽 𝑁 subject to 𝜃 𝑇 𝐒𝑌 𝜃 = 1 • この解がLDAに一致することを見るために，まずはCCAとの関係を見る 2022/8/7 統計サマーセミナー2022 10

OSとCCAの関係（分類問題に限らず成立） • OSの目的関数は制約下で次のように書き換えられる： 1 𝐘𝜃 − 𝐗𝛽 2 = 𝜃 𝑇 𝐒𝑌 𝜃 + 𝛽𝑇 𝐒𝑋 𝛽 − 2𝜃 𝑇 𝐒𝑌𝑋 𝛽 = 1 − 𝛽𝑇 𝐒𝑋 𝛽 − 2𝜃 𝑇 𝐒𝑌𝑋 𝛽 𝑁 • 対して，CCAの問題は次で定義される： max 𝜃 𝑇 𝐒𝑌𝑋 𝛽 subject to 𝜃 𝑇 𝐒𝑌 𝜃 = 1, 𝛽𝑇 𝐒𝑋 𝛽 = 1 𝜃,𝛽 • ゆえに， OSとCCAは同一の 𝜃 を与える． • 𝛽はスケールだけが異なり，それぞれ𝛽𝑐𝑐𝑎 , 𝛽𝑜𝑠 とすると以下を得る． 𝛽𝑜𝑠 𝛽𝑐𝑐𝑎 = 𝑇 𝛽𝑜𝑠 𝐒𝑋 𝛽𝑜𝑠 2022/8/7 統計サマーセミナー2022 11

分類問題におけるCCAとLDA • はじめにNotationを用意しておく（すこし記号をabuse） • 𝐘 ∈ 0,1 𝑁×𝐾 は中心化せず，𝐗 ∈ ℝ𝑁×𝑝 は中心化済とする 1 𝑇 𝐒𝑋 = 𝐗 𝐗, 𝑁 1 𝑇 𝑁𝑘 𝐒𝑌 = 𝐘 𝐘 = diag , 𝑁 𝑁 1 𝑇 𝑇 𝐒𝑌𝑋 = 𝒀 𝐗 = 𝐒𝑋𝑌 𝑁 −1 𝑇 𝑻 =𝐘 𝐘 𝐘 𝐘 𝑃𝑌 𝐌 = 𝐒Y−1 𝐒𝑌𝑋 = 𝑋ത1 , … , 𝑋ത𝐾 𝑇 1 𝑺𝑏𝑒𝑡 = 𝑃𝑌 𝐗 𝑇 𝑃𝑌 𝐗 = 𝐒𝑋𝑌 𝐒Y−1 𝐒𝑌𝑋 = 𝐌𝑇 𝐒𝑌 𝐌 𝑁 1 𝑺𝑤𝑖𝑡ℎ = 𝐼 − 𝑃𝑌 𝐗 𝑇 𝐼 − 𝑃𝑌 𝐗 = 𝐒𝑋 − 𝑺𝑏𝑒𝑡 𝑁 • LDAは max 𝛽𝑇 𝑺𝑏𝑒𝑡 𝛽 subject to 𝛽𝑇 𝐒𝑤𝑖𝑡ℎ 𝛽 = 1 𝛽 2022/8/7 統計サマーセミナー2022 12

CCAとLDA • −1/2 ෨ −1/2 ෨ CCAにおいて，パラメータを 𝜃 = 𝐒𝑌 𝜃, 𝛽 = 𝐒𝑋 𝛽 と変換すると −1/2 −1/2 max 𝜃෨ 𝑇 𝐒𝑌 𝐒𝑌𝑋 𝐒𝑋 𝛽෨ subject to 𝜃෨ = 𝛽෨ = 1 𝜃,𝛽 • これは通常の特異値分解であるから， ෩ 𝐁 ෩ とすると，次の特異値分解を得る． 変換後CCAの解を集めた行列を𝚯, −1/2 −1/2 ෩ 𝛼𝐁 ෩𝑇 𝐒𝑌 𝐒𝑌𝑋 𝐒𝑋 = 𝚯𝐷 • パラメータを元に戻すと， 以下を得る． 𝐒𝑌−1 𝐒𝑌𝑋 𝐒𝑋−1 = 𝚯𝐃𝛼 𝐁𝑇 , 𝚯𝑇 𝐒𝑌 𝚯 = 𝐁𝑇 𝐒𝑋 𝐁 = 𝐈𝐿 (𝐿 = min 𝑁, 𝐾 ) 2022/8/7 統計サマーセミナー2022 13

CCAとLDA 𝐒𝑌−1 𝐒𝑌𝑋 𝐒𝑋−1 = 𝚯𝐃𝛼 𝐁𝑇 , 𝚯𝑇 𝐒𝑌 𝚯 = 𝐁𝑇 𝐒𝑋 𝐁 = 𝐈𝐿 (𝐿 = min 𝑁, 𝐾 ) • 上式を使って変形すると，次の関係を容易に示せる： 𝐒𝑌−1 𝐒𝑌𝑋 𝐁 = 𝚯𝐃𝛼 𝐒𝑋−1 𝐒𝑋𝑌 𝚯 = 𝐁𝐃𝛼 𝚯𝑇 𝐒𝑌𝑋 𝐁 = 𝐃𝛼 𝚯𝑇 𝐒𝑌𝑋 𝐒𝑋−1 𝐒𝑋𝑌 𝚯 = 𝐃𝛼2 𝐁𝑇 𝐒𝑋𝑌 𝐒𝑌−1 𝐒𝑌𝑋 𝐁 = 𝐃𝛼2 𝑇 • CCAとOSの関係 𝛽𝑐𝑐𝑎 = 𝛽𝑜𝑠 / 𝛽𝑜𝑠 𝐒𝑋 𝛽𝑜𝑠 を踏まえると， 𝐁𝑜𝑠 = 𝐁𝐃𝛼 2022/8/7 統計サマーセミナー2022 14

CCAとLDA • さらに， 𝑺𝑏𝑒𝑡 = 1 𝑁 𝑃𝑌 𝐗 𝑇 𝑃𝑌 𝐗 = 𝐒𝑋𝑌 𝐒Y−1 𝐒𝑌𝑋 を使えば， 𝐁𝑇 𝑺𝑏𝑒𝑡 𝐁 = 𝐁𝑇 𝐒𝑋𝑌 𝐒𝑌−1 𝐒𝑌𝑋 𝐁 = 𝐃𝛼2 • また， 𝑺𝑤𝑖𝑡ℎ = 𝐒𝑋 − 𝑺𝑏𝑒𝑡 より， 𝐁𝑇 𝑺𝑤𝑖𝑡ℎ 𝐁 = 𝐈𝐿 − 𝐃𝛼2 • ゆえに， 𝐁𝐿𝐷𝐴 = 𝐁𝐃 1−𝛼2 −1/2 とすれば， 𝐁𝐿𝐷𝐴 の各列はLDAの解となる． 𝑇 𝐁𝐿𝐷𝐴 𝑺𝑤𝑖𝑡ℎ 𝐁𝐿𝐷𝐴 = 𝐈𝐿 , 𝑇 𝐁𝐿𝐷𝐴 𝑺𝑏𝑒𝑡 𝐁𝐿𝐷𝐴 = 𝐃𝛼2 / 1−𝛼 2 • したがって，Optimal Scoringの解はLDAの解とスケール以外一致する 𝐁𝑜𝑠 = 𝐁𝐃𝛼 = 𝐁𝐿𝐷𝐴 𝐃𝛼 1−𝛼 2 2022/8/7 統計サマーセミナー2022 15

OSでLDAを構成すると何が嬉しいのか？ • Optimal Scoringは 𝜃 と 𝛽 の交互更新で解ける． • 𝜃𝑘 の更新は以下の繰り返し（逆行列計算や特異値分解なし） 1. 推定済みの 𝜃𝑙 (𝑙 < 𝑘) に対する直交化 2. 目的関数の最小化とNormamalization • 𝛽𝑘 の更新は 𝐘𝜃 をターゲットとする（正則化）最小二乗法 • 二乗誤差最小化問題であるから，さまざまな正則化がそのまま使える 実装もglmnetほか，スパース推定パッケージを組み込むだけ． 例：Elastic Net型 1 min 𝐘𝜃 − 𝐗𝛽 2 + 𝜆 1 − 𝛼 𝛽 1 + 𝛼 𝛽 22 𝜃,𝛽 𝑁 subject to 𝜃 𝑇 𝐒𝑌 𝜃 = 1 2022/8/7 統計サマーセミナー2022 16

別の形のLeast Square型LDA • 実は，2クラスならもっと簡単にできる（DSDA; Mai, Zou, and Yuan, 2012, Biometrika） • クラス1に対して 𝑦𝑖 = −𝑛/𝑛1 ， 𝑦𝑖 = 𝑛/𝑛2 として，以下の回帰問題を解く 𝑛 𝛽መ 𝐷𝑆𝐷𝐴 𝜆 = argmin ෍ 𝑦𝑖 − 𝛽0 − 𝑋𝑖 𝛽 2 + 𝜆 𝛽 1 𝛽 𝑖=1 • 多クラス拡張は Mai, Yang, and Zou (2019, Statistica Sinica) だが， こちらはもう少し煩雑 • Sparse Discriminant Analysis (Clemmensen+, 2011) 2022/8/7 統計サマーセミナー2022 17

Experiment • OS + glmnetによるLDAを実装したので，簡単な設定だが分析結果を紹介． • 𝑁 = 50 • 𝑝 = 100 • 𝐾=3 • 𝑋は多変量正規分布 • X1,X2,X3,X4だけ分類に寄与 • X1とX2は強く相関 2022/8/7 統計サマーセミナー2022 18

Experiment • LDAの結果 2022/8/7 統計サマーセミナー2022 19

Experiment • Sparse LDA （Elastic net型： 𝜆 = 0.1, 𝛼 = 0.5）の結果 • スパースな解を得ている 2022/8/7 統計サマーセミナー2022 20

References 1. Frank, I. E. and Friedman, J. H. (1993). A Statistical View of Some Chemometrics Regression Tools. Technometrics 35 109–35. 2. Hastie, T., Buja, A. and Tibshirani, R. (1995). Penalized Discriminant Analysis. AoS 23 73–102. 3. Barker, M. and Rayens, W. (2003). Partial least squares for discrimination. J. Chemom. 17 166–73. 4. Sun, L., Ji, S., Yu, S. and Ye, J. (2009). On the Equivalence between Canonical Correlation Analysis and Orthonormalized Partial Least Squares. In IJCAI. 5. Clemmensen, L., Hastie, T., Witten, D. and Ersbøll, B. (2011). Sparse Discriminant Analysis. Technometrics 53 406–13. 6. Kawano, S., Fujisawa, H., Takada, T. and Shiroishi, T. (2015). Sparse principal component regression with adaptive loading. Comput. Stat. Data Anal. 89 192–203. 7. Mai, Q., Yang, Y. and Zou, H. (2019). MULTICLASS SPARSE DISCRIMINANT ANALYSIS. Stat. Sin. 29 97–111. 2022/8/7 統計サマーセミナー2022 21

References 8. Mai, Q. and Zou, H. (2013). A Note On the Connection and Equivalence of Three Sparse Linear Discriminant Analysis Methods. Technometrics 55 243–6. 9. Mai, Q., Zou, H. and Yuan, M. (2012). A direct approach to sparse discriminant analysis in ultra-high dimensions. Biometrika 99 29–42. 2022/8/7 統計サマーセミナー2022 22