測試在nlme中的線性混合模型中的相互作用的顯着性

我在nlmeR包中使用lme函數來測試因子items的水平是否與因子condition的水平具有顯着的相互作用。因子condition有兩個等級：Control和Treatment，因子items有3個等級：E1,...,E3。我使用以下代碼：測試在nlme中的線性混合模型中的相互作用的顯着性

f.lme = lme(response ~ 0 + factor(condition) * factor(items), random = ~1|subject)

其中subject是隨機效應。這樣一來，當我運行：

summary(f.lme)$tTable

我會得到下面的輸出：

factor(condition)Control 
factor(condition)Treatment 
factor(items)E2 
factor(items)E3 
factor(condition)Treatment:factor(items)E2 
factor(condition)Treatment:factor(items)E3

與Value, Std.Error, DF, t-value, p-value列在一起。我有兩個問題：

如果我想比較Control與Treatment，我將只使用estimable()功能gmodels做出的(-1,1,0,0,0,0)對比？
我很感興趣，是否items水平，即E1, E2, E3是跨condition不同的，所以我很感興趣的交互項是否顯著（僅通過檢查p-value列??）：

factor(condition)Treatment:factor(items)E2 factor(condition)Treatment:factor(items)E3

但是，如何知道factor(condition)Treatment:factor(items)E1是否顯着？它沒有顯示在彙總輸出中，我認爲它與R中使用的對比度有關...非常感謝！

來源

2013-07-22 alittleboy

我恭敬地@斯文海恩斯坦

就R不同意，對於範疇變量的默認編碼是治療的對比編碼。在治療對比中，第一個水平是參考水平。所有其餘的因子水平與參考水平進行比較。

首先，固定效應在這裏指定爲零截距，... ~ 0 + ...。這意味着condition編碼不再是contr.treatment。如果我沒有記錯的Control和Treatment的主要作用是，現在解釋作爲從組各自偏差的意思是......

在你的模型，該因素的項目有三個層次：E1，E2和E3。這兩種對比測試了（a）E2和E1以及（b）E3和E1之間的差異。這些對比度的主要影響是在因子條件的水平控制下估算的，因爲這是該因子的參考類別。

...當items的值在其參考水平E1！因此：

主要作用Control = Control:E1意見從項目E1的平均多少偏離。
主效應Treatment =有多少Treatment:E1觀測值偏離項目E1的平均值。
主效應E2 =有多少Control:E2觀測值偏離項目E2的平均值。
主效應E3 =有多少Control觀測值偏離項目E3的平均值。
互動Treatment:E2 = Treatment:E2意見從項目E2
互動Treatment:E3 = Treatment:E3意見從項目E3的平均多少偏離的均值多少偏離。

感謝您指向estimable，我以前沒有嘗試過。對於自定義對比，我已經從multcomp包得到(ab)using glht。

來源

2013-07-25 20:59:45 f1r3br4nd

感謝您的建議！看來，產出因素（條件）治療：因素（項目）E2'等都與參考相比較，對嗎？那麼，我怎麼能知道'item'這個因素的水平在'condition'中是不同的？也就是說，我仍然不知道如何獲得測試的p值，說E1在控制方面與在處理方面不同... – alittleboy

你也許可以嘗試省略'...〜0 +。 ..'來自模型。那麼'因素（條件）治療'的主要作用就是它與因素（項目）E2'上的控制之間的差異（如果它與零有顯着差異，那麼治療和控制是不同的）。麻煩的是，那麼相互作用是多少或多或少的治療與對照不同的是E1的差異。你可以嘗試multcomp的glht - http://stackoverflow.com/a/17867984/945039 – f1r3br4nd

好的答案。我忽略了失蹤的攔截。事實上，當沒有攔截時，對比的解釋會改變。 –

測試相互作用是否顯着的常用方法是進行似然比檢驗（例如，參見discussion on R-Sig-ME）。

要做到這一點，你必須也沒有估計的交互模型，你也將不得不使用method="ML"：

f0 = lme(response ~ 0 + factor(condition) * factor(items), 
    random = ~1|subject, method="ML") 
f1 = lme(response ~ 0 + factor(condition) + factor(items), 
    random = ~1|subject, method="ML")

可以使用anova然後比較：

anova(f0, f1)

另見this blog post

來源

2013-07-31 09:50:06

您需要先解決有關交互的第二個問題。您可以像Jan van der Laan的答案那樣設定似然比檢驗。您也可以直接在合適的lme物體上使用anova。有關更多信息，請參閱anova.lme的幫助頁面。

在解釋你的係數方面，我經常發現它有時需要製作一個羣組平均值的彙總表，以便適當地找出模型中哪些線性係數組合代表每個羣組。在你的問題中，我將展示一個攔截刪除的例子，儘管我發現這很少幫助我在模型中有兩個因素後找出我的係數。這裏是我用Orthodont數據集（我決定做平衡）的一個例子：

require(nlme) 

# Make dataset balanced 
Orthodont2 = Orthodont[-c(45:64),] 
# Factor age 
Orthodont2$fage = factor(Orthodont2$age) 

# Create a model with an interaction using lme; remove the intercept 
fit1 = lme(distance ~ Sex*fage - 1, random = ~1|Subject, data = Orthodont2) 
summary(fit1)

這裏是估計的固定效應。但是每個這些係數代表什麼？

Fixed effects: distance ~ Sex * fage - 1 
        Value Std.Error DF t-value p-value 
SexMale   23.636364 0.7108225 20 33.25213 0.0000 
SexFemale  21.181818 0.7108225 20 29.79903 0.0000 
fage10   0.136364 0.5283622 61 0.25809 0.7972 
fage12   2.409091 0.5283622 61 4.55954 0.0000 
fage14   3.727273 0.5283622 61 7.05439 0.0000 
SexFemale:fage10 0.909091 0.7472171 61 1.21664 0.2284 
SexFemale:fage12 -0.500000 0.7472171 61 -0.66915 0.5059 
SexFemale:fage14 -0.818182 0.7472171 61 -1.09497 0.2778

總結一下小組意味着幫助解決這個問題。

require(plyr) 
ddply(Orthodont2, .(Sex, age), summarise, dist = mean(distance)) 

    Sex fage  dist 
1 Male 8 23.63636 
2 Male 10 23.77273 
3 Male 12 26.04545 
4 Male 14 27.36364 
5 Female 8 21.18182 
6 Female 10 22.22727 
7 Female 12 23.09091 
8 Female 14 24.09091

注意，第一固定效應係數，稱爲SexMale，對於年齡8名男性的平均距離。固定效應SexFemale是8歲女性平均距離。這些是最容易看到的（我總是從簡單的開始），但其餘的都不難看出來。 10歲男性的平均距離是第一個係數加上第三個係數（fage10）。 10歲女性的平均距離是係數SexFemale,fage10和SexFemale:fage10的總和。其餘部分沿着相同的路線。

一旦您知道如何爲組平均值創建係數的線性組合，您可以使用estimable來計算任何感興趣的比較。當然，這裏有一些關於主要效果的注意事項，交互的統計證據，離開交互的理論原因等等。這是由您來決定的。但是如果我離開模型中的相互作用（請注意，沒有統計證據表明相互作用，請參見anova(fit1)），但想要比較Male與Female的整體均值，我會寫出以下線性係數組合：

# male/age group means 
male8 = c(1, 0, 0, 0, 0, 0, 0, 0) 
male10 = c(1, 0, 1, 0, 0, 0, 0, 0) 
male12 = c(1, 0, 0, 1, 0, 0, 0, 0) 
male14 = c(1, 0, 0, 0, 1, 0, 0, 0) 

# female/age group means 
female8 = c(0, 1, 0, 0, 0, 0, 0, 0) 
female10 = c(0, 1, 1, 0, 0, 1, 0, 0) 
female12 = c(0, 1, 0, 1, 0, 0, 1, 0) 
female14 = c(0, 1, 0, 0, 1, 0, 0, 1) 

# overall male group mean 
male = (male8 + male10 + male12 +male14)/4 
# overall female group mean 
female = (female8 + female10 + female12 + female14)/4 

require(gmodels) 
estimable(fit1, rbind(male - female))

您可以檢查您的整體組意味着確保您正確地使係數的線性組合。

ddply(Orthodont2, .(Sex), summarise, dist = mean(distance)) 

    Sex  dist 
1 Male 25.20455 
2 Female 22.64773

來源

2013-07-31 16:38:45 aosmith

測試在nlme中的線性混合模型中的相互作用的顯着性

回答

相關問題