Intervalle de confiance pour le nombre de patients atteints de coronavirus (calcul de la mortalité)

Un argument populaire pour une publication virale sur le coronavirus - comment obtenir des statistiques pour trois cas? Vous ne pouvez pas tirer de conclusions sur de si petits échantillons! Tous ceux qui ont étudié les sciences sociales ont absorbé cette histoire sur la taille des échantillons avec du lait maternel d'Alma. Et cela est vrai dans les situations que nous traitons habituellement - avec des statistiques sélectives.

Pour le cas des trois défunts, ces statistiques sont très indirectement liées. Dans les années où j'enseignais encore les méthodes mathématiques aux psychologues à l'université, j'ai toujours essayé de m'arrêter à cet endroit - quelque chose dont tout ce cours traite n'est pas lié à des données réelles. Seulement au problème, quand nous devons tirer une conclusion sur la population générale à partir d'un échantillon aléatoire.

Et nous avons ici le nombre 3. Trois morts, pas une sorte de vecteur, pas une table ou un échantillon. C'est un fait. Trois des morts sont venus complètement par hasard. Ils sont morts.

Ainsi, nous considérerons l'une des méthodes les plus simples pour déterminer le nombre de cas - par le taux de mortalité et le nombre de décès. Supposons que nous connaissions la mortalité et qu'elle soit de 1%. Dans cette situation, il serait logique et correct de considérer que le nombre de personnes récupérées sera de 297. Mais quelle est la fiabilité de ce jugement? Pouvons-nous simplement écarter le fait que nous avons trois morts, déclarant que trois ne sont pas des statistiques?

image

Une distribution binomiale négative et son prophète répondront à cette question - Wikipedia. Il y a beaucoup de lettres grecques, si vous, comme moi, en avez peur, alors je vais vous dire ce qui se passe. Cette distribution répond juste à la question du nombre de fois qu'il faut lancer le dé pour que le six tombe cinq fois. J'utilise le langage de programmation R pour les calculs, dans lequel il existe une fonction prête à l'emploi qui vous permet d'évaluer l'intervalle de confiance.

qnbinom(p=c(.025,.975),size=3, prob=0.01)

Ici, p est de 2,5% par le bas et 2,5% par le haut, entre lesquels se situe la plage souhaitée.

Le résultat est un intervalle de confiance de 60 à 717. Pas si mal! Il est probable que les trois morts ne signifient pas du tout 297 récupérés, mais seulement soixante! Mais peut-être pour sept cents. :-(

Pour les plus méfiants, qui ne croient pas à une distribution binomiale négative, je peux proposer une modélisation numérique. En général, si vous ne savez pas calculer par formules et distributions, modèle! Dans toute situation incompréhensible, modèle, Monte Carlo vous attend.

Nous écrirons la fonction random_infected, qui simule une situation de maladie et de mort.

random_infected <- function(deaths, fatality_rate)
{
  dead = 0
  all = 1
  while (dead < deaths) {
    if (runif(1) < fatality_rate) {
      all = all + 1
      dead = dead + 1
    } else
      all = all + 1
  }
  return(all)
}

Cette fonction fait ce qui suit - lance un cube «n-face» (en utilisant une distribution uniforme). Si l'un est lâché, cela augmente le nombre de morts et le nombre d'un. Et si ce n'était pas le cas, alors seulement le nombre tout. Chaque lancer de ce dé est une personne malade qui peut mourir ou récupérer. Dès que nous avons le nombre de décès spécifié par le paramètre de décès, nous nous arrêtons et signalons le nombre de fois que nous lançons le dé (le nombre total). La probabilité que l'un tombe sur notre cube imaginaire est la mortalité, dans notre cas le paramètre fatality_rate.

infected_sizes<-replicate(100000,random_infected(deaths=3,fatality_rate=0.01))

Et maintenant, calculons ce nombre 100 mille fois. J'ai un vieil ordinateur portable, donc j'hésite à attendre qu'un million soit compté.

Après cela, vous pouvez calculer la moyenne arithmétique des nombres obtenus. J'ai obtenu 301,2 - très similaire au nombre attendu de 300. Voici à quoi ressemble la distribution du nombre de rouleaux de notre cube de la mort:

library(ggplot2)
theme_set(theme_classic())

g <- ggplot(data.frame(infected_sizes=infected_sizes), aes(infected_sizes))
g + geom_density(alpha=0.8,fill="plum")

image

La voici - une distribution binomiale négative, veuillez aimer et favoriser. Sur la base de ces données, il est possible de donner des réponses approximatives aux questions - quelle est la probabilité que le nombre total de cas soit inférieur à cinquante (1,2%) ou supérieur à 1000 (0,3%).

Bien sûr, ce ne sont que des estimations. Ils sont basés sur des données qui peuvent être incorrectes. Nous ne connaissons pas la véritable mortalité du coronavirus. Mais plus le taux de mortalité est faible, plus il y a de cas de maladie par personne décédée et plus on estime l'ampleur de la pandémie.

Permettez-moi de vous rappeler que nous lançons ce dé instantanément. Pour le modèle de calcul de la mortalité, qui a été utilisé dans l'article sensationnel de Thomas Pueyo, j'ai une petite plainte. Là, nous supposons que, sur la base des 3 décès le jour X, un taux de mortalité de 1% et la connaissance que le délai moyen entre l'infection et le décès est de 17 jours, 300 personnes ont été infectées le jour X-17. Cependant, un tel calcul n'est valable que si le nombre de personnes malades est le même chaque jour. Étant donné que 17 jours n'est pas un nombre strict, il comporte également des intervalles de confiance et des erreurs. Si nous constatons une augmentation rapide du nombre de cas, alors parmi ceux qui sont décédés le jour X, nous avons un certain nombre de personnes infectées non pas il y a 17 jours, mais 16 ou 15 jours, et peut-être 10 jours. Il y en a peut-être encore plus que ceux qui ont été infectés il y a 17 jours. De cette façon,dans une situation d'augmentation rapide du nombre de cas, un tel calcul de retour peut conduire à des estimations surestimées de la prévalence de la maladie. En général, tout est compliqué.

PS Merci à Gregory Demin pour un indice sur le type de distribution.

All Articles