Régression logistique: variables groupées et non groupées (en utilisant R)

9

Je lis A. Agresti (2007), An Introduction to Categorical Data Analysis , 2nd. édition, et je ne sais pas si je comprends bien ce paragraphe (p.106, 4.2.1) (bien que cela devrait être facile):

Dans le tableau 3.1 sur le ronflement et les maladies cardiaques du chapitre précédent, 254 sujets ont signalé des ronflements chaque nuit, dont 30 avaient des maladies cardiaques. Si le fichier de données contient des données binaires groupées, une ligne du fichier de données signale ces données comme 30 cas de maladie cardiaque sur un échantillon de 254. Si le fichier de données contient des données binaires non groupées, chaque ligne du fichier de données fait référence à un sujet séparé, donc 30 lignes contiennent un 1 pour les maladies cardiaques et 224 lignes contiennent un 0 pour les maladies cardiaques. Les estimations ML et SE sont les mêmes pour les deux types de fichiers de données.

Transformer un ensemble de données non groupées (1 dépendant, 1 indépendant) prendrait plus qu'une "ligne" pour inclure toutes les informations!?

Dans l'exemple suivant, un ensemble de données simple (irréaliste!) Est créé et un modèle de régression logistique est créé.

À quoi ressembleraient réellement les données groupées (onglet variable?)? Comment créer le même modèle à l'aide de données groupées?

> dat = data.frame(y=c(0,1,0,1,0), x=c(1,1,0,0,0))
> dat
  y x
1 0 1
2 1 1
3 0 0
4 1 0
5 0 0
> tab=table(dat)
> tab
   x
y   0 1
  0 2 1
  1 1 1
> mod1=glm(y~x, data=dat, family=binomial())
Banquise
la source

Réponses:

11

Le tableau 3.1 est reproduit ci-dessous:

entrez la description de l'image ici

Agresti a considéré les scores numériques suivants pour le niveau de ronflement: {0,2,4,5}.

Il existe deux façons d'adapter un GLM à R: soit votre résultat est fourni sous la forme d'un vecteur de 0/1 ou d'un facteur à deux niveaux, avec les prédicteurs sur les RH de votre formule; ou vous pouvez donner une matrice avec deux colonnes de décomptes de réussite / échec comme les lhs de la formule. Cette dernière correspond à ce que Agresti appelle des données «groupées». Une troisième méthode, qui s'applique également aux milieux groupés, consisterait à utiliser l' weights=argument pour indiquer le nombre de résultats positifs et négatifs observés pour chaque catégorie du tableau de classification.

Les données en vue matricielle se liraient comme suit:

snoring <- matrix(c(24,35,21,30,1355,603,192,224), nc=2)

À partir de cela, nous pouvons générer un data.frameformat long (2484 lignes = sum(snoring)observations) comme suit:

snoring.df <- data.frame(snoring=gl(4, 1, labels=c("Never", "Occasional",
                                                   "Nearly every night", 
                                                   "Every night")),
                         disease=gl(2, 4, labels=c("Yes", "No")),
                         counts=as.vector(snoring))
snoring.df <- snoring.df[rep(seq_len(nrow(snoring.df)), snoring.df$counts), 1:2]

Et les deux modèles suivants donneront des résultats identiques:

levels(snoring.df$snoring) <- c(0, 2, 4, 5)
y <- abs(as.numeric(snoring.df$disease)-2)
x <- as.numeric(as.character(snoring.df$snoring))
fit.glm1 <- glm(y ~ x, family=binomial)

fit.glm2 <- glm(snoring ~ c(0, 2, 4, 5), family=binomial)

Autrement dit, , en utilisant la notation d'Agresti.logit[π^(X)]=-3,87+0,40X

La deuxième notation est fréquemment utilisée sur une table agrégée avec une instruction comme cbind(a, b), où aet bsont des colonnes de décomptes pour un événement binaire (voir par exemple, Modèles linéaires généralisés ). Il semble que cela fonctionnerait également lors de l'utilisation de table au lieu de matrice (comme dans votre exemple), par exemple

glm(as.table(snoring) ~ c(0, 2, 4, 5), family=binomial)
chl
la source
Merci beaucoup! Une réponse parfaite! J'ai un ajout simple: au lieu de as.table (ronflement), je suggérerais un tableau (x, y, dnn = c ('ronflement', 'maladie')) comme équivalent de mon exemple, car la transformation de non groupé en groupé les données étaient également intéressantes.
FloE
1
@FloE Vous avez raison. Vous aurez encore besoin de construire votre rhs ad hoc . Par exemple, quelque chose comme tab <- table(x,y, dnn=c('snoring','disease')); glm(tab ~ as.numeric(rownames(tab)), family=binomial)cela fonctionnerait (inversion du signe moins pour les coefficients parce que "Oui" est codé 0 au lieu de 1).
chl