Manière appropriée de gérer un tableau d'urgence à 3 niveaux

J'ai un tableau de contingence à trois niveaux, avec des données de comptage pour plusieurs espèces, la plante hôte à partir de laquelle elles ont été collectées et si cette collecte a eu lieu un jour de pluie (cela importe vraiment!). En utilisant R, de fausses données pourraient ressembler à ceci:

count    <- rpois(8, 10)
species  <- rep(c("a", "b"), 4)
host     <- rep(c("c","c", "d", "d"), 2)
rain     <- c(rep(0,4), rep(1,4))
my.table <- xtabs(count ~ host + species + rain)


, , rain = 0

    species
host  a  b
   c 12 15
   d 10 13

, , rain = 1

    species
host  a  b
   c 11 12
   d 12  7

Maintenant, je veux savoir deux choses: les espèces sont-elles associées aux plantes hôtes? "Pluie ou pas" affecte-t-elle cette association? Je loglm()de MASScette:

 # Are species independent to host plants, given the effect of rain?
loglm(~species + host + rain + species*rain + host*rain, data=my.table)

 # Given any relationship between host plants and species, does rain change it?
loglm(~species + host + rain + species*host)

C'est un peu en dehors de mon niveau de confort, et je voulais vérifier que j'avais bien configuré les modèles et que c'était la meilleure façon d'aborder ces questions.

Il y a deux façons d'interpréter votre première question, qui se reflètent dans les deux façons que vous lui avez posées: «Les espèces sont-elles associées aux plantes hôtes?» »et« Les espèces sont-elles indépendantes des plantes hôtes, étant donné l'effet de la pluie? »

La première interprétation correspond à un modèle d' indépendance conjointe , qui stipule que les espèces et les hôtes sont dépendants, mais conjointement indépendamment de la pluie:

$\quad p_{shr} = p_{sh} p_r$

où est la probabilité qu'une observation tombe dans le cellule où espèces d'index, type d'hôte, et la valeur de la pluie, est la probabilité marginale de la cellule où l'on s'effondre sur la variable pluie, et est la probabilité marginale de pluie.$p_{shr}$$(s,h,r)$$s$$h$$r$$p_{sh}$$(s,h,\cdot)$$p_r$

La deuxième interprétation correspond à un modèle d' indépendance conditionnelle , qui stipule que les espèces et les hôtes sont indépendants selon qu'il a plu:

$\quad p_{sh|r} = p_{s|r}p_{h|r}$ ou$p_{shr} = p_{sr}p_{hr} / p_r$

où est la probabilité conditionnelle de la cellule , étant donné une valeur de . ( s , h , r ) $p_{sh|r}$$(s,h,r)$$r$

Vous pouvez tester ces modèles en R ( loglincela fonctionnerait bien aussi mais je suis plus familier glm):

count <- c(12,15,10,13,11,12,12,7)
species <- rep(c("a", "b"), 4)
host <- rep(c("c","c", "d", "d"), 2)
rain <- c(rep(0,4), rep(1,4))
my.table <- xtabs(count ~ host + species + rain)
my.data <- as.data.frame.table(my.table)
mod0 <- glm(Freq ~ species + host + rain, data=my.data, family=poisson())
mod1 <- glm(Freq ~ species * host + rain, data=my.data, family=poisson())
mod2 <- glm(Freq ~ (species + host) * rain, data=my.data, family=poisson())
anova(mod0, mod1, test="Chi") #Test of joint independence
anova(mod0, mod2, test="Chi") #Test of conditional independence

Ci-dessus, mod1correspond à l'indépendance conjointe et mod2correspond à l'indépendance conditionnelle, alors que mod0correspond à un modèle d' indépendance mutuelle . Vous pouvez voir les estimations des paramètres à l'aide de , etc. Comme d'habitude, vous devez vérifier si les hypothèses du modèle sont respectées. Dans les données que vous avez fournies, le modèle nul correspond en fait de manière adéquate.$p_{shr} = p_s p_h p_r$summary(mod2)

Une autre façon d'aborder votre première question serait d'effectuer le test exact de Fischer ( fisher.test(xtabs(count ~ host + species))) sur la table 2x2 repliée (première interprétation) ou le test de Mantel-Haenszel ( mantelhaen.test(xtabs(count ~ host + species + rain))) pour les tables 2x2 stratifiées ou d'écrire un test de permutation qui respecte la stratification (deuxième interprétation).

Pour paraphraser votre deuxième question, la relation entre l'espèce et l'hôte dépend-elle de la pluie?

mod3 <- glm(Freq ~ species*host*rain - species:host:rain, data=my.data, family=poisson())
mod4 <- glm(Freq ~ species*host*rain, data=my.data, family=poisson())
anova(mod3, mod4, test=”Chi”)
pchisq(deviance(mod3), df.residual(mod3), lower=F)

Le modèle complet mod4est saturé, mais vous pouvez tester l'effet en question en regardant la déviance mod3comme je l'ai fait ci-dessus.

La régression logistique semble appropriée à votre problème. La variable que vous essayez de prédire est la probabilité qu'une observation (qui est soit l'espèce A soit l'espèce B) soit l'espèce A. Les covariables sont l' , la et éventuellement la .r a i n h o s t ∗ r a i $host$$rain$$host*rain$

La commande R serait:

glm (formule = espèce ~ hôte + pluie, famille = binôme (logit), poids = comptes)

et vous serez intéressé par les valeurs des pistes. Gardez à l'esprit que vous testez plusieurs hypothèses, cependant.$p$

Au début, j'ai suggéré d'essayer l'une des techniques d'ordination contraintes du veganpackage, mais après une seconde réflexion, je doute que cela soit utile, car vous avez en fait 2 tables de contingence. J'espère que la deuxième partie de cet exemple [PDF: Démonstration R - Analyse catégorique] sera utile.