Simulations de la décroissance de l’incidence

Données utilisées

Le temps de génération du SARS-CoV-2 est d’environ 5.3.

L’incidence selon SPF au 28 octobre : 399/100000

Le 1er confinement avait duré du 17 mars au 11 mai, soit : 14+30+11=55 jours.

Par ailleurs, 56 jours (1+30+25) séparent le 30 octobre, date du début du second confinement, et le 25 décembre.

Dans ce qui suit nous nous demandons quelle serait l’efficacité d’un confinement pendant 56 jours, selon 2 différents \(R_t\) : un \(R_t\) à 0.7, cohérent à ce qui a été estimé pendant le premier confinement, et un \(R_t\) à 0.9.

Fonctions de calcul de l’évolution du nombre de personnes infectieuses

## Nouvelles infections par jour 
add_n <- function (rt) {
  return ( rt/5.3 )
}

## Personnes qui cessent d'être infectieuses par jour 
remove_n <- function () {
  return ( -1/5.3 )
}


## Changement dans le nombre d'individus infectés par jour.
dn <- function (rt) {
  return ( add_n(rt) + remove_n() )
}

## Simulation du nombre de nouveaux individus infectés, sur num_days jours
nouv_infec <- function (rt, num_days, initial_num_infec) {
  current_num_infec = initial_num_infec
  nouveaux_infec = c()
  for (i in 1:num_days) {
    nouveaux_infec = c(nouveaux_infec, current_num_infec*add_n(rt))
    current_num_infec = current_num_infec * (1 + dn(rt) )
  }
  return (nouveaux_infec)
}

## Simulation du nombre d'individus infectés, sur num_days jours
num_infec <- function (rt, num_days, initial_num_infec) {
  current_num_infec = initial_num_infec
    infecs=c(current_num_infec)
  for (i in 1:num_days) {
    current_num_infec = current_num_infec * (1 + dn(rt) )
    infecs = c(infecs, current_num_infec)
  }
  return (infecs)
}

Simulation de la décroissance du taux d’incidence avec \(R_t\) = 0.7

Un Rt à 0.7 est ce qui a été estimé par plusieurs études, par exemple : Duchemin et col., Sofonea et col., Salje et col..

num_infecs = num_infec(0.7, 56, 399)

plot(1:length(num_infecs), num_infecs, main="Décroissance avec Rt=0.7", ylab="Incidence", xlab="Jours", t="l", ylim=c(0, 500))
abline(v=57, col="red", lwd=2)
abline(h=num_infecs[57], col="red", lwd=2, lty=2)
legend("topleft", col="red", lwd=2, legend="Jour de Noël")

Avec un \(R_t\) à 0.7, le confinement réussirait à diviser le taux d’incidence par :

399 / num_infecs[57]
## [1] 26.12934

soit environ 26.

J’évalue ci-dessous cette diminution du taux d’incidence par un facteur 26 au regard de ce qui a été observé lors du premier confinement (voir ci-dessous pour le raisonnement).

Quand le seuil des 5000 nouvelles infections par jour serait-il atteint, si le second confinement est aussi efficace que le premier ?

Il y a 67 081 000 habitants en France. On peut utiliser les formules ci-dessous pour étudier le nombre de nouvelles infections par jour, et évaluer quand ce nombre passera en-dessous de 5000, le seuil qui a été choisi par le gouvernement pour arrêter le confinement.

num_infecs = nouv_infec(0.7, 56, 399*67081000/100000)

plot(1:length(num_infecs), num_infecs, main="Décroissance avec Rt=0.7", ylab="Incidence", xlab="Jours", t="l", ylim=c(0, 50000))
abline(v=57, col="red", lwd=2)
abline(h=5000, col="red", lwd=2, lty=2)
legend("topleft", col="red", lwd=2, legend="Jour de Noël")

Avec un Rt similaire à celui qu’on a eu lors du premier confinement, la levée du confinement, selon le critère des 5000 nouvelles infections par jour, tomberait après environ 35 jours, soit début décembre.

Simulation de la décroissance du taux d’incidence avec \(R_t\) = 0.9

Le second confinement est moins contraignant que le premier confinement. Par exemple, les écoles restent ouvertes, et les citoyens sont invités à continuer le travail, y compris en présentiel quand le travail à distance n’est pas possible. Pour cette raison, on peut s’attendre à ce que le second confinement conduise à un \(R_t\) moins faible que le premier. Pour cette raison nous testons un \(R_t\)=0.9. Pour rappel, si le \(R_t\) était supérieur à 1.0, alors le taux d’incidence ne diminuerait pas mais continuerait de monter.

num_infecs = num_infec(0.9, 56, 399)
plot(1:length(num_infecs), num_infecs, main="Décroissance avec Rt=0.9", ylab="Incidence", xlab="Jours", t="l", ylim=c(0, 500))
abline(v=57, col="red", lwd=2)
abline(h=num_infecs[57], col="red", lwd=2, lty=2)
legend("topleft", col="red", lwd=2, legend="Jour de Noël")

Avec un \(R_t\) à 0.9, le confinement réussirait à diviser le taux d’incidence par :

399 / num_infecs[57]
## [1] 2.905771

Seulement environ 3.

Quand le seuil des 5000 nouvelles infections par jour serait-il atteint, si le second confinement conduit à un Rt de 0.9 ?

num_infecs = nouv_infec(0.9, 56, 399*67081000/100000)

plot(1:length(num_infecs), num_infecs, main="Décroissance avec Rt=0.9", ylab="Incidence", xlab="Jours", t="l", ylim=c(0, 50000))
abline(v=57, col="red", lwd=2)
abline(h=5000, col="red", lwd=2, lty=2)
legend("topleft", col="red", lwd=2, legend="Jour de Noël")

Si le Rt est à 0.9, le modèle prédit que nous serons encore loin du seuil de levée du confinement au moment de Noël.

Date de levée du confinement en fonction du Rt

On a vu que la date de levée du confinement est très sensible à la valeur du Rt. Observons ses valeurs lorsque Rt prend des valeurs entre 0.7 et 0.9.

## Simulation du nombre de nouveaux individus infectés, sur num_days jours
jour_seuil <- function (rt, num_days, initial_num_infec) {
  current_num_infec = initial_num_infec
  for (i in 1:num_days) {
    nouveaux_infec = current_num_infec*add_n(rt)
    if (nouveaux_infec <= 5000) {
      return (i)
    }
    current_num_infec = current_num_infec * (1 + dn(rt) )
  }
  return (num_days)
}


jours_seuil = c()
for (i in 70:90) {
  jours_seuil = c(jours_seuil, jour_seuil(i/100, 150, 399*67081000/100000))
}
par(mar=c(4, 5, 4, 2) + 0.1)
plot(70:90/100, jours_seuil, ylab="Jour d'atteinte du seuil de\n5000 nouvelles infections", xlab="Valeur de Rt pendant le second confinement")
abline(h=56, col="red", lwd=2)
legend("topleft", col="red", lwd=2, legend="Jour de Noël")

Si le Rt est inférieur à 0.8 lors de ce second confinement, alors le seuil des 5000 nouvelles infections par jour pourrait être atteint avant Noël. S’il est de 0.9, alors il faudrait 120 jours.

Réalisme du modèle de simulation

Dans cette section on cherche à évaluer le réalisme du modèle de simulation. Pour ce faire, on compare la décroissance du taux d’incidence simulée avec ce modèle à la décroissance estimée sur les données, réelles, lors du premier confinement. On a vu que la simulation avec le modèle le plus simple, avec un \(R_t\) de 0.7 sur 56 jours donnait une diminution par 26. Sur 55 jours, la durée du premier confinement, la diminution serait de :

num_infecs = num_infec(0.7, 55, 399)

plot(1:length(num_infecs), num_infecs, main="Simulation du 1er confinement avec Rt=0.7", ylab="Incidence", xlab="Jours", t="l", ylim=c(0, 500))

Avec un \(R_t\) à 0.7, selon le modèle, le premier confinement aurait réussi à diviser le taux d’incidence par :

399 / num_infecs[56]
## [1] 24.65032

Environ 25.

Afin d’évaluer le réalisme du modèle de simulation, nous voulons comparer ce chiffre de 25 au chiffre observé empiriquement lors du premier confinement. Malheureusement, en mars, pour le premier confinement, SPF ne fournit pas de données d’incidences. A la place, nous nous basons sur les données de mortalité. Le raisonnement est le suivant : la mortalité à un temps t est grosso modo proportionnelle à l’incidence 23 jours plus tôt (la médiane de la distribution utilisée par Flaxman et al. est de 23 jours).

d<-read.table("donnees-hospitalieres-nouveaux-covid19-2020-10-31-19h00.csv", h=T, sep=";")
dd = aggregate(. ~ jour, data = d, sum)
summary(dd)
##          jour          dep         incid_hosp       incid_rea    
##  2020-03-19:  1   Min.   :5151   Min.   :  14.0   Min.   :  1.0  
##  2020-03-20:  1   1st Qu.:5151   1st Qu.: 134.0   1st Qu.: 17.0  
##  2020-03-21:  1   Median :5151   Median : 317.0   Median : 47.0  
##  2020-03-22:  1   Mean   :5151   Mean   : 775.9   Mean   :125.1  
##  2020-03-23:  1   3rd Qu.:5151   3rd Qu.: 993.0   3rd Qu.:154.0  
##  2020-03-24:  1   Max.   :5151   Max.   :4281.0   Max.   :771.0  
##  (Other)   :221                                                  
##     incid_dc       incid_rad     
##  Min.   :  0.0   Min.   :   5.0  
##  1st Qu.: 17.5   1st Qu.: 197.0  
##  Median : 47.0   Median : 345.0  
##  Mean   :112.9   Mean   : 540.3  
##  3rd Qu.:143.0   3rd Qu.: 748.0  
##  Max.   :613.0   Max.   :2224.0  
## 
plot(dd$incid_dc, t="l", ylab="Mortalité en France", xlab="Jours depuis le 19 mars 2020")
abline(h=35, col="green", lwd=2)
legend("topright", col="green", lwd=2, legend="Mortalité correspondant\nau seuil de 5000 infections")

23 jours après le début du confinement, le 9 avril, combien de décès y avait-il ? Nous faisons une somme sur 5 jours pour diminuer l’effet de la stochasticité hebdomadaire des rapports de morts.

initial_deaths = dd$incid_dc[which(dd$jour=="2020-04-07")] + dd$incid_dc[which(dd$jour=="2020-04-08")] + dd$incid_dc[which(dd$jour=="2020-04-09")] + dd$incid_dc[which(dd$jour=="2020-04-10")] + dd$incid_dc[which(dd$jour=="2020-04-11")]  
initial_deaths
## [1] 2524

23 jours après la fin du confinement, le 3 juin, combien de décès y avait-il ?

final_deaths = dd$incid_dc[which(dd$jour=="2020-06-01")] + dd$incid_dc[which(dd$jour=="2020-06-02")] + dd$incid_dc[which(dd$jour=="2020-06-03")] + dd$incid_dc[which(dd$jour=="2020-06-04")] + dd$incid_dc[which(dd$jour=="2020-06-05")]  
final_deaths
## [1] 294

En considérant qu’il faut environ trois semaines entre l’infection et le décès, et que le pourcentage des personnes infectées en mars qui mouraient est le même que celui des personnes infectées en mai qui mouraient, on conclut que le premier confinement, avec un \(R_t\) à environ 0.7 a divisé l’incidence par :

initial_deaths / final_deaths 
## [1] 8.585034

soit environ 9.

Le modèle, qui prédisait 25, semble donc peu réaliste : il prédit une décroissance beaucoup plus vive qu’observée.

Salje et al. estiment que 0.7% des personnes infectées meurent. Donc on peut estimer combien de personnes ont été infectées sur la base du nombre de personnes qui meurent 23 jours plus tard. Donc, à la fin du 1er confinement, il y avait :

final_deaths/0.007/5
## [1] 8400

8400 nouvelles infections par jour.

Quand a-t-on atteint le seuil de 5000 infections par jour ? On regarde la situation tout au long du 1er confinement, du 17 mars au 11 mai. 23 jours après le 17 mars, on était le 9 avril, le jour 22 dans la table dd ; 23 jours après la fin, on était au jour 77.

infecs = c()
for (i in 22:77) {
  infecs=c(infecs, sum(dd$incid_dc[(i-2):(i+2)])/0.007/5)
}
plot((22:77)-22, infecs, ylab="Infections par jour", xlab="Jours de confinement")
abline(h=5000, col="red", lwd=2)

Jamais le nombre de nouvelles infections par jour n’était < 5000.