关于R语言中混合线性模型summary()结果中交互作用beta值的含义
本文以2*2的实验设计为例,利用lmerTest包在R中进行混合线性模型分析,采用sum的因子编码方式,简单介绍一下在summary的结果中,交互作用的beta值的含义。
数据准备:
library(tidyverse);
library(lmertest)
DF = read_csv('https://raw.githubusercontent.com/usplos/Eye-movement-related/master/DataforShiny.csv')
DF[c('A','B')] = lapply(DF[c('A','B')], factor)
数据中,A和B是自变量(first level),Y是因变量,Sub和Item是群组变量(second level)
采用sum编码方式,有两种操作方式,下面分别介绍:
在option()中直接设置
options(contrasts = c('contr.sum','contr.poly')) # 这里需要同时依次设置无序因子和有序因子的编码方式
查看因子A和因子B的编码方式:
> contrasts(DF$A) [,1]
A1 1
A2 -1
############################
> contrasts(DF$B) [,1]
B1 1
B2 -1
即A因素中以A2为基线,B因素中以B2为基线。
建立模型:
Model = lmer(data = DF, Y~A*B+(1|Sub)+(1|Item)) # 这里采用最简单的模型为例
获取summary中固定效应的结果
> summary(Model) %>% coef(.) %>% round(x = .,digits = 3)
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 248.476 10.116 9.671 24.563 0.000
A1 -3.218 2.551 1092.286 -1.262 0.207
B1 -3.906 2.550 1093.941 -1.531 0.126
A1:B1 -2.382 3.028 58.248 -0.787 0.435
即A因素beta值为3.218,B因素为3.906,AB交互作用为2.382。
注意,这里A的beta值不是A1和A2均值差距,而是 (A1-A2)/2
进行A因素的事后检验,验证上面的结论:
> emmeans::emmeans(Model, pairwise~A) %>% .$contr
NOTE: Results may be misleading due to involvement in interactions
contrast estimate SE df t.ratio p.value
A1 - A2 -6.44 5.12 1091 -1.257 0.2091 Results are averaged over the levels of: B
可以看到:6.44 = 2 * 3.22.
同理,B的beta值为(B1-B2)/2.(均值差的一半)。
那么交互作用呢?beta值为2.382,代表什么意思?这里进行简单效应分析:
> emmeans::emmeans(Model, pairwise~A|B) %>% .$contr
B = B1:
contrast estimate SE df t.ratio p.value
A1 - A2 -11.20 8.00 146 -1.400 0.1637 B = B2:
contrast estimate SE df t.ratio p.value
A1 - A2 -1.67 7.86 136 -0.213 0.8318
对交互作用了解比较深刻的人似乎可以察觉到 2.382 和 11.2、1.67这两个数值之间有一些关系。再细心一点,就可以看到下面的关系:
> (11.2-1.67)/4
[1] 2.3825
诶?11.2–1.67 代表的是条件A在条件B的不同水平上效应的差异,也就是交互作用在数值上的体现。而summary结果中的2.382是这一数值的1/4,为什么是1/4是想必大家也懂了——既然A和B的beta值是其效应量的一半,那么交互作用肯定是其效应量的1/4了。
那么有什么办法可以直接从summary中获取真实的效应量,而不是减半呢?有!
在建模命令中自定义设置为sum编码
重新建立模型:
Model = lmer(data = Df,
Y~A*B+(1|Sub)+(1|Item),
contrasts = list(A = contr.sum(2)/2,
B = contr.sum(2)/2))
再来看一下A和B的编码方式:
> contr.sum(2)/2
[,1]
1 0.5
2 -0.5
即从原来的1和-1,变成了0.5和-0.5.
获取summary中固定效应的结果:
> summary(Model) %>% coef(.) %>% round(x = .,digits = 3)
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 248.476 10.116 9.671 24.563 0.000
A1 -6.437 5.102 1092.286 -1.262 0.207
B1 -7.811 5.101 1093.941 -1.531 0.126
A1:B1 -9.530 12.110 58.248 -0.787 0.435
再来做一下A和B的事后检验:
> emmeans::emmeans(Model, pairwise~A)$contr contrast estimate SE df t.ratio p.value A1 - A2 -6.44 5.12 1091 -1.257 0.2091
>
> emmeans::emmeans(Model, pairwise~B)$contr contrast estimate SE df t.ratio p.value B1 - B2 -7.81 5.12 1092 -1.526 0.1274
可以看到,此时的beta值等于对应的效应量。
再做一下简单效应:
> emmeans::emmeans(Model, pairwise~A|B)$contr
B = B1:
contrast estimate SE df t.ratio p.value
A1 - A2 -11.20 8.00 146 -1.400 0.1637 B = B2:
contrast estimate SE df t.ratio p.value
A1 - A2 -1.67 7.86 136 -0.213 0.8318
验证一下:
此时交互作用的beta值也等于交互作用的效应量。
总结一下:
交互作用其实是看因素A在因素B不同水平上的效应量的差异;
sum编码下,summary中固定效应的效应量不一定是真实的效应,应是编码方式而定;
有两种sum编码的方式:第一种是在环境中声明——option()命令中设置,第二种是建模中自定义设定(通过设置contrasts参数,contr.sum(n)/n);
第一种设置下,主效应的 beta值 为真实效应量的1/2,交互作用的 beta 值为真实效应量的1/4;
第二种设置下,主效应和交互作用的beta值等于真实效应量;
编码类型的不同影响anova和summary中p值的一致性;
sum编码方式下,矩阵的数值影响summary中固定效应beta值和真实效应量的一致性;
Shiny Performer中默认为第一种设置方式(因为第一种typing codes起来方便,我懒…)。
小红书,“黄”了
苹果谷歌双双被曝,你的手机正在窃听你的生活
遇事不决赖毛子,美国这次打算封杀变脸APP
关于R语言中混合线性模型summary()结果中交互作用beta值的含义相关推荐
- R语言对数线性模型loglm函数_使用R语言进行混合线性模型(mixed linear model) 分析代码及详解...
1.混合线性模型简介 混合线性模型,又名多层线性模型(Hierarchical linear model).它比较适合处理嵌套设计(nested)的实验和调查研究数据.此外,它还特别适合处理带有被试内 ...
- R语言ggplot2可视化移除数据中的NA值再可视化实战:消除图形中非常突出的NA柱状图、使用subset函数、使用drop_na函数
R语言ggplot2可视化移除数据中的NA值再可视化实战:消除图形中非常突出的NA柱状图.使用subset函数.使用drop_na函数 目录
- R语言使用edit函数在Rsudio中生成数据编辑器(在windows中生成编辑器)、在编辑器中输出需要的数据生成最终的dataframe
R语言使用edit函数在Rsudio中生成数据编辑器(在windows中生成编辑器).在编辑器中输出需要的数据生成最终的dataframe 目录
- R语言ggplot2可视化:计算dataframe中每个数据列缺失值的个数、使用堆叠的条形图(Stacked Barplot)可视化每个数据列的缺失值的情况(自定义堆叠条形图的形式)
R语言ggplot2可视化:计算dataframe中每个数据列缺失值的个数.使用堆叠的条形图(Stacked Barplot)可视化每个数据列的缺失值的情况(自定义堆叠条形图的形式) 目录
- R语言ggplot2可视化:为图像中的均值竖线、中位数竖线、 geom_vline添加图例(legend)
R语言ggplot2可视化:为图像中的均值竖线.中位数竖线. geom_vline添加图例(legend) 目录
- R语言ggplot2可视化在可视化图像中添加上限线条、下限线条、添加上下限图例实战
R语言ggplot2可视化在可视化图像中添加上限线条.下限线条.添加上下限图例实战 目录
- R语言ggplot2可视化抑制可视化网格中的竖线输出、抑制可视化网格中的横线线输出、抑制背景网格输出实战
R语言ggplot2可视化抑制可视化网格中的竖线输出.抑制可视化网格中的横线线输出.抑制背景网格输出实战 目录
- R语言ggplot2可视化在箱图中为箱图添加均值的标签及对应数值实战
R语言ggplot2可视化在箱图中为箱图添加均值的标签及对应数值实战 目录 R语言ggplot2可视化在箱图中为箱图添加均值的标签及对应数值实战
- R语言ggplot2可视化移除图例中的a字符实战
R语言ggplot2可视化移除图例中的a字符实战 目录 R语言ggplot2可视化移除图例中的a字符实战
最新文章
- python使用fpdf创建pdf文件包含:页眉、页脚并嵌入logo图片、设置使用中文字体
- SQLBulkCopy 性能统计
- stm32之iap实现应用(基于串口,上位机,详细源码)
- Go内置库模块 flag
- centos安装ES(elasticsearch)
- java并发之线程池
- 在mysql支持关系模型中_MySQL支持关系模型中、和三种不同的完整性约束
- 思科Smart Software Manager高权限登录凭证遭暴露
- (备忘)Java web项目迁移到Centos7中验证码无法显示
- 全国行政区划代码(json版)
- 终端代码重复率检测实践
- CEC2018:动态多目标测试函数DF6~DF9的PS及PF(提供Matlab代码)
- HP-UX之MP管理
- 微商分销系统哪家好,要怎么做?
- js和python前景比较好_Python,Java和JavaScript这3个编程语言未来哪个更有前景?
- C语言程序设计ncre,NCRE二级C语言程序设计辅导
- 企业上云计划:上云前应该考虑哪些因素
- 使用程序自动调用ANSYS并运行命令流文件
- easyexcel结合zip 导出压缩文件(包含多个excel)
- TOP 1比不加TOP慢的疑惑