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

Rを通じて統計学を学ぶ備忘録ブログ

SPSSからRに移行したい私のような人向けのR解説ブログ兼学習用備忘録。

第76回 重回帰分析と交互作用〜説明変数が2つ以上の時の作図〜

久しぶりの更新です。このブログのアクセス解析ページを見てみると、重回帰分析と交互作用のページへのアクセスが多いようなので、それに関連するトピックで記事を書いていきます。過去の関連記事は以下から参照ください。

第31回 回帰分析と重回帰分析 - SPSS→R備忘録ブログ

第34回 重回帰分析 - SPSS→R備忘録ブログ

第44回 階層的重回帰分析 - SPSS→R備忘録ブログ

第45回 重回帰分析と交互作用 - SPSS→R備忘録ブログ

第46回の補足の補足 重回帰分析×交互作用項(質的×質的) - SPSS→R備忘録ブログ

補足の回 重回帰分析の交互作用の検討 R実践 - SPSS→R備忘録ブログ

 

重回帰分析と交互作用〜説明変数が2つ以上の時の作図〜

重回帰分析を行った後に、他の変数を平均値で固定し、注目している1つの変数を動かして、y への影響の変化を観察したい時があります。例えば、個人のエンターテインメント傾向を数値化し、各種情報メディアの使用有無との交互作用変数を作成して個人の知識量の変化を観察する、というような感じです。(→この仮説は私の修論の一部です。笑 プリンストン大学の先生の論文に似たような研究があり、それの日本版として分析しました。)

毎度ざっくりな説明で恐縮ですが、違う例を考えると、CMに接触したかどうかで、違うメディア、リスティングなどの広告の接触量を増加させると、ブランドの想起度や認知度が変わるかどうか。みたいな感じでしょうか。仮に高まるのであれば、CMに予算を追加してクリエイティブを重視して、リスティングで接触回数を増やして、ブランド力を高めていく、みたいな戦略も考えれるかもしれません。(もっとちゃんと検討する必要がありますが。あくまでも例なのでご勘弁ください。)

このような仮説を考えた時に、変化量を観察するためには、図示するのが手っ取り早いです。今回は、x1x2classx1*classの4変数を用いた重回帰分析を実施し、x1を動かした時のyの振る舞いを観察しています。

d <- read.csv("inter.data.csv", header = TRUE)

# 重回帰、x1*classは交互作用
model <- lm( y ~ x1 + x2 + class + x1*class, d)
const <- aggregate(x2 ~ class, d, mean)  # classとx2代入用
min(d$x1) #xの範囲を知るため
max(d$x1) #xの範囲を知るため
xpara <- seq(1, 8, 0.1) #x1を動かす用

#x2を平均で固定し、x1を動かした時のclass_aの信頼区間
a_conf <- predict(model, data.frame(x1 = xpara, x2 = const[1, 2], class = const[1,1]), interval="confidence")               
#x2を平均で固定し、x1を動かした時のclass_bの信頼区間
b_conf <- predict(model, data.frame(x1 = xpara, x2 = const[2, 2], class = const[2,1]), interval="confidence")  
  
#塗りつぶし。border=NAは縁をつけない場合
plot(0, 0, xlim=c(2, 7), ylim=c(2, 12), type="n", xlab = "x1", ylab = "y", main = "Multi regression with interaction")  
lines(xpara, a_conf[,1], col="2", lwd=2)       
lines(xpara, b_conf[,1], col="4", lwd=2)       
polygon( c(xpara, rev(xpara)),  c(a_conf[,2], rev(a_conf[,3]) ),col="#FF000020", border=NA )
polygon( c(xpara, rev(xpara)),  c(b_conf[,2], rev(b_conf[,3]) ),col="#0000FF20", border=NA )
 
#点線かつfor文で楽する際のコード。
plot(0, 0, xlim=c(2, 7), ylim=c(2, 12), type="n", ann=F)  
for(i in 1:3) lines(xpara, a_conf[,i], col = 2, lty = 1)
for(i in 1:3) lines(xpara, b_conf[,i], col = 4, lty = 1)
 
#予測値は実線がいい場合はline()関数で先に書く。
lines(xpara, a_conf[,1], col="2", lwd=2)       
lines(xpara, a_conf[,1], col="4", lwd=2)  

 

f:id:teruaki-sugiura:20160302234602p:plain

以上で今回はおしまい。