交互作用效应(p for Interaction)在SCI文章中可以算是一个必杀技,几乎在高分的SCI中必出现,因为把人群分为亚组后再进行统计可以增强文章结果的可靠性,进行可视化后可以清晰的表明变量之间的关系。不仅如此,交互作用还可以使用来进行数据挖掘。在既往文章中,我们已经介绍了怎么使用R语言对交互作用进行可视化分析
今天咱们来介绍一下sjPlot包,可以使用它可以轻松绘制各种回归模型中交互项的边际效应。它接受许多模型对象,例如 lm、glm、lme, lmerMod等。下面咱们来演示一下,先导入数据和R包,数据使用sjPlot包自带的efc数据
library(sjPlot)
library(sjmisc)
library(sjlabelled)
library(ggplot2)
data(efc)
这是一个有关于老年人护理的数据,我介绍一下等会我要用到的变量,neg_c_7:7个项目的负面影响,c12hour:每周平均护理时数,BARTHTOT:总分 BARTHEL INDEX,c161sex :照顾者的性别,
如果想了解更详细的数据情况,可以
get_label(efc)
分析前先设定下绘图背景
theme_set(theme_sjplot())
分类变量转成因子
efc$c161sex <- to_factor(efc$c161sex)
假设咱们想了解不同性别种BARTHEL INDEX分数和负面影响的关系
设立模型
fit <- lm(neg_c_7 ~ c12hour + barthtot * c161sex, data = efc)
绘图,如果type = “pred”, terms要制定两个矩阵参数
plot_model(fit, type = "pred", terms = c("barthtot", "c161sex"))
上图表明了BARTHEL INDEX分数不断升高,男性和女性的负面影响评分下降,男性下降更快。
模型默认交互项的第一个作为绘图类型,咱们可以把交互项换一下位置
fit <- lm(neg_c_7 ~ c12hour + c161sex * barthtot, data = efc)
plot_model(fit, type = "int")
这又是另一种风格来表示。如果咱们使用mdrt.values = “meansd”,它会默认按均值分成3组。
plot_model(fit, type = "int", mdrt.values = "meansd")
除了2项交互,terms- 参数还可以接受三个模型项,因此您还可以计算三向交互的边际效应。
假设护理时间,barthtot分数和性别存在交互
fit <- lm(neg_c_7 ~ c12hour * barthtot * c161sex, data = efc)
绘图,这里把barthtot 分数分层3段,相当于分类变量了
plot_model(fit, type = "pred", terms = c("c12hour", "barthtot [30,50,70]", "c161sex"))
得出不同性别,在不同barthtot 分数段中,护理时间和负面影响的关系,我们可以看到barthtot 分数70分这段,护理时间越长,负面影响越高,男女都是一样,儿30分这段随着护理时间延长没有什么变化。
本期结束啦,内容少了点,主要是最近实在太忙了。