Creating Database only with summer months

villes_s<-list()

for(i in 1:length(villes)) {
  villes_s[[i]]<-villes[[i]][villes[[i]]$saison=="summer",]
  villes_s[[i]]$temps <- ave(villes_s[[i]]$time,villes_s[[i]]$annee,  FUN = seq_along)
  }

names(villes_s)<-cities

table of N obs per months to check to have only summer month

for (i in 1:length(villes)){
print(table(villes_s[[i]]$mois))
}
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0 
## 
##   1   2   3   4   5   6   7   8   9  10  11  12 
##   0   0   0   0   0 540 558 558 540   0   0   0

Creating Heat Wave variable

#seuil_IBMax<-c(33,35,34,34,33,33,34,35,35,34,34,31,31,34,33,34,36)
#seuil_IBMin<-c(18,21,19,19,18,18,20,24,22,18,20,24,21,19,19,19,21)

## imputation NA pour Marseille
villes_s[[6]]$tempmax[which(is.na(villes_s[[6]]$tempmax))] <- villes_s[[6]]$tempmaxmoy7j[which(is.na(villes_s[[6]]$tempmax))]
villes_s[[6]]$tempmin[which(is.na(villes_s[[6]]$tempmin))]<-mean(villes_s[[6]]$tempmin)

## imputation NA pour rouen
villes_s[[13]]$tempmax[which(is.na(villes_s[[13]]$tempmax))] <- villes_s[[13]]$tempmaxmoy7j[which(is.na(villes_s[[13]]$tempmax))]


trshld975<-c()

# for (i in 1:length(villes_s)){
#   trshld975[i]<-quantile(villes_s[[i]]$tempmax, probs=c(0.975),na.rm=TRUE)
#   villes_s[[i]]$heat_wave<-0
#   villes_s[[i]]$heat_wave[which(villes_s[[i]]$tempmax > trshld975[i])]<- 1
# }

for (i in 1:length(villes_s)){
  trshld975[i]<-quantile(villes_s[[i]]$tempmax, probs=c(0.975),na.rm=TRUE)
  for (r in 2:nrow(villes_s[[i]])){
  villes_s[[i]]$heat_wave[r]<-  if (villes_s[[i]]$tempmax[r] >= trshld975[i] & villes_s[[i]]$tempmax[r-1] >= trshld975[i]) 1 else 0
  }}

 

Number of Heat Wave for each city:

Heat wave=0 Heat wave=1
bordeaux 2174 21
clermont 2168 27
grenoble 2165 30
lehavre 2181 14
lyon 2168 27
marseille 2180 15
montpellier 2181 14
nancy 2169 26
nantes 2178 17
nice 2181 14
paris 2170 25
rennes 2177 18
rouen 2173 22
strasbourg 2172 23
toulouse 2174 21

Analysis On Non-Accidental Mortality

BOOTSTRAP METHOD

Code:

# some data management
 for (i in 1:length(villes_s)){
#   villes_s[[i]]<-villes_s[[i]][villes_s[[i]]$Dates>="2007-06-01",] ## create dataframe without missing values for SpF data
   villes_s[[i]]$time<-1:nrow(villes_s[[i]]) # Create a new variable for time 
 } 
   

    
  # start bootstrap
     
  nreps=1000 #number of bootstrap reps
  results<-matrix(NA,nrow = 15,ncol=10)
  colnames(results)<-c("city","NDE","NDE_2.5","NDE_97.5","NIE","NIE_2.5","NIE_97.5","PMM","PMM_2.5","PMM_97.5")
  
for (i in (1:length(villes_s))){
  boot.nde<- rep(NA,nreps)
  boot.nie <- rep(NA,nreps)
  boot.pmm <- rep(NA,nreps)
  for (rep in 1:nreps)
      
    {
      dat_resample <- sample(1:nrow(villes_s[[i]]),nrow(villes_s[[i]]),replace=T)
      dat_boot <- villes_s[[i]][dat_resample,]
        
      # outcome model - Poisson Regression: Mortality, HW and O3
      m.outcome<-glm(nocc_tot~heat_wave+o3+ns(time,df=round(3*length(time)/122))+Jours+Vacances+hol,data=dat_boot,family=quasipoisson)
      
      # mediator model - multinomial linear regression
      m.mediator<-lm(o3~heat_wave+no2moy+ns(time,df=round(3*length(time)/122))+Jours+Vacances+hol,data=dat_boot) 
      
      # save coefficients
      theta1<-coef(m.outcome)[2]
      
      theta2<-coef(m.outcome)[3]
     
      beta1<-coef(m.mediator)[2]
      
      # estimate CDE and NIE
      
      boot.nde[rep] <- exp(theta1) #NDE
      
      boot.nie[rep] <-exp(theta2*beta1) #NIE 
      
      boot.pmm[rep] <-boot.nde[rep]*(boot.nie[rep]-1)/(boot.nde[rep]*boot.nie[rep]-1) #NIE 
                   
      
    } #end bootstrap

  results[i,1]<-cities[[i]]
  results[i,2]<-median(boot.nde, na.rm=T)
  results[i,3]<-quantile(boot.nde, 0.025, na.rm=T)
  results[i,4]<-quantile(boot.nde, 0.975, na.rm=T)
  results[i,5]<-median(boot.nie, na.rm=T)
  results[i,6]<-quantile(boot.nie, 0.025, na.rm=T)
  results[i,7]<-quantile(boot.nie, 0.975, na.rm=T)
  results[i,8]<-median(boot.pmm, na.rm=T)
  results[i,9]<-quantile(boot.pmm, 0.025, na.rm=T)
  results[i,10]<-quantile(boot.pmm, 0.975, na.rm=T)
       }
  
  # output CDE and NIE and confidence intervals

Results:

res<-data.frame(results)

for (i in (2:10)){
res[,i]<-as.numeric(as.character(res[,i]))
}

res[,2:10]<-round(res[,2:10],digits = 4)

kable(res)%>%kable_styling()%>%
  column_spec(1, bold = T, border_right = T)%>%
  column_spec(4, border_right = T)%>%
  column_spec(7, border_right = T)
city NDE NDE_2.5 NDE_97.5 NIE NIE_2.5 NIE_97.5 PMM PMM_2.5 PMM_97.5
bordeaux 1.6423 1.4345 1.9452 1.0317 1.0057 1.0650 0.0746 0.0130 0.1597
clermont 1.3533 1.1510 1.5998 1.0216 0.9912 1.0585 0.0774 -0.0337 0.2363
grenoble 1.1653 1.0043 1.3387 1.0354 1.0092 1.0648 0.2045 0.0502 0.7252
lehavre 1.0596 0.7650 1.3645 1.0531 1.0054 1.1146 0.2805 -2.9903 4.0148
lyon 1.4171 1.1683 1.7226 1.0370 1.0189 1.0598 0.1130 0.0534 0.2189
marseille 1.2165 1.0679 1.4095 0.9974 0.9849 1.0105 -0.0142 -0.1594 0.0797
montpellier 0.9103 0.7069 1.1173 1.0021 0.9967 1.0133 -0.0124 -0.3775 0.3520
nancy 1.1946 0.9854 1.4830 1.0002 0.9743 1.0300 0.0022 -0.4441 0.5494
nantes 1.1626 0.9345 1.4308 1.0877 1.0378 1.1475 0.3809 0.1302 2.1579
nice 1.0937 0.9114 1.3107 1.0124 1.0011 1.0289 0.1079 -0.4646 0.8282
paris 1.6807 1.2694 2.3071 1.0637 1.0349 1.1103 0.1355 0.0701 0.3082
rennes 0.9503 0.7251 1.2206 1.1188 1.0384 1.2113 0.6087 -9.6659 10.5483
rouen 1.4853 1.2538 1.7620 1.0745 1.0346 1.1243 0.1886 0.0946 0.3153
strasbourg 1.4436 1.2196 1.7045 1.0287 1.0062 1.0560 0.0865 0.0178 0.1978
toulouse 1.1703 1.0047 1.3506 1.0404 1.0169 1.0675 0.2210 0.0808 0.7578 
write_xlsx(res,"C:/Users/Anna/Dropbox/Heat wave and O3/97.5 perc/mortality/nonaccidental/results_2.xlsx")

Forestplot

Natural Direct Effect

Natural Indirect Effect

Analysis On Cardiovascular Mortality

Results:

city NDE NDE_2.5 NDE_97.5 NIE NIE_2.5 NIE_97.5 PMM PMM_2.5 PMM_97.5
bordeaux 1.4349 1.1510 1.7796 1.0408 0.9950 1.0945 0.1225 -0.0161 0.3530
clermont 1.5732 1.2023 2.0007 1.0011 0.9406 1.0575 0.0031 -0.1952 0.1914
grenoble 1.2044 0.8852 1.5681 1.0662 1.0177 1.1216 0.2646 -0.9674 1.7940
lehavre 1.0982 0.5253 1.9316 1.1392 1.0372 1.2761 0.3452 -3.5196 4.0003
lyon 1.7250 1.3093 2.2538 1.0339 1.0021 1.0698 0.0742 0.0045 0.1712
marseille 1.1062 0.8454 1.4669 0.9974 0.9759 1.0142 -0.0139 -0.8427 0.4577
montpellier 0.5888 0.1782 1.0852 1.0025 0.9926 1.0202 -0.0030 -0.0893 0.0375
nancy 1.3394 0.9088 1.9159 0.9962 0.9401 1.0496 -0.0104 -0.8139 0.6846
nantes 1.1866 0.8158 1.6112 1.1023 1.0079 1.2068 0.3692 -2.5805 2.3688
nice 1.1155 0.8155 1.5073 1.0245 1.0036 1.0603 0.1270 -1.5480 1.5838
paris 1.7906 1.2907 2.4632 1.0630 1.0250 1.1166 0.1238 0.0517 0.2864
rennes 0.9656 0.5762 1.4224 1.0914 0.9511 1.2497 0.1255 -8.4460 5.9450
rouen 1.4187 1.0228 1.8597 1.0855 1.0032 1.1804 0.2286 -0.0014 0.6608
strasbourg 1.4149 1.0336 1.8623 0.9945 0.9519 1.0395 -0.0196 -0.2732 0.1993
toulouse 1.2221 0.9398 1.5563 1.0617 1.0235 1.1113 0.2486 -0.2930 1.4578

Forestplot

Natural Direct Effect

Natural Indirect Effect

Analysis On Respiratory Mortality

## Warning: glm.fit: algorithm did not converge

## Warning: glm.fit: algorithm did not converge

## Warning: glm.fit: algorithm did not converge

Results:

city NDE NDE_2.5 NDE_97.5 NIE NIE_2.5 NIE_97.5 PMM PMM_2.5 PMM_97.5
bordeaux 2.3166 1.4840 3.6128 1.1268 1.0088 1.2718 0.1859 0.0119 0.4112
clermont 1.3745 0.5454 2.7220 0.9452 0.8196 1.0908 -0.1133 -2.4135 1.8459
grenoble 1.3195 0.6285 2.3352 0.9734 0.8663 1.0987 -0.0570 -2.3743 1.3128
lehavre 1.8509 0.3982 4.2844 0.8658 0.6788 1.0489 -0.2633 -3.7908 2.7884
lyon 1.8923 1.2851 2.6105 1.0982 1.0316 1.1876 0.1757 0.0593 0.3457
marseille 1.2680 0.8049 1.8163 0.9994 0.9736 1.0165 -0.0013 -0.2923 0.2396
montpellier 0.6348 0.0000 1.8056 1.0150 0.9856 1.0681 -0.0038 -0.6232 0.5420
nancy 1.1421 0.5403 2.3685 0.9711 0.8631 1.0924 -0.0069 -2.0232 1.9896
nantes 1.1217 0.5515 1.9655 1.3395 1.1077 1.6891 0.6823 -3.7880 8.9663
nice 1.0004 0.4996 1.7214 1.0257 0.9806 1.0844 0.0042 -1.2248 1.2727
paris 2.2047 1.6029 2.9815 1.0883 1.0432 1.1448 0.1392 0.0693 0.2422
rennes 2.0912 0.6487 4.5626 1.1604 0.8876 1.5504 0.2323 -0.4553 1.1497
rouen 2.2582 1.2678 3.7399 1.0068 0.8331 1.1795 0.0115 -0.4331 0.3532
strasbourg 2.2008 1.2799 3.6728 1.0239 0.9343 1.1290 0.0442 -0.1402 0.2661
toulouse 1.5673 1.0141 2.3270 1.0152 0.9218 1.1123 0.0400 -0.3443 0.5420

Forestplot

Natural Direct Effect

Natural Indirect Effect