文書の過去の版を表示しています。

目次

概要

今週の実験の内容は

  • データマイニング：保険のデータのV86の予測モデルの構築
    • 回帰分析
    • ロジスティック回帰
    • 決定木

である。3週間の流れがデータ分析(第1週)を行って、回帰分析等によるモデル構築(第2週)、そして第3週の他の手法も組み合わせたモデル構築と結論へと繋がるため、今週は多変量解析との内容の重複を避けずにおいた。また同科目の履修を前提としていないため、回帰分析の学習も自習内容として含めてある。

  • 先週に続いてデータの把握、特にV86を中心に。
  • Rコマンダーを用いた回帰分析の2つの課題(自習、練習に相当)
  • 解析データを用いた1つの課題(本番)

に取り組んで貰う。

実験の流れ

  1. 配付資料とRの実行結果を見比べながら、lm関数を用いた回帰分析で出力される情報のどれが、配付資料のどれに対応するのかを把握する
    • 回帰係数の推定値: Estimate
      • 切片: Intercept
    • 寄与率: R-Squared
    • 自由度調整済み寄与率: Adjusted R-squared
    • F値: F-statistic
    • 自由度: DF, degrees of freedom
    • t値: t value
    • P値: Pr(>|t|), p-value
    • 標準誤差: Std. Error, standard error
    • 残差: Residuals
    • (回帰)係数: Coefficients
    • てこ比: Leverage
    • 標準化残差: Studentized Residuals
    • 変数増減法
  2. rpart関数を用いた決定木分析の出力を見ながら、決定木分析はどのようなモデルを構築するのかを検討する
  3. 保険データの回帰分析、に取り組む (保険データの回帰分析)
    1. 回帰係数の推定
    2. 変数の加工(今回初)
    3. 分散分析によるモデルの有意性の検討や回帰係数の有意性の検討
    4. てこ比や標準化残差などの検討
    5. 変数の増減 (stepAIC関数を使う?)
    6. 以上を繰り返す
  4. 保険のデータの決定木分析、に取り組む (時間外課題)

データの説明

TIC2000

tic.data.txtからの要約。

  • CoIL 2000 Challengeで用いられた保険会社の顧客に関するデータ。86個の変数は、契約状況(V44-V85)と社会人口統計学的な変数(V1-V43)を含んでいる。この調査は “Can you predict who would be interested in buying a caravan insurance policy and give an explanation why?” という問いに答えるように集められた。
  • このデータはオランダのデータマイニング会社Sentinent Machine Researchから提供され、現実のビジネスの問題に基づいている。学習用データ(ticdata2000.txt)は5000レコードでcaravan insurance policyの契約の有無(V86)を含んでおり、検証用データ(ticeval2000.txt)は4000レコードで契約の有無(V86)は含んでいない。検証用データの正解は、CoIL 2000 Challengeの開催時には公開されていなかったが、現在はテストデータ(tictest2000.txt)として公開されている。
  • V1-V43のうち、コード化が指定されていない変数はすべて、郵便番号の一桁目のエリアを指している。たとえばV30が9ならばその顧客は郵便番号が9で始まるエリアに家を借りていることを、V31が5ならば郵便番号が5のエリアに持ち家があることを意味する。職業、社会層などもすべて、該当するエリアの箇所が郵便番号の一桁目で埋まっている。
変数

dictionary.txtからの抜粋と要約、の日本語版。

変数分類メモ
V1顧客分類２L0でコード化されている、数字の大きさに意味なし
V2住居数大きいほど住む箇所が多い
V3世帯構成員数の平均人数
V4世帯構成員の平均年齢L1でコード化されている、年齢
V5顧客分類１L2でコード化されている、数字の大きさに意味なし
V6-V9宗教L3でコード化されている、V6+V7+V8+V9は9から12の間。それぞれの宗教を信じる割合？
V10-V13結婚場所を表す変数, 例えばV10が0ならば無し？
V14-V15世帯の大きさL3でコード化されている、なぜかV14+V15は10以下。割合？
V16-V18教育水準L3でコード化されている、なぜかV16+V17+V18はほぼ10、それぞれの年数？割合？
V19-V24職業L3でコード化されている、なぜかV19+V20+V21+V22+V23+V24は9から13の間
V25-V29社会層L3でコード化されている、なぜかV25+V26+V27+V28+V29は9から12の間
V30-V31住居L3でコード化されている、なぜかV30+V31は9か10
V32-V34自動車L3でコード化されている、なぜかV32+V33+V34は9から11の間
V35-V36健康保険L3でコード化されている、なぜかV35+V36は9か10
V37-V41収入L3でコード化されている、なぜかV37+V38+V39+V40+V41は9から13の間
V42平均収入L3でコード化されている
V43購買力L3でコード化されている、1から8の間。
V44-V64各種保険支払い額L4でコード化
V65-V85各種保険契約件数件数

メモの確認用のコード。 

table((tic.learn$V16+tic.learn$V17+tic.learn$V18))
table((tic.learn$V19+tic.learn$V20+tic.learn$V21+tic.learn$V22+tic.learn$V23+tic.learn$V24))
table((tic.learn$V25+tic.learn$V26+tic.learn$V27+tic.learn$V28+tic.learn$V29))
table(tic.learn$V30+tic.learn$V31)
table(tic.learn$V32+tic.learn$V33+tic.learn$V34)
table(tic.learn$V35+tic.learn$V36)
table(tic.learn$V37+tic.learn$V38+tic.learn$V39+tic.learn$V40+tic.learn$V41)
各変数のコーディング

L0:分類を表す数字なので、大小関係に意味がなく、名義尺度である。そのままでは説明変数にならない。

ValueLabel
1High Income, expensive child
2Very Important Provincials
3High status seniors
4Affluent senior apartments
5Mixed seniors
6Career and childcare
7Dinki's (double income no kids)
8Middle class families
9Modern, complete families
10Stable family
11Family starters
12Affluent young families
13Young all american family
14Junior cosmopolitan
15Senior cosmopolitans
16Students in apartments
17Fresh masters in the city
18Single youth
19Suburban youth
20Etnically diverse
21Young urban have-nots
22Mixed apartment dwellers
23Young and rising
24Young, low educated
25Young seniors in the city
26Own home elderly
27Seniors in apartments
28Residential elderly
29Porchless seniors: no front yard
30Religious elderly singles
31Low income catholics
32Mixed seniors
33Lower class large families
34Large family, employed child
35Village families
36Couples with teens 'Married with children'
37Mixed small town dwellers
38Traditional families
39Large religous families
40Large family farms
41Mixed rurals

L1:大きさが年齢の順なので、そのまま説明変数に使える。

120-30 years
230-40 years
340-50 years
450-60 years
560-70 years
670-80 years

L2:数字は分類を表すだけなので、連続尺度でも順序尺度でもなく、名義尺度。そのままでは説明変数にならない。

1Successful hedonists
2Driven Growers
3Average Family
4Career Loners
5Living well
6Cruising Seniors
7Retired and Religeous
8Family with grown ups
9Conservative families
10Farmers

L3:順序尺度。このまま連続尺度の説明変数として用いる。

00%
11 - 10%
211 - 23%
324 - 36%
437 - 49%
550 - 62%
663 - 75%
776 - 88%
889 - 99%
9100%

L4: 順序尺度。今回はこのまま連続尺度の変数として用いる。

0f 0
1f 1 - 49
2f 50 - 99
3f 100 - 199
4f 200 - 499
5f 500 - 999
6f 1000 - 4999
7f 5000 - 9999
8f 10.000 - 19.999
9f 20.000 - ?
参考

kernlabパッケージに、加工済みのデータが入っていて、それを使うこともできる。 

install.packages(c("kernlab"), dependencies=TRUE)
tic.learn <- ticdata[1:5822,]
tic.eval <- ticdata[5823:9822,]

今回の課題

概要
  • 回帰分析について学ぶ。
  • TICデータは、V65からV86までが22種類の保険商品の契約件数である。今回はV86の契約に漕ぎ着けるためのモデル構築、を目的とする。
  • 重回帰分析を行い、分析結果を考察せよ。
  • 重回帰分析の結果を用いて、指定された保険商品の契約に関する予測式を構築せよ。
  • 線形モデルと重回帰分析を用いて、このデータをデータマイニングすることについて、考察せよ。

時間内課題１：回帰分析の自習

単回帰分析

配付資料の表4.1のデータの入力は次の通り。 

x <- c(2.2,4.1,5.5,1.9,3.4,2.6,4.2,3.7,4.9,3.2)
y <- c(71,81,86,72,77,73,80,81,85,74)
table.4.1 <- data.frame(x=x,y=y)

回帰分析の実行にはlm関数を用いる。 

lm(y~x, data=table.4.1)

「y~x」は、変数yを変数xで説明する回帰分析を行うための表現。これは <jsmath> Y = \beta_0+\beta_1 X+\epsilon </jsmath> という回帰式を推定せよ、という意味である。

回帰分析の詳細を示すには、lm関数の実行結果を別の変数lm.4.1などに代入してから、詳細を出力させるsummary関数を用いる。 

lm.4.1 <- lm(y~x, data=table.4.1)
print(lm.4.1)
summary(lm.4.1)

問１：このsummary関数の出力結果を、配付資料と対比させて理解せよ。

回帰分析の結果を診断するために、幾つかのグラフを出力する機能がある。lm関数の実行結果をplot関数にかければ、そのようなグラフが4枚、順次出力される。 

plot(test.lm)

問２：このplot関数の出力結果を、配付資料と対比させて理解せよ。

重回帰分析

表5.1のデータは次のように入力する。 

x.1 <- c(51,38,57,51,53,77,63,69,72,73)
x.2 <- c(16,4,16,11,4,22,5,5,2,1)
y <- c(3.0,3.2,3.3,3.9,4.4,4.5,4.5,5.4,5.4,6.0)
table.5.1 <- data.frame(x.1=x.1,x.2=x.2,y=y)

回帰分析の実行には、説明変数が複数あっても(重回帰でも)lm関数を用いる。 

lm(y~x.1+x.2, data=table.5.1)

今度は説明変数が2つあるので、「y~x.1+x.2」となる。これは <jsmath> Y = \beta_0+\beta_1 X_1 + \beta_2 X_2+\epsilon </jsmath> という回帰式を推定せよ、という意味である。

回帰分析の詳細を示すには、lm関数の実行結果を別の変数lm.5.1などに代入してから、詳細を出力させるsummary関数を用いる。 

lm.5.1 <- lm(y~x.1+x.2, data=table.5.1)
print(lm.5.1)
summary(lm.5.1)

問４：このplot関数の出力結果を、配付資料と対比させて理解せよ。

重回帰応用(水準変数が説明変数に含まれる場合)

上のデータを、広さを「広め(w)」「狭め(n)」とし、築年数も「新しめ(new)」「古め(old)」にする 

x.1 <- c("n","n","n","n","n","w","w","w","w","w")
x.2 <- c("old","new","old","old","new","old","new","new","new","new")
y <- c(3.0,3.2,3.3,3.9,4.4,4.5,4.5,5.4,5.4,6.0)
table.5.1.c <- data.frame(x.1=x.1,x.2=x.2,y=y)

問５：このデータで、回帰分析を行い、グラフも作成して、先のモデルと何が異なるか、検討せよ。 

lm.5.1.c <- lm(y~x.1+x.2, data=table.5.1.c)
print(lm.5.1.c)
summary(lm.5.1.c)
plot(lm.5.1.c)

時間内課題２：決定木分析の自習

まず、次の一行を実行しておく。 

library(mvpart)

上の実行例で「lm」とあるところをすべて「rpart」で置き換える。 

rpart.4.1 <- rpart(y~x, data=table.4.1)
print(rpart.4.1)
summary(rpart.4.1)
rpart.5.1 <- rpart(y~x.1+x.2, data=table.5.1)
print(rpart.5.1)
summary(rpart.5.1)
rpart.5.1.c <- rpart(y~x.1+x.2, data=table.5.1.c)
print(rpart.5.1.c)
summary(rpart.5.1.c)

問６：これらがどのようなモデルか、データと出力に照らして検討せよ。 特にlm.5.1、lm.5.1.c、rpart.5.1、rpart.5.1.cの間の違いは考察せよ。

レポート課題

この課題ではMASSライブラリのみ、使う可能性がある。 

library(MASS)

1つ目のデータは、Rに次の命令を実行させておく。 

x <- c(2.2,4.1,5.5,1.9,3.4,2.6,4.2,3.7,4.9,3.2)
y <- c(71,81,86,72,77,73,80,81,85,74)
data.1 <- data.frame(x=x,y=y)
rm(x,y)

２つ目のデータは、Rに次の命令を実行させておく。 

x1 <- c(51,38,57,51,53,77,63,69,72,73)
x2 <- c(16,4,16,11,4,22,5,5,2,1)
y <- c(3.0,3.2,3.3,3.9,4.4,4.5,4.5,5.4,5.4,6.0)
data.2 <- data.frame(x1=x1,x2=x2,y=y)
rm(x1,x2,y)

演習で用いる保険データは、Rに次の命令を実行させておく。 

Sys.setenv("http_proxy"="http://130.153.8.66:8080/")
tic.learn <- read.table("http://kdd.ics.uci.edu/databases/tic/ticdata2000.txt")
tic.eval <- read.table("http://kdd.ics.uci.edu/databases/tic/ticeval2000.txt")
tic.test <- read.table("http://kdd.ics.uci.edu/databases/tic/tictgts2000.txt")
tic.eval <- cbind(tic.eval, tic.test)
colnames(tic.eval)[86] <- "V86"
rm(tic.test)

あとは、上の回帰分析と同様に行えばよいのだが、V1とV5が特殊なので、これらを例に「分類変数のグループ化」について、少し手順を解説する。

分類変数のグループ化

  • 分類尺度の変数の、分類数を減らす (回帰分析でも決定木分析でも数千のサンプルで40以上の分類は細かすぎる)

V1とV5を「分類尺度」の変数に変換し、それぞれV1fおよびV5fと名付ける。 

tic.learn <- read.table("http://kdd.ics.uci.edu/databases/tic/ticdata2000.txt")
tic.learn$V1f <- as.factor(tic.learn$V1)
tic.learn$V5f <- as.factor(tic.learn$V5)

以下、幾つかのグラフを描いているので、1行ずつ貼って実行していくことを進める。

V1やV5のヒストグラムは、hist関数を用いて 

hist(tic.learn$V1)
hist(tic.learn$V5)

などで描ける。一方、V1fやV5fのヒストグラムは 

hist(tic.learn$V1f)
hist(tic.learn$V5f)

を実行するとエラーになる。これを得るには、table関数で集計してから、barplot関数で棒グラフを描く。 

barplot(table(tic.learn$V1f))
barplot(table(tic.learn$V5f)

更にカテゴリごとの契約の有無の割合を見るには 

barplot(table(tic.learn$V86,tic.learn$V1f))
barplot(table(tic.learn$V86,tic.learn$V5f))

とする。色の濃い部分が契約無、色の薄い部分が契約有である。 

library(mvpart)
table(tic.learn$V1f, tic.learn$V86)
table(tic.learn$V5f, tic.learn$V86)
rpart(V86~V1f, data=tic.learn, control=rpart.control(cp=0))
<code>

 V1の値を幾つかのグループにまとめる

  * 契約にあまり寄与しないカテゴリには、コード0を割り当ててしまおう。
  * 契約にしそうなカテゴリは残そう。

決定木(実際には回帰木)を用いて、V86の増加に寄与しそうなカテゴリのみを残す。

<code>
library(mvpart)
rpart(V86~V1f, data=tic.learn, control=rpart.control(cp=0))

剪定せずにそのまま出力した回帰樹は次の通り。 

n= 5822 

node), split, n, deviance, yval
      * denotes terminal node

  1) root 5822 327.1989000 0.05977327  
    2) V1f=2,4,5,7,9,10,11,13,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41 4880 220.0654000 0.04733607  
      4) V1f=4,5,9,15,16,17,18,19,21,22,23,24,25,26,27,28,29,30,31,35,40,41 2075  56.3788000 0.02795181  
        8) V1f=15,16,17,18,19,21,23,25,26,27,28,29,40,41 885  14.7457600 0.01694915  
         16) V1f=15,16,17,18,19,21,28,40 163   0.0000000 0.00000000 *
         17) V1f=23,25,26,27,29,41 722  14.6883700 0.02077562  
           34) V1f=23 251   3.9362550 0.01593625 *
           35) V1f=25,26,27,29,41 471  10.7431000 0.02335456  
             70) V1f=26,27 98   1.9591840 0.02040816  
              140) V1f=27 50   0.9800000 0.02000000 *
              141) V1f=26 48   0.9791667 0.02083333 *
             71) V1f=25,29,41 373   8.7828420 0.02412869  
              142) V1f=29 86   1.9534880 0.02325581 *
              143) V1f=25,41 287   6.8292680 0.02439024 *
        9) V1f=4,5,9,22,24,30,31,35 1190  41.4462200 0.03613445  
         18) V1f=24,30,31 503  14.5526800 0.02982107  
           36) V1f=24,31 385  10.6857100 0.02857143  
             72) V1f=24 180   4.8611110 0.02777778 *
             73) V1f=31 205   5.8243900 0.02926829 *
           37) V1f=30 118   3.8644070 0.03389831 *
         19) V1f=4,5,9,22,35 687  26.8588100 0.04075691  
           38) V1f=4,35 266   9.6240600 0.03759398  
             76) V1f=35 214   7.7009350 0.03738318 *
             77) V1f=4 52   1.9230770 0.03846154 *
           39) V1f=5,9,22 421  17.2304000 0.04275534  
             78) V1f=22 98   3.8367350 0.04081633 *
             79) V1f=5,9 323  13.3931900 0.04334365  
              158) V1f=9 278  11.4820100 0.04316547 *
              159) V1f=5 45   1.9111110 0.04444444 *
      5) V1f=2,7,10,11,13,20,32,33,34,36,37,38,39 2805 162.3301000 0.06167558  
       10) V1f=10,11,32,33,34,39 1779  94.3788600 0.05621135  
         20) V1f=34 182   8.5549450 0.04945055 *
         21) V1f=10,11,32,33,39 1597  85.8146500 0.05698184  
           42) V1f=10 165   8.5090910 0.05454545 *
           43) V1f=11,32,33,39 1432  77.3044700 0.05726257  
             86) V1f=32,33 951  50.9337500 0.05678233  
              172) V1f=32 141   7.5460990 0.05673759 *
              173) V1f=33 810  43.3876500 0.05679012 *
             87) V1f=11,39 481  26.3700600 0.05821206  
              174) V1f=39 328  17.8993900 0.05792683 *
              175) V1f=11 153   8.4705880 0.05882353 *
       11) V1f=2,7,13,20,36,37,38 1026  67.8060400 0.07115010  
         22) V1f=7,38 383  24.2349900 0.06788512  
           44) V1f=38 339  21.4395300 0.06784661 *
           45) V1f=7 44   2.7954550 0.06818182 *
         23) V1f=2,13,20,36,37 643  43.5645400 0.07309487  
           46) V1f=2,13,36 486  32.4794200 0.07201646  
             92) V1f=36 225  14.8622200 0.07111111 *
             93) V1f=2,13 261  17.6168600 0.07279693  
              186) V1f=13 179  12.0558700 0.07262570 *
              187) V1f=2 82   5.5609760 0.07317073 *
           47) V1f=20,37 157  11.0828000 0.07643312  
             94) V1f=37 132   9.2424240 0.07575758 *
             95) V1f=20 25   1.8400000 0.08000000 *
    3) V1f=1,3,6,8,12 942 102.4682000 0.12420380  
      6) V1f=1,3,6 492  44.9187000 0.10162600  
       12) V1f=3,6 368  33.2798900 0.10054350  
         24) V1f=3 249  22.4899600 0.10040160 *
         25) V1f=6 119  10.7899200 0.10084030 *
       13) V1f=1 124  11.6371000 0.10483870 *
      7) V1f=8,12 450  57.0244400 0.14888890  
       14) V1f=12 111  13.6936900 0.14414410 *
       15) V1f=8 339  43.3274300 0.15044250 *

下記は、V1の値をV1grに代入するだけのコード。これで幾つかを0にしたり、幾つかの変数に同じ値を割り振るなどして、「V1grの値の種類を減らす」工夫をする。 

tic.learn$V1gr <- 0
tic.learn$V1gr[tic.learn$V1==1] <- 1
tic.learn$V1gr[tic.learn$V1==2] <- 2
tic.learn$V1gr[tic.learn$V1==3] <- 3
tic.learn$V1gr[tic.learn$V1==4] <- 4
tic.learn$V1gr[tic.learn$V1==5] <- 5
tic.learn$V1gr[tic.learn$V1==6] <- 6
tic.learn$V1gr[tic.learn$V1==7] <- 7
tic.learn$V1gr[tic.learn$V1==8] <- 8
tic.learn$V1gr[tic.learn$V1==9] <- 9
tic.learn$V1gr[tic.learn$V1==10] <- 10
tic.learn$V1gr[tic.learn$V1==11] <- 11
tic.learn$V1gr[tic.learn$V1==12] <- 12
tic.learn$V1gr[tic.learn$V1==13] <- 13
#tic.learn$V1gr[tic.learn$V1==14] <- 14
tic.learn$V1gr[tic.learn$V1==15] <- 15
tic.learn$V1gr[tic.learn$V1==16] <- 16
tic.learn$V1gr[tic.learn$V1==17] <- 17
tic.learn$V1gr[tic.learn$V1==18] <- 18
tic.learn$V1gr[tic.learn$V1==19] <- 19
tic.learn$V1gr[tic.learn$V1==20] <- 10
tic.learn$V1gr[tic.learn$V1==21] <- 21
tic.learn$V1gr[tic.learn$V1==22] <- 22
tic.learn$V1gr[tic.learn$V1==23] <- 23
tic.learn$V1gr[tic.learn$V1==24] <- 24
tic.learn$V1gr[tic.learn$V1==25] <- 25
tic.learn$V1gr[tic.learn$V1==26] <- 26
tic.learn$V1gr[tic.learn$V1==27] <- 27
tic.learn$V1gr[tic.learn$V1==28] <- 28
tic.learn$V1gr[tic.learn$V1==29] <- 29
tic.learn$V1gr[tic.learn$V1==30] <- 30
tic.learn$V1gr[tic.learn$V1==31] <- 31
tic.learn$V1gr[tic.learn$V1==32] <- 32
tic.learn$V1gr[tic.learn$V1==33] <- 33
tic.learn$V1gr[tic.learn$V1==34] <- 34
tic.learn$V1gr[tic.learn$V1==35] <- 35
tic.learn$V1gr[tic.learn$V1==36] <- 36
tic.learn$V1gr[tic.learn$V1==37] <- 37
tic.learn$V1gr[tic.learn$V1==38] <- 38
tic.learn$V1gr[tic.learn$V1==39] <- 39
tic.learn$V1gr[tic.learn$V1==40] <- 40
tic.learn$V1gr[tic.learn$V1==41] <- 41
tic.learn$V1gr <- as.factor(tic.learn$V1gr)

V5についても同様。 

rpart(V86~V5f, data=tic.learn, control=rpart.control(cp=0))

V5の値も幾つかのグループにまとめる。 

tic.learn$V5gr <- 0
tic.learn$V5gr[tic.learn$V5==1] <- 1
tic.learn$V5gr[tic.learn$V5==2] <- 2
tic.learn$V5gr[tic.learn$V5==3] <- 3
tic.learn$V5gr[tic.learn$V5==4] <- 4
tic.learn$V5gr[tic.learn$V5==5] <- 5
tic.learn$V5gr[tic.learn$V5==6] <- 6
tic.learn$V5gr[tic.learn$V5==7] <- 7
tic.learn$V5gr[tic.learn$V5==8] <- 8
tic.learn$V5gr[tic.learn$V5==9] <- 9
tic.learn$V5gr[tic.learn$V5==10] <- 10
tic.learn$V5gr <- as.factor(tic.learn$V5gr)

参考までに、V1とV5で回帰樹を作ってみると・・・。 

rpart(V86~V1f+V5f, data=tic.learn, control=rpart.control(cp=0))

レポート

レポート提出要領:下記「XXXXXXX」は各自の学籍番号(半角文字)で置き換えること

項目指定
提出期限実験実施の翌週の月曜日の午後7時0分まで
提出方法電子メールに添付 (宛先は配付資料に記載)
ファイル形式Wordファイル (LaTeXで作成する場合は、dvipdfmxでPDFに変換すること)
メールの件名統計工学実験2レポート提出(XXXXXXX)
レポートファイルの名称統計工学実験2_XXXXXXX.doc あるいは 統計工学実験2_XXXXXXX.docx
提出部数レポートは各自1通ずつ。レポートの表紙に、共同実験者の学籍番号と氏名を記すこと。

参考文献

サポート欄

  • data.2のyがひとつ足りなかったのを、追加しました。(1122a)
  • tic.learnというデータ名をtic.leaanとミスタイプしていたのを修正しました。(1122a)
  • インターネットに繋がらないパソコンを使っている人は、TAさんから次の2つのファイルを貰ってください。(1122a)
    • このWikiページのPDFファイル
    • ticdata2000.txt
  • 回帰分析の結果から標準化残差とテコ比の散布図を描くとき、配布資料では残差を標準化するのに、「残差の平方和を残差の自由度で割ったもの」を誤差分散の推定値としますが、Rでは「残差の標本分散」を誤差分散の推定値としています。第5章の例題では、それぞれ「残差平方和/7」と「残差平方和/9」ですので、Rが描く標準化残差のグラフはすべて、配布資料よりも9/7だけ原点から拡大されることになります。(0211p)
  • V1は使わないのがおすすめ。番号の順序に意味がなく、各コードごとの頻度を集計させると、次のようになるため。(0237p) 
    > table(tic.learn$V1)
  1   2   3   4   5   6   7   8   9  10  11  12  13  15  16  17  18  19  20 
124  82 249  52  45 119  44 339 278 165 153 111 179   5  16   9  19   3  25 
 21  22  23  24  25  26  27  28  29  30  31  32  33  34  35  36  37  38  39 
 15  98 251 180  82  48  50  25  86 118 205 141 810 182 214 225 132 339 328 
 40  41 
 71 205
  • 保険商品ごとに難易度が異なります。V86が一番簡単。
V0123456
V75542638214
V7655291731001181
V77579131
V78578438
V795799194
V80266630171267321
V8158193
V825789312
V835675111342
V845777441

参考

少し加工する

以下の6行は、実行しない方がいい場合もある。 

tic.learn$V1 <- as.factor(tic.learn$V1)
tic.learn$V5 <- as.factor(tic.learn$V5)
tic.learn$V86 <- as.factor(tic.learn$V86)
tic.eval$V1 <- as.factor(tic.eval$V1)
tic.eval$V5 <- as.factor(tic.eval$V5)
tic.eval$V86 <- as.factor(tic.eval$V86)

あとはそのまま。

考えたルールに基づく対象限定

各変数に閾値を設けてルールを生成したとする。 たとえば、「V47が5.5以上かつV44が1未満」または「V47が5.5以上かつV1が{1,3,6,8,12,20}のどれか」、というルールは 次のように記す。 

(tic.learn$V47>5.5 & tic.learn$V44<1) | (tic.learn$V47>5.5 & (tic.learn$V1==1 |tic.learn$V1==3 | tic.learn$V1==6 | tic.learn$V1==8 | tic.learn$V1==12 | tic.learn$V1==20) )

「&」が「かつ(AND)」、「|」が「または(OR)」である。

このルールを検証用データに適用するには、 

tic.learn.visit <- (tic.learn$V47>5.5 & tic.learn$V44<1) | (tic.learn$V47>5.5 & (tic.learn$V1==1 |tic.learn$V1==3 | tic.learn$V1==6 | tic.learn$V1==8 | tic.learn$V1==12 | tic.learn$V1==20) )

と、訪問するか否かを二値(TRUE, FALSE)で表すオブジェクトを生成する。 このモデルに予測に基づいた訪問の成果を検証するには、訪問対象のリストtic.visitと検証用データの正解V86のクロス集計を行えばよい。 

table(tic.learn.visit)

FALSE  TRUE 
 3029   971 

table(tic.learn.visit, tic.learn$V86)

tic.learn.visit    0    1
         FALSE 2878  151
         TRUE   884   87

ここでは、訪問対象に884+87=971人を選定し、そのうちの87人が実際に契約してくれる人だったことになる。 契約率は87/971=8.96%。また誤判別率は 

(884+151)/4000

で25.9%となる。

モデルに基づく対象限定

学習したモデルに基づいて、訪問対象を狭めるには、predict()という関数を用いて、訪問対象か否かというリストを作成する。 まず、設定まで調整したモデルを、学習用データ(tic.learn)から得る。 

tic.rpart <- rpart(V86~., data=tic.learn, control=c(cp=0.005))

次に、このモデル(ここではtic.rpart)を検証用データ(tic.eval)に適用して、契約してくれるか否かの予測を行う。 この際、0.05という閾値も調整の必要がある。 

tic.eval.visit <- predict(tic.rpart, newdata=tic.eval)[,2]>0.05

このモデルに予測に基づいた訪問の成果を検証するには、訪問対象のリストtic.visitと検証用データの正解V86のクロス集計を行えばよい。 V86の保険商品が分析対象の場合は、 

table(tic.eval.visit)

tic.eval.visit 
FALSE  TRUE
 2389  1611

table(tic.eval.visit, tic.eval$V86)

tic.eval.visit    0    1
         FALSE 2310   79
         TRUE  1452  159

ここでは、訪問対象に1452+159=1611人を選定し、そのうちの159人が実際に契約してくれる人だったことになる。契約率は159/1452=11.0%。 また誤判別率は 

(79+1452)/4000

で38.275%となる。

想定される困難

次の1行を実行すると、かなり時間がかかってエラーになる。 

tic.glm.step <- step(glm(V86~., family="binomial", data=tic.learn)

次の４行、いずれもエラーになる。変数間の関係が悪すぎるよう。変数の意味を考えて、追加しないといけないかも。 

tic.glm <- glm(V86~V1+V2+V3+V4+V5+V6+V7+V8+   V10+
V11+V12+    V14+    V16+V17+    V19+V20+
V21+V22+V23+    V25+V26+V27+V28+    V30+
        V33+V34+V35+    V37+V38+V39+V40+
    V42+V43+V44+V45+V46+V47+V48+V49+V50+
V51+V52+V53+V54+V55+V56+V57+V58+V59+V60+
V61+V62+V63+V64+V65+V66+V67+V68+V69+V70+
V71+V72+V73+V74+V75+V76+V77+V78+V79+V80+
V81+V82+V83+V84+V85, family="binomial", data=tic.learn)
table(predict(tic.glm, newdata=tic.eval)>0.5)
tic.glm <- glm(V86~     V2+V3+V4+    V6+V7+V8+   V10+
V11+V12+    V14+    V16+V17+    V19+V20+
V21+V22+V23+    V25+V26+V27+V28+    V30+
        V33+V34+V35+    V37+V38+V39+V40+
    V42+V43+V44+V45+V46+V47+V48+V49+V50+
V51+V52+V53+V54+V55+V56+V57+V58+V59+V60+
V61+V62+V63+V64+V65+V66+V67+V68+V69+V70+
V71+V72+V73+V74+V75+V76+V77+V78+V79+V80+
V81+V82+V83+V84+V85, family="binomial", data=tic.learn)
tic.glm <- glm(V86~V44+V45+V46+V47+V48+V49+V50+
V51+V52+V53+V54+V55+V56+V57+V58+V59+V60+
V61+V62+V63+V64+V65+V66+V67+V68+V69+V70+
V71+V72+V73+V74+V75+V76+V77+V78+V79+V80+
V81+V82+V83+V84+V85, family="binomial", data=tic.learn)
tic.glm <- glm(V86~V44+V45+V46+V47+V48+V49+V50+
V51+V52+V53+V54+V55+V56+V57+V58+V59+V60+
V61+V62+V63+V64, family="binomial", data=tic.learn)