語法
給個(gè)鏈接到之前的文檔：如何定義兩物種之間基因表達(dá)量的保守性（內(nèi)有l(wèi)me4包使用說明）
基本的表達(dá)式為：resp ~ FEexpr + (REexpr1 | factor1) + (REexpr2 | factor2) + ...
固定效應(yīng)的截距和斜率不會由于隨機(jī)效應(yīng)項(xiàng)的變化而變化
這里介紹以下幾種常用的方法：
fm1 <- lmer(Reaction ~ Days + ( Days | Subject), dat)
fm1 <- lmer(Reaction ~ Days + (1 | Subject), dat)
fm1 <- lmer(Reaction ~ Days + (Subject | Days), dat)
fm1 <- lmer(Reaction ~ Days + (1 | Days), dat)
fm1 <- lmer(Reaction ~ Subject + (Days | Subject), dat)
fm1 <- lmer(Reaction ~ Subject + (Subject | Days), dat)
fm1 <- lmer(Reaction ~ Subject + (1 | Subject), dat)
fm1 <- lmer(Reaction ~ Subject + (1 | Days), dat)
要點(diǎn)：
- 固定效應(yīng)部分寫 1 代表只考慮固定效應(yīng)截距項(xiàng)而不考慮固定效應(yīng)的斜率；若賦予變量則代表同時(shí)考慮固定效應(yīng)截距項(xiàng)和固定效應(yīng)的斜率
- 隨機(jī)效應(yīng)部分 (隨機(jī)效應(yīng)項(xiàng) | 隨機(jī)分組因子)，隨機(jī)效應(yīng)項(xiàng)寫 1 代表只考慮固定效應(yīng)截距項(xiàng)而不考慮固定效應(yīng)的斜率尸闸；若賦予變量則代表同時(shí)考慮隨機(jī)效應(yīng)截距項(xiàng)和隨機(jī)效應(yīng)的斜率益缎；隨機(jī)分組因子則是對隨機(jī)效應(yīng)項(xiàng)進(jìn)行進(jìn)一步分組与境，比方隨機(jī)分組因子有A块饺，B，C三個(gè)水平分組恩袱，那么線性模型將會分A夫晌，B雕薪，C作擬合
fm1 <- lmer(Reaction ~ Subject + (Subject | Subject), dat)等價(jià)于fm1 <- lmer(Reaction ~ Subject + (1 | Subject), dat)；fm1 <- lmer(Reaction ~ Days+ (Days | Days), dat)等價(jià)于fm1 <- lmer(Reaction ~ Days+ (1 | Days), dat)
例子
其中 Days為連續(xù)變量晓淀，Subject為因子變量
對應(yīng)第一個(gè)表達(dá)式
library(lme4)
data(sleepstudy)
dat = sleepstudy
dat$Subject = as.factor(c(rep('A',60),rep('B',60),rep('C',60)))
fm1 <- lmer(Reaction ~ Days + (Days | Subject), dat)
summary(fm1)
coef(fm1)
> coef(fm1)
$Subject
(Intercept) Days
A 249.2140 8.502512
B 252.3907 11.351088
C 252.6106 11.548257
attr(,"class")
[1] "coef.mer"
表達(dá)式fm1 <- lmer(Reaction ~ Days + ( Days | Subject), dat)的意義是以Days為固定效應(yīng)（考慮固定效應(yīng)的截距和斜率）所袁，還是以Days為隨機(jī)效應(yīng)，即考慮隨機(jī)效應(yīng)的截距和斜率凶掰，分組的隨機(jī)因子為Subject燥爷；因此這種情況相當(dāng)于沒有隨機(jī)效應(yīng)項(xiàng)蜈亩，只有分組，即分組的固定效應(yīng)回歸模型前翎，沒有隨機(jī)效應(yīng)的斜率和截距稚配，全是固定效應(yīng)的斜率和截距
- 當(dāng)Days = 0，Subject = A 時(shí)港华，Reaction 值為 8.502512 × 0 + 249.2140
- 當(dāng)Days = 1道川，Subject = A 時(shí)，Reaction 值為 8.502512 × 1 + 249.2140
- 當(dāng)Days = 2立宜，Subject = A 時(shí)冒萄，Reaction 值為 8.502512 × 2 + 249.2140
- 當(dāng)Days = 0，Subject = B 時(shí)赘理，Reaction 值為 11.351088 × 0 + 252.3907
- 當(dāng)Days = 1宦言，Subject = B 時(shí)，Reaction 值為 11.351088 × 1 + 252.3907
- 當(dāng)Days = 2商模，Subject = B 時(shí)，Reaction 值為 11.351088 × 2 + 252.3907
- 當(dāng)Days = 0蜘澜，Subject = C 時(shí)施流，Reaction 值為 11.548257 × 0 + 252.6106
- 當(dāng)Days = 1，Subject = C 時(shí)鄙信，Reaction 值為 11.548257 × 1 + 252.6106
- 當(dāng)Days = 2瞪醋，Subject = C 時(shí)，Reaction 值為 11.548257 × 2 + 252.6106
因此它的線性關(guān)系需要分組來看：
對應(yīng)第二個(gè)表達(dá)式
library(lme4)
data(sleepstudy)
dat = sleepstudy
dat$Subject = as.factor(c(rep('A',60),rep('B',60),rep('C',60)))
fm1 <- lmer(Reaction ~ Days + (1 | Subject), dat)
summary(fm1)
coef(fm1)
> coef(fm1)
$Subject
(Intercept) Days
A 241.1497 10.46729
B 255.8238 10.46729
C 257.2418 10.46729
attr(,"class")
[1] "coef.mer"
表達(dá)式fm1 <- lmer(Reaction ~ Days + ( 1 | Subject), dat)的意義是以Days為固定效應(yīng)（考慮固定效應(yīng)的截距和斜率）装诡，隨機(jī)效應(yīng)項(xiàng)寫 1 代表不考慮隨機(jī)效應(yīng)的斜率項(xiàng)银受，考慮其截距項(xiàng)（即將隨機(jī)因子的影響轉(zhuǎn)換為截距加入到方程中），分組的隨機(jī)因子為Subject鸦采，例如 A 中的截距 241.1497宾巍，B中的截距為255.8238 受不同Subject分組而不同，并且它們的截距為固定效應(yīng)截距與隨機(jī)效應(yīng)截距之和渔伯，而斜率為固定效應(yīng)斜率
因此它的線性關(guān)系需要分組來看：
對應(yīng)第三個(gè)表達(dá)式
library(lme4)
data(sleepstudy)
dat = sleepstudy
dat$Subject = as.factor(c(rep('A',60),rep('B',60),rep('C',60)))
fm1 <- lmer(Reaction ~ Days + (Subject | Days), dat)
summary(fm1)
coef(fm1)
> coef(fm1)
$Days
SubjectB SubjectC (Intercept) Days
0 2.399840e-08 2.920663e-08 251.4051 10.46729
1 4.067253e-08 4.970416e-08 251.4051 10.46729
2 -8.451885e-09 -1.019267e-08 251.4051 10.46729
3 9.618661e-09 1.181750e-08 251.4051 10.46729
4 -9.571652e-09 -1.162255e-08 251.4051 10.46729
5 4.944233e-08 6.047095e-08 251.4051 10.46729
6 -6.091995e-09 -7.603941e-09 251.4051 10.46729
7 2.267801e-08 2.766016e-08 251.4051 10.46729
8 5.934868e-08 7.253005e-08 251.4051 10.46729
9 7.305345e-08 8.927133e-08 251.4051 10.46729
attr(,"class")
[1] "coef.mer"
表達(dá)式fm1 <- lmer(Reaction ~ Days + ( Subject | Days), dat)的意義是以Days為固定效應(yīng)顶霞，隨機(jī)效應(yīng)項(xiàng)為因子變量Subject（即將隨機(jī)效應(yīng)項(xiàng)的斜率作為影響距加入到方程中，忽略隨機(jī)效應(yīng)截距）锣吼，分組的隨機(jī)因子為Days选浑，按照Days進(jìn)行分組計(jì)算；由于此時(shí)模型將Subject當(dāng)作因子變量玄叠，所以默認(rèn)SubjectA為參考物古徒，SubjectB和SubjectC分別和SubjectA對比（此時(shí)因子變量忽略截距，即忽略隨機(jī)效應(yīng)截距）读恃；由于Subject為因子變量（此時(shí)因子變量忽略截距）隧膘，截距=251.4051為固定效應(yīng)的截距崎苗，斜率=10.46729為固定效應(yīng)的斜率，此時(shí)隨機(jī)效應(yīng)的影響用斜率體現(xiàn)
因此它的線性關(guān)系需要分組來看：
- 當(dāng)Days = 0舀寓，Subject = A 時(shí)胆数，Reaction 值為 10.46729 × 0 + 251.4051
- 當(dāng)Days = 0，Subject = B 時(shí)互墓，Reaction 值為 10.46729 × 0 + 251.4051 + 2.399840e-08
- 當(dāng)Days = 0必尼，Subject = C 時(shí)，Reaction 值為 10.46729 × 0 + 251.4051 + 2.920663e-08
- 當(dāng)Days = 1篡撵，Subject = A 時(shí)判莉，Reaction 值為 10.46729 × 1 + 251.4051
- 當(dāng)Days = 1，Subject = B 時(shí)育谬，Reaction 值為 10.46729 × 1 + 251.4051 + 2.399840e-08
- 當(dāng)Days = 1券盅，Subject = C 時(shí)，Reaction 值為 10.46729 × 1 + 251.4051 + 2.920663e-08
- 當(dāng)Days = 2膛檀，Subject = A 時(shí)锰镀，Reaction 值為 10.46729 × 2 + 251.4051
- 當(dāng)Days = 2，Subject = B 時(shí)咖刃，Reaction 值為 10.46729 × 2 + 251.4051 + 2.399840e-08
- 當(dāng)Days = 2泳炉，Subject = C 時(shí)，Reaction 值為 10.46729 × 2 + 251.4051 + 2.920663e-08
對應(yīng)第四個(gè)表達(dá)式
library(lme4)
data(sleepstudy)
dat = sleepstudy
dat$Subject = as.factor(c(rep('A',60),rep('B',60),rep('C',60)))
fm1 <- lmer(Reaction ~ Days + (1 | Days), dat)
summary(fm1)
coef(fm1)
> coef(fm1)
$Days
(Intercept) Days
0 251.4051 10.46729
1 251.4051 10.46729
2 251.4051 10.46729
3 251.4051 10.46729
4 251.4051 10.46729
5 251.4051 10.46729
6 251.4051 10.46729
7 251.4051 10.46729
8 251.4051 10.46729
9 251.4051 10.46729
attr(,"class")
[1] "coef.mer"
表達(dá)式fm1 <- lmer(Reaction ~ Days + ( 1 | Days), dat)的意義是以Days為固定效應(yīng)嚎杨，隨機(jī)效應(yīng)項(xiàng)寫 1 代表不考慮隨機(jī)效應(yīng)的斜率花鹅，而考慮隨機(jī)效應(yīng)的截距；此時(shí)隨機(jī)效應(yīng)的截距等同固定效應(yīng)的截距枫浙，因此此時(shí)只有固定效應(yīng)的斜率和截距刨肃，并且隨機(jī)分組因子為Days，相當(dāng)于完全不考慮隨機(jī)效應(yīng)項(xiàng)箩帚，因此結(jié)果完全為固定效應(yīng)的系數(shù)（斜率和截距）
- 當(dāng)Days = 0真友，Reaction 值為 251.4051 + 10.46729 × 0
- 當(dāng)Days = 1，Reaction 值為 251.4051 + 10.46729 × 1
- 當(dāng)Days = 2膏潮，Reaction 值為 251.4051 + 10.46729 × 2
- 當(dāng)Days = 3锻狗，Reaction 值為 251.4051 + 10.46729 × 3
- 當(dāng)Days = 4，Reaction 值為 251.4051 + 10.46729 × 4
- 當(dāng)Days = 5焕参，Reaction 值為 251.4051 + 10.46729 × 5
- 當(dāng)Days = 6轻纪，Reaction 值為 251.4051 + 10.46729 × 6
- 當(dāng)Days = 7，Reaction 值為 251.4051 + 10.46729 × 7
- 當(dāng)Days = 8叠纷，Reaction 值為 251.4051 + 10.46729 × 8
- 當(dāng)Days = 9刻帚，Reaction 值為 251.4051 + 10.46729 × 9
對應(yīng)第五個(gè)表達(dá)式
library(lme4)
data(sleepstudy)
dat = sleepstudy
dat$Subject = as.factor(c(rep('A',60),rep('B',60),rep('C',60)))
fm1 <- lmer(Reaction ~ Subject + (Days | Subject), dat)
summary(fm1)
coef(fm1)
> coef(fm1)
$Subject
Days (Intercept) SubjectB SubjectC
A 7.680847 250.1575 1.854289 4.273369
B 11.291513 251.7820 1.854289 4.273369
C 11.187903 251.7354 1.854289 4.273369
attr(,"class")
[1] "coef.mer"
表達(dá)式fm1 <- lmer(Reaction ~ Subject + (Days | Subject), dat)的意義是以Subject 為固定效應(yīng)（由于Subject為因子變量，所以不考慮固定效應(yīng)的截距）涩嚣，隨機(jī)效應(yīng)項(xiàng)為Days（即考慮隨機(jī)效應(yīng)的斜率和截距）崇众，分組的隨機(jī)因子為Days掂僵，按照Days進(jìn)行分組計(jì)算；由于此時(shí)模型將Subject當(dāng)作因子變量顷歌，所以默認(rèn)SubjectA為參考物锰蓬，SubjectB和SubjectC分別和SubjectA對比（此時(shí)因子變量忽略截距，忽略隨機(jī)效應(yīng)截距）眯漩；因此固定效應(yīng)只有統(tǒng)一的斜率項(xiàng)：斜率B = 1.854289芹扭，斜率C = 4.273369
- 當(dāng)Days = 0，Subject = A 時(shí)赦抖，Reaction 值為 7.680847 × 0 + 250.1575
- 當(dāng)Days = 0舱卡，Subject = B 時(shí)，Reaction 值為 11.291513 × 0 + 251.7820 + 1.854289
- 當(dāng)Days = 0队萤，Subject = C 時(shí)轮锥，Reaction 值為 11.187903 × 0 + 251.7354 + 4.273369
- 當(dāng)Days = 1，Subject = A 時(shí)要尔，Reaction 值為 7.680847 × 1 + 250.1575
- 當(dāng)Days = 1舍杜，Subject = B 時(shí)，Reaction 值為 11.291513 × 1 + 251.7820 + 1.854289
- 當(dāng)Days = 1盈电，Subject = C 時(shí)蝴簇，Reaction 值為 11.187903 × 1 + 251.7354 + 4.273369
- 當(dāng)Days = 2，Subject = A 時(shí)匆帚，Reaction 值為 7.680847 × 2 + 250.1575
- 當(dāng)Days = 2，Subject = B 時(shí)旁钧，Reaction 值為 1.291513 × 2 + 251.7820 + 1.854289
- 當(dāng)Days = 2吸重，Subject = C 時(shí)，Reaction 值為 11.187903 × 2 + 251.7354 + 4.273369
對應(yīng)第六個(gè)表達(dá)式
library(lme4)
data(sleepstudy)
dat = sleepstudy
dat$Subject = as.factor(c(rep('A',60),rep('B',60),rep('C',60)))
fm1 <- lmer(Reaction ~ Subject + (Subject | Days), dat)
summary(fm1)
coef(fm1)
> coef(fm1)
$Days
(Intercept) SubjectB SubjectC
0 257.1725 5.629783 8.06793
1 263.1822 8.705003 11.02711
2 262.9781 8.600570 10.92662
3 273.9903 14.235679 16.34910
4 277.5391 16.051596 18.09649
5 291.3354 23.111346 24.88986
6 292.6354 23.776574 25.52998
7 298.5781 26.817544 28.45621
8 310.3491 32.840869 34.25225
9 319.4526 37.499269 38.73487
attr(,"class")
[1] "coef.mer"
表達(dá)式fm1 <- lmer(Reaction ~ Subject + (Subject | Days), dat)的意義是以Subject 為固定效應(yīng)歪今，還是以Subject 為隨機(jī)效應(yīng)嚎幸，即考慮隨機(jī)效應(yīng)的截距和斜率，分組的隨機(jī)因子為Days寄猩，按照Days分組計(jì)算嫉晶；因此這種情況相當(dāng)于沒有隨機(jī)效應(yīng)項(xiàng)，只有分組田篇，即分組的固定效應(yīng)回歸模型替废，沒有隨機(jī)效應(yīng)的斜率和截距，全是分組后的固定效應(yīng)的斜率和截距
- 當(dāng) Subject = A 時(shí)泊柬，Days = 0椎镣，Reaction = 257.1725
- 當(dāng) Subject = B 時(shí)，Days = 0兽赁，Reaction = 257.1725 + 5.629783
- 當(dāng) Subject = C 時(shí)状答，Days = 0冷守，Reaction = 257.1725 + 8.06793
- 當(dāng) Subject = A 時(shí)，Days = 1惊科，Reaction = 263.1822
- 當(dāng) Subject = B 時(shí)拍摇，Days = 1，Reaction = 263.1822 + 8.705003
- 當(dāng) Subject = C 時(shí)馆截，Days = 1充活，Reaction = 263.1822 + 11.02711
- 當(dāng) Subject = A 時(shí)，Days = 2孙咪，Reaction = 262.9781
- 當(dāng) Subject = B 時(shí)堪唐，Days = 2，Reaction = 262.9781 + 8.600570
- 當(dāng) Subject = C 時(shí)翎蹈，Days = 2淮菠，Reaction = 262.9781 + 10.92662
對應(yīng)第七個(gè)表達(dá)式
當(dāng)然如果將只考慮Subject的固定效應(yīng)項(xiàng)
library(lme4)
data(sleepstudy)
dat = sleepstudy
dat$Subject = as.factor(c(rep('A',60),rep('B',60),rep('C',60)))
fm1 <- lmer(Reaction ~ Subject + (1 | Subject), dat)
summary(fm1)
coef(fm1)
> coef(fm1)
$Subject
(Intercept) SubjectB SubjectC
A 284.7213 19.72682 21.63304
B 284.7213 19.72682 21.63304
C 284.7213 19.72682 21.63304
attr(,"class")
[1] "coef.mer"
表達(dá)式fm1 <- lmer(Reaction ~ Subject + ( 1 | Subject), dat)的意義是以Subject 為固定效應(yīng)，隨機(jī)效應(yīng)項(xiàng)寫 1 代表不考慮隨機(jī)效應(yīng)的斜率荤堪，而考慮隨機(jī)效應(yīng)的截距合陵；此時(shí)隨機(jī)效應(yīng)的截距等同固定效應(yīng)的截距，因此此時(shí)只有固定效應(yīng)的斜率和截距澄阳，并且隨機(jī)分組因子為 Subject拥知，相當(dāng)于完全不考慮隨機(jī)效應(yīng)項(xiàng)，因此結(jié)果完全為固定效應(yīng)的系數(shù)（斜率和截距）碎赢；由于Subject為因子變量低剔，所以默認(rèn) Subject A 為參考物，Subject B和Subject C分別與Subject A作比較
- 當(dāng) Subject = A 時(shí)肮塞，Reaction = 284.7213
- 當(dāng) Subject = B 時(shí)襟齿，Reaction = 284.7213 + 19.72682
- 當(dāng) Subject = C 時(shí)，Reaction = 284.7213 + 21.63304
對應(yīng)第八個(gè)表達(dá)式
library(lme4)
data(sleepstudy)
dat = sleepstudy
dat$Subject = as.factor(c(rep('A',60),rep('B',60),rep('C',60)))
fm1 <- lmer(Reaction ~ Subject + (1 | Days), dat)
summary(fm1)
coef(fm1)
> coef(fm1)
$Days
(Intercept) SubjectB SubjectC
0 248.0516 19.72682 21.63304
1 254.9236 19.72682 21.63304
2 255.6824 19.72682 21.63304
3 271.1280 19.72682 21.63304
4 276.0844 19.72682 21.63304
5 293.4914 19.72682 21.63304
6 296.6977 19.72682 21.63304
7 302.4557 19.72682 21.63304
8 318.1192 19.72682 21.63304
9 330.5787 19.72682 21.63304
attr(,"class")
[1] "coef.mer"
表達(dá)式fm1 <- lmer(Reaction ~ Subject + ( 1 | Days), dat)的意義是以Subject 為固定效應(yīng)枕赵，隨機(jī)效應(yīng)項(xiàng)寫 1 代表不考慮隨機(jī)效應(yīng)的斜率猜欺，而考慮隨機(jī)效應(yīng)的截距（即將隨機(jī)因子的影響轉(zhuǎn)換為截距加入到方程中），并且隨機(jī)分組因子為 Days拷窜，按Days進(jìn)行分組計(jì)算开皿；由于Subject為因子變量，所以默認(rèn) Subject A 為參考物篮昧，Subject B和Subject C分別與Subject A作比較赋荆；而Days作為連續(xù)變量帶入對應(yīng)數(shù)值運(yùn)算；因此結(jié)果的截距為固定效應(yīng)的截距和隨機(jī)效應(yīng)的截距之和恋谭，斜率為固定效應(yīng)的斜率：斜率B=19.72682糠睡，斜率C=21.63304
- 當(dāng) Subject = A 時(shí)，Days =0疚颊，Reaction = 248.0516
- 當(dāng) Subject = A 時(shí)狈孔，Days =1信认，Reaction = 254.9236
- 當(dāng) Subject = A 時(shí)，Days =2均抽，Reaction = 255.6824
- 當(dāng) Subject = B 時(shí)嫁赏，Days =0，Reaction = 248.0516 + 19.72682
- 當(dāng) Subject = B 時(shí)油挥，Days =1潦蝇，Reaction = 254.9236 + 19.72682
- 當(dāng) Subject = B 時(shí)，Days =2深寥，Reaction = 255.6824 + 19.72682
- 當(dāng) Subject = C 時(shí)攘乒，Days =0，Reaction = 248.0516 + 21.63304
- 當(dāng) Subject = C 時(shí)惋鹅，Days =1则酝，Reaction = 254.9236 + 21.63304
- 當(dāng) Subject = C 時(shí)，Days =2闰集，Reaction = 255.6824 + 21.63304
復(fù)合情形
library(lme4)
data(sleepstudy)
dat = sleepstudy
dat$Subject = as.factor(c(rep('A',60),rep('B',60),rep('C',60)))
fm1 <- lmer(Reaction ~ Subject + (1 | Days) + (1 | Subject), dat)
summary(fm1)
coef(fm1)
> coef(fm1)
$Days
(Intercept) SubjectB SubjectC
0 248.0518 19.72682 21.63304
1 254.9238 19.72682 21.63304
2 255.6826 19.72682 21.63304
3 271.1280 19.72682 21.63304
4 276.0844 19.72682 21.63304
5 293.4914 19.72682 21.63304
6 296.6977 19.72682 21.63304
7 302.4556 19.72682 21.63304
8 318.1190 19.72682 21.63304
9 330.5784 19.72682 21.63304
$Subject
(Intercept) SubjectB SubjectC
A 284.7213 19.72682 21.63304
B 284.7213 19.72682 21.63304
C 284.7213 19.72682 21.63304
attr(,"class")
[1] "coef.mer"
這種隨機(jī)效應(yīng)相加的情形相當(dāng)于將這兩項(xiàng)隨機(jī)效應(yīng)分開來看沽讹，詳見第七個(gè)表達(dá)式和第八個(gè)表達(dá)式
