Ich habe ein Experiment durchgeführt, bei dem ich verschiedene Familien aus zwei verschiedenen Bevölkerungsgruppen großgezogen habe. Jede Familie erhielt eine von zwei Behandlungen. Nach dem Experiment habe ich mehrere Merkmale an jedem Individuum gemessen. Um die Wirkung einer Behandlung oder einer Quelle sowie deren Wechselwirkung zu testen, verwendete ich ein lineares Mischeffektmodell mit der Familie als Zufallsfaktor, d. H
lme(fixed=Trait~Treatment*Source,random=~1|Family,method="ML")
Soweit so gut. Jetzt muss ich die relativen Varianzkomponenten berechnen, dh den Prozentsatz der Variation, der durch die Behandlung oder die Quelle sowie die Wechselwirkung erklärt wird.
Ohne einen zufälligen Effekt könnte ich leicht die Quadratsummen (SS) verwenden, um die Varianz zu berechnen, die durch jeden Faktor erklärt wird. Da es für ein gemischtes Modell (mit ML-Schätzung) keine SS gibt, dachte ich, ich könnte Treatment und Source auch als Zufallseffekte verwenden, um die Varianz zu schätzen, d. H
lme(fixed=Trait~1,random=~(Treatment*Source)|Family, method="REML")
In einigen Fällen konvergiert lme jedoch nicht. Daher habe ich lmer aus dem Paket lme4 verwendet:
lmer(Trait~1+(Treatment*Source|Family),data=DATA)
Wo ich die Varianzen aus dem Modell mit der Zusammenfassungsfunktion extrahiere:
model<-lmer(Trait~1+(Treatment*Source|Family),data=regrexpdat)
results<-VarCorr(model)
variances<-results[,3]
Ich erhalte die gleichen Werte wie bei der VarCorr-Funktion. Ich benutze dann diese Werte, um den tatsächlichen Prozentsatz der Variation zu berechnen, wobei die Summe als Gesamtvariation genommen wird.
Ich habe Probleme mit der Interpretation der Ergebnisse aus dem ersten Modell (mit Behandlung und Quelle als festen Effekten) und dem Zufallsmodell zur Schätzung der Varianzkomponenten (mit Behandlung und Quelle als zufälligen Effekten). Ich finde in den meisten Fällen, dass der Prozentsatz der durch jeden Faktor erklärten Varianz nicht der Signifikanz des festen Effekts entspricht.
Zum Beispiel für das Merkmal HD deutet der erste Abschnitt auf eine Tendenz zur Wechselwirkung sowie auf eine Bedeutung für die Behandlung hin. Ich stelle fest, dass die Behandlung nach einem Rückwärtsverfahren eine beinahe signifikante Tendenz aufweist. Bei der Schätzung der Varianzkomponenten stelle ich jedoch fest, dass die Quelle die höchste Varianz aufweist und bis zu 26,7% der Gesamtvarianz ausmacht.
Die lme:
anova(lme(fixed=HD~as.factor(Treatment)*as.factor(Source),random=~1|as.factor(Family),method="ML",data=test),type="m")
numDF denDF F-value p-value
(Intercept) 1 426 0.044523 0.8330
as.factor(Treatment) 1 426 5.935189 0.0153
as.factor(Source) 1 11 0.042662 0.8401
as.factor(Treatment):as.factor(Source) 1 426 3.754112 0.0533
Und der lmer:
summary(lmer(HD~1+(as.factor(Treatment)*as.factor(Source)|Family),data=regrexpdat))
Linear mixed model fit by REML
Formula: HD ~ 1 + (as.factor(Treatment) * as.factor(Source) | Family)
Data: regrexpdat
AIC BIC logLik deviance REMLdev
-103.5 -54.43 63.75 -132.5 -127.5
Random effects:
Groups Name Variance Std.Dev. Corr
Family (Intercept) 0.0113276 0.106431
as.factor(Treatment) 0.0063710 0.079819 0.405
as.factor(Source) 0.0235294 0.153393 -0.134 -0.157
as.factor(Treatment)L:as.factor(Source) 0.0076353 0.087380 -0.578 -0.589 -0.585
Residual 0.0394610 0.198648
Number of obs: 441, groups: Family, 13
Fixed effects:
Estimate Std. Error t value
(Intercept) -0.02740 0.03237 -0.846
Daher ist meine Frage, ist es richtig, was ich tue? Oder sollte ich eine andere Methode verwenden, um den Betrag der Varianz zu schätzen, der durch jeden Faktor erklärt wird (dh Behandlung, Quelle und deren Wechselwirkung). Wären die Effektgrößen beispielsweise ein geeigneterer Weg?
Antworten:
Eine gebräuchliche Methode, um den relativen Beitrag jedes Faktors zu einem Modell zu bestimmen, besteht darin, den Faktor zu entfernen und die relative Wahrscheinlichkeit mit einem Chi-Quadrat-Test zu vergleichen:
Da die Art und Weise, wie Wahrscheinlichkeiten zwischen Funktionen berechnet werden, leicht unterschiedlich sein kann, vergleiche ich sie normalerweise nur zwischen derselben Methode.
quelle
1-pchisq(logLik(model1) - logLik(model2), 1)
?