Calcul de la taille de l'échantillon pour la régression logistique univariée

Comment calcule-t-on la taille de l'échantillon nécessaire pour une étude dans laquelle une cohorte de sujets aura une seule variable continue mesurée au moment d'une intervention chirurgicale puis deux ans plus tard, ils seront classés comme résultat fonctionnel ou résultat altéré.

Nous aimerions voir si cette mesure aurait pu prédire le mauvais résultat. À un moment donné, nous pouvons vouloir dériver un point de coupure dans la variable continue au-dessus de laquelle nous tenterions d'intervenir pour diminuer la probabilité de l'issue altérée.

Des idées? Toute implémentation R.

Les calculs de taille d'échantillon pour la régression logistique sont complexes. Je n'essaierai pas de le résumer ici. Des solutions raisonnablement accessibles à ce problème se trouvent dans:

Hsieh FY. Exemples de tableaux de tailles pour la régression logistique. Statistiques en médecine. 1989 juil; 8 (7): 795-802.

Hsieh FY et al. Une méthode simple de calcul de la taille de l'échantillon pour la régression linéaire et logistique. Statistiques en médecine. 30 juil.1998; 17 (14): 1623-34.

Une discussion accessible des problèmes avec des exemples de calculs peut être trouvée dans le dernier chapitre (Section 8.5 pp 339-347) de Hosmer & Lemeshow's Applied Logistic Regression .

Je trouve généralement plus facile et plus rapide d'exécuter une simulation. Les articles mettent beaucoup de temps à lire, à comprendre et à conclure finalement qu'ils ne s'appliquent pas dans le cas particulier qui nous intéresse.

Par conséquent, je voudrais simplement choisir un certain nombre de sujets, simuler la covariable qui vous intéresse (distribuée comme vous le pensez), simuler de bons / mauvais résultats en fonction de la forme fonctionnelle que vous posez (effets de seuil de la covariable? Non-linéarité?) avec la taille d'effet minimale (cliniquement) significative que vous souhaitez détecter, exécutez le résultat à travers votre analyse et voyez si l'effet se trouve à votre alpha. Réexécutez cela 10 000 fois et regardez si vous avez trouvé l'effet dans 80% des simulations (ou quelle que soit la puissance dont vous avez besoin). Ajustez le nombre de sujets, répétez jusqu'à ce que vous ayez un pouvoir qui vous convient.

Cela a l'avantage d'être très général, donc vous n'êtes pas confiné à une forme fonctionnelle spécifique ou à un nombre ou une distribution spécifique de covariables. Vous pouvez inclure les abandons, voir le commentaire de chl ci-dessus, soit au hasard, soit influencé par la covariable ou le résultat. Vous codez essentiellement l'analyse que vous allez faire sur l'échantillon final au préalable, ce qui aide parfois à concentrer ma réflexion sur la conception de l'étude. Et cela se fait facilement en R (vectoriser!).

Dans le prolongement du post de Stephan Kolassa (je ne peux pas ajouter ceci en tant que commentaire), j'ai un code alternatif pour une simulation. Cela utilise la même structure de base, mais est un peu plus explosé, donc c'est peut-être un peu plus facile à lire. Il est également basé sur le code de Kleinman et Horton pour simuler la régression logistique.

nn est le nombre dans l'échantillon. La covariable doit être distribuée normalement en continu et normalisée pour signifier 0 et sd 1. Nous utilisons rnorm (nn) pour générer cela. Nous sélectionnons un rapport de cotes et le stockons dans odds.ratio. Nous choisissons également un numéro pour l'interception. Le choix de ce nombre détermine quelle proportion de l'échantillon a vécu "l'événement" (par exemple 0,1, 0,4, 0,5). Vous devez jouer avec ce nombre jusqu'à ce que vous obteniez la bonne proportion. Le code suivant vous donne une proportion de 0,1 avec une taille d'échantillon de 950 et un OR de 1,5:

nn <- 950
runs <- 10000
intercept <- log(9)
odds.ratio <- 1.5
beta <- log(odds.ratio)
proportion  <-  replicate(
              n = runs,
              expr = {
                  xtest <- rnorm(nn)
                  linpred <- intercept + (xtest * beta)
                  prob <- exp(linpred)/(1 + exp(linpred))
                  runis <- runif(length(xtest),0,1)
                  ytest <- ifelse(runis < prob,1,0)
                  prop <- length(which(ytest <= 0.5))/length(ytest)
                  }
            )
summary(proportion)

le résumé (proportion) confirme que la proportion est ~ 0,1

Puis en utilisant les mêmes variables, la puissance est calculée sur 10000 runs:

result <-  replicate(
              n = runs,
              expr = {
                  xtest <- rnorm(nn)
                  linpred <- intercept + (xtest * beta)
                  prob <- exp(linpred)/(1 + exp(linpred))
                  runis <- runif(length(xtest),0,1)
                  ytest <- ifelse(runis < prob,1,0)
                  summary(model <- glm(ytest ~ xtest,  family = "binomial"))$coefficients[2,4] < .05
                  }
            )
print(sum(result)/runs)

Je pense que ce code est correct - je l'ai vérifié par rapport aux exemples donnés dans Hsieh, 1998 (tableau 2), et il semble être d'accord avec les trois exemples qui y sont donnés. Je l'ai également testé contre l'exemple des pages 342 à 343 de Hosmer et Lemeshow, où il a trouvé une puissance de 0,75 (contre 0,8 à Hosmer et Lemeshow). Il se peut donc que dans certaines circonstances, cette approche sous-estime le pouvoir. Cependant, lorsque j'ai exécuté le même exemple dans cette calculatrice en ligne , j'ai constaté qu'il était d'accord avec moi et non le résultat dans Hosmer et Lemeshow.

Si quelqu'un peut nous dire pourquoi c'est le cas, je serais intéressé de savoir.

$\theta = 1$$_0: \theta = 0$

en fait, il semble que votre étude sera menée de manière séquentielle. dans ce cas, il peut être avantageux d'en faire une partie explicite de l'expérience. l'échantillonnage séquentiel peut souvent être plus efficace qu'une expérience à taille d'échantillon fixe [moins d'observations nécessaires, en moyenne].

farrel: j'ajoute ceci en réponse à votre commentaire.

pour obtenir une taille d'échantillon, on spécifie généralement une sorte de critère de précision pour une estimation [telle que la longueur d'un CI] OU la puissance à une alternative spécifiée d'un test à effectuer sur les données. vous semblez avoir mentionné ces deux critères. il n'y a rien de mal à cela, en principe: il vous suffit ensuite de faire deux calculs de taille d'échantillon - un pour atteindre la précision d'estimation souhaitée - et un autre pour obtenir la puissance souhaitée à l'alternative indiquée. alors la plus grande des deux tailles d'échantillon est ce qui est requis. [btw - à part dire 80% de puissance - vous ne semblez pas avoir mentionné quel test vous prévoyez d'effectuer - ou l'alternative à laquelle vous voulez la puissance de 80%.]

quant à l'utilisation de l'analyse séquentielle: si les sujets sont tous inscrits à l'étude en même temps, alors une taille d'échantillon fixe est logique. mais si les sujets sont peu nombreux, il peut prendre un an ou deux [ou plus] pour obtenir le nombre requis inscrit. le procès pouvait donc durer trois ou quatre ans [ou plus]. dans ce cas, un schéma séquentiel offre la possibilité de s'arrêter plus tôt que cela - si le ou les effets que vous recherchez deviennent statistiquement significatifs plus tôt dans le procès.