混合线性模型如何检测固定因子和随机因子的显著性以及计算R2_混合效应模型怎么计算r2-程序员宅基地
技术标签： 混合线性模型
很多朋友写信问我， 像要知道固定因子的显著性和随机因子的显著性如何计算，他们使用的是lme4这个R包, 但是这个包使用anova时没有P值,还要手动计算, 随机因子也需要自己计算loglikehood值, 然后使用LRT的卡方检验进行显著性检验, 其实lme4包有扩展的包可以非常友好的做这件事情.
1. 载入数据和软件包
###载入软件包和数据
library(lme4)
library(lmerTest)
library(sjstats)
library(learnasreml)
data(fm)
2. 软件包介绍
lme4
- R语言中最流行的混合线性包
- 结果不太友好, 所以才有下面两个包作为辅助
- 安装方法
install.packages("lme4")
lmerTest
- 主要是用于检测lme4对象的固定因子和随机因子,它有两个函数:
lmerTest::anova.lmerModLmerTest用于检测固定因子的显著性, 方差分析表采用III平方和的形式.lmerTest::ranova用于检测随机因子的显著性, 使用的是LRT检验, 给出的是卡方结果.- 安装方法
install.packages("lmerTest")
sjstats
- 可以计算R2
- 可以提取方差组分
- 安装方法
install.packages("lmerTest")
3. 使用lme4进行混合线性分析
模型介绍
- 固定因子: Spacing + Rep
- 随机因子: Fam
建模
### 固定因子: Spacing+Rep, 随机因子: Fam
fm1 <- lmer(h1 ~Spacing + Rep + (1|Fam), fm)
固定因子检验
anova(fm1) # 固定因子显著性检验
| Sum Sq | Mean Sq | NumDF | DenDF | F value | Pr(>F) | |
|---|---|---|---|---|---|---|
| Spacing | 1775.888 | 1775.888 | 1 | 786.2359 | 34.98368 | 4.957481e-09 |
| Rep | 5339.632 | 1334.908 | 4 | 781.5386 | 26.29670 | 1.762292e-20 |
可以看到Spacing 和Rep都达到极显著
随机因子显著性检验
ranova(fm1) # 随机因子显著性检验,LRT
| npar | logLik | AIC | LRT | Df | Pr(>Chisq) | |
|---|---|---|---|---|---|---|
| <none> | 8 | -2831.504 | 5679.008 | NA | NA | NA |
| (1 | Fam) | 7 | -2876.161 | 5766.323 | 89.31466 | 1 | 3.367543e-21 |
可以看到Fam达到极显著
计算R2
r2(fm1) # 计算R2
R-Squared for Generalized Linear Mixed Model
[34mFamily : gaussian (identity)
Formula: h1 ~ Spacing + Rep + (1 | Fam)
[39m Marginal R2: 0.116
Conditional R2: 0.277
计算固定因子每个水平的P值
p_value(fm1) # 计算每个水平的显著性
| term | p.value | std.error |
|---|---|---|
| (Intercept) | 1.535094e-127 | 0.7915991 |
| Spacing3 | 4.957481e-09 | 0.5463546 |
| Rep2 | 2.886600e-01 | 0.8082299 |
| Rep3 | 7.443430e-08 | 0.8218056 |
| Rep4 | 1.720753e-10 | 0.7995633 |
| Rep5 | 4.635631e-01 | 0.7663026 |
提取方差组分
re_var(fm1) # 计算方差组分
-
_sigma_2
- 50.7633345854535 Fam_tau.00
- 11.3168413429073
4. 使用asreml进行对照
建模
library(asreml)
fm2 = asreml(h1 ~ Spacing + Rep, random = ~ Fam, data=fm,trace=F)
固定因子检验
anova(fm2) # 固定因子显著性检验, 这里anova 是anova.asreml
| Df | Sum of Sq | Wald statistic | Pr(Chisq) | |
|---|---|---|---|---|
| (Intercept) | 1 | 253059.46646 | 4985.08358 | 0.000000e+00 |
| Spacing | 1 | 1164.71720 | 22.94406 | 1.667844e-06 |
| Rep | 4 | 5339.63197 | 105.18678 | 0.000000e+00 |
| residual (MS) | NA | 50.76333 | NA | NA |
随机因子显著性检验
这里首先构建一个空模型, 然后使用LRT检验
fm_Null = asreml(h1 ~ Spacing + Rep, data=fm,trace=F)
lrt.asreml(fm2,fm_Null) # 随机因子显著性检验LRT
| Df | LR statistic | Pr(Chisq) | |
|---|---|---|---|
| fm2/fm_Null | 1 | 89.31466 | 0 |
summary(fm2)$varcomp[,1:2] # 方差组分
| gamma | component | |
|---|---|---|
| Fam!Fam.var | 0.2229334 | 11.31684 |
| R!variance | 1.0000000 | 50.76333 |
5. 关于混合线性模型计算R2
还有一个包叫MuMIn,也可以计算R2
library(MuMIn)
r.squaredLR(fm1)#计算R2
0.217233511687581
6. 完整代码分享
# 混合线性模型, 如何检测固定因子和随机因子
###载入数据
library(lme4)
library(lmerTest)
library(sjstats)
library(learnasreml)
data(fm)
str(fm)
### 固定因子: Spacing+Rep, 随机因子: Fam
fm1 <- lmer(h1 ~Spacing + Rep + (1|Fam), fm)
summary(fm1)
anova(fm1) # 固定因子显著性检验
ranova(fm1) # 随机因子显著性检验,LRT
r2(fm1) # 计算R2
p_value(fm1) # 计算每个水平的显著性
re_var(fm1) # 计算方差组分
### 对比asreml
fm2 = asreml(h1 ~ Spacing + Rep, random = ~ Fam, data=fm)
anova(fm2) # 固定因子显著性检验, 这里anova 是anova.asreml
fm_Null = asreml(h1 ~ Spacing + Rep, data=fm)
lrt.asreml(fm2,fm_Null) # 随机因子显著性检验LRT
summary(fm2)$varcomp[,1:2] # 方差组分
library(MuMIn)
r.squaredLR(fm1)#计算R2
智能推荐
JWT(Json Web Token)实现无状态登录_无状态token登录-程序员宅基地
文章浏览阅读685次。1.1.什么是有状态？有状态服务，即服务端需要记录每次会话的客户端信息，从而识别客户端身份，根据用户身份进行请求的处理，典型的设计如tomcat中的session。例如登录：用户登录后，我们把登录者的信息保存在服务端session中，并且给用户一个cookie值，记录对应的session。然后下次请求，用户携带cookie值来，我们就能识别到对应session，从而找到用户的信息。缺点是什么？服务端保存大量数据，增加服务端压力 服务端保存用户状态，无法进行水平扩展 客户端请求依赖服务.._无状态token登录
SDUT OJ逆置正整数-程序员宅基地
文章浏览阅读293次。SDUT OnlineJudge#include<iostream>using namespace std;int main(){int a,b,c,d;cin>>a;b=a%10;c=a/10%10;d=a/100%10;int key[3];key[0]=b;key[1]=c;key[2]=d;for(int i = 0;i<3;i++){ if(key[i]!=0) { cout<<key[i.
年终奖盲区_年终奖盲区表-程序员宅基地
文章浏览阅读2.2k次。年终奖采用的平均每月的收入来评定缴税级数的，速算扣除数也按照月份计算出来，但是最终减去的也是一个月的速算扣除数。为什么这么做呢，这样的收的税更多啊，年终也是一个月的收入，凭什么减去12*速算扣除数了？这个霸道（不要脸）的说法，我们只能合理避免的这些跨级的区域了，那具体是那些区域呢？可以参考下面的表格：年终奖一列标红的一对便是盲区的上下线，发放年终奖的数额一定一定要避免这个区域，不然公司多花了钱..._年终奖盲区表
matlab 提取struct结构体中某个字段所有变量的值_matlab读取struct类型数据中的值-程序员宅基地
文章浏览阅读7.5k次，点赞5次，收藏19次。matlab结构体struct字段变量值提取_matlab读取struct类型数据中的值
Android fragment的用法_android reader fragment-程序员宅基地
文章浏览阅读4.8k次。1，什么情况下使用fragment通常用来作为一个activity的用户界面的一部分例如, 一个新闻应用可以在屏幕左侧使用一个fragment来展示一个文章的列表,然后在屏幕右侧使用另一个fragment来展示一篇文章 – 2个fragment并排显示在相同的一个activity中,并且每一个fragment拥有它自己的一套生命周期回调方法,并且处理它们自己的用户输_android reader fragment
FFT of waveIn audio signals-程序员宅基地
文章浏览阅读2.8k次。FFT of waveIn audio signalsBy Aqiruse An article on using the Fast Fourier Transform on audio signals. IntroductionThe Fast Fourier Transform (FFT) allows users to view the spectrum content of _fft of wavein audio signals
随便推点
Awesome Mac：收集的非常全面好用的Mac应用程序、软件以及工具_awesomemac-程序员宅基地
文章浏览阅读5.9k次。https://jaywcjlove.github.io/awesome-mac/ 这个仓库主要是收集非常好用的Mac应用程序、软件以及工具，主要面向开发者和设计师。有这个想法是因为我最近发了一篇较为火爆的涨粉儿微信公众号文章《工具武装的前端开发工程师》，于是建了这么一个仓库，持续更新作为补充，搜集更多好用的软件工具。请Star、Pull Request或者使劲搓它 issu_awesomemac
java前端技术---jquery基础详解_简介java中jquery技术-程序员宅基地
文章浏览阅读616次。一.jquery简介 jQuery是一个快速的，简洁的javaScript库，使用户能更方便地处理HTML documents、events、实现动画效果，并且方便地为网站提供AJAX交互 jQuery 的功能概括1、html 的元素选取2、html的元素操作3、html dom遍历和修改4、js特效和动画效果5、css操作6、html事件操作7、ajax_简介java中jquery技术
Ant Design Table换滚动条的样式_ant design ::-webkit-scrollbar-corner-程序员宅基地
文章浏览阅读1.6w次，点赞5次，收藏19次。我修改的是表格的固定列滚动而产生的滚动条引用Table的组件的css文件中加入下面的样式：.ant-table-body{ &amp;::-webkit-scrollbar { height: 5px; } &amp;::-webkit-scrollbar-thumb { border-radius: 5px; -webkit-box..._ant design ::-webkit-scrollbar-corner
javaWeb毕设分享 健身俱乐部会员管理系统【源码+论文】-程序员宅基地
文章浏览阅读269次。基于JSP的健身俱乐部会员管理系统项目分享:见文末!
论文开题报告怎么写？_开题报告研究难点-程序员宅基地
文章浏览阅读1.8k次，点赞2次，收藏15次。同学们，是不是又到了一年一度写开题报告的时候呀？是不是还在为不知道论文的开题报告怎么写而苦恼?Take it easy!我带着倾尽我所有开题报告写作经验总结出来的最强保姆级开题报告解说来啦，一定让你脱胎换骨，顺利拿下开题报告这个高塔，你确定还不赶快点赞收藏学起来吗？_开题报告研究难点
原生JS 与 VUE获取父级、子级、兄弟节点的方法 及一些DOM对象的获取_获取子节点的路径 vue-程序员宅基地
文章浏览阅读6k次，点赞4次，收藏17次。原生先获取对象var a = document.getElementById("dom");vue先添加ref <div class="" ref="divBox">获取对象let a = this.$refs.divBox获取父、子、兄弟节点方法var b = a.childNodes; 获取a的全部子节点 var c = a.parentNode; 获取a的父节点var d = a.nextSbiling; 获取a的下一个兄弟节点 var e = a.previ_获取子节点的路径 vue