Comment dois-je modéliser les interactions entre les variables explicatives lorsque l'une d'entre elles peut avoir des termes quadratiques et cubiques?

10

J'espère sincèrement avoir formulé cette question de manière à ce qu'elle puisse être résolue définitivement - sinon, faites-le moi savoir et je vais réessayer! Je devrais également, je suppose, noter que j'utiliserai R pour ces analyses.

J'ai plusieurs mesures plant performance (Ys)dont je soupçonne qu'elles ont été influencées par quatre traitements que j'ai imposés - flower thinning (X1), fertilization (X2), leaf clipping (X3), et biased flower thinning (X4). Pour tous les Y possibles, N est au moins 242, donc mes tailles d'échantillon étaient grandes. Toutes les parcelles ont été éclaircies ou non, mais chaque parcelle a également été soumise à l'un (et un seul) des trois autres traitements (ou non - il y avait également des parcelles témoins). L'idée de cette conception était de tester si les trois autres traitements étaient capables de "masquer" ou "d'améliorer" les effets de l'amincissement. Ainsi, par conception, les trois derniers traitements (X2-X4) ne pouvaient pas interagir les uns avec les autres car ils n'étaient pas croisés, mais ils peuvent chacun interagir avec l'éclaircissement des fleurs - et ils le font probablement.

Mes hypothèses explicites sont que 1) l'éclaircie de la floraison sera significative et que 2) les termes d'interaction X1*X2, X1*X3, and X1*X4,entre l'éclaircie des fleurs et les trois autres traitements seront également significatifs. C'est-à-dire que l'éclaircissage des fleurs devrait avoir de l'importance, mais les façons dont il importe devraient être considérablement modifiées par ce que les trois autres traitements ont fait.

Je voudrais inclure toutes ces informations dans un modèle mixte:

Y ~ X0 + X1 + X2 + X3 + X4 + X1*X2 + X1*X3 + X1*X4 + (Up to three random effects)

Mais il y a un problème: j'ai de bonnes raisons de croire que les effets de l'amincissement sur Y ne sont pas linéaires. Ils sont probablement quadratiques mais peut-être même cubiques dans certains cas. En effet, les effets de l'amincissement sur les performances sont très susceptibles d'augmenter plus rapidement à des niveaux d'amincissement plus élevés. Si j'essaie de modéliser cette relation non linéaire via l'équation ci-dessus en ajoutant des termes quadratiques et cubiques pour X1, je ne sais pas comment modéliser les termes d'interaction - suis-je censé inclure toutes les combinaisons possibles de X1, (X1) ^ 2 et (X1) ^ 3 * X2, X3 et X4? Parce que cela semble être un grand nombre de paramètres à essayer d'estimer, même avec le nombre de points de données que j'ai, et je ne sais pas comment interpréter les résultats que j'obtiendrais. Cela dit, je n'ai aucune raison biologique de penser que ce serait une façon imprudente de modéliser la situation.

J'ai donc trois réflexions sur la façon de résoudre ce problème:

  1. Ajustez d'abord un modèle plus petit, par exemple Y ~ X1 + X1^2 + X^3 + Random effects, dans le seul but de déterminer si la relation entre l'amincissement et Y est linéaire, quadratique ou cubique, puis transformez l'amincissement via une racine carrée ou cubique pour linéariser la relation de manière appropriée. De là, les termes d'interaction peuvent être modélisés comme ci-dessus avec la variable transformée.
  2. Supposons que les interactions significatives, si elles se produisent, n'affectent qu'un seul des termes X1 (c'est-à-dire uniquement le terme linéaire, quadratique ou cubique), et modélisez les interactions en conséquence. Je ne sais même pas si cette approche a du sens.
  3. Il suffit d'adapter le «modèle complet» avec tous les termes d'interaction possibles entre les termes d'amincissement et les autres traitements, comme discuté ci-dessus. Ensuite, élaguez les termes d'interaction insignifiants et utilisez des graphiques et d'autres techniques pour interpréter les résultats.

Laquelle de ces approches, le cas échéant, a le plus de sens et pourquoi, étant donné que je suis intéressé par les tests d'hypothèses et non par la sélection de modèles? En particulier, si # 1 ci - dessus n'a pas de sens, pourquoi? J'ai lu cet article et cet article et j'ai essayé de comprendre ce qu'ils pourraient signifier pour moi, mais toutes les sources pour une lecture plus approfondie seraient également très appréciées!

Bajcz
la source

Réponses:

7

Aucune de ces approches ne fonctionnera correctement. L'approche 3. s'est rapprochée, mais vous avez dit que vous élagueriez les termes insignifiants. Cela est problématique parce que les colinéarités ne permettent pas de trouver les termes à supprimer, et parce que cela vous donnerait les mauvais degrés de liberté dans les tests d'hypothèse si vous souhaitez conserver l'erreur de type I.

En fonction de la taille effective de l'échantillon et du rapport signal / bruit dans votre problème, je suggère d'ajuster un modèle avec tous les termes de produit et d'effet principal, et d'interpréter le modèle à l'aide de tracés et de "tests de blocs" (plusieurs tests df de termes connexes, c'est-à-dire, un test d'interaction globale, test d'interaction non linéaire, test d'effet global incluant effet principal + interaction, etc.). Le rmspackage R facilite cette opération pour les modèles univariés standard et pour les modèles longitudinaux lorsque est normal multivarié. Exemple:Y

# Fit a model with splines in x1 and x2 and tensor spline interaction surface
# for the two.  Model is additive and linear in x3.
# Note that splines typically fit better than ordinary polynomials
f <- ols(y ~ rcs(x1, 4) * rcs(x2, 4) + x3)
anova(f)   # get all meaningful hypothesis tests that can be inferred
           # from the model formula
bplot(Predict(f, x1, x2))    # show joint effects
plot(Predict(f, x1, x2=3))   # vary x1 and hold x2 constant

Lorsque vous voyez le anovatableau, vous verrez des lignes étiquetées All Interactionsqui, pour l'ensemble du modèle, testent l'influence combinée de tous les termes d'interaction. Pour un prédicteur individuel, cela n'est utile que lorsque le prédicteur interagit avec plusieurs variables. Il y a une option dans la printméthode pour anova.rmsmontrer par chaque ligne du tableau exactement quels paramètres sont testés par rapport à zéro. Tout cela fonctionne avec des mélanges de prédicteurs catégoriques et continus.

Si vous souhaitez utiliser des polynômes ordinaires, utilisez polplutôt que rcs.

Malheureusement, je n'ai pas implémenté de modèles à effets mixtes.

Frank Harrell
la source
1
Merci pour cette réponse. Je n'ai jamais utilisé de splines auparavant, mais je pense que je comprends votre exemple. J'ai quelques questions de suivi, si ça va? 1. Lorsque vous regardez les résultats anova de ols, comme dans votre exemple, que signifie "Toutes les interactions" sous un facteur? Autrement dit, toutes les interactions avec quoi? 2. Une approche similaire sera-t-elle autorisée dans une approche de modélisation mixte? Je pense que je suis coincé avec besoin de facteurs aléatoires. Votre exemple est-il compatible avec, par exemple, lme4? 3. Est-ce que cela fonctionnera si certains des traitements en interaction sont catégoriques? Par exemple, que se passerait-il si X2 était un facteur à 2 niveaux?
Bajcz
2

Je suis un fan de l'utilisation des régressions de lissage non paramétriques pour évaluer les formes fonctionnelles des relations entre les variables dépendantes et les prédicteurs, même lorsque je vais ensuite estimer des modèles de régression paramétrique. Bien que j'aie très souvent trouvé des relations non linéaires, je n'ai jamais trouvé de terme d'interaction d'interaction non linéaire, même lorsque les principaux effets sont fortement non linéaires. Mon retour à la maison: les effets d'interaction n'ont pas besoin d'être composés des mêmes formes fonctionnelles que les prédicteurs dont ils sont constitués.

Alexis
la source
Donc, pour clarifier, votre retour à la maison est que si je choisis l'option # 2, je peux en toute sécurité simplement inclure des termes d'interaction avec le terme linéaire X1 et ne pas me soucier des "termes d'interaction d'ordre supérieur", par exemple X1 ^ 2 * X3 et ainsi de suite?
Bajcz
1
@Bajcz Eh bien ... je suppose que je dis deux choses: (1) j'ai réussi à me débrouiller dans les ensembles de données que j'ai rencontrés avec des interactions linéaires uniquement, mais aussi (2) j'aime regarder (en utilisant des régressions non paramétriques) et laissez les données me dire si je dois ou non envisager des alternatives non linéaires. [L'adoption d'une approche d'ajustement de modèle ou de test d'hypothèse à des termes non linéaires n'est pas la bonne façon de procéder de l'OMI, car cela implique, par exemple, une inférence basée, par exemple, sur un ensemble arbitraire de termes polynomiaux, plutôt que sur les données elles-mêmes.]
Alexis
3
Il n'y a pas de grande raison de croire que les interactions sont plus susceptibles d'être linéaires. J'ai rencontré de grands exemples d'interactions non linéaires. L'idée de «regarder» et de «laisser les données vous dire» est confrontée à des problèmes d'inférence, y compris des problèmes de couverture de l'intervalle de confiance.
Frank Harrell
1
@FrankHarrell Merci! Votre première phrase est exactement le point que j'essayais de faire passer dans mon (2) dans le commentaire ci-dessus (mon expérience passée peut varier considérablement à l'avenir). OTOH: ne pas laisser les données parler est une excellente stratégie pour imputer des inférences sur les artefacts des hypothèses de modélisation sur des inférences sur les données réelles.
Alexis