R數據分析:雙連續變量交互作用的簡單斜率圖作圖及解釋

Interaction: When the effect of one independent variable differs based on the level or magnitude of another independent variable

對於交互作用相信很多人都不陌生，論文中也會經常出現，大家自己也會想要看看感興趣變量之間的交互，交互作用的意思就是在一個自變量的不同水平，另一個自變量的效應大小會有差別，這個和調節作用統計上基本上是一樣的，兩者只有理論意義上的差異。

今天就帶大家用一個實際例子做一個交互作用

y = A + B + A*B

上面的這個公式就是最簡單的，有交互作用的回歸，對於此方面的詳盡解釋，大家可以參考文獻：Jaccard & Turrisi 2003 Interaction Effects in Multiple Regression

今天重點給大家寫兩個連續變量的簡單斜率圖和解釋：

兩連續變量的交互

先模擬出我們今天的數據集：

library(car)

#此部分為數據的模擬
n <- 250
#模擬兩個正態分布的連續變量
X <- rnorm(n, 2.75, .75)
Z <- rnorm(n, 15, 15)
#模擬因變量
Y <- .7*X + .3*Z + 2.5*X*Z + rnorm(n, sd = 5)
#因變量轉化
Y = (Y - min(Y)) / (max(Y) - min(Y))*4
#生成我們的數據
GPA.Data <- data.frame(GPA=Y, Work.Ethic=X, IQ=Z)

上面的代碼生成一個模擬的數據框，我們假設應變量是學生的GPA，兩個自變量分別是學習態度和智商，我們文章要探討的就是學習態度和智商在影響GPA時的交互作用。

此時智商和態度都是正態分布的連續變量。

首先擬合模型：

GPA.Data$IQ.C <- scale(GPA.Data$IQ, center = TRUE, scale = FALSE)
GPA.Data$Work.Ethic.C <- scale(GPA.Data$Work.Ethic, center = TRUE, scale = FALSE)


GPA.Model.1 <- lm(GPA~IQ.C+Work.Ethic.C, GPA.Data)
GPA.Model.2 <- lm (GPA~IQ.C*Work.Ethic.C, GPA.Data)

library(stargazer)
stargazer(GPA.Model.1, GPA.Model.2,type="html",
column.labels = c("Main Effects", "Interaction"),
intercept.bottom = FALSE,
single.row=FALSE,
notes.append = FALSE,
header=FALSE,
out="test.html",
out.header=TRUE)

上面的代碼進行了有交互和沒有交互時的模型擬合，同時我還用了stargazer來輸出模型結果：

可以看到交互作用是有的，下面準備畫簡單斜率圖：

簡單斜率圖

兩個連續變量放在回歸方程中很好理解，在做交互的簡單斜率圖的時候我們就得指定水平，比如AB兩個連續變量有交互作用，我們如果以A為x軸做簡單斜率圖，我們需要表達的是A的效應隨著B的不同水平的不同而不同。

然而B是一個連續變量，所以我們此時得給B指定作圖的水平。

指定水平時一般有三種方法：hand picking, quantiles, standard deviation。

我們先來看hand picking：

library(effects)
Inter.HandPick <- effect('IQ.C*Work.Ethic.C', GPA.Model.2,
xlevels=list(IQ.C = c(-15, 0, 15),
Work.Ethic.C = c(-1.1, 0, 1.1)),
se=TRUE, confidence.level=.95, typical=mean)

Inter.HandPick <- as.data.frame(Inter.HandPick)
head(Inter.HandPick)

可以看到在我們指定的不同水平都有擬合係數均值：

有了上面的數據我們就可以做簡單斜率圖了：

Inter.HandPick$IQ <- factor(Inter.HandPick$IQ.C,
levels=c(-15, 0, 15),
labels=c("1 SD Below Population Mean", "Population Mean", "1 SD Above Population Mean"))


Inter.HandPick$Work.Ethic <- factor(Inter.HandPick$Work.Ethic.C,
levels=c(-1.1, 0, 1.1),
labels=c("Poor Worker", "Average Worker", "Hard Worker"))

library(ggplot2)
Plot.HandPick<-ggplot(data=Inter.HandPick, aes(x=Work.Ethic, y=fit, group=IQ))+
geom_line(size=2, aes(color=IQ))+
ylim(0,4)+
ylab("GPA")+
xlab("Work Ethic")+
ggtitle("Hand Picked Plot")


Plot.HandPick

運行以上代碼即得到簡單斜率圖：

因為我們選擇的B的水平是levels=c(-15, 0, 15)，而我們模擬的B也就是IQ的分布是一個以15為均值15為標準差的正態分布，而後進行了中心化，所以我們畫簡單斜率圖選擇的這個水平levels=c(-15, 0, 15)就是均值和加減一個標準差的水平。這個是我們自己選的，所以叫做hand picking。

對於這個簡單斜率圖的解釋如下：

對於IQ均值在總體均值一個標準差以上的這些人，他們的學習態度越好那麼GPA也越好，在普通IQ的學生中也有這麼一種關係，但是比較弱一點，但是對於那些IQ低於人群一個標準差的同學，他們的學習態度再好，GPA好像也不增加。

我們接著看另外一種劃分水平的方法---quantiles

首先我們將我們的B的水平化出來：

IQ.Quantile <- quantile(GPA.Data$IQ.C, probs=c(0,.25,.50,.75,1))
IQ.Quantile <- round(IQ.Quantile, 2)

其餘的步驟和基本就一樣了，依然還是先跑我們的B的不同水平的係數均值：

library(effects)

Inter.Quantile <- effect('IQ.C*Work.Ethic.C', GPA.Model.2,
xlevels=list(IQ.C = c(-35.44, -9.78, -0.04, 9.89, 41.90),
Work.Ethic.C = c(-1.1, 0, 1.1)),
se=TRUE, confidence.level=.95, typical=mean)

Inter.Quantile <- as.data.frame(Inter.Quantile)


Inter.Quantile$IQ<-factor(Inter.Quantile$IQ.C,
levels=c(-35.44, -9.78, -0.04, 9.89, 41.90),
labels=c("0%", "25%", "50%", "75%", "100%"))

Inter.Quantile$Work.Ethic<-factor(Inter.Quantile$Work.Ethic.C,
levels=c(-1.1, 0, 1.1),
labels=c("Poor Worker", "Average Worker", "Hard Worker"))

然後再畫圖：

library(ggplot2)
Plot.Quantile<-ggplot(data=Inter.Quantile, aes(x=Work.Ethic, y=fit, group=IQ))+
geom_line(size=2, aes(color=IQ))+
ylab("GPA")+
xlab("Work Ethic")+
scale_color_manual(values=c("#42c5f4","#54f284","#f45dcc",
"#ff9d35","#d7afff"))+
theme_bw()+
theme(text = element_text(family="Impact", size=14, color="black"))+ #可以在這換字體
ggtitle("Quantile Plot")

Plot.Quantile

對於上面圖的解釋相信大家都會了，這兒不多羅嗦，直接繼續看第三種劃分水平的方法

第三種劃分水平的方法叫做Standard Deviation

其實我們第一種方法就是按照標準差劃分的，所以這個做出來的圖和第一種方法基本沒有區別，首先還是水平劃分：

IQ.SD <- c(mean(GPA.Data$IQ.C)-sd(GPA.Data$IQ.C),
mean(GPA.Data$IQ.C),
mean(GPA.Data$IQ.C)+sd(GPA.Data$IQ.C))

IQ.SD <- round(IQ.SD, 2)

然後做交互，得出係數均值：

Inter.SD <- effect(c("IQ.C*Work.Ethic.C"), GPA.Model.2,
xlevels=list(IQ.C=c(-14.75, 0, 14.75),
Work.Ethic.C=c(-1.1, 0, 1.1)))

Inter.SD <- as.data.frame(Inter.SD)


Inter.SD$IQ<-factor(Inter.SD$IQ.C,
levels=c(-14.75, 0, 14.75),
labels=c("1 SD Below Mean", "Mean", "1 SD Above Mean"))

Inter.SD$Work.Ethic<-factor(Inter.SD$Work.Ethic.C,
levels=c(-1.1, 0, 1.1),
labels=c("Poor Worker", "Average Worker", "Hard Worker"))

然後再出圖：

Plot.SD<-ggplot(data=Inter.SD, aes(x=Work.Ethic, y=fit, group=IQ))+
geom_line(size=1, aes(color=IQ))+
geom_point(aes(colour = IQ), size=2)+
geom_ribbon(aes(ymin=fit-se, ymax=fit+se),fill="gray",alpha=.6)+
ylim(0,4)+
ylab("GPA")+
xlab("Work Ethic")+
ggtitle("Standard Deviation Plot")+
theme_bw()+ #Removes the gray background
theme(panel.grid.major=element_blank(),
panel.grid.minor=element_blank(),
legend.key = element_blank())+ #Removes the lines
scale_fill_grey()
Plot.SD

此圖的解釋和第一種方法一模一樣哈。

小結

今天給大家寫了回歸中雙連續變量交互作用的簡單斜率圖的畫法和解釋，之後會給大家寫有分類變量的交互，感謝大家耐心看完。發表這些東西的主要目的就是督促自己，希望大家關注評論指出不足，一起進步。內容我都會寫的很細，用到的數據集也會在原文中給出連結，你只要按照文章中的代碼自己也可以做出一樣的結果，一個目的就是零基礎也能懂，因為自己就是什麼編程基礎沒有從零學Python和R的，加油。數據分析問題諮詢，代處理請私信。

（站外連結發不了，請關注後私信回復「數據連結」獲取本頭條號所有使用數據）

往期內容：

R數據分析：雙因素方差分析與交互作用檢驗

python機器學習：機器學習模型評價-交叉驗證與留一驗證