반복 측정을위한 varIdent가있는 lme

Question

반복 측정과 varIdent를 사용하여 lme를 사용할 때 이상한 결과가 나타납니다. 이 모든 도움은 매우 감사하겠습니다!

나는 시계열에 따른 잎의 13C 신호가 2 종 (A와 B)간에 다른지 여부를 테스트하고 있습니다. 나는 기본적으로 특정 시점이 아닌 종간의 전반적인 차이에 관심이있다.

Block Species time X13C 
1 B 2 0.775040865 
2 B 2 0.343913792 
3 B 2 0.381053614 
1 A 2 0.427101597 
2 A 2 0.097743662 
3 A 2 0.748345826 
1 B 24 0.416700446 
2 B 24 0.230773558 
3 B 24 0.681386484 
1 A 24 0.334026511 
2 A 24 0.866426406 
3 A 24 0.606346215 
1 B 48 0.263085491 
2 B 48 0.083323709 
3 B 48 0.534697801 
1 A 48 0.30594443 
2 A 48 0.024555489 
3 A 48 0.790670392 
1 B 96 0.158090804 
2 B 96 0.254880689 
3 B 96 0.082666799 
1 A 96 0.139189281 
2 A 96 0.300340119 
3 A 96 0.233149535 
1 B 192 0.055421148 
2 B 192 0.082582155 
3 B 192 0.136636735 
1 A 192 0.03641637 
2 A 192 0.06082544 
3 A 192 0.126029308 

나는 다음과 같은 모델을 적용하고있다 : 여기

내 데이터 세트입니다 시간의 잔차의 이질성이 존재하므로

bulk<-lme(X13Catex ~ Species*time, random = ~1|Block/Species, method='REML', na.action=na.exclude, data=VacL, corAR1())

을, 나는 모델을 개선 varIdent이 적용 적합 (AIC). 정규화 된 잔차 그림 또한 좋게 보입니다.

bulk.var<-lme(X13Catex ~ Species*time, random = ~1|Block/Species, method='REML', na.action=na.exclude, data=VacL, corAR1(), weights=varIdent(form=~1|time))

문제는이 코드로 내가 종의 중요한 p- 값을 얻을 수 있지만,이 종은 전혀 다르다는 것을하지 않는 것 내 데이터보고 ... 나는 그것이 매우 이상한 생각이다 오차 막대가 각 시점에서 중첩되고 어떤 시점 A에서 B가 B보다 크고 다른 어떤 점에서는 반올림하는 등 낮은 p 값을 얻습니다. 문제는 각 샘플링 시간에서 각 종족에 대한 낮은 복제 될 수 있다면 나는 다른 유사한 변수를 분석 할 때 다시 일어

> anova(bulk.var) 
         numDF denDF F-value p-value 
(Intercept)     1 15 13.25772 0.0024 
SpeciesCode     1  2 67.08281 0.0146 
SamplingTime     4 15 4.42320 0.0147 
SpeciesCode:SamplingTime  4 15 1.27659 0.3227 

...

가 궁금 (N = 3). 복제가 적은 varIdent 및 "비교적 복잡한"모델을 적용하면 중요한 p 값을 발견 할 수 있습니까? 이것에 대처하는 방법에 대한 제안?

감사합니다.

Answer 1

확인하겠습니다.

우선, 귀하의 상관 관계 구조가 나에게 맞지 않는 것 같습니다. 거기에 시간 공변량을 사용해야합니다.

fit0 <- lme(X13C ~ Species*time, random = ~1|Block/Species, method='REML', 
      na.action=na.exclude, data=VacL, 
      corAR1(0.9, form = ~ time | Block/Species)) 
summary(fit0) 

그러면 중첩 된 임의 효과의 분산은 다소 작아 보입니다. 이 매개 변수를 제거해 봅시다.

fit1 <- lme(X13C ~ Species*time, random = ~1|Block, method='REML', 
      na.action=na.exclude, data=VacL, 
      corAR1(0.9, form = ~ time | Block/Species)) 
summary(fit1) 

anova(fit0, fit1) 
#  Model df  AIC  BIC logLik Test  L.Ratio p-value 
#fit0  1 8 35.47003 45.5348 -9.735014       
#fit1  2 7 33.47003 42.2767 -9.735014 1 vs 2 8.37192e-09 0.9999 

plot(fit1) 

줄거리

은 참으로 강력한 이질성을 보여줍니다. 이 시점에서 나는 GLMM을 진지하게 고려할 것이다. 델타 13C는 (변형 된) 분수 13C/12C임을 기억하십시오. 정상적인 가정은 약간의 모호성이있는 것처럼 보입니다. (가끔 델타 값으로 사용하기도합니다.) 그러나 필자는 적합 값의 의존성을 모델링 할 수있는 것처럼 보입니다.

anova(fit1, fit2) 
#  Model df  AIC  BIC logLik Test L.Ratio p-value 
#fit1  1 7 33.47003 42.27670 -9.735014       
#fit2  2 8 11.34319 21.40796 2.328405 1 vs 2 24.12684 <.0001 

하지 나쁘지

fit2 <- lme(X13C ~ Species*time, random = ~1|Block, method='REML', 
      na.action=na.exclude, data=VacL, 
      corAR1(0.9, form = ~ time | Block/Species), 
      weights = varPower()) 

plot(fit2, resid(., type = "normalized") ~ fitted(.)) 

. p- 값을 확인해 봅시다.

coef(summary(fit2)) 
#      Value Std.Error DF t-value  p-value 
#(Intercept) 0.3906798322 0.0640391495 24 6.1006405 2.661703e-06 
#SpeciesB  -0.0303078937 0.0777180616 24 -0.3899723 6.999965e-01 
#time   -0.0016541839 0.0003059863 24 -5.4060727 1.491893e-05 
#SpeciesB:time 0.0002578782 0.0004048368 24 0.6369930 5.301592e-01 

절편이나 슬로프는 크게 다릅니다. 이제 ANOVA를 사용해 봅시다. 분산 구조가없는 비교를 위해

anova(fit2) 
      numDF denDF F-value p-value 
(Intercept)  1 24 9.0061 0.0062 
Species   1 24 525.7457 <.0001 
time    1 24 56.5135 <.0001 
Species:time  1 24 0.4058 0.5302 

:

anova(fit1) 
      numDF denDF F-value p-value 
(Intercept)  1 24 29.536428 <.0001 
Species   1 24 0.319802 0.5770 
time    1 24 17.602173 0.0003 
Species:time  1 24 0.041482 0.8403 

그래서, 당신은 더 적은 네 개의 매개 변수를 사용하는 모델과 같은 문제를 얻을. 그리고 현재 분산 매개 변수가 포함되어있는 경우 해당 모델 매개 변수가 중요하지 않지만 왜 종 효과가 순차적 ANOVA에서 중요한지는 알 수 없습니다. Pinheiro & Bates 2000을 공부하고 자신을 찾거나 Cross Validated에 문의하십시오.