Model efektów mieszanych z zagnieżdżaniem

34

Mam dane zebrane z eksperymentu zorganizowane w następujący sposób:

Dwa miejsca, każde z 30 drzewami. 15 jest leczonych, 15 kontroluje w każdym miejscu. Z każdego drzewa pobieramy próbki trzech kawałków łodygi i trzech kawałków korzeni, a więc 6 poziomów 1 próbki na drzewo, które jest reprezentowane przez jeden z dwóch poziomów czynników (korzeń, łodyga). Następnie z tych próbek łodygi / korzenia pobieramy dwie próbki, dzieląc różne tkanki w próbce, co jest reprezentowane przez jeden z dwóch poziomów czynników dla typu tkanki (typ tkanki A, typ tkanki B). Próbki te są mierzone jako zmienna ciągła. Całkowita liczba obserwacji wynosi 720; 2 stanowiska * 30 drzew * (trzy próbki łodygi + trzy próbki korzenia) * (jedna próbka tkanki A + jedna próbka tkanki B). Dane wyglądają tak ...

        ï..Site Tree Treatment Organ Sample Tissue Total_Length
    1        L  LT1         T     R      1 Phloem           30
    2        L  LT1         T     R      1  Xylem           28
    3        L  LT1         T     R      2 Phloem           46
    4        L  LT1         T     R      2  Xylem           38
    5        L  LT1         T     R      3 Phloem          103
    6        L  LT1         T     R      3  Xylem           53
    7        L  LT1         T     S      1 Phloem           29
    8        L  LT1         T     S      1  Xylem           21
    9        L  LT1         T     S      2 Phloem           56
    10       L  LT1         T     S      2  Xylem           49
    11       L  LT1         T     S      3 Phloem           41
    12       L  LT1         T     S      3  Xylem           30

Próbuję dopasować model efektów mieszanych przy użyciu R i Lme4, ale jestem nowy w modelach mieszanych. Chciałbym modelować odpowiedź jako Czynnik + Czynnik Poziomu 1 (łodyga, korzeń) + Czynnik Poziomu 2 (tkanka A, tkanka B), z losowymi efektami dla określonych próbek zagnieżdżonych na dwóch poziomach.

W R robię to za pomocą lmer, jak następuje

fit <- lmer(Response ~ Treatment + Organ + Tissue + (1|Tree/Organ/Sample)) 

Z mojego zrozumienia (... co nie jest pewne i dlaczego publikuję!) Termin:

(1|Tree/Organ/Sample)

Określa, że ​​„Próbka” jest zagnieżdżona w próbkach organów, które są zagnieżdżone w drzewie. Czy tego rodzaju zagnieżdżanie jest ważne / ważne? Przepraszamy, jeśli to pytanie nie jest jasne, jeśli tak, proszę podać, gdzie mogę je rozwinąć.

Erik
źródło

Odpowiedzi:

33

Myślę, że to prawda.

  • (1|Tree/Organ/Sample)rozwija się do / jest równoważne (1|Tree)+(1|Tree:Organ)+(1|Tree:Organ:Sample)(gdzie :oznacza interakcję).
  • Ustalone czynniki Treatment, Organi Tissueautomatycznie się obchodzić na właściwym poziomie.
  • Prawdopodobnie powinieneś uwzględnić Sitejako ustalony efekt (koncepcyjnie jest to efekt losowy, ale nie jest praktyczne, aby spróbować oszacować wariancję między witrynami tylko z dwoma witrynami); to nieznacznie zmniejszy wariancję między drzewami.
  • Prawdopodobnie powinieneś zawrzeć wszystkie dane w ramce danych i przekazać je jawnie lmerza pomocą data=my.data.frameargumentu.

Pomocne może być najczęściej zadawane pytania glmm (dotyczy GLMM, ale zawiera także rzeczy istotne dla liniowych modeli mieszanych).

Ben Bolker
źródło
Co jeśli Erik chciałby określić strukturę kowariancji dla tych przechwyceń? Tzn. Można oczekiwać, że próbka z dodatnim przechwyceniem drzewa również będzie miała pozytywne przechwycenie narządu. Czy zagnieżdżanie rozwiązuje ten problem automatycznie? Jeśli nie, to jak określić taką strukturę?
Sheridan Grant
Myślę, że jeśli spróbujesz zapisać równania dla tej sprawy, przekonasz się, że została rozwiązana.
Ben Bolker,
13

Przeczytałem to pytanie i odpowiedź doktora Bolkera i próbowałem powielić dane (szczerze mówiąc, nie dbając o to, co „długość” reprezentuje w kategoriach biologicznych lub jednostkach, a następnie dopasowałem model jak powyżej. Zamieszczam wyniki tutaj do dzielenia się i szukania informacji zwrotnej na temat prawdopodobnych nieporozumień.

Kod, którego użyłem do wygenerowania tych fikcyjnych danych, można znaleźć tutaj , a zestaw danych ma wewnętrzną strukturę OP:

     site     tree treatment organ sample tissue    length
1    L       LT01         T  root      1  phloem  108.21230
2    L       LT01         T  root      1  xylem   138.54267
3    L       LT01         T  root      2  phloem   68.88804
4    L       LT01         T  root      2  xylem   107.91239
5    L       LT01         T  root      3  phloem   96.78523
6    L       LT01         T  root      3  xylem    88.93194
7    L       LT01         T  stem      1  phloem  101.84103
8    L       LT01         T  stem      1  xylem   118.30319

Struktura jest następująca:

 'data.frame':  360 obs. of  7 variables:
     $ site     : Factor w/ 2 levels "L","R": 1 1 1 1 1 1 1 1 1 1 ...
 $ tree     : Factor w/ 30 levels "LT01","LT02",..: 1 1 1 1 1 1 1 1 1 1 ...
     $ treatment: Factor w/ 2 levels "C","T": 2 2 2 2 2 2 2 2 2 2 ...
 $ organ    : Factor w/ 2 levels "root","stem": 1 1 1 1 1 1 2 2 2 2 ...
     $ sample   : num  1 1 2 2 3 3 1 1 2 2 ...
 $ tissue   : Factor w/ 2 levels "phloem","xylem": 1 2 1 2 1 2 1 2 1 2 ...
     $ length   : num  108.2 138.5 68.9 107.9 96.8 ...

Zestaw danych został „sfałszowany” (mile widziane tutaj opinie):

  1. Ponieważ treatmentistnieje ustalony efekt z dwoma wyraźnymi punktami przechwytywania dla leczenia w porównaniu do kontroli (w 100porównaniu 70) i brak efektów losowych.
  2. Ustawiam wartości tissuez widocznymi stałymi efektami z bardzo różnymi punktami przechwytywania dla phloemkontra xylem( 3kontra 6) i losowymi efektami z sd = 3.
  3. organN.(0,3))sd = 36rootstem
  4. Ponieważ treemamy po prostu losowe efekty z sd = 7.
  5. Dla samplePróbowałem założyć tylko z efektów losowych sd = 5.
  6. Dla siterównież po prostu losowych efektów z sd = 3.

Nie ustalono nachyleń ze względu na kategoryczny charakter zmiennych.

Wyniki modelu mieszanych efektów:

fit <- lmer(length ~ treatment + organ + tissue + (1|tree/organ/sample), data = trees) 

byli:

 Random effects:
 Groups              Name        Variance  Std.Dev. 
 sample:(organ:tree) (Intercept) 9.534e-14 3.088e-07
 organ:tree          (Intercept) 0.000e+00 0.000e+00
 tree                (Intercept) 4.939e+01 7.027e+00
 Residual                        3.603e+02 1.898e+01
Number of obs: 360, groups:  sample:(organ:tree), 180; organ:tree, 60; tree, 30

Fixed effects:
            Estimate Std. Error       df t value Pr(>|t|)    
(Intercept)  79.8623     2.7011  52.5000  29.567  < 2e-16 ***
treatmentT   21.4368     3.2539  28.0000   6.588 3.82e-07 ***
organstem     0.1856     2.0008 328.0000   0.093    0.926    
tissuexylem   3.1820     2.0008 328.0000   1.590    0.113    
---
Signif. codes:  0***0.001**0.01*0.05 ‘.’ 0.1 ‘ ’ 1

Jak to zadziałało:

  1. Dla treatmentprzechwytywania bez leczenia było 79.8623(ustawiłem średnią 70), a przy leczeniu było 79.8623 + 21.4368 = 101.2991(ustawiliśmy średnią) 100.
  2. Dla tissuenastąpił 3.1820wpis do uprzejmości osią xylem, a ja ustawić różnicę pomiędzy phloemi xylemod 3. Losowe efekty nie były częścią modelu.
  3. Dla organ, próbek z stemzwiększyło przechwycenia przez 0.1856- Miałem założyć żadnej różnicy w trwałych efektów pomiędzy stemi root. Standardowe odchylenie tego, co chciałem działać jako efekty losowe, nie zostało odzwierciedlone.
  4. Te treelosowe efekty z SD 7powierzchniowe jak ładnie 7.027.
  5. Jeśli chodzi o sample, inicjał sdnie 5został podkreślony jako 3.088.
  6. site nie był częścią modelu.

Ogólnie rzecz biorąc, wydaje się, że model pasuje do struktury danych.

Antoni Parellada
źródło