Gérer l'ajustement singulier dans les modèles mixtes

16

Disons que nous avons un modèle

mod <- Y ~ X*Condition + (X*Condition|subject)

# Y = logit variable  
# X = continuous variable  
# Condition = values A and B, dummy coded; the design is repeated 
#             so all participants go through both Conditions  
# subject = random effects for different subjects 

summary(model)
Random effects:
 Groups  Name             Variance Std.Dev. Corr             
 subject (Intercept)      0.85052  0.9222                    
         X                0.08427  0.2903   -1.00            
         ConditionB       0.54367  0.7373   -0.37  0.37      
         X:ConditionB     0.14812  0.3849    0.26 -0.26 -0.56
Number of obs: 39401, groups:  subject, 219

Fixed effects:
                 Estimate Std. Error z value Pr(>|z|)    
(Intercept)       2.49686    0.06909   36.14  < 2e-16 ***
X                -1.03854    0.03812  -27.24  < 2e-16 ***
ConditionB       -0.19707    0.06382   -3.09  0.00202 ** 
X:ConditionB      0.22809    0.05356    4.26 2.06e-05 ***

Nous observons ici un ajustement singulier, car la corrélation entre l'ordonnée à l'origine et les effets aléatoires x est -1. Maintenant, selon ce lien utile, une façon de traiter ce modèle est de supprimer les effets aléatoires d'ordre supérieur (par exemple, X: ConditionB) et de voir si cela fait une différence lors du test de singularité. L'autre consiste à utiliser l'approche bayésienne, par exemple le blmepackage pour éviter la singularité.

Quelle est la méthode préférée et pourquoi?

Je pose cette question car l'utilisation du premier ou du second conduit à des résultats différents - dans le premier cas, je supprimerai l'effet aléatoire X: ConditionB et je ne pourrai pas estimer la corrélation entre les effets aléatoires X et X: ConditionB. En revanche, l'utilisation blmeme permet de conserver X: ConditionB et d'estimer la corrélation donnée. Je ne vois aucune raison pour laquelle je devrais même utiliser les estimations non bayésiennes et supprimer les effets aléatoires lorsque des ajustements singuliers se produisent alors que je peux tout estimer avec l'approche bayésienne.

Quelqu'un peut-il m'expliquer les avantages et les problèmes de l'une ou l'autre méthode pour gérer les crises singulières?

Je vous remercie.

User33268
la source
Qu'est-ce qui vous inquiète à ce sujet corr = -1? C'est une corrélation entre des effets aléatoires.
user158565
Donc, chaque sujet vous donne deux mesures de Y, une sous condition A et une sous condition B? Si cela est vrai, pouvez-vous également nous dire si la valeur de la variable continue X change pour un sujet donné entre les conditions A et B?
Isabella Ghement
Pourquoi mettez-vous Condition dans les effets aléatoires? Avez-vous testé si cela est nécessaire?
Dimitris Rizopoulos
@ user158565 oui mais cela indique la singularité ...
User33268
@IsabellaGhement En effet. Oui, x change pour un sujet donné entre A et B. En outre, il y a une raison théorique de supposer que la régression de Y sur X est différente pour chaque sujet
User33268

Réponses:

21

Lorsque vous obtenez un ajustement singulier, cela indique souvent que le modèle est sur-ajusté - c'est-à-dire que la structure des effets aléatoires est trop complexe pour être prise en charge par les données, ce qui conduit naturellement à conseiller de supprimer la partie la plus complexe des effets aléatoires. structure (pentes généralement aléatoires). L'avantage de cette approche est qu'elle conduit à un modèle plus parcimonieux qui n'est pas trop adapté.

Cependant, avant de faire quoi que ce soit, avez-vous une bonne raison de vouloir X, Conditionet leur interaction, varier en fonction du sujet en premier lieu? La théorie de la façon dont les données sont générées le suggère-t-elle?

Si vous souhaitez ajuster le modèle avec la structure d'effets aléatoires maximale et obtenir lme4un ajustement singulier, alors l'ajustement du même modèle dans un cadre bayésien pourrait très bien vous expliquer pourquoi vous avez lme4rencontré des problèmes, en inspectant les tracés de trace et la convergence des différentes estimations de paramètres. . L'avantage de l'approche bayésienne est que ce faisant, vous pouvez découvrir un problème avec le modèle d'origine, c'est-à-dire. la raison pour laquelle la structure d'effets aléatoires maximale n'est pas prise en charge par les données) ou elle pourrait découvrir pourquoi lme4est incapable d'adapter le modèle. J'ai rencontré des situations où un modèle bayésien ne converge pas bien, à moins que des priorités informatives ne soient utilisées - qui peuvent ou non être OK.

Bref, les deux approches ont du mérite.

Cependant, je partirais toujours d'un endroit où le modèle initial est parcimonieux et informé par des connaissances de domaine expertes pour déterminer la structure d'effets aléatoires la plus appropriée. La spécification des variables de regroupement est relativement facile, mais les pentes aléatoires ne doivent généralement pas être incluses. Incluez-les uniquement si elles ont un sens théorique solide ET si elles sont étayées par les données.

Edit: Il est mentionné dans les commentaires qu'il existe de bonnes raisons théoriques d'adapter la structure maximale des effets aléatoires. Ainsi, un moyen relativement facile de procéder avec un modèle bayésien équivalent consiste à échanger l'appel vers glmeravec stan_glmerdu rstanarmpackage - il est conçu pour être plug and play. Il a des priorités par défaut, vous pouvez donc rapidement installer un modèle. Le paquet contient également de nombreux outils pour évaluer la convergence. Si vous trouvez que tous les paramètres convergent vers des valeurs plausibles, alors vous êtes tous bons. Cependant, il peut y avoir un certain nombre de problèmes - par exemple, une variance estimée à zéro ou en dessous, ou une estimation qui continue de dériver. Le site mc-stan.org possède une mine d'informations et un forum d'utilisateurs.

Robert Long
la source
1
Oui, j'ai de bonnes raisons théoriques de supposer que la régression de Y sur X devrait différer d'un sujet à l'autre pour les conditions A et B. La régression implique un style de traitement. Pouvez-vous me donner plus d'informations sur la façon d'interpréter les tracés de trace comme un outil de diagnostic pour les causes de singularité?
User33268
11

Ceci est un fil très intéressant, avec des réponses et des commentaires intéressants! Comme cela n'a pas encore été soulevé, je voulais souligner que nous avons très peu de données pour chaque sujet (si je comprends bien). En effet, chaque sujet n'a que deux valeurs pour chacune de la variable de réponse Y, la variable catégorielle Condition et la variable continue X. En particulier, nous savons que les deux valeurs de Condition sont A et B.

Si nous devions poursuivre la modélisation de régression en deux étapes au lieu de la modélisation à effets mixtes, nous ne pourrions même pas adapter un modèle de régression linéaire aux données d'un sujet spécifique, comme illustré dans l'exemple de jouet ci-dessous pour l'un des sujets:

y <- c(4, 7)
condition <- c("A", "B")
condition <- factor(condition)
x <- c(0.2, 0.4)

m <- lm(y ~ condition*x)
summary(m)

Le résultat de ce modèle thématique serait:

Call:
lm(formula = y ~ condition * x)

Residuals:
ALL 2 residuals are 0: no residual degrees of freedom!

Coefficients: (2 not defined because of singularities)
         Estimate Std. Error t value Pr(>|t|)
(Intercept)         4         NA      NA       NA
conditionB          3         NA      NA       NA
x                  NA         NA      NA       NA
conditionB:x       NA         NA      NA       NA

Residual standard error: NaN on 0 degrees of freedom
Multiple R-squared:      1,     Adjusted R-squared:    NaN 
F-statistic:   NaN on 1 and 0 DF,  p-value: NA

Notez que l'ajustement du modèle souffre de singularités, car nous essayons d'estimer 4 coefficients de régression plus l'écart-type d'erreur en utilisant seulement 2 observations.

Les singularités persisteraient même si nous observions ce sujet deux fois - plutôt qu'une fois - sous chaque condition. Cependant, si nous observions le sujet 3 fois dans chaque condition, nous nous débarrasserions des singularités:

y <- c(4, 7, 3, 5, 1, 2)
condition <- c("A", "B", "A","B","A","B")
condition <- factor(condition)
x <- c(0.2, 0.4, 0.1, 0.3, 0.3, 0.5)

m2 <- lm(y ~ condition*x)
summary(m2)

Voici la sortie R correspondante pour ce deuxième exemple, dont les singularités ont disparu:

>     summary(m2)

Call:
lm(formula = y ~ condition * x)

Residuals:
    1       2       3       4       5       6 
1.3333  2.3333 -0.6667 -1.1667 -0.6667 -1.1667 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)
(Intercept)     4.667      3.555   1.313    0.320
conditionB      6.000      7.601   0.789    0.513
x             -10.000     16.457  -0.608    0.605
conditionB:x   -5.000     23.274  -0.215    0.850

Residual standard error: 2.327 on 2 degrees of freedom
Multiple R-squared:  0.5357,    Adjusted R-squared:  -0.1607 
F-statistic: 0.7692 on 3 and 2 DF,  p-value: 0.6079

Bien sûr, le modèle à effets mixtes ne correspond pas à des modèles de régression linéaire indépendants et non liés pour chaque sujet - il correspond à des modèles "apparentés" dont les intersections et / ou les pentes s'écartent de manière aléatoire autour d'une interception et / ou d'une pente typique, de sorte que les écarts aléatoires par rapport à la l'ordonnée à l'origine typique et / ou la pente typique suivent une distribution normale avec un zéro moyen et un écart type inconnu.

Malgré cela, mon intuition suggère que le modèle à effets mixtes a du mal avec le petit nombre d'observations - seulement 2 - disponibles pour chaque sujet. Plus le modèle est chargé de pentes aléatoires, plus il a probablement du mal. Je soupçonne que, si chaque sujet apportait 6 observations au lieu de 2 (c'est-à-dire 3 par condition), il n'aurait plus de mal à s'adapter à toutes les pentes aléatoires.

Il me semble que cela pourrait être (?) Un cas où la conception actuelle de l'étude ne prend pas en charge les ambitions de modélisation complexes - pour soutenir ces ambitions, plus d'observations seraient nécessaires dans chaque condition pour chaque sujet (ou au moins pour certains des sujets?). Ceci est juste mon intuition donc j'espère que d'autres pourront ajouter leurs idées à mes observations ci-dessus. Merci d'avance!

Isabella Ghement
la source
Je dois vous corriger - chaque participant a 30 observations pour X et Y, dans les conditions A et B!
User33268
2
Oh, cela n'était pas indiqué dans votre réponse initiale, il m'aurait donc été impossible de deviner combien d'observations vous avez réellement par sujet et par condition. Il se passe autre chose alors. Avez-vous essayé de standardiser votre variable X? Est-ce que cela aide le lme à s'adapter? De plus, avez-vous examiné les graphiques de Y contre X (ou X normalisé) pour Condition = A contre Condition = B séparément pour chaque sujet? Cela pourrait vous donner des indices supplémentaires sur ce qui se passe.
Isabella Ghement
Je n'ai pas standardisé x car ce sont des données de temps de réaction et c'est important pour l'interprétation du coefficient de régression. Cependant, les données ont été centralisées. Je vais regarder dans les parcelles individuelles, et voir ...
User33268
1
@ User33268 Je suis un peu en retard à la fête, mais vous pouvez toujours interpréter les coefficients standardisés, il vous suffit de stocker les valeurs utilisées pour la mise à l'échelle, puis de retransformer après avoir exécuté le modèle.
Frans Rodenburg