読者です 読者をやめる 読者になる 読者になる

東京に棲む日々

データ分析、統計、ITを勉強中。未だ世に出ず。

多項ロジット(Multinomial Logit), R - mlogit 使用メモ

R Statistics

Rのmlogitパッケージで多項ロジット(Multinomial Logit)を使用する際のメモ。


まず、用語の整理。参考文献(A) p.8より。
-------------------------
A model with only individual specific variables is sometimes called a multinomial logit (多項ロジット)model, one with only alternative specific variables a conditional logit(条件付ロジット) model and one with both kind of variables a mixed logit(混合ロジット) model.
This is seriously misleading: conditional logit model is also a logit model for longitudinal data in the statistical literature and mixed logit is one of the names of a logit model with random parameters. Therefore, in what follow, we'll use the name multinomial logit model for the model we've just described whatever the nature of the explanatory variables included in the model.
-------------------------
とのことなので、多項ロジット(Multinomial Logit)と呼ぶことにする。


指定する説明変数がややこしいので、整理しておく。

以下のように効用関数が定義されるとする。

 

f:id:High_School_Student:20150301111544j:plain

i : 個人の添え字
j : 選択肢の添え字

αj は、選択肢jの基本的効用と考えられる。

β : 個人iと選択肢jによって異なる変数Xijに対して、共通の係数
γj : 個人iに固有の変数Ziに対して、選択肢j別に異なる係数
σj : 個人iと選択肢jによって異なる変数Wijに対して、選択肢j別に異なる係数

 

library(mlogit)


data(Fishing)

head(Fishing, 10)

      mode price.beach price.pier price.boat price.charter catch.beach
1  charter     157.930    157.930    157.930       182.930      0.0678
2  charter      15.114     15.114     10.534        34.534      0.1049
3     boat     161.874    161.874     24.334        59.334      0.5333
4     pier      15.134     15.134     55.930        84.930      0.0678
5     boat     106.930    106.930     41.514        71.014      0.0678
6  charter     192.474    192.474     28.934        63.934      0.5333
7    beach      51.934     51.934    191.930       220.930      0.0678
8  charter      15.134     15.134     21.714        56.714      0.0678
9     boat      34.914     34.914     34.914        53.414      0.2537
10    boat      28.314     28.314     28.314        46.814      0.2537
   catch.pier catch.boat catch.charter   income
1      0.0503     0.2601        0.5391 7083.332
2      0.0451     0.1574        0.4671 1250.000
3      0.4522     0.2413        1.0266 3750.000
4      0.0789     0.1643        0.5391 2083.333
5      0.0503     0.1082        0.3240 4583.332
6      0.4522     0.1665        0.3975 4583.332
7      0.0789     0.1643        0.5391 8750.001
8      0.0789     0.0102        0.0209 2083.333
9      0.1498     0.0233        0.0219 3750.000
10     0.1498     0.0233        0.0219 2916.667

一列が1個人となっている”wide型”と呼ばれるデータ。

str(Fishing)

'data.frame':   1182 obs. of  10 variables:
 $ mode         : Factor w/ 4 levels "beach","pier",..: 4 4 3 2 3 4 1 4 3 3 ...
 $ price.beach  : num  157.9 15.1 161.9 15.1 106.9 ...
 $ price.pier   : num  157.9 15.1 161.9 15.1 106.9 ...
 $ price.boat   : num  157.9 10.5 24.3 55.9 41.5 ...
 $ price.charter: num  182.9 34.5 59.3 84.9 71 ...
 $ catch.beach  : num  0.0678 0.1049 0.5333 0.0678 0.0678 ...
 $ catch.pier   : num  0.0503 0.0451 0.4522 0.0789 0.0503 ...
 $ catch.boat   : num  0.26 0.157 0.241 0.164 0.108 ...
 $ catch.charter: num  0.539 0.467 1.027 0.539 0.324 ...
 $ income       : num  7083 1250 3750 2083 4583 ...

 

Mlogitによってモデルをあてはめるには、mlogit.data型に変換する必要がある。

Fish <- mlogit.data(Fishing, varying = c(2:9), shape = "wide", choice ="mode")

Fishingはサンプルデータで、実データで実施すると上手く行かないことが多いが、とりあえず続ける。

ここら辺の変換はhelpを参考にすればよい。

head(Fish, 12)

           mode   income     alt   price  catch chid
1.beach   FALSE 7083.332   beach 157.930 0.0678    1
1.boat    FALSE 7083.332    boat 157.930 0.2601    1
1.charter  TRUE 7083.332 charter 182.930 0.5391    1
1.pier    FALSE 7083.332    pier 157.930 0.0503    1
2.beach   FALSE 1250.000   beach  15.114 0.1049    2
2.boat    FALSE 1250.000    boat  10.534 0.1574    2
2.charter  TRUE 1250.000 charter  34.534 0.4671    2
2.pier    FALSE 1250.000    pier  15.114 0.0451    2
3.beach   FALSE 3750.000   beach 161.874 0.5333    3
3.boat     TRUE 3750.000    boat  24.334 0.2413    3
3.charter FALSE 3750.000 charter  59.334 1.0266    3
3.pier    FALSE 3750.000    pier 161.874 0.4522    3

4行(選択肢の総数)で1個人のデータとなる。

str(Fish)

Classes ‘mlogit.data’ and 'data.frame':       4728 obs. of  6 variables:
 $ mode  : logi  FALSE FALSE TRUE FALSE FALSE FALSE ...
 $ income: num  7083 7083 7083 7083 1250 ...
 $ alt   : chr  "beach" "boat" "charter" "pier" ...
 $ price : num  157.9 157.9 182.9 157.9 15.1 ...
 $ catch : num  0.0678 0.2601 0.5391 0.0503 0.1049 ...
 $ chid  : int  1 1 1 1 2 2 2 2 3 3 ...
 - attr(*, "reshapeLong")=List of 4
  ..$ varying:List of 2
  .. ..$ price: chr  "price.beach" "price.pier" "price.boat" "price.charter"
  .. ..$ catch: chr  "catch.beach" "catch.pier" "catch.boat" "catch.charter"
  .. ..- attr(*, "v.names")= chr  "price" "catch"
  .. ..- attr(*, "times")= chr  "beach" "pier" "boat" "charter"
  ..$ v.names: chr  "price" "catch"
  ..$ idvar  : chr "chid"
  ..$ timevar: chr "alt"
 - attr(*, "index")='data.frame':       4728 obs. of  2 variables:
  ..$ chid: Factor w/ 1182 levels "1","2","3","4",..: 1 1 1 1 2 2 2 2 3 3 ...
  ..$ alt : Factor w/ 4 levels "beach","boat",..: 1 2 3 4 1 2 3 4 1 2 ...
 - attr(*, "choice")= chr "mode"

目的変数となるmodeはlogical型になる。

 

 

もう一つ、サンプルデータ。

library(AER)
data(TravelMode)

str(TravelMode)

'data.frame':   840 obs. of  9 variables:
 $ individual: Factor w/ 210 levels "1","2","3","4",..: 1 1 1 1 2 2 2 2 3 3 ...
 $ mode      : Factor w/ 4 levels "air","train",..: 1 2 3 4 1 2 3 4 1 2 ...
 $ choice    : Factor w/ 2 levels "no","yes": 1 1 1 2 1 1 1 2 1 1 ...
 $ wait      : int  69 34 35 0 64 44 53 0 69 34 ...
 $ vcost     : int  59 31 25 10 58 31 25 11 115 98 ...
 $ travel    : int  100 372 417 180 68 354 399 255 125 892 ...
 $ gcost     : int  70 71 70 30 68 84 85 50 129 195 ...
 $ income    : int  35 35 35 35 30 30 30 30 40 40 ...
 $ size      : int  1 1 1 1 2 2 2 2 1 1 ...

head(TravelMode, 10)

   individual  mode choice wait vcost travel gcost income size
1           1   air     no   69    59    100    70     35    1
2           1 train     no   34    31    372    71     35    1
3           1   bus     no   35    25    417    70     35    1
4           1   car    yes    0    10    180    30     35    1
5           2   air     no   64    58     68    68     30    2
6           2 train     no   44    31    354    84     30    2
7           2   bus     no   53    25    399    85     30    2
8           2   car    yes    0    11    255    50     30    2
9           3   air     no   69   115    125   129     40    1
10          3 train     no   34    98    892   195     40    1

“long型”と呼ばれるデータ。

 

mlogit.data型に変換。

TM <- mlogit.data(TravelMode, choice = "choice", shape = "long", id.var = "individual", alt.levels = c("air", "train", "bus", "car"))

str(TM)

Classes ‘mlogit.data’ and 'data.frame':       840 obs. of  9 variables:
 $ individual: Factor w/ 210 levels "1","2","3","4",..: 1 1 1 1 2 2 2 2 3 3 ...
 $ mode      : Factor w/ 4 levels "air","train",..: 1 2 3 4 1 2 3 4 1 2 ...
 $ choice    : logi  FALSE FALSE FALSE TRUE FALSE FALSE ...
 $ wait      : int  69 34 35 0 64 44 53 0 69 34 ...
 $ vcost     : int  59 31 25 10 58 31 25 11 115 98 ...
 $ travel    : int  100 372 417 180 68 354 399 255 125 892 ...
 $ gcost     : int  70 71 70 30 68 84 85 50 129 195 ...
 $ income    : int  35 35 35 35 30 30 30 30 40 40 ...
 $ size      : int  1 1 1 1 2 2 2 2 1 1 ...
 - attr(*, "index")='data.frame':       840 obs. of  3 variables:
  ..$ chid: Factor w/ 210 levels "1","2","3","4",..: 1 1 1 1 2 2 2 2 3 3 ...
  ..$ alt : Factor w/ 4 levels "air","train",..: 1 2 3 4 1 2 3 4 1 2 ...
  ..$ id  : Factor w/ 210 levels "1","2","3","4",..: 1 1 1 1 2 2 2 2 3 3 ...
 - attr(*, "choice")= chr "choice"

head(TM, 12)

        individual  mode choice wait vcost travel gcost income size
1.air            1   air  FALSE   69    59    100    70     35    1
1.train          1 train  FALSE   34    31    372    71     35    1
1.bus            1   bus  FALSE   35    25    417    70     35    1
1.car            1   car   TRUE    0    10    180    30     35    1
2.air            2   air  FALSE   64    58     68    68     30    2
2.train          2 train  FALSE   44    31    354    84     30    2
2.bus            2   bus  FALSE   53    25    399    85     30    2
2.car            2   car   TRUE    0    11    255    50     30    2
3.air            3   air  FALSE   69   115    125   129     40    1
3.train          3 train  FALSE   34    98    892   195     40    1
3.bus            3   bus  FALSE   35    53    882   149     40    1
3.car            3   car   TRUE    0    23    720   101     40    1

4行(選択肢の総数)で1個人。

Helpにデータ項目の簡単な説明がある。
210人のデータである。

 

選択に関するクロス集計。
( actChoice <- table(TM$mode, TM$choice) )

        FALSE TRUE
  air     152   58
  train   147   63
  bus     180   30
  car     151   59

割合で出す場合。
( actChoiceRatio <- actChoice[,2]/210 )      # 210 = #individuals

      air     train       bus       car 
0.2761905 0.3000000 0.1428571 0.2809524 

 

 

では、多項ロジットモデルのあてはめを行う。

FitTM <- mlogit(choice ~ vcost | income | travel, data=TM)

formulaは、“|”が仕切りとなって分かれている。
これは、上で説明したように、Xij | Zi | Wij とそれぞれの変数の役割に従って指定する。

vcostは各個人iが提示された各選択肢jの値段(各乗り物の値段)である。個人iによって同じ選択肢jでも値段が異なっていることに注意。もしi間で共通の変数であれば、選択肢と全く同じとなってしまうので、モデルに含めることはできない。
vcost変数の、選択肢jによって異ならない共通の効用への影響度を推定する(β)。

incomeは収入で、選択肢jに依存しない、個人i依存の変数である。
incom変数の、効用への影響度を選択肢j別に推定する(γj)。

travelは各個人iが選択肢jを選択した場合に経験するトラベルタイムである。
travel変数の、効用への影響度を選択肢j別に推定する(σj)。トラベルタイムに対して感じる効用が、選択肢(乗り物)によって異なると仮定する場合はここに変数を入れるのだと思われる。


結果の確認。
summary( FitTM )

Call:
mlogit(formula = choice ~ vcost | income | travel, data = TM, 
    method = "nr", print.level = 0)

Frequencies of alternatives:
    air   train     bus     car 
0.27619 0.30000 0.14286 0.28095 

nr method
5 iterations, 0h:0m:0s 
g'(-H)^-1g = 7.87E-05 
successive function values within tolerance limits 

Coefficients :
                    Estimate Std. Error t-value  Pr(>|t|)    
train:(intercept)  1.6866632  0.8477430  1.9896  0.046636 *  
bus:(intercept)    0.1449168  0.9201414  0.1575  0.874856    
car:(intercept)   -0.7699463  0.8161644 -0.9434  0.345491    
vcost             -0.0056710  0.0066967 -0.8468  0.397089    
train:income      -0.0598253  0.0129413 -4.6228 3.786e-06 ***
bus:income        -0.0414433  0.0136923 -3.0268  0.002472 ** 
car:income        -0.0129009  0.0111799 -1.1539  0.248527    
air:travel        -0.0367774  0.0067680 -5.4340 5.509e-08 ***
train:travel      -0.0074593  0.0011859 -6.2898 3.180e-10 ***
bus:travel        -0.0068219  0.0013397 -5.0922 3.539e-07 ***
car:travel        -0.0066324  0.0011331 -5.8531 4.825e-09 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Log-Likelihood: -232.48
McFadden R^2:  0.18071 
Likelihood ratio test : chisq = 102.56 (p.value = < 2.22e-16)

ここで出てくる“Frequencies of alternatives:”は上で計算したとおり、実測選択の割合である。予測選択の割合ではない。


切片(intercept)とγ、σは一つの水準が0に固定される。
0にする水準を入れ替える場合は、引数reflevelを指定する。
summary( mlogit(choice ~ vcost | income | travel, reflevel=4, data=TM) )

Call:
mlogit(formula = choice ~ vcost | income | travel, data = TM, 
    reflevel = 4, method = "nr", print.level = 0)

Frequencies of alternatives:
    car     air   train     bus 
0.28095 0.27619 0.30000 0.14286 

nr method
5 iterations, 0h:0m:0s 
g'(-H)^-1g = 7.87E-05 
successive function values within tolerance limits 

Coefficients :
                    Estimate Std. Error t-value  Pr(>|t|)    
air:(intercept)    0.7699463  0.8161644  0.9434 0.3454908    
train:(intercept)  2.4566096  0.6244004  3.9343 8.342e-05 ***
bus:(intercept)    0.9148631  0.7662845  1.1939 0.2325191    
vcost             -0.0056710  0.0066967 -0.8468 0.3970894    
air:income         0.0129009  0.0111799  1.1539 0.2485269    
train:income      -0.0469243  0.0121817 -3.8520 0.0001171 ***
bus:income        -0.0285424  0.0129991 -2.1957 0.0281118 *  
car:travel        -0.0066324  0.0011331 -5.8531 4.825e-09 ***
air:travel        -0.0367774  0.0067680 -5.4340 5.509e-08 ***
train:travel      -0.0074593  0.0011859 -6.2898 3.180e-10 ***
bus:travel        -0.0068219  0.0013397 -5.0922 3.539e-07 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Log-Likelihood: -232.48
McFadden R^2:  0.18071 
Likelihood ratio test : chisq = 102.56 (p.value = < 2.22e-16)

 

適合度の指標に関して。

“Log-Likelihood:”は対数尤度。
“McFadden R^2:”は参考文献(B) p.124では、McFadden ρ^2 と書かれている。
“Likelihood ratio test :”は切片モデルに対する検定。

“McFadden R^2:”を計算してみる。

( ll <- summary( FitTM )$logLik )

'log Lik.' -232.4804 (df=11)

 

切片モデル。
FitTM0 <- mlogit(choice ~ 1, data=TM)

summary( FitTM0 )

Call:
mlogit(formula = choice ~ 1, data = TM, method = "nr", print.level = 0)

Frequencies of alternatives:
    air   train     bus     car 
0.27619 0.30000 0.14286 0.28095 

nr method
4 iterations, 0h:0m:0s 
g'(-H)^-1g = 0.000622 
successive function values within tolerance limits 

Coefficients :
                   Estimate Std. Error t-value Pr(>|t|)   
train:(intercept)  0.082692   0.181974  0.4544 0.649529   
bus:(intercept)   -0.659237   0.224888 -2.9314 0.003374 **
car:(intercept)    0.017094   0.184907  0.0924 0.926341   
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Log-Likelihood: -283.76
McFadden R^2:  0 
Likelihood ratio test : chisq = 0 (p.value = 1)

 

( ll0 <- summary( FitTM0 )$logLik )

 'log Lik.' -283.7588 (df=3)

 

as.numeric( 1 - (ll / ll0) )

[1] 0.1807111

となる。

対数尤度を利用した計算からも分かるとおり、モデルに変数を加えれば加えるほど値が高くなっていくことに注意。


この指標の読み方であるが、参考文献(B)が(C)を参照している。

 

f:id:High_School_Student:20150301114657j:plain

 (参考文献(C) p.124より)

 横軸が“McFadden R^2:”で、縦軸が通常のR^2である。

どうやって導き出したかは不明だが、経験則的にこのように判断できるとなっている。0.5程度のモデルが作れれば、十分合格点ということになるだろう。

 

各個人に対する、各選択肢の選択予測値は以下の操作で出力を行う。
( fitRes <- FitTM$probabilities )
# fitted(FitTM, outcome=FALSE)     # これでも同じ

               air      train         bus         car
  [1,] 0.119817364 0.23053702 0.090576734 0.559068885
  [2,] 0.317031379 0.28838285 0.102175448 0.292410320
  [3,] 0.651515719 0.04527835 0.049321451 0.253884481
  [4,] 0.492368062 0.03999761 0.029395135 0.438239197
  [5,] 0.117289099 0.49291692 0.237541845 0.152252133
  [6,] 0.087258878 0.53618734 0.131474081 0.245079699
  [7,] 0.829788104 0.02232870 0.027580156 0.120303036
  [8,] 0.251498768 0.34596333 0.176161979 0.226375922
  [9,] 0.154158841 0.21992794 0.094726675 0.531186546
 [10,] 0.526192120 0.04275839 0.031663300 0.399386190
 [11,] 0.179378500 0.42617592 0.112891226 0.281554353
 [12,] 0.114347018 0.29971102 0.117604796 0.468337167
 [13,] 0.144166110 0.17389664 0.089848946 0.592088309
 [14,] 0.223915354 0.27859546 0.119081403 0.378407786
 [15,] 0.498162161 0.25977304 0.200121678 0.041943118
..... skip .....
[207,] 0.218547728 0.24899471 0.211234675 0.321222889
[208,] 0.056502916 0.56022179 0.276306236 0.106969061
[209,] 0.170163004 0.37815127 0.183175726 0.268510000
[210,] 0.405440673 0.04138912 0.057691327 0.495478876


各個人の実際の選択に対する予測値のみを返したい場合。
FitTM$fitted.values
# fitted( FitTM, outcome=TRUE )      # これでも同じ

      1.air       2.air       3.air       4.air       5.air       6.air 
0.559068885 0.292410320 0.253884481 0.438239197 0.152252133 0.536187342 
      7.air       8.air       9.air      10.air      11.air      12.air 
0.829788104 0.226375922 0.531186546 0.399386190 0.281554353 0.468337167 
     13.air      14.air      15.air      16.air      17.air      18.air 
0.592088309 0.378407786 0.041943118 0.625623550 0.558350920 0.738326313 
..... skip .....
    199.air     200.air     201.air     202.air     203.air     204.air 
0.510837490 0.361912483 0.258362415 0.315766795 0.202403945 0.267393122 
    205.air     206.air     207.air     208.air     209.air     210.air 
0.166773948 0.128082868 0.218547728 0.276306236 0.268510000 0.495478876 


混同行列(Confusion Matrix)を作成してみる。
標準出力にあるのかもしれないが、分からないので、選択予測確率(fitResに格納している)から作成してみる。

このfitResの各行において、最大の確率に対応する列名を取り出す関数を作成。

## fitResの各行において、maxとなる列名をまとめたvectorを返す関数
fn_findPred1 <- function(df0) {
     df1 <- as.data.frame(df0)
     maxLevel <- apply(df1, 1, max)
     nc <- ncol(df1)

     temp_return <- NULL
     for(i in 1:nc) {
          #print(i)
          temp_col <- ifelse( df1[i]==maxLevel, names(df1[i]), "" )
          #print(head(temp_col))

          if(i==1){
          temp_return <- temp_col
          }else{
          temp_return <- paste(temp_return, temp_col, sep="")
          }
          #print(head(temp_return))
     }

     return( temp_return )
}

fitResを代入して関数を実行。
( predChoiceCol <- fn_findPred1(fitRes) )

  [1] "car"   "air"   "air"   "air"   "train" "train" "air"   "train" "car"  
 [10] "air"   "train" "car"   "car"   "car"   "air"   "train" "train" "train"
 [19] "train" "air"   "train" "train" "air"   "air"   "train" "air"   "air"  
 [28] "train" "train" "train" "train" "train" "air"   "train" "train" "train"
..... skip .....
[181] "train" "train" "train" "train" "train" "air"   "train" "train" "car"  
[190] "air"   "train" "train" "car"   "bus"   "train" "train" "train" "air"  
[199] "air"   "air"   "train" "air"   "train" "train" "air"   "air"   "car"  
[208] "train" "train" "car"  

各個人の予測選択肢が出力される。

この予測選択による選択肢の割合を計算するには以下。
table(predChoiceCol)/210

predChoiceCol
       air        bus        car      train 
0.30000000 0.04761905 0.25238095 0.40000000 

 

各個人の実際の選択をデータ(TM)から作成。

( actChoiceCol <- as.character(TM$mode[TM$choice==TRUE]) )

  [1] "car"   "car"   "car"   "car"   "car"   "train" "air"   "car"   "car"  
 [10] "car"   "car"   "car"   "car"   "car"   "car"   "train" "train" "train"
 [19] "train" "train" "train" "train" "air"   "air"   "air"   "air"   "air"  
 [28] "air"   "train" "train" "train" "train" "train" "train" "train" "train"
..... skip .....
[181] "car"   "car"   "train" "train" "train" "bus"   "air"   "air"   "car"  
[190] "car"   "car"   "air"   "air"   "bus"   "train" "train" "car"   "air"  
[199] "air"   "air"   "bus"   "car"   "bus"   "bus"   "bus"   "car"   "air"  
[208] "bus"   "car"   "car"  

 

predChoiceColとactChoiceColでクロス集計。

( ConfMatrix <- table(actChoiceCol, predChoiceCol) )

            predChoiceCol
actChoiceCol air bus car train
       air    25   0  20    13
       bus     9  10   3     8
       car    17   0  30    12
       train  12   0   0    51

と混同行列(Confusion Matrix)が作成できる。


誤分類率を求める場合は、混同行列の対角要素を利用し、以下のように計算できる。
1 - sum(diag(ConfMatrix))/210      # 210 = #individuals

[1] 0.447619

 

 

では、変数を増やしたモデルも作成しておく。


FitTM2 <- mlogit(choice ~ wait+vcost | income+size | travel, data=TM)

summary( FitTM2 )

Call:
mlogit(formula = choice ~ wait + vcost | income + size | travel, 
    data = TM, method = "nr", print.level = 0)

Frequencies of alternatives:
    air   train     bus     car 
0.27619 0.30000 0.14286 0.28095 

nr method
5 iterations, 0h:0m:0s 
g'(-H)^-1g = 7.32E-07 
gradient close to zero 

Coefficients :
                    Estimate Std. Error t-value  Pr(>|t|)    
train:(intercept) -2.1240906  1.2159544 -1.7469 0.0806633 .  
bus:(intercept)   -3.3311281  1.3317591 -2.5013 0.0123739 *  
car:(intercept)   -7.5393851  1.2585583 -5.9905 2.092e-09 ***
wait              -0.0973882  0.0113005 -8.6180 < 2.2e-16 ***
vcost             -0.0058938  0.0086924 -0.6780 0.4977476    
train:income      -0.0727496  0.0172146 -4.2260 2.378e-05 ***
bus:income        -0.0333524  0.0176288 -1.8919 0.0585002 .  
car:income        -0.0163247  0.0141088 -1.1571 0.2472481    
train:size         1.1771968  0.3094744  3.8039 0.0001425 ***
bus:size           0.8055721  0.3962361  2.0331 0.0420464 *  
car:size           1.0019958  0.2754798  3.6373 0.0002755 ***
air:travel        -0.0329873  0.0072772 -4.5330 5.816e-06 ***
train:travel      -0.0069413  0.0013858 -5.0090 5.472e-07 ***
bus:travel        -0.0067683  0.0016560 -4.0872 4.367e-05 ***
car:travel        -0.0070640  0.0012872 -5.4879 4.068e-08 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Log-Likelihood: -162.79
McFadden R^2:  0.42631 
Likelihood ratio test : chisq = 241.94 (p.value = < 2.22e-16)

McFrdden R^2など大きく上昇している。

 

各個人における、選択肢の予測確率。
fitRes2 <- fitted(FitTM2, outcome=FALSE)

上で定義した関数を適用。
predChoiceCol2 <- fn_findPred1(fitRes2) 

 

混同行列。

( ConfMatrix2 <- table(actChoiceCol, predChoiceCol2) )

            predChoiceCol2
actChoiceCol air bus car train
       air    40   1  14     3
       bus     5  22   2     1
       car     5   0  46     8
       train   5   1   9    48

 

誤分類率。
1 - sum(diag(ConfMatrix2))/210

 [1] 0.2571429

 

 

実データで扱う場合は、それはそれで大変なのだが、今回は使い方メモということで…。

 

---------- 参考文献 ----------

(A) Estimation of multinomial logit models in R : The mlogit Packages
Yves Croissant
http://cran.r-project.org/web/packages/mlogit/vignettes/mlogit.pdf


(B) マーケティングのデータ分析―分析手法と適用事例
岡太彬訓, 守口剛
http://www.asakura.co.jp/books/isbn/978-4-254-12822-2/


(C) Urban Travel Demand: A Behavioral Analysis
Tom Domencich, Daniel L. McFadden
http://eml.berkeley.edu/~mcfadden/travel.html
Ch.5
http://eml.berkeley.edu/~mcfadden/travel/ch5.pdf