交互作用モデルと、最初から変数として用意した重回帰分析の違いを調べてみた

まだ試してませんが、仕事でどうも変数同士が交互作用しているような気がする変数があります。例えば、中古の携帯電話の価格の推定をするとした場合、説明変数に発売年とカメラの画素数を加えるとすると、去年発売されたものの画素数が200万ぐらいなら価値がかなり下がるが、10年ぐらい前だとむしろ価値があがるといった具合です。でも、これを普通の回帰分析でやってしまうと、画素数の係数は一つの値しかでないので、発売年が去年でも10年前でも同じ値だけプラスすることになります。
で、これを解決してくれそうなのが交互作用モデルを用いた重回帰分析というわけです。
いや、実は言うと交互作用モデルについてはよく分かってない部分があるので、本当にこの解釈があっているかどうか自身を持ってはいえないのですが、そんなに間違ってないとは思います。多分。

ちなみに、交互作用モデルというのを簡単にいうと、上の携帯電話の価格の例でいうと、普通の重回帰分析では下記のような式を考えます。
『予測価格=a*発売年+b*カメラの画素数』
普通の重回帰分析ではこのaとbを求めるわけですが、交互作用モデルを用いた重回帰分析では下記のような式を考えます。
『予測価格=a*発売年+b*カメラの画素数+c*発売年*カメラの画素数』
この、aとbとcを求めるわけです。

ただ、ここで思ったわけです。最初から、『発売年*カメラの画素数』というような変数があった場合は、結果が異なるのかどうなのかと。例えば、『発売年*カメラの画素数』を含んだ『年画素』という変数が会った場合、
『予測価格=a’*発売年+b’*カメラの画素数+c’*年画素』
という式になるaとa’、bとb’、cとc’は違うのかどうなのかと。カテゴリー変数得られた係数の数値をそれぞれのカテゴリーに割り当てて重回帰分析した結果ですら同じになったので、多分同じになるんじゃないかと思いながらもやってみました。

利用したファイルは下記。
nico4
先ほど、ニコニコ動画の24時間の総合ランキングから取得したデータです。普段使っているニコ動のそれぞれのカテゴリーの合計データを使おうと思ったのですが、オーバーフローしてしまったので新しくとってきました。
今回試したのは、コメント数とマイリスト数を用いて再生数を推定するというもの。
結果は以下。

> nico <- read.table("http://am-yu.net/wp-content/uploads/2013/12/nico4.txt", header=T, sep="\t")
> head(nico)
  rank                                                                 title    dougaid   view
1    1  島風&金剛にハッピーシンセサイザ歌わせてみた【人力VOCALOID&艦これ】 sm22400980 179468
2    2                              完全版【神回】ずんこ 民主党をフルボッコ sm22393110 304734
3    3 【Minecraft】 方向音痴のマインクラフト Season5 Part7 【ゆっくり実況】 sm22408930 136729
4    4                                                  ひぐらし煌なのです。 sm18719109  50161
5    5                   最強のバッターをつくろう!【パワプロ2012実況】part6 sm22404904 177542
6    6                  【実況】行ってはいけない――“約束の場所”【Part10】 sm22410955  82262
  comment mylist
1    1918  18084
2   19212   9353
3    4560   3623
4     474   8249
5    5806   3501
6   11201   2338
> nico.lm <- lm(view~comment+mylist+comment*mylist, data=nico) #交互作用モデル
> summary(nico.lm)

Call:
lm(formula = view ~ comment + mylist + comment * mylist, data = nico)

Residuals:
    Min      1Q  Median      3Q     Max 
-135968  -30149  -11546   20103  140435 

Coefficients:
                    Estimate    Std. Error t value   Pr(>|t|)    
(Intercept)    30609.0821597  7651.1796297   4.001   0.000124 ***
comment            5.9690678     1.2651654   4.718 0.00000809 ***
mylist             8.0245064     1.7196667   4.666 0.00000993 ***
comment:mylist     0.0006418     0.0002093   3.067   0.002813 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 48030 on 96 degrees of freedom
Multiple R-squared:  0.7452,	Adjusted R-squared:  0.7373 
F-statistic:  93.6 on 3 and 96 DF,  p-value: < 0.00000000000000022

> nico$commyl <- nico$comment * nico$mylist
> nico.lm <- lm(view~comment+mylist+commyl, data=nico) #普通の重回帰分析
> summary(nico.lm)

Call:
lm(formula = view ~ comment + mylist + commyl, data = nico)

Residuals:
    Min      1Q  Median      3Q     Max 
-135968  -30149  -11546   20103  140435 

Coefficients:
                 Estimate    Std. Error t value   Pr(>|t|)    
(Intercept) 30609.0821597  7651.1796297   4.001   0.000124 ***
comment         5.9690678     1.2651654   4.718 0.00000809 ***
mylist          8.0245064     1.7196667   4.666 0.00000993 ***
commyl          0.0006418     0.0002093   3.067   0.002813 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 48030 on 96 degrees of freedom
Multiple R-squared:  0.7452,	Adjusted R-squared:  0.7373 
F-statistic:  93.6 on 3 and 96 DF,  p-value: < 0.00000000000000022

最初の『view~comment+mylist+comment*mylist』は『view~comment*mylist』でもOKです。
というわけで、結論。
全く同じなようです。自由度調整済み決定係数ぐらい異なるかもしれないと思っていましたが、それすら同じでした。

最初から変数同士の掛け算を入れた変数を用意できるなら、そうやるほうがいいような気がします。これだと、Rで、コメント数とコメント数とマイリスト数を掛けた値のみの重回帰分析というのも簡単にできるので(自分がやり方を知らないだけだとは思いますが)。

> nico.lm <- lm(view~comment+comment*mylist, data=nico)
> summary(nico.lm)

Call:
lm(formula = view ~ comment + comment * mylist, data = nico)

Residuals:
    Min      1Q  Median      3Q     Max 
-135968  -30149  -11546   20103  140435 

Coefficients:
                    Estimate    Std. Error t value   Pr(>|t|)    
(Intercept)    30609.0821597  7651.1796297   4.001   0.000124 ***
comment            5.9690678     1.2651654   4.718 0.00000809 ***
mylist             8.0245064     1.7196667   4.666 0.00000993 ***
comment:mylist     0.0006418     0.0002093   3.067   0.002813 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 48030 on 96 degrees of freedom
Multiple R-squared:  0.7452,	Adjusted R-squared:  0.7373 
F-statistic:  93.6 on 3 and 96 DF,  p-value: < 0.00000000000000022
> nico.lm <- lm(view~comment+commyl, data=nico)
> summary(nico.lm)

Call:
lm(formula = view ~ comment + commyl, data = nico)

Residuals:
    Min      1Q  Median      3Q     Max 
-205341  -32452  -12019   21851  174492 

Coefficients:
                 Estimate    Std. Error t value          Pr(>|t|)    
(Intercept) 52851.0342302  6594.9156759   8.014 0.000000000002507 ***
comment         3.4083338     1.2561082   2.713           0.00788 ** 
commyl          0.0013430     0.0001606   8.365 0.000000000000448 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 52920 on 97 degrees of freedom
Multiple R-squared:  0.6874,	Adjusted R-squared:  0.681 
F-statistic: 106.7 on 2 and 97 DF,  p-value: < 0.00000000000000022

コメント

タイトルとURLをコピーしました