Heat wave and Ozone, Mediation analysis on mortality (95 perc threshold)

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))]


trshld95<-c()

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

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

Number of Heat Wave for each city:

	Heat wave=0	Heat wave=1
bordeaux	2154	41
clermont	2140	55
grenoble	2134	61
lehavre	2159	36
lyon	2133	62
marseille	2151	44
montpellier	2153	42
nancy	2142	53
nantes	2150	45
nice	2154	41
paris	2147	48
rennes	2153	42
rouen	2150	45
strasbourg	2140	55
toulouse	2157	38

Analysis On Non-Accidental Mortality

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.3704	1.1933	1.5554	1.0280	1.0028	1.0557	0.0938	0.0109	0.2106
clermont	1.1949	1.0348	1.3764	1.0213	0.9916	1.0506	0.1176	-0.0556	0.5118
grenoble	1.0730	0.9660	1.1760	1.0369	1.0104	1.0661	0.3394	-0.0131	2.2385
lehavre	1.2024	1.0026	1.4557	1.0344	0.9863	1.0834	0.1717	-0.0845	0.7863
lyon	1.2077	1.0666	1.3470	1.0361	1.0193	1.0535	0.1780	0.0886	0.3802
marseille	1.0777	1.0001	1.1671	1.0046	0.9989	1.0135	0.0587	-0.0436	0.3187
montpellier	1.0614	0.9197	1.2028	1.0037	0.9960	1.0143	0.0366	-0.6995	0.6776
nancy	1.0869	0.9129	1.2703	1.0012	0.9723	1.0342	0.0087	-1.2408	1.9490
nantes	1.2276	1.0499	1.4013	1.0624	1.0246	1.1064	0.2489	0.0905	0.6078
nice	1.0813	0.9644	1.2067	1.0082	1.0000	1.0209	0.0870	-0.5860	0.8907
paris	1.5114	1.2564	1.7863	1.0527	1.0340	1.0816	0.1384	0.0929	0.2041
rennes	1.1014	0.9071	1.3224	1.0838	1.0200	1.1532	0.4550	-0.6685	4.0600
rouen	1.3379	1.1876	1.5080	1.0620	1.0228	1.1053	0.1978	0.0736	0.3443
strasbourg	1.3257	1.2078	1.4702	1.0233	1.0007	1.0476	0.0886	0.0022	0.1857
toulouse	1.1164	0.9908	1.2611	1.0407	1.0180	1.0658	0.2827	0.0959	0.9706

# write_xlsx(res,"C:/Users/Anna/Dropbox/Heat wave and O3/95 perc/mortality/nonaccidental/results.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.1887	0.9647	1.4202	1.0360	1.0015	1.0801	0.1803	-0.0933	0.9203
clermont	1.2518	1.0120	1.5720	1.0015	0.9506	1.0590	0.0067	-0.4353	0.4264
grenoble	1.1785	0.9622	1.4192	1.0610	1.0131	1.1100	0.2770	0.0231	1.3567
lehavre	1.0062	0.6822	1.3913	1.1326	1.0345	1.2471	0.5320	-6.0701	7.4871
lyon	1.3089	1.1064	1.5434	1.0343	1.0045	1.0667	0.1323	0.0190	0.3264
marseille	1.0493	0.9154	1.1992	1.0057	0.9987	1.0213	0.0628	-1.2717	1.2252
montpellier	0.8091	0.5776	1.0926	1.0056	0.9912	1.0286	-0.0198	-0.5020	0.2949
nancy	1.0306	0.7875	1.3409	1.0035	0.9505	1.0594	0.0074	-2.9880	2.7303
nantes	1.2720	1.0060	1.6023	1.0713	0.9946	1.1530	0.2548	-0.0343	0.8758
nice	1.0404	0.8608	1.2508	1.0146	1.0004	1.0385	0.1068	-1.7605	1.8335
paris	1.5620	1.2794	1.9182	1.0517	1.0249	1.0915	0.1286	0.0716	0.2033
rennes	1.2784	0.8409	1.8359	1.0423	0.9281	1.1680	0.1473	-1.1432	1.4096
rouen	1.1532	0.8815	1.4220	1.0871	1.0122	1.1788	0.3828	-2.0174	2.4647
strasbourg	1.3254	1.1010	1.5754	0.9887	0.9411	1.0311	-0.0476	-0.3458	0.1734
toulouse	1.1229	0.9248	1.3339	1.0630	1.0222	1.1097	0.3495	-0.6798	2.3967

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	1.7288	1.1921	2.4353	1.1065	1.0119	1.2207	0.2065	0.0246	0.4931
clermont	0.9481	0.4493	1.6214	0.9592	0.8473	1.0834	0.0532	-3.8926	3.4337
grenoble	1.5048	0.9724	2.2577	0.9624	0.8578	1.0852	-0.1152	-1.2745	1.1527
lehavre	1.3573	0.6126	2.4887	0.8804	0.7031	1.0607	-0.2788	-4.7697	4.8184
lyon	1.8258	1.4269	2.3207	1.0739	1.0077	1.1500	0.1408	0.0158	0.2851
marseille	1.2179	0.8832	1.5909	1.0026	0.9890	1.0245	0.0114	-0.1905	0.2441
montpellier	1.0698	0.5423	1.7642	1.0264	0.9919	1.0886	0.0508	-1.2063	1.2224
nancy	0.9223	0.5271	1.5150	0.9771	0.8676	1.1017	0.0435	-4.1058	2.8290
nantes	1.1747	0.7149	1.7574	1.2762	1.0952	1.5225	0.6245	-1.0326	4.6067
nice	1.1478	0.7602	1.6254	1.0140	0.9882	1.0517	0.0494	-1.0601	0.9653
paris	1.6928	1.3631	2.0995	1.0798	1.0430	1.1283	0.1657	0.0882	0.2801
rennes	1.8629	0.8844	3.2681	1.1058	0.8711	1.4344	0.1848	-0.4361	1.1337
rouen	1.6465	1.0182	2.4145	0.9967	0.8415	1.1641	-0.0074	-0.6162	0.6609
strasbourg	1.5605	1.0529	2.2303	1.0212	0.9232	1.1355	0.0596	-0.2905	0.5404
toulouse	1.4791	1.0241	2.0780	1.0082	0.9165	1.1062	0.0258	-0.4622	0.4503

Forestplot

Natural Direct Effect

Natural Indirect Effect

Heat wave and Ozone, Mediation analysis on mortality (95 perc threshold)

Anna Alari

04/11/2020

Creating Database only with summer months

Creating Heat Wave variable

Analysis On Non-Accidental Mortality

Forestplot

Analysis On Cardiovascular Mortality

Forestplot

Analysis On Respiratory Mortality

Forestplot