基本的な考え方
(Press ?
for help, n
and p
for next and previous slide)
村田 昇
ID | V1 | V2 | V3 | V4 | V5 | V6 | V7 | V8 | V9 | class |
---|---|---|---|---|---|---|---|---|---|---|
1000025 | 5 | 1 | 1 | 1 | 2 | 1 | 3 | 1 | 1 | benign |
1002945 | 5 | 4 | 4 | 5 | 7 | 10 | 3 | 2 | 1 | benign |
1015425 | 3 | 1 | 1 | 1 | 2 | 2 | 3 | 1 | 1 | benign |
1016277 | 6 | 8 | 8 | 1 | 3 | 4 | 3 | 7 | 1 | benign |
1017023 | 4 | 1 | 1 | 3 | 2 | 1 | 3 | 1 | 1 | benign |
1017122 | 8 | 10 | 10 | 8 | 7 | 10 | 9 | 7 | 1 | malignant |
1018099 | 1 | 1 | 1 | 1 | 2 | 10 | 3 | 1 | 1 | benign |
1018561 | 2 | 1 | 2 | 1 | 2 | 1 | 3 | 1 | 1 | benign |
1033078 | 2 | 1 | 1 | 1 | 2 | 1 | 1 | 1 | 5 | benign |
1033078 | 4 | 2 | 1 | 1 | 2 | 1 | 2 | 1 | 1 | benign |
1035283 | 1 | 1 | 1 | 1 | 1 | 1 | 3 | 1 | 1 | benign |
1036172 | 2 | 1 | 1 | 1 | 2 | 1 | 2 | 1 | 1 | benign |
1041801 | 5 | 3 | 3 | 3 | 2 | 3 | 4 | 4 | 1 | malignant |
1043999 | 1 | 1 | 1 | 1 | 2 | 3 | 3 | 1 | 1 | benign |
1044572 | 8 | 7 | 5 | 10 | 7 | 9 | 5 | 5 | 4 | malignant |
Figure 1: 9種類の生検の散布図
Figure 2: 主成分分析を用いて2次元に縮約した生検と良性・悪性の関係
Figure 3: 線形判別関数による判別関数値の分布
判別分析 (discriminant analysis) の目的
個体の特徴量からその個体の属する クラス を予測する関係式を構成する方法
\(X=\boldsymbol{x}\) の下で \(Y=k\) となる 条件付確率 を計算
\begin{equation} p_{k}(\boldsymbol{x})=P(Y=k|X=\boldsymbol{x}) \end{equation}
観測データ : \(n\) 個の \((Y,X_{1},\dots,X_{q})\) の組
\begin{equation} \{(y_{i},x_{i1},\dots,x_{iq})\}_{i=1}^n \end{equation}
事象 \(X=\boldsymbol{x}\) が起きたという条件の下で 事象 \(Y=k\) が起きる条件付確率
\begin{equation} p_{k}(\boldsymbol{x}) = P(Y=k|X=\boldsymbol{x}) = \frac{P(Y=k,X=\boldsymbol{x})}{P(X=\boldsymbol{x})} \end{equation}
\(Y=k\) の下での \(X\) の条件付き確率質量関数
\begin{equation} f_{k}(\boldsymbol{x}) = P(X=\boldsymbol{x}|Y=k)=\frac{P(X=\boldsymbol{x},Y=k)}{P(Y=k)} \end{equation}
のモデル化を通じて \(p_{k}(\boldsymbol{x})\) をモデル化する
\(f_{k}(\boldsymbol{x})\) から \(p_{k}(\boldsymbol{x})\) を得る数学的原理
原因 \(X=\boldsymbol{x}\) から結果 \(Y=k\) が生じる確率 を 結果 \(Y=k\) が生じる原因が \(X=\boldsymbol{x}\) である確率 から計算する方法
Bayes の公式 (Bayes’ formula)
\begin{equation} p_{k}(\boldsymbol{x}) = P(Y=k|X=\boldsymbol{x}) = \frac{f_{k}(\boldsymbol{x})P(Y=k)}{\sum_{l=1}^{K}f_l(\boldsymbol{x})P(Y=l)} \end{equation}
定義より
\begin{equation} f_{k}(\boldsymbol{x}) = P(X=\boldsymbol{x}|Y=k) = \frac{P(X=\boldsymbol{x},Y=k)}{P(Y=k)} \end{equation}
求める条件付確率
\begin{equation} p_{k}(\boldsymbol{x}) = P(Y=k|X=\boldsymbol{x}) = \frac{f_{k}(\boldsymbol{x})P(Y=k)}{P(X=\boldsymbol{x})} \end{equation}
分母の展開
\begin{align} P(X=\boldsymbol{x}) &= \sum_{l=1}^{K}P(X=\boldsymbol{x},Y=l)\\ &= \sum_{l=1}^{K}f_l(\boldsymbol{x})P(Y=l) \end{align}
Bayes の公式による書き換え
\begin{equation} p_{k}(\boldsymbol{x}) = \frac{f_{k}(\boldsymbol{x})\pi_{k}}{\sum_{l=1}^{K}f_l(\boldsymbol{x})\pi_l} = \frac{f_{k}(\boldsymbol{x})}{\sum_{l=1}^{K}f_l(\boldsymbol{x})\pi_l} \cdot\pi_{k} \end{equation}
事前に特別な情報がない場合
データから自然に決まる確率
\begin{equation} \pi_{k} = \frac{\text{\(Y=k\)のサンプル数}}{\text{全サンプル数}} \end{equation}
事前に情報がある場合
食事・運動・飲酒・ストレスなどの生活の特徴から生活習慣病か否かを判別
- 健常者の食事・運動・飲酒・ストレスなどの特徴量を収集
- 罹患者の食事・運動・飲酒・ストレスなどの特徴量を収集
- 事前確率は 別の調査の日本人の罹患率 を利用
判別関数 : \(\delta_{k}(\boldsymbol{x})\) (\(k=1,\dots,K\))
\begin{equation} p_{k}(\boldsymbol{x}) < p_l(\boldsymbol{x}) \Leftrightarrow \delta_{k}(\boldsymbol{x}) < \delta_l(\boldsymbol{x}) \end{equation}
共分散行列 \(\Sigma\) : すべてのクラスで共通
\begin{equation} f_{k}(\boldsymbol{x}) = \frac{1}{(2\pi)^{q/2}\sqrt{\det\Sigma}} \exp\left(-\frac{1}{2}(\boldsymbol{x}-\boldsymbol{\mu}_{k})^{\mathsf{T}} \Sigma^{-1}(\boldsymbol{x}-\boldsymbol{\mu}_{k})\right) \end{equation}
線形判別関数 : \(\boldsymbol{x}\) の1次式
\begin{equation} \delta_{k}(\boldsymbol{x}) = \boldsymbol{x}^{\mathsf{T}}\Sigma^{-1}\boldsymbol{\mu}_{k} -\frac{1}{2}\boldsymbol{\mu}_{k}^{\mathsf{T}}\Sigma^{-1}\boldsymbol{\mu}_{k} +\log\pi_{k} \end{equation}
事後確率と判別関数の関係
\begin{align} &p_{k}(\boldsymbol{x}) < p_{l}(\boldsymbol{x})\\ &\Leftrightarrow f_{k}(\boldsymbol{x})\pi_{k} < f_l(\boldsymbol{x})\pi_l \quad\text{(分母は共通)}\\ &\Leftrightarrow \log f_{k}(\boldsymbol{x})+\log\pi_{k} < \log f_l(\boldsymbol{x})+\log\pi_l\\ &\Leftrightarrow -\frac{1}{2}(\boldsymbol{x}-\boldsymbol{\mu}_{k})^{\mathsf{T}} \Sigma^{-1}(\boldsymbol{x}-\boldsymbol{\mu}_{k})+\log\pi_{k}\\ &\phantom{\Leftrightarrow}\quad < -\frac{1}{2}(\boldsymbol{x}-\boldsymbol{\mu}_l)^{\mathsf{T}} \Sigma^{-1}(\boldsymbol{x}-\boldsymbol{\mu}_l)+\log\pi_l\\ &\Leftrightarrow \delta_{k}(\boldsymbol{x}) < \delta_l(\boldsymbol{x}) \quad\text{(2次の項は右辺と左辺で共通)} \end{align}
平均の推定 (クラスごとに行う)
\begin{equation} \hat{\boldsymbol{\mu}}_{k} = \frac{1}{n_{k}}\sum_{i:y_{i}=k}\boldsymbol{x}_{i} \end{equation}
分散の推定 (まとめて行う)
\begin{equation} \hat{\Sigma} = \frac{1}{n{-}K}\sum_{k=1}^{K}\sum_{i:y_{i}=k} (\boldsymbol{x}_{i}-\hat{\boldsymbol{\mu}}_{k}) (\boldsymbol{x}_{i}-\hat{\boldsymbol{\mu}}_{k})^{\mathsf{T}} \end{equation}
関数 MASS::lda()
library(MASS) # パッケージの読み込み lda(formula, data, ..., subset, na.action) #' formula: モデル式 (ラベル ~ 判別に用いる変数) #' data: 必要な情報を含むデータフレーム #' 詳細は '?MASS::lda' を参照
判別結果取得のための関数
predict(object, newdata, prior = object$prior, dimen, method = c("plug-in", "predictive", "debiased"), ...) #' object: ldaの返すオブジェクト #' newdata: 予測の対象とするデータフレーム #' prior: ラベルの事前分布 #' 返値はリスト型で以下の項目がある #' class: 判別関数による予測ラベル #' posterior: ラベルの事後確率 #' x: 判別関数値 #' 詳細は '?MASS::predict.lda' を参照
9月と10月の気温と湿度のデータを抽出する
tw_data <- read_csv("data/tokyo_weather.csv") tw_subset <- tw_data |> filter(month %in% c(9,10)) |> select(temp, humid, month) |> mutate(month = as_factor(month)) # 月を因子化
半分のデータを用いて線形判別関数を構成し,残りのデータを用いて判別を行う
library(MASS) idx <- seq(2, 60, by = 2) tw_train <- tw_subset[ idx,] # 訓練データ tw_test <- tw_subset[-idx,] # 試験データ tw_lda <- lda(month ~ temp + humid, data = tw_train) # 線形判別関数の構成 tw_lda_fitted <- predict(tw_lda) # 判別関数によるクラス分類結果の取得 tw_lda_predict <- predict(tw_lda, newdata = tw_test) # 新しいデータの予測
共分散行列 \(\Sigma_{k}\) : クラスごとに異なる
\begin{equation} f_{k}(\boldsymbol{x}) = \frac{1}{(2\pi)^{q/2}\sqrt{\det\Sigma_{k}}} \exp\left(-\frac{1}{2}(\boldsymbol{x}-\boldsymbol{\mu}_{k})^{\mathsf{T}} \Sigma_{k}^{-1}(\boldsymbol{x}-\boldsymbol{\mu}_{k})\right) \end{equation}
2次判別関数 : \(\boldsymbol{x}\) の2次式
\begin{equation} \delta_{k}(\boldsymbol{x}) = -\frac{1}{2}\det\Sigma_{k} -\frac{1}{2}(\boldsymbol{x}-\boldsymbol{\mu}_{k})^{\mathsf{T}} \Sigma_{k}^{-1}(\boldsymbol{x}-\boldsymbol{\mu}_{k}) +\log\pi_{k} \end{equation}
事後確率と判別関数の関係
\begin{align} &p_{k}(\boldsymbol{x}) < p_l(\boldsymbol{x})\\ &\Leftrightarrow f_{k}(\boldsymbol{x})\pi_{k} < f_l(\boldsymbol{x})\pi_l \quad\text{(分母は共通)}\\ &\Leftrightarrow \log f_{k}(\boldsymbol{x})+\log\pi_{k} < \log f_l(\boldsymbol{x})+\log\pi_l\\ &\Leftrightarrow -\frac{1}{2}\det\Sigma_{k} -\frac{1}{2}(\boldsymbol{x}-\boldsymbol{\mu}_{k})^{\mathsf{T}} \Sigma_{k}^{-1}(\boldsymbol{x}-\boldsymbol{\mu}_{k}) +\log\pi_{k}\\ &\phantom{\Leftrightarrow}\quad < -\frac{1}{2}\det\Sigma_l -\frac{1}{2}(\boldsymbol{x}-\boldsymbol{\mu}_l)^{\mathsf{T}} \Sigma_l^{-1}(\boldsymbol{x}-\boldsymbol{\mu}_l) +\log\pi_l\\ &\Leftrightarrow \delta_{k}(\boldsymbol{x}) < \delta_l(\boldsymbol{x}) \end{align}
平均の推定 (クラスごとに行う)
\begin{equation} \hat{\boldsymbol{\mu}}_{k} = \frac{1}{n_{k}}\sum_{i:y_{i}=k}\boldsymbol{x}_{i} \end{equation}
分散の推定 (クラスごとに行う)
\begin{equation} \hat{\Sigma}_{k} = \frac{1}{n_{k}-1}\sum_{i:y_{i}=k} (\boldsymbol{x}_{i}-\hat{\boldsymbol{\mu}}_{k}) (\boldsymbol{x}_{i}-\hat{\boldsymbol{\mu}}_{k})^{\mathsf{T}} \end{equation}
関数 MASS::qda()
library(MASS) # 既に読み込んでいる場合は不要 qda(formula, data, ..., subset, na.action) #' formula: モデル式 (ラベル ~ 判別に用いる変数) #' data: 必要な情報を含むデータフレーム #' 詳細は '?MASS::qda' を参照
判別結果取得のための関数
predict(object, newdata, prior = object$prior, method = c("plug-in", "predictive", "debiased", "looCV"), ...) #' object: qda の返すオブジェクト #' newdata: 予測対象のデータフレーム.省略すると推定に用いたデータフレーム #' 返値はリスト型で以下の項目がある #' class: 判別関数による予測ラベル #' posterior: ラベルの事後確率 #' 詳細は '?MASS::predict.qda' を参照
前問と同様な設定で2次判別を行いなさい
tw_qda <- qda(month ~ temp + humid, data = tw_train) # 2次判別関数の構成 tw_qda_fitted <- predict(tw_qda) # 判別関数によるクラス分類結果の取得 tw_qda_predict <- predict(tw_qda, newdata = tw_test) # 新しいデータの予測
変動の関係
\begin{equation} \text{(全変動)} = \text{(群内変動)} + \text{(群間変動)} \end{equation}\begin{equation} A = W + B \end{equation}
Fisherの基準
\begin{equation} \text{maximize}\quad \boldsymbol{\alpha}^{\mathsf{T}} B\boldsymbol{\alpha} \quad\text{s.t.}\quad \boldsymbol{\alpha}^{\mathsf{T}} W\boldsymbol{\alpha}=\text{const.} \end{equation}
\(K=2\) の場合 : 最大固有値を用いる (線形判別と一致)
\begin{equation} \boldsymbol{\alpha}\propto W^{-1}(\boldsymbol{\mu}_{1}-\boldsymbol{\mu}_2) =\Sigma^{-1}(\boldsymbol{\mu}_{1}-\boldsymbol{\mu}_2) \end{equation}
9月,10月,11月の気温と湿度のデータを用いて判別関数を作成しなさい.
tw_subset <- tw_data |> filter(month %in% c(9,10,11)) |> select(temp, humid, month) |> mutate(month = as_factor(month))
別の月や変数を用いて判別分析を行いなさい
#' 雨の有無を識別する例 tw_subset2 <- tw_data |> mutate(rain = factor(rain > 0), month = as_factor(month)) |> # 雨の有無でラベル化する select(rain, temp, solar, wind, month) tw_lda2 <- lda(rain ~ ., data = tw_subset2) # 'rain' をそれ以外で判別