2016-05-09 7 views
5

3つのサイトでの前(2007年)および後(2011年)の測定値を持つBACI実験があります。コントラストを設定する際に問題があるので、3つのサイトのそれぞれで、治療のために年の効果が対照と異なるかどうかを確認します(コントロール参照)。マルチコムを使用した3者対話の対照設定

例をvignetteに拡張しようとしましたが、成功しませんでした。 lsmeansパッケージも試してみました。

例dfと以下のコードと、私がテストしたい明示的な対比。

library(nlme) 
library(multcomp) 
require(lsmeans) 

# I have these versions: 
# R version 3.2.4 Revised (2016-03-16 r70336) 
# Platform: x86_64-w64-mingw32/x64 (64-bit) 
# Running under: Windows 7 x64 (build 7601) Service Pack 1 
# nlme_3.1-126 
# multcomp_1.4-4 



exdf <- expand.grid(Site = c("Z", "Y", "X"), 
        Plot = 1:30, 
        Year = c("2007", "2011"), 
        Treatment = c("Control", "A", "B", "A+B")) 
exdf$Site <- factor(exdf$Site) 
exdf$Plot <- factor(exdf$Plot) 
exdf$Year <- factor(exdf$Year) 
exdf$Treatment <- factor(exdf$Treatment) 
exdf$Treatment <- relevel(exdf$Treatment, ref = "Control") 

exdf$Response <- rnorm(mean = 1150, sd = 130, n = nrow(exdf)) 

mod1 <- lme(data = exdf, 
      fixed = Response ~ Year * Treatment * Site, 
      random = ~1|Plot) 


## would like to test (obviously that's not actual code) 

## How do the treatments A, B, A+B differ in their interaction with 
## year, compared to the control treatment (ie do the slopes differ) 

SiteX:Year:TreatmentControl - SiteX:Year:TreatmentA == 0 
SiteX:Year:TreatmentControl - SiteX:Year:TreatmentB == 0 
SiteX:Year:TreatmentControl - SiteX:Year:TreatmentA+B == 0 

SiteY:Year:TreatmentControl - SiteY:Year:TreatmentA == 0 
SiteY:Year:TreatmentControl - SiteY:Year:TreatmentB == 0 
SiteY:Year:TreatmentControl - SiteY:Year:TreatmentA+B == 0 

SiteZ:Year:TreatmentControl - SiteZ:Year:TreatmentA == 0 
SiteZ:Year:TreatmentControl - SiteZ:Year:TreatmentB == 0 
SiteZ:Year:TreatmentControl - SiteZ:Year:TreatmentA+B == 0 

## Also potentially - does the interaction between site and year 
## differ across sites? 
SiteX:Year:TreatmentA - SiteY:Year:TreatmentA == 0 
SiteZ:Year:TreatmentA - SiteY:Year:TreatmentA == 0 
SiteZ:Year:TreatmentA - SiteX:Year:TreatmentA == 0 

SiteX:Year:TreatmentB - SiteY:Year:TreatmentB == 0 
SiteZ:Year:TreatmentB - SiteY:Year:TreatmentB == 0 
SiteZ:Year:TreatmentB - SiteX:Year:TreatmentB == 0 

SiteX:Year:TreatmentA+B - SiteY:Year:TreatmentA+B == 0 
SiteZ:Year:TreatmentA+B - SiteY:Year:TreatmentA+B == 0 
SiteZ:Year:TreatmentA+B - SiteX:Year:TreatmentA+B == 0 


## Started doing the example from 
#https://cran.r-project.org/web/packages/multcomp/vignettes/multcomp-examples.pdf 
# Got a bit stuck 

tmp <- expand.grid(Year = levels(exdf$Year), 
        Treatment = levels(exdf$Treatment), 
        Site = levels(exdf$Site)) 
X <- model.matrix(~ Year * Treatment * Site, data = tmp) 
glht(mod1, linfct = X) 


## Gives all the interactions not in the form I want.. 
tt <- lsmeans(mod1, specs = ~Year:Treatment|Site) 
comps <- pairs(tt, interaction = TRUE) 
summary(comps, adjust = "holm") 
+0

「ではないについてあなたはどんな形をしたいですか?おそらくあなたは〜年|治療*サイト、そしておそらくその結果のペア()のような何かが欲しいのだろうか? – rvl

+0

ありがとうございます。 「フォーム」とは、コード内に設定された比較をテストすることを意味します。私は上記の比較をテストしないあなたの提案(私のコードの最後にあるものに似ています)を試しました。少し異なる式 'ls < - lsmeans(mod1、specs =〜Year *治療|サイト);ペア(tt、interaction = FALSE)は適切な構造を提供しますが、年と治療のあらゆる組み合わせをテストしています。 – Walden

答えて

3

ステップアップしましょう。まず、他の要素の組み合わせごとに2年間の違いを知りたいとします。

pairs(dd, by = "Site") 
pairs(dd, by = "Treatment") 

を...またはわずかにそれらを比較:これは簡単に

tt = lsmeans(mod1, specs = ~ Year | Treatment:Site) 
dd = pairs(tt, reverse = TRUE) # (2011) - (2008) for each T*S comb 
summary(dd, by = NULL) 

を使用してLS平均で行われているあなたは、あなたが好きなようになりましたこれらの違いを比較することができ

pairs(lsmeans(dd, "Treatment")) 
pairs(lsmeans(dd, "Site")) 
+0

非常に明確な答えは、ありがとう! – Walden

関連する問題