5.6 判別分析 (discriminant analysis)
5.6.2 線形判別分析による予測 predict.lda()
v 6.692857 3.021429 5.603571 1.9964286
Coefficients of linear discriminants: # 線形判別子の係数
LD1 LD2
Sepal.L. 0.9055751 0.3902900 Sepal.W. 1.5851015 -2.5279530 Petal.L. -2.1554253 0.1300322 Petal.W. -3.0802797 -1.4416391 Proportion of trace:
LD1 LD2 0.9912 0.0088
x
最大dimen
個の判別変数に関するテスト例のスコア# データセットiris3を使用する
> library(MASS) # MASSパッケージ読み込み
> tr <- sample(1:50, 25) # 訓練用データ．25例の添字をランダムに選ぶ
# 訓練用データ(3種類のアヤメからそれぞれ添字trに該当するものを抜き出す)
> train <- rbind(iris3[tr, , 1], iris3[tr, , 2], iris3[tr, , 3])
# テスト用データ(残り，3種類のアヤメ25例の4種類の性質)
> test <- rbind(iris3[-tr, , 1], iris3[-tr, , 2], iris3[-tr, , 3])
> cl <- factor(c(rep("s", 25), rep("c", 25), rep("v", 25))) # 種類を表す因子
> z <- lda(train, cl) # 線形判別解析実行
# 結果の線形判別子を用いてテスト用データを判別(つまり4種類の性質から種類を当てる)
> predict(z, test)$class
[1] s s s s s s s s s s s s s s s s s s s s s s s s s c c c c c c c c c c c [37] c c c c v c c c c c c c c c v v v v v v v v v v v v v v v v v v v v v v [73] v v v
Levels: c s v
> (cl == predict(z, test)$class) # 判別結果は正しかったか?(正答率74/75) [1] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
[13] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [25] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [37] TRUE TRUE TRUE TRUE FALSE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [49] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [61] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [73] TRUE TRUE TRUE
5.6.3 2 次判別関数 qda()
qda()
は2
次関数(Maharanobis
距離による楕円状の判別領域を用いる)
による判別分 析(quadratic discriminant analysis)
を行う．この関数は標準パッケージには無く，パッ ケージMASS
中にある．QR
分解を使い，もしグループ内共分散行列が特定のグループに 対し特異なら，エラーメッセージを与える．書式：
qda(x, ...)
#
クラス"formula"
用のS3
メソッドqda(formula, data, ..., subset, na.action)
#
既定のS3
メソッドqda(x, grouping, prior = proportions, method, CV = FALSE, nu, ...)
#
クラス"data.frame"
用の既定メソッドqda(x, ...)
#
クラス"matrix"
用の既定メソッドqda(x, grouping, ..., subset, na.action)
引数：
formula groups ~ x1+x2+...
の形のモデル式．つまり目的変数はグルーピング因 子で，右辺は(
因子でない)
判別子を指定するdata formula
中で指定された説明変数が優先的に選ばれるデータフレームx
行列，データフレーム，もしくは説明変数を含む行列(
もし如何なるモデル式も主要引き数として与えられないとき必要
)
grouping
各観測値に対するクラスを指示する因子(
もし如何なるモデル式も主要引き数として与えられないとき必要
)
prior
クラスへの所属可能性に対する事前確率．もし指定されないと訓練データに対するクラス比率が使われる．もし与えられると，確率は因子水準の順番で並べる 必要がある
subset
訓練用サンプルとして使われる観測値例を指定する添字ベクトル(
もし与えられるなら，この引き数は名前付きでなければならない
)
na.action
もし欠損値NA
が見つかったときの処理法を指定する関数．既定の処理は手続きを失敗させる．もう一つの選択肢は
na.omit
で，必要とされる変数に欠 損値が含まれる例があれば，それを除外する(
もし与えられるなら，この引き数 は名前付きでなければならない)
method
平 均・分 散 を 標 準 的 な 推 定 す る な ら"moment"
，最 尤 推 定 な ら"mle"
，cov.mve()
を使えば"mve"
，t
分布に基づく頑健推定量を使えば"t"
CV
もし真なら，「一時に一つ取り除く」クロスバリデーションに対する結果(
クラス と事後確率)
を返す．もし事前分布を推定するなら，全データセット中の比率が 使われることを注意するnu method="t"
に対する自由度...
他のメソッドへ(
から)
引き渡される引き数 返り値： クラス"qda"
のリストで，次の成分を含む：prior
使用された事前確率means
グループ平均scaling
各グループi
に対して，scaling[,,i]
はグループ内共分散行列が球形に なるように観測値を変換する配列であるldet
偏差行列の行列式の対数値の半分のベクトルlev
グループ因子の水準terms (
もしformula
がモデル式なら)
モデル式を要約するモード表現とクラス項のオブジェクト
call
マッチした関数呼び出しclass MAP(
最大事後確率)
分類(
因子)
posterior
クラスに対する事後確率# あやめデータiris3を使用
> library(MASS) # MASSパッケージを読み込む
> tr <- sample(1:50, 25) # 訓練用データ25例の添字をランダムに選ぶ
# 訓練用データ(3種類のアヤメからそれぞれ添字trに該当するものを抜き出す)
> train <- rbind(iris3[tr, , 1], iris3[tr, , 2], iris3[tr, , 3])
# テスト用データ(残り，3種類のアヤメ25例の4種類の性質)
> test <- rbind(iris3[-tr, , 1], iris3[-tr, , 2], iris3[-tr, , 3])
> cl <- factor(c(rep("s", 25), rep("c", 25), rep("v", 25))) # 種類を表す因子
> zq <- qda(train, cl) # 2次判別解析実行
> str(zq) # 返り値の持つ全情報の簡易要約
List of 8
$ prior : Named num [1:3] 0.333 0.333 0.333 # 事前確率 ..- attr(*, "names")= chr [1:3] "c" "s" "v" # 対応する名前
$ counts : Named int [1:3] 25 25 25 # 種類毎の度数 ..- attr(*, "names")= chr [1:3] "c" "s" "v" # 種類名
$ means : num [1:3, 1:4] 6.00 5.04 6.47 2.79 3.47 ... # グループ平均 ..- attr(*, "dimnames")=List of 2 # 次元属性は行列 .. ..$ : chr [1:3] "c" "s" "v" # 因子名
.. ..$ : chr [1:4] "Sepal L." "Sepal W." "Petal L." "Petal W." # 変数名
# 観測値変換配列
$ scaling: num [1:4, 1:4, 1:3] -1.79 0.00 0.00 0.00 1.01 ...
..- attr(*, "dimnames")=List of 3 # 次元属性は3次元配列
# 第一(変数)次元名 .. ..$ : chr [1:4] "Sepal L." "Sepal W." "Petal L." "Petal W."
.. ..$ : chr [1:4] "1" "2" "3" "4" # 第二次元名 .. ..$ : chr [1:3] "c" "s" "v" # 第三(種類)次元名
$ ldet : num [1:3] -10.70 -13.01 -8.84
$ lev : chr [1:3] "c" "s" "v" # 種類因子の3水準
$ N : int 75 # 訓練用の観測値総数
$ call : language qda(x=train, grouping=cl) # 呼び出し式 - attr(*, "class")= chr "qda" # 返り値の属性
# 結果の2次判別子を用いてテスト用データを判別(つまり4種類の性質から種類を当てる)
> predict(zq, test)$class
[1] s s s s s s s s s s s s s s s s s s s s s s s s s c c c c c c c c v c c [37] c c c c v c c c c c c c c c v v v v v v v v v v v v v v v v v v v v v v [73] v v v
Levels: c s v
5.6.4 2 次判別関数による予測 predict.qda()
この関数はクラス
"qda"
に対する総称的関数predict()
の一つのメソッドである．これは適当なクラスのオブジェクト
x
に対して呼び出しpredict()
を行うか，または 直接にx
のクラスに無関係に呼び出しpredict.qda(x)
を行うことにより起動される．newdata
中の欠損値は，線形判別子が計算できなければNA
を返す．もしnewdata
が省略され，欠損値処理
na.action
が欠損値の無視を指示していれば，これらは予測に対 して無視される．2
次判別分析関数qda()
による判別子を用いて多変量観測値の判別を 行う．書式：
#
クラス"qda"
に対するS3
メソッドpredict(object, newdata, prior = object$prior,
method = c("plug-in","predictive","debiased","looCV"), ...)
引数：
object
クラス”qda”
のオブジェクトnewdata
判別すべき例のデータフレーム，もしくは，object
がモデル式を持てば，使われた変数と同じ列名を持つデータフレーム．ベクトルは行ベクトルと見なされ
る．もし
newdata
が無ければ，引数qda
のオブジェクトの当てはめに使われたデータを探そうと試みる
prior
クラスに対する事前確率．既定では訓練データにおける比率か，引数lda
の呼び出しで設定された値
method
これはパラメータの推定をどうするか決める．既定の”plug-in”
では通常の不偏推定量が使われ，正確な値と仮定される．
”debiased”
では対数事後確率の不偏 推定量が使われる．”predictive”
では曖昧な事前分布を用い，パラメータ推定量 が「integrate out
」される...
他のメソッドから(
へ)
引き渡される追加引数返り値： 次の成分を持つリスト：
class MAP(
最大事後確率)
判別(
因子)
posterior
クラスに対する事後確率x
最大dimen
個の判別変数に関するテスト例のスコア# データiris3を使用
> library(MASS) # MASSパッケージを読み込む
> tr <- sample(1:50, 25) # 訓練用データ25例の添字をランダムに選ぶ
# 訓練用データ(3種類のアヤメからそれぞれ添字 tr に該当するものを抜き出す)
> train <- rbind(iris3[tr, , 1], iris3[tr, , 2], iris3[tr, , 3])
# テスト用データ(残り，3種類のアヤメ25例の4種類の性質)
> test <- rbind(iris3[-tr, , 1], iris3[-tr, , 2], iris3[-tr, , 3])
> cl <- factor(c(rep("s", 25), rep("c", 25), rep("v", 25))) # 種類を表す因子
> zq <- qda(train, cl) # 2次判別解析実行
# 結果の2次判別子を用いてテスト用データを判別(つまり4種類の性質から種類を当てる)
> predict(zq, test)$class
[1] s s s s s s s s s s s s s s s s s s s s s s s s s c c c c c c (以下省略) [39] c c c c c c c c c c c c v v v v v v v v v v v v v v v v v v v (以下省略) Levels: c s v
# 判別結果は正しかったか? (正答率 72/75)
# 注意：訓練用データをランダムに選ぶので毎回結果が異なる可能性がある
> (cl == predict(z, test)$class)
[1] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [13] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [25] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE FALSE FALSE TRUE [37] TRUE TRUE TRUE TRUE TRUE FALSE TRUE TRUE TRUE TRUE TRUE TRUE [49] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE [61] TRUE TRUE TRUE TRUE TRUE TRUE FALSE TRUE TRUE TRUE TRUE TRUE [73] TRUE TRUE TRUE
第 6 章
線形回帰モデル
R
は回帰分析関連の関数を多数持ち，統計解析機能の中心的位置を占める．以下では，線形