差分

このページの2つのバージョン間の差分を表示します。

この比較画面へのリンク

次のリビジョン		前のリビジョン
mselab:2016:stat:week2 [2016/11/08 09:21] – created watalumselab:2016:stat:week2 [2016/11/29 11:57] (現在) – watalu
行 1: 行 1:
-=== 第2週 ===+==== 第2週 ====
  
 今週の{{:mselab:2013:stat:week3:note-statistics-3-20131029.pdf|配付資料}}。グラフと表のみによるデータ解析と比較的単純な学習機械によるデータ解析に続いて、今回は次の2つの課題に取り組んでもらう。 今週の{{:mselab:2013:stat:week3:note-statistics-3-20131029.pdf|配付資料}}。グラフと表のみによるデータ解析と比較的単純な学習機械によるデータ解析に続いて、今回は次の2つの課題に取り組んでもらう。
  
-  - より柔軟な学習機械と単純な学習機械との比較+PPDACの5つの要素のうち、Pが変わって、Dにデータが追加されて、Aも変更になる。 
 + 
 +  - P(roblem)は「V86の契約の有無と他の変数との関係の分析」から「V86の契約の有無の他の変数からの予測モデルの構築」に変わる 
 +    - 先週は「学習機械の当てはめ」によるデータの分析、今週は「学習機械による予測結果」を評価してどれを用いるかを検討 
 +    - 各学習機械には様々な「パラメータ」があり、それを変えることで学習結果も予測精度も変わる 
 +    - 線形判別分析(lda)、二次判別分析(qda)、決定木(rpart)だけでなく、より柔軟な学習機械も用いる候補に含める 
 +  - D(ata)に先週に解析してもらったtic.learnに加えて、tic.evalを追加する 
 +    - tic.learnは学習用、tic.evalは検証用 
 +  - A(nalysis)に「予測精度の評価」が加わる 
 +    - 前回は学習しただけ、今回は学習したモデルを予測に用いる(仮想的な現場投入) 
 +    - tic.learnで学習したモデルを、tic.evalの予測に用いて、予測精度を評価する 
 +  - C(onclusion)は「どのモデルが予測精度が良いか」とその考察
   - キャンペーンの提案書の作成   - キャンペーンの提案書の作成
 +    - 一番上手に予測できた「制御パラメータを調整済み」の学習機械に基づいて、予測精度を自慢しつつ、効率が良いと思われるキャンペーンを提案する提案書を起草せよ。
 +
 +１回目と２回目とで、各手法に対する印象が少し異なってくると嬉しい。
 +
 +=== 学習と予測 ===
 +== 学習 ==
 +
 +例えばあるデータに回帰分析を行う。
 +
 +<code>
 +tic.lm <- lm(V86~.,data=tic.learn)
 +</code>
 +
 +学習理論の言葉では、これは「回帰分析」という学習機械を「最小二乗法」という学習アルゴリズムで学習させたことに相当する。
 +そして学習の精度は、「実際の個々の観測値(tic.learn$V86)」と、モデルを推定した結果による個別の平均の推定値である「当てはめ値(fitted(tic.lm))」とを比較して検討する。
 +
 +まずは横軸を観測値、縦軸を当てはめ値にした散布図。
 +<code>
 +plot(tic.learn$V86, fitted(tic.lm),
 +     xlab="V86",
 +     ylab="Fitted Value")
 +</code>
 +
 +あまりよくわからないので、箱ひげ図に直してみる。
 +<code>
 +plot(fitted(tic.lm)~as.factor(tic.learn$V86),
 +     xlab="V86",
 +     ylab="Fitted Value")
 +</code>
 +
 +ヒストグラムを描いてもいい。
 +<code>
 +par(mfrow=c(2,1))
 +hist(fitted(tic.lm)[tic.learn$V86==0],
 +     main="V86=0",
 +     xlab="Fitted Value")
 +hist(fitted(tic.lm)[tic.learn$V86==1],
 +     main="V86=1",
 +     xlab="Fitted Value")
 +</code>
 +
 +あまり変わらないように見える。
 +平均の差のt検定(Welchの検定)ぐらいをやってみると、この差に統計的な有意差はあることが分かる。
 +
 +<code>
 +t.test(fitted(tic.lm)[tic.learn$V86==0],fitted(tic.lm)[tic.learn$V86==1],var.equal=F)
 +</code>
 +さりとて、当てはめ値が0付近なのはあまりよい回帰分析とは言えない。
 +
 +当てはまりの良さ(悪さ)は、残差という量でも検討できる。残差の定義は、「実際の個々の観測値(tic.learn$V86)-当てはめ値(fitted(tic.lm))」であり、データの個々の観測値の変動のうちのモデルで説明しきれなかった部分である。
 +関数residuals()を用いて計算させることができる。
 +<code>
 +plot(residuals(tic.lm))
 +</code>
 +単に、レコード番号順に残差を計算させるとこうなる。
 +
 +残差の大きさや傾向を検討をするときには、観測値と残差の傾向をみる散布図
 +<code>
 +plot(tic.learn$V86,residuals(tic.lm),
 +     xlab="V86",
 +     ylab="Residuals")
 +</code>
 +や、当てはめ値と残差の散布図
 +<code>
 +plot(fitted(tic.lm),residuals(tic.lm),
 +     xlab="Fitted Values",
 +     ylab="Residuals")
 +</code>
 +を見て、何か説明できていない傾向が残っていないかどうかを検討する。
 +他に、てこ比を見たり、半正規プロットなどを検討することもある。
 +
 +重回帰分析に限って、上で紹介したグラフの多くは
 +<code>
 +plot(tic.lm)
 +</code>
 +で自動的に描画される。
 +
 +他にも要約の出力
 +<code>
 +summary(tic.lm)
 +</code>
 +を見て、学習精度を検討する。
 +
 +== 予測 ==
 +
 +予測とは、学習した結果(学習済みの学習機械)を別のデータに基づく予測に用いることを言う。
 +そのためには、学習用データと同じ構造を持つ、別のデータを用意する。
 +予測精度を評価するので、ここでは評価用データと呼ぶ。
 +今年の実験ではtic.learnが学習用データ、tic.evalが評価用データであり、前回のデータ読み込みの手順で両方ともダウンロードされる。
 +
 +予測値の算出には関数predictを用いる。学習済み機械tic.lmに、tic.evalの中の各レコードのV1からV85に基づいてV86を予測させるにはpredict()を
 +<code>
 +predict(tic.lm, newdata=tic.eval)
 +</code>
 +と用いる。学習にtic.evalを用いていないことから、この結果は予測値として扱える。
 +
 +lmを用いて予測したので、
 +<code>
 +v86.lm <- predict(tic.lm, newdata=tic.eval)
 +</code>
 +と名付けておく。
 +
 +<code>
 +plot(tic.eval$V86, v86.lm,
 +     xlab="V86",
 +     ylab="Predicted Value")
 +</code>
 +
 +ヒストグラムを描いてもいい。
 +<code>
 +par(mfrow=c(2,1))
 +hist(v86.lm[tic.eval$V86==0],
 +     main="V86=0",
 +     xlab="Predicted Value")
 +hist(v86.lm[tic.eval$V86==1],
 +     main="V86=1",
 +     xlab="Predicted Value")
 +</code>
 +なんかダメそうに見える。
 +
 +さて残差は「当てはめ値-観測値」だったが、予測の良さを測るには予測誤差を用いる。
 +予測誤差の定義は「予測値-正解」である。
 +<code>
 +plot(tic.eval$V86, tic.eval$V86-v86.lm,
 +     xlab="V86",
 +     ylab="Prediction Error")
 +</code>
 +予測誤差にも傾向があり、とてもではないが、重回帰分析はこのままではよい予測を与えてくれる気がしない。
 +箱ひげ図を描くまでもないが、念のため。
 +
 +<code>
 +plot((tic.eval$V86-v86.lm)~as.factor(tic.eval$V86),
 +     xlab="V86",
 +     ylab="Prediction Error")
 +</code>
 +
 +あまり変わらないように見える。
 +平均の差のt検定(Welchの検定)ぐらいをやってみると、この差に統計的な有意差はあることが分かる。
 +
 +<code>
 +t.test((tic.eval$V86-v86.lm)[tic.eval$V86==0],
 +       (tic.eval$V86-v86.lm)[tic.eval$V86==1],
 +       var.equal=F)
 +</code>
 +
 +まあだめそう。。。
  
-== 学習機械の予測精度の評価 ==+== 線形判別分析(lda)の場合 ==
  
 学習機械は予測精度が高くて初めて、有用性が認められる。そのためには「学習に用いなかったデータ」に関して予測を行い、精度を評価する必要がある。 学習機械は予測精度が高くて初めて、有用性が認められる。そのためには「学習に用いなかったデータ」に関して予測を行い、精度を評価する必要がある。
行 77: 行 234:
 </code> </code>
 のようにクロス集計を実行して、予測結果(v86.lda$class)と実際(tic.eval$V86)を比較する。 のようにクロス集計を実行して、予測結果(v86.lda$class)と実際(tic.eval$V86)を比較する。
 +
 +== 二次判別分析(qda)の場合 ==
  
 他の学習機械についても同様に 他の学習機械についても同様に
行 95: 行 254:
 などと、ひとつひとつ増やしていってみるか、例えば などと、ひとつひとつ増やしていってみるか、例えば
 <code> <code>
-tic.qda <- qda(V86~V1,data=tic.learn[,c(1:30,86)])+tic.qda <- qda(V86~.,data=tic.learn[,c(1:30,86)])
 </code> </code>
 のように最初の30変数を使って86番目を予測する、とするなどの工夫が必要になる。 のように最初の30変数を使って86番目を予測する、とするなどの工夫が必要になる。
 +
 +== 決定木(rpart)の場合 ==
  
 再帰的分割(rpart)を用いる場合は、次のようなデータの変換が必要になる。 再帰的分割(rpart)を用いる場合は、次のようなデータの変換が必要になる。
行 107: 行 268:
 <code> <code>
 tic.rpart <- rpart(V86~.,data=tic.learn.2) tic.rpart <- rpart(V86~.,data=tic.learn.2)
-v86.rpart <- predict(tic.rpart,nedata=tic.eval,type="class") +v86.rpart <- predict(tic.rpart,newdata=tic.eval,type="class") 
-table(v86.rpart$class, tic.eval$V86)+table(v86.rpart, tic.eval$V86) 
 +</code> 
 + 
 +rpartの制御はたとえば 
 +<code> 
 +tic.rpart <- rpart(V86~.,data=tic.learn.2, control=rpart.control(cp=0.00001, minsplit=5))
 </code> </code>
 +のように指定する。
  
-謝判別率(誤った予測をした割合)などを算出してみよ。+誤判別率(誤った予測をした割合)などを算出しつつ、各学習機械の制御パラメータを調整してみよ。
  
-== 柔軟な学習機械を使ってみる ==+=== 柔軟な学習機械を使ってみる ===
  
 ひとつ目の課題は、より柔軟な学習機械を適用してもらい、先週までの方法との違いを検討してもらう。こちらについては今回も、同志社大学の金(じん)先生が公開されてらっしゃる ひとつ目の課題は、より柔軟な学習機械を適用してもらい、先週までの方法との違いを検討してもらう。こちらについては今回も、同志社大学の金(じん)先生が公開されてらっしゃる
行 124: 行 291:
 サポートベクターマシン(ksvm)はkernlabパッケージに、バギング(bagging)はipredパッケージに、ブースト(ada)はadaパッケージに、ランダムフォレスト(randomForest)はrandomForestパッケージに、それぞれ含まれているので、まずは必要なパッケージをインストールする。 サポートベクターマシン(ksvm)はkernlabパッケージに、バギング(bagging)はipredパッケージに、ブースト(ada)はadaパッケージに、ランダムフォレスト(randomForest)はrandomForestパッケージに、それぞれ含まれているので、まずは必要なパッケージをインストールする。
 <code> <code>
-bank <- read.table("c:¥¥Users¥¥Student¥¥bank-full.csv",  
-                   header=TRUE,  
-                   sep=";") 
 Sys.setenv("http_proxy"="http://130.153.8.19:8080/") Sys.setenv("http_proxy"="http://130.153.8.19:8080/")
 install.packages("kernlab", dependencies=TRUE) install.packages("kernlab", dependencies=TRUE)
行 140: 行 304:
 サポートベクターマシンに学習させる実行例： サポートベクターマシンに学習させる実行例：
 <code> <code>
-bank.svm <- ksvm(y~., data=bank, +tic.svm <- ksvm(V86~., data=tic.learn.2, 
                  kernel="rbfdot",                  kernel="rbfdot",
                  kpar=list(sigma=0.01))                  kpar=list(sigma=0.01))
-bank.svm.p<-predict(bank.svm)  +v86.svm <- predict(tic.svm, newdata=tic.eval)  
-table(bank.svm.p, bank$y)+table(v86.svm,tic.eval$V86)
 </code> </code>
 +kparの中身を変更しないと、全く予測精度が稼げない感じ。
  
 バギングに学習させる実行例： バギングに学習させる実行例：
 <code> <code>
-bank.bag<-bagging(y~., data=bank,+tic.bagging <-bagging(V86~., data=tic.learn.2,
                   nbagg=40)                    nbagg=40) 
-bank.bag.p <- predict(bank.bag,type= "class")  +v86.bagging <- predict(tic.bagging,newdata=tic.eval,type= "class")  
-table(bank.bag.p, bank$y)+table(v86.bagging, tic.eval$V86)
 </code> </code>
 +これもバッグの数が40でいいか分からない微妙な感じ。
  
 Adaブーストに学習させる実行例： Adaブーストに学習させる実行例：
 <code> <code>
-bank.ada <- ada(y~., data=bank,+tic.ada <- ada(V86~., data=tic.learn.2,
                 iter=20)                  iter=20) 
-bank.ada.p <- predict(bank.ada,type= "vector")  +v86.ada <- predict(tic.ada, newdata=tic.eval, type= "vector")  
-table(bank.ada.p, bank$y)+table(v86.ada, tic.eval$V86)
 </code> </code>
 +これも繰り返し数(iter)が20ではぜんぜんだめな感じ。
  
 ランダムフォレストに学習させる実行例： ランダムフォレストに学習させる実行例：
 <code> <code>
-bank.rf <- randomForest(y~., data=bank)  +tic.randomForest <- randomForest(V86~., data=tic.learn.2)  
-bank.rf.p <- predict(bank.rf)  +v86.randomForest <- predict(tic.randomForest, newdata=tic.eval)  
-table(bank.rf.p, bank$y)+table(v86.randomForest, tic.eval$V86)
 </code> </code>
  
-以上は「弱い学習機械」の指定と、個々の学習機械の設定パラメータが、デフォルトのままであることとは注意しておく。パラメータのチューニングが必要なはずである。+以上は「弱い学習機械」の指定と、個々の学習機械の設定パラメータが、デフォルトのままであることとは注意しておく。それぞれの予測結果を見れば、パラメータのチューニングが必要なことは明らか。
  
 なお、この課題では「過学習」について言及していないが、過学習は大事な問題であるので、各自、少し調べて気にすることを進める。 なお、この課題では「過学習」について言及していないが、過学習は大事な問題であるので、各自、少し調べて気にすることを進める。
  
-== 提案書 ==+=== 提案書 ===
  
 2つ目の課題は、前回までのレポートの考察として、どのような顧客層に重点的にキャンペーンを展開するのがよいか、提案書を起草してもらう。「提案書」という書類の形式については、成書を参照するのがよいが、参考までに幾つかのウェブサイトへのリンクを掲げておく。 2つ目の課題は、前回までのレポートの考察として、どのような顧客層に重点的にキャンペーンを展開するのがよいか、提案書を起草してもらう。「提案書」という書類の形式については、成書を参照するのがよいが、参考までに幾つかのウェブサイトへのリンクを掲げておく。