Wie kann ich testen, ob ein zufälliger Effekt signifikant ist?

34

Ich versuche zu verstehen, wann ein zufälliger Effekt verwendet werden soll und wann er nicht erforderlich ist. Ich habe eine Faustregel erhalten, wenn Sie 4 oder mehr Gruppen / Individuen haben, die ich tue (15 einzelne Elche). Einige dieser Elche wurden zwei- oder dreimal für insgesamt 29 Versuche experimentiert. Ich möchte wissen, ob sie sich anders verhalten, wenn sie sich in Risikolandschaften befinden, als nicht. Also dachte ich, ich würde das Individuum als zufälligen Effekt einstellen. Jetzt wird mir jedoch gesagt, dass es nicht erforderlich ist, die Person als zufälligen Effekt einzubeziehen, da ihre Reaktion nicht sehr unterschiedlich ist. Was ich nicht herausfinden kann, ist, wie man testet, ob wirklich etwas berücksichtigt wird, wenn man Individuum als zufälligen Effekt einstellt. Vielleicht ist eine erste Frage: Mit welchem ​​Test / welcher Diagnose kann ich herausfinden, ob Individual eine gute erklärende Variable ist und ob es sich um einen festen Effekt handelt - qq-Diagramme? Histogramme? Streudiagramme? Und wonach würde ich in diesen Mustern suchen?

Ich habe das Modell mit dem Individuum als zufälligen Effekt und ohne ausgeführt, aber dann habe ich http://glmm.wikidot.com/faq gelesen, wo es heißt:

Vergleichen Sie lmer Modelle nicht mit den entsprechenden lm-Passungen oder glmer / glm; Die Log-Wahrscheinlichkeiten sind nicht angemessen (dh sie enthalten verschiedene additive Ausdrücke)

Und hier gehe ich davon aus, dass dies bedeutet, dass Sie kein Modell mit oder ohne Zufallseffekt vergleichen können. Aber ich würde sowieso nicht wirklich wissen, was ich zwischen ihnen vergleichen sollte.

In meinem Modell mit dem Zufallseffekt habe ich auch versucht, anhand der Ausgabe festzustellen, welche Art von Beweisen oder Bedeutung die RE hat

lmer(Velocity ~ D.CPC.min + FD.CPC + (1|ID), REML = FALSE, family = gaussian, data = tv)

Linear mixed model fit by maximum likelihood 
Formula: Velocity ~ D.CPC.min + FD.CPC + (1 | ID) 
   Data: tv 
    AIC    BIC logLik deviance REMLdev
 -13.92 -7.087  11.96   -23.92   15.39
Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.00000  0.00000 
 Residual             0.02566  0.16019 
Number of obs: 29, groups: ID, 15

Fixed effects:
              Estimate Std. Error t value
(Intercept)  3.287e-01  5.070e-02   6.483
D.CPC.min   -1.539e-03  3.546e-04  -4.341
FD.CPC       1.153e-04  1.789e-05   6.446

Correlation of Fixed Effects:
          (Intr) D.CPC.
D.CPC.min -0.010       
FD.CPC    -0.724 -0.437

Sie sehen, dass meine Varianz und SD aus der individuellen ID als Zufallseffekt = 0 ist. Wie ist das möglich? Was bedeutet 0? Ist das richtig? Dann ist mein Freund, der sagte, "da es keine Variation gibt, die ID als zufälligen Effekt verwendet, unnötig", richtig? Würde ich es dann als festen Effekt verwenden? Aber würde die Tatsache, dass es so wenig Variationen gibt, nicht bedeuten, dass es uns sowieso nicht viel sagen wird?

Kerry
quelle
Informationen zum Erzielen einer Abweichung eines zufälligen Effekts von 0 finden Sie unter stats.stackexchange.com/questions/115090 .
Amöbe sagt Reinstate Monica

Antworten:

21

Die geschätzte IDVarianz = 0 gibt an, dass das Niveau der Variabilität zwischen Gruppen nicht ausreicht, um die Einbeziehung zufälliger Effekte in das Modell zu rechtfertigen. dh Ihr Modell ist entartet.

Wie Sie sich richtig identifizieren: höchstwahrscheinlich ja; IDda ein zufälliger effekt nicht nötig ist. Es fallen mir nur wenige Dinge ein, um diese Annahme zu testen:

  1. Sie können REML = Fden AIC (oder Ihren Lieblings-IC im Allgemeinen) zwischen einem Modell mit und ohne zufällige Effekte vergleichen ( immer verwenden) und sehen, wie dies funktioniert.
  2. Sie würden sich die anova()Ausgabe der beiden Modelle ansehen .
  3. Sie können einen parametrischen Bootstrap mit der von Ihrem Originalmodell definierten posterioren Verteilung durchführen.

Bedenken Sie, dass die Entscheidungen 1 und 2 ein Problem haben: Sie überprüfen, ob sich etwas an den Grenzen des Parameterbereichs befindet, sodass diese technisch nicht stimmen. Trotzdem glaube ich nicht, dass Sie von ihnen falsche Erkenntnisse erhalten, und viele Leute verwenden sie (z. B. Douglas Bates, einer der Entwickler von lme4, verwendet sie in seinem Buch , weist jedoch eindeutig darauf hin, dass Parameterwerte getestet werden an der Grenze der Menge möglicher Werte). Wahl 3 ist die langweiligste der 3, gibt Ihnen aber tatsächlich die beste Vorstellung davon, was wirklich los ist. Einige Leute sind versucht, auch nicht-parametrische Bootstraps zu verwenden, aber ich denke, angesichts der Tatsache, dass Sie zu Beginn parametrische Annahmen treffen, könnten Sie diese auch verwenden.

usεr11852 sagt Reinstate Monic
quelle
6
Das RLRsim-Paket ist eine sehr praktische Möglichkeit, zufällige Effekte mithilfe simulationsbasierter Likelihood-Ratio-Tests zu testen.
Atrichornis
@atrichornis: +1. Interessantes Paket; Ich habe nichts davon gewusst. Ich habe mir nur den Code angesehen, ganz einfach, könnte ich sagen. Ich wünschte, sie würden es (oder so etwas) lme4besonders jetzt integrieren, wo mcmcsamp()es kaputt ist und die Leute nur noch ihre eigenen Ad-hoc-Bootstrap-Implementierungen haben, um ein paar anständige p-Werte usw. herauszubekommen.
usεr11852 sagt Reinstate Monic
Richtig, gemischte Modelle sind in R nicht einfach. Viele Annäherungen und Problemumgehungen ... Obwohl ich SAS usw. sammle, beschönige ich nur einige der gleichen Unsicherheiten? Ben Bolker ist Mitverfasser beider Pakete, er hat möglicherweise Gründe, warum er sie nicht einbezieht. Wahrscheinlich Zeit!
Atrichornis
4
Der Bootstrap an der Grenze des Parameterraums hat seine eigenen Probleme und Probleme, die zu Inkonsistenzen führen . Der Bootstrap ist kein Allheilmittel und sollte nicht leichtfertig in die Tasche geworfen werden, da sonst alles gelöst wird.
StasK
2
Schauen Sie doch mal rein, das Argument ist sehr subtil. Soweit ich mich erinnern kann, läuft es darauf hinaus, dass Sie den Bootstrap von einer Distribution ausführen, die sich von der Null unterscheidet. Angesichts der an der Grenze ermittelten Nicht-Standardverteilungen werden die Regelmäßigkeitsbedingungen verletzt, und die Bootstrap-Verteilung konvergiert nicht zum Ziel. Ich denke, dass nicht-parametrisches Bootstrap hier immer noch konstruiert werden kann, indem die Gruppenmittel der Residuen herausgenommen werden. Bei Verletzung der Unabhängigkeit von Beobachtungen zwischen Gruppen kann jedoch eine weitere Komplikationsschicht auftreten.
StasK
3

Ich bin mir nicht sicher, ob der Ansatz, den ich vorschlagen werde, vernünftig ist, daher korrigieren mich diejenigen, die mehr über dieses Thema wissen, wenn ich mich irre.

Mein Vorschlag ist, eine zusätzliche Spalte in Ihren Daten zu erstellen, die einen konstanten Wert von 1 hat:

IDconst <- factor(rep(1, each = length(tv$Velocity)))

Anschließend können Sie ein Modell erstellen, das diese Spalte als Zufallseffekt verwendet:

fm1 <- lmer(Velocity ~ D.CPC.min + FD.CPC + (1|IDconst), 
  REML = FALSE, family = gaussian, data = tv)

Zu diesem Zeitpunkt könnten Sie Ihr ursprüngliches Modell mit dem Zufallseffekt ID(nennen wir es fm0) mit dem neuen Modell vergleichen (AIC) , das dies nicht berücksichtigt, IDda IDconstes für alle Ihre Daten gleich ist.

anova(fm0,fm1)

Aktualisieren

user11852 hat nach einem Beispiel gefragt, da der obige Ansatz seiner Meinung nach nicht einmal ausgeführt werden kann. Im Gegenteil, ich kann zeigen, dass der Ansatz funktioniert (zumindest mit dem lme4_0.999999-0, den ich gerade benutze).

set.seed(101)
dataset <- expand.grid(id = factor(seq_len(10)), fac1 = factor(c("A", "B"),
  levels = c("A", "B")), trial = seq_len(10))
dataset$value <- rnorm(nrow(dataset), sd = 0.5) +
      with(dataset, rnorm(length(levels(id)), sd = 0.5)[id] +
      ifelse(fac1 == "B", 1.0, 0)) + rnorm(1,.5)
    dataset$idconst <- factor(rep(1, each = length(dataset$value)))

library(lme4)
fm0 <- lmer(value~fac1+(1|id), data = dataset)
fm1 <- lmer(value~fac1+(1|idconst), data = dataset)

anova(fm1,fm0)

Ausgabe:

  Data: dataset
  Models:
  fm1: value ~ fac1 + (1 | idconst)
  fm0: value ~ fac1 + (1 | id)

      Df    AIC    BIC  logLik  Chisq Chi Df Pr(>Chisq)
  fm1  4 370.72 383.92 -181.36                      
  fm0  4 309.79 322.98 -150.89 60.936      0  < 2.2e-16 ***

Nach diesem letzten Test sollten wir den Zufallseffekt beibehalten, da das fm0Modell sowohl den niedrigsten AIC als auch den niedrigsten BIC aufweist.

Update 2

Der gleiche Ansatz wird übrigens von NW Galwey in "Einführung in die gemischte Modellierung: Jenseits von Regression und Varianzanalyse" auf den Seiten 213-214 vorgeschlagen.

VLC
quelle
Hast du deine Idee getestet? Bitte beweise, dass ich falsch liege, aber ich denke, deine Idee wird nicht einmal umgesetzt. Wenn das IDconstfür alle Ihre Daten gleich ist, haben Sie keine Gruppierung. Sie benötigen einen Gruppierungsfaktor, um mindestens eine Stichprobenebene zu haben, und die Art und Weise, wie Sie das Modell einrichten, hat keine. Ich könnte vielleicht glauben, dass es sinnvoll ist, eine "zufällige Gruppierung" zu verwenden, aber das ist alles in allem ein anderes Ballspiel. Testen Sie Ihren Ansatz mit einigen Dummy-Daten. Ich bin der festen Überzeugung, dass mit Ihrem vorgeschlagenen Setup lmer()nicht ausgeführt werden kann. (Ich benutze lme4_0.99999911-1)
usεr11852 sagt Reinstate Monic
@ user11852 Bitte sehen Sie sich mein Update an und lassen Sie uns wissen, ob dieser Ansatz auch mit funktioniert lme4_0.99999911-1.
VLC
Z
3
Ja, ich habe getan, was Sie vorgeschlagen haben. es wird nicht funktionieren / rechnen. Error in lFormula(formula = value ~ fac1 + (1 | idconst), data = dataset) : grouping factors must have at least 1 sampled level. Und wie gesagt, konzeptionell ist es falsch. Es geht nicht darum, die Software dazu zu bringen, Zahlen preiszugeben, es geht darum, ob Sie es für sinnvoll halten. Sie müssen kein zweites gemischtes Modell vergleichen, wenn in diesem Modell der Zufallseffekt konstruktionsbedingt eine Konstante ist. Sie können es auch ausschließen und stattdessen ein lineares Modell verwenden.
usεr11852 sagt Reinstate Monic
1
Aktualisieren Sie die Abstimmung, indem Sie eine einzelne Gruppen-Zufallsvariable in definieren lme4. Dies kann geschehen, wenn Sie die Option festgelegt: control=lmerControl(check.nlev.gtr.1="ignore"). Ben Bolker erwähnt es hier: github.com/lme4/lme4/issues/411 .
Robin Beaumont
1

Ich würde gerne die anfänglichere Frage beantworten.

Wenn Sie aufgrund einiger Faktoren eine Heterogenität in der Varianz zwischen den abhängigen Variablen vermuten, sollten Sie fortfahren und die Daten mithilfe von Streu- und Boxplots zeichnen. Nach einigen gebräuchlichen Mustern zu suchen, habe ich diese Liste aus verschiedenen Quellen im Internet unten abgelegt.

Heteroskedastizitätsmuster

Zeichnen Sie außerdem Ihre abhängige Variable nach Faktor / Behandlungsgruppe, um festzustellen, ob eine konstante Varianz vorliegt. Wenn nicht, möchten Sie möglicherweise zufällige Effekte oder gewichtete Regressionen untersuchen. Zum Beispiel. Diese Grafik unten ist ein Beispiel für eine trichterförmige Varianz in meinen Behandlungsgruppen. Also entscheide ich mich für zufällige Effekte und teste die Auswirkungen auf Steigung und Abschnitte.

Boxplot zur Überprüfung der Heteroskedastizität

Von hier aus beantworten die obigen Antworten Ihre Hauptfrage. Es gibt auch Tests, die auf Heteroskedastizität prüfen, einer davon ist hier - https://dergipark.org.tr/download/article-file/94971 . Ich bin mir jedoch nicht sicher, ob Tests existieren, um eine Heteroskedastizität auf Gruppenebene festzustellen.

Arun
quelle
Bitte verwenden Sie nur das Feld "Ihre Antwort", um die Frage des OP zu beantworten. CV ist eine strenge Q & A-Site, kein Diskussionsforum. Der letztere, fettgedruckte Teil Ihres Beitrags ist eine neue Frage, keine Antwort auf diese Frage. Wenn Sie eine neue Frage haben, klicken Sie oben auf das graue ASK QUESTIONSymbol und stellen Sie sie dort. Da Sie neu hier sind, möchten Sie vielleicht an unserer Tour teilnehmen , die Informationen für neue Benutzer enthält.
gung - Wiedereinsetzung von Monica