diff --git a/Codes/Conditional_incidence_and_reconstruction.R b/Codes/Conditional_incidence_and_reconstruction.R
index f3f0f9dec96cdcf935bba458c8ad61b733439c3f..c2699e25e7268f26e8931daf91e13d41940c4f29 100644
--- a/Codes/Conditional_incidence_and_reconstruction.R
+++ b/Codes/Conditional_incidence_and_reconstruction.R
@@ -118,7 +118,7 @@ Pars_cond_prev_quantiles<- c(variable_x,variable_y,variable_z,"_Quantiles.csv",b
Conditional_incidence_quantiles_info<-Conditional_incidence_Quantiles(Env_laboratory_data,Env_Pathogen_data,Pars_cond_prev_quantiles)
## plot Conditional incidence for bins divided in quantiles
-Pars_plot_cond_prev_quantiles <- c(variable_x,variable_y,variable_z)
+Pars_plot_cond_prev_quantiles <- c(variable_x,variable_y,variable_z,"_quantile")
source("Function_Plot_Conditional_Incidence_Quantiles.R", chdir=TRUE)
Plot_Conditional_incidence_quantiles(Conditional_incidence_quantiles_info,Pars_plot_cond_prev_quantiles)
diff --git a/Codes/Conditional_incidence_and_reconstruction_for_constant_situation.R b/Codes/Conditional_incidence_and_reconstruction_for_constant_situation.R
new file mode 100644
index 0000000000000000000000000000000000000000..5e77509021321ec869bb768bfcfb263cbec0bcd1
--- /dev/null
+++ b/Codes/Conditional_incidence_and_reconstruction_for_constant_situation.R
@@ -0,0 +1,223 @@
+# The code estimate the incidence of the diseases conditional to local environmental variables. It perform a series of diagnostic and
+# then reconstruct the time series of cases given the local environmental factors in England and Wales
+
+### Require Relevant R Packages ###########
+rm(list=ls(all=TRUE))
+
+library(data.table)
+library(readr)
+library(dplyr)
+library(lubridate)
+library(plyr)
+library(dplR)
+library(MALDIquant)
+library(ggplot2)
+library(slider)
+library(WaveletComp)
+library(gdata)
+library(wesanderson)
+library(zoo)
+
+# Global Parameters
+
+# Set length of the desired time-lag
+
+time_lag <<-14
+time_lag_char <<-paste(time_lag)
+
+# Set the years you are focusing on
+
+year_first<<-1990
+year_last<<-2009
+first_day<<-as.Date(paste(as.character(year_first),"-01-01",sep=""))
+last_day<<-as.Date(paste(as.character(year_last),"-12-31",sep=""))
+dates_s<<-seq(first_day,last_day,by="1 day")
+
+## Gather Information on laboratory Postcodes
+
+variable_df_1<-read.csv("../Data_Base/Laboratories_Post_Codes.csv")
+
+All_PC_s<-variable_df_1[,1]
+dates<-rep(dates_s,times=length(All_PC_s))
+All_PC<<-rep(All_PC_s,each=length(dates_s))
+
+
+# Selection of the variables of interest.
+#Change HERE depending on the variable object of study, by copy-pasting from the list below
+
+variable_x<-"Maximum_air_temperature"
+variable_y<-"Relative_humidity"
+#variable_y<-"Mean_Precipitation"
+#variable_x<-"Mean_wind_speed"
+variable_z<-"daylength"
+#variable_z<-"Relative_humidity"
+
+### Set the size of the bins for Uniform conditional incidence and resolution in rounding the values of the bins.
+### For example:
+## if the size of the bin for humidity is 5 and we want the interval of he bins like ...70,75,80... we choose delta_variable<-5, round<-0
+## if the size of the bin for precipitation is 0.1 and we want the interval of he bins like ...0,0.1.0.2 ..., we choose delta_variable<-0.1, round_x<-1
+## to distinguishes one decimal point
+
+delta_variable_x<-1
+round_x<-0
+delta_variable_y<-3
+round_y<-0
+delta_variable_z<-0.5
+round_z<-1
+## partitions of the quantiles
+by_z<-0.25
+by_y<-0.25
+by_x<-0.1
+
+# Potential Variables to use:
+#"Maximum_air_temperature",
+#"Minimum_air_temperature",
+#"Mean_wind_speed",
+#"Mean_Precipitation"
+#"Relative_humidity"
+#"daylength"
+
+
+
+
+#const_variable<-"const_temp_20_hum_76"
+#const_variable<-"const_hum_76_daylength_15"
+const_variable<-"const_temp_20_daylength_15"
+
+
+if (const_variable=="const_hum_76_daylength_15"){
+ Relative_humidity_value<-76
+ daylength_value <-15
+}
+if (const_variable=="const_temp_20_hum_76"){
+ Maximum_air_temperature_value <-20
+ Relative_humidity_value<-76
+}
+if (const_variable=="const_temp_20_daylength_15"){
+ Maximum_air_temperature_value <-20
+ daylength_value <-15
+}
+
+
+
+
+######## Input of data ####
+source("Function_Input_files.R", chdir=TRUE)
+
+
+######## Perform Wavelet analysis for the selected 17 years to check for periodicity in the data ####
+###### helpful but not necessary. So this is commented
+source("Function_Wavelet_analysis.R", chdir=TRUE)
+
+Pars_wavelet_analysis<-c("../Graphs/Power_Spectrum.tiff","../Graphs/Power_Average.tiff","No")
+Env_Pathogen_data$Cases_computed<-Env_Pathogen_data$Cases
+campylobacter_data_national<-wavelet_analysis(Env_Pathogen_data, Pars_wavelet_analysis)
+
+
+###### Calculate Conditional incidence Uniform ######
+# The code does look at how the risk of the disease in humans depends on environmental variables
+# Analysis was done following an uniform division of the range of the environmental variables, independently of the number of observations
+# This is used for the reconstruction
+
+#
+source("Function_Conditional_Incidence_Uniform.R", chdir=TRUE)
+#Conditional incidence, Uniform bins
+Pars_cond_incidence <- c(variable_x,variable_y,variable_z,"_Uniform.csv",delta_variable_x,round_x,delta_variable_y,round_y,delta_variable_z,round_z)
+Conditional_incidence_unif_info<-Conditional_incidence_Uniform(Env_laboratory_data,Env_Pathogen_data,Pars_cond_incidence)
+
+## plot Conditional incidence for uniform, for diagnostic only
+#Pars_plot_cond_prev_uniform <- c(variable_x,variable_y,variable_z,time_lag_char)
+#source("Function_Plot_Conditional_incidence_Uniform.R", chdir=TRUE)
+#Plot_Conditional_incidence_Uniform(Conditional_incidence_unif_info,Pars_plot_cond_prev_uniform)
+
+#
+#source("Function_Conditional_Incidence_Quantiles.R", chdir=TRUE)
+#Conditional incidence, bins divided in quantiles
+#Pars_cond_prev_quantiles<- c(variable_x,variable_y,variable_z,"_Quantiles.csv",by_x,by_y,by_z)
+#Conditional_incidence_quantiles_info<-Conditional_incidence_Quantiles(Env_laboratory_data,Env_Pathogen_data,Pars_cond_prev_quantiles)
+
+if (const_variable=="const_hum_76_daylength_15"){
+
+
+
+
+ ###########
+ wt<-which(abs((Conditional_incidence_unif_info[[1]]$daylength) - daylength_value) ==
+ min(abs((Conditional_incidence_unif_info[[1]]$daylength) - daylength_value)))
+ Conditional_incidence_unif_info_subset1<-Conditional_incidence_unif_info[[1]][wt,]
+
+ ww<-which(abs((Conditional_incidence_unif_info_subset1$Relative_humidity) - Relative_humidity_value) ==
+ min(abs((Conditional_incidence_unif_info_subset1$Relative_humidity) - Relative_humidity_value)))
+ Conditional_incidence_unif_info_subset<-Conditional_incidence_unif_info_subset1[ww,]
+
+ ## plot Conditional incidence for bins divided in unif
+ Pars_plot_cond_prev_unif <- c(variable_x,paste("_",const_variable))
+ source("Function_Plot_Conditional_Incidence_Uniform.R", chdir=TRUE)
+Plot_Conditional_incidence_unif_constant(Conditional_incidence_unif_info_subset,Pars_plot_cond_prev_unif)
+
+}
+
+if (const_variable=="const_temp_20_hum_76"){
+
+ wt<-which(abs((Conditional_incidence_unif_info[[1]]$Maximum_air_temperature) - Maximum_air_temperature_value) == min(abs((Conditional_incidence_unif_info[[1]]$Maximum_air_temperature) - Maximum_air_temperature_value)))
+ Conditional_incidence_unif_info_subset1<-Conditional_incidence_unif_info[[1]][wt,]
+
+ ww<-which(abs((Conditional_incidence_unif_info_subset1$Relative_humidity) - Relative_humidity_value) ==
+ min(abs((Conditional_incidence_unif_info_subset1$Relative_humidity) - Relative_humidity_value)))
+ Conditional_incidence_unif_info_subset<-Conditional_incidence_unif_info_subset1[ww,]
+
+ ## plot Conditional incidence for bins divided in quantiles
+ Pars_plot_cond_prev_unif <- c(variable_z,paste("_",const_variable))
+ source("Function_Plot_Conditional_Incidence_Uniform.R", chdir=TRUE)
+ Plot_Conditional_incidence_unif_constant(Conditional_incidence_unif_info_subset,Pars_plot_cond_prev_unif)
+}
+if (const_variable=="const_temp_20_daylength_15"){
+ wt<-which(abs((Conditional_incidence_unif_info[[1]]$Maximum_air_temperature) - Maximum_air_temperature_value) == min(abs((Conditional_incidence_unif_info[[1]]$Maximum_air_temperature) - Maximum_air_temperature_value)))
+ Conditional_incidence_unif_info_subset1<-Conditional_incidence_unif_info[[1]][wt,]
+
+ ww<-which(abs((Conditional_incidence_unif_info_subset1$daylength) - daylength_value) == min(abs((Conditional_incidence_unif_info_subset1$daylength) - daylength_value)))
+ Conditional_incidence_unif_info_subset<-Conditional_incidence_unif_info_subset1[ww,]
+
+ ## plot Conditional incidence for bins divided in quantiles
+ Pars_plot_cond_prev_unif <- c(variable_y,paste("_",const_variable))
+ source("Function_Plot_Conditional_Incidence_Uniform.R", chdir=TRUE)
+ Plot_Conditional_incidence_unif_constant(Conditional_incidence_unif_info_subset,Pars_plot_cond_prev_unif)
+
+}
+
+
+
+
+
+
+###### Reconstruction ######
+
+
+
+if (const_variable=="const_hum_76_daylength_15"){
+ Env_laboratory_data$Relative_humidity<-Relative_humidity_value
+ Env_laboratory_data$daylength <-daylength_value
+ }
+if (const_variable=="const_temp_20_hum_76"){
+ Env_laboratory_data$Maximum_air_temperature <-Maximum_air_temperature_value
+ Env_laboratory_data$Relative_humidity<-Relative_humidity_value
+}
+if (const_variable=="const_temp_20_daylength_15"){
+ Env_laboratory_data$Maximum_air_temperature <-Maximum_air_temperature_value
+ Env_laboratory_data$daylength <-daylength_value
+}
+
+
+source("Function_Reconstruction.R")
+Pars_reconstruction <- c(variable_x,variable_y,variable_z,"_reconstructed.csv")
+Reconstruction_time_series<-Reconstruction_time_series(Conditional_incidence_unif_info,Env_laboratory_data,Pars_reconstruction)
+
+############## End reconstructions ####################
+
+###### Plot Reconstruction ######
+source("Function_Plot_Reconstruction.R")
+Pars_Plot_Reconstruction <- c(paste("../Graphs/",const_variable,"_",time_lag_char,"_time_series.tiff",sep=""),paste("../Graphs/",const_variable,"_",time_lag_char,"_yearly_average.tiff",sep=""),"Campylobacteriosis")
+Plot_Reconstruction_time_series(Reconstruction_time_series,campylobacter_data_national,Pars_Plot_Reconstruction)
+#
+
+############## End Plot Reconstruction ####################
diff --git a/Codes/Conditional_incidence_and_reconstruction_two_factors.R b/Codes/Conditional_incidence_and_reconstruction_two_factors.R
index 17050175ba69b17a850160542ed47d546ebd9919..127e7fb508f3218f9ff5ed98a3d9496069ec4a96 100644
--- a/Codes/Conditional_incidence_and_reconstruction_two_factors.R
+++ b/Codes/Conditional_incidence_and_reconstruction_two_factors.R
@@ -45,7 +45,7 @@ All_PC<<-rep(All_PC_s,each=length(dates_s))
# Selection of the variables of interest.
#Change HERE depending on the variable object of study, by copy-pasting from the list below
-variable_x<-"Difference_temperature"
+#variable_x<-"Difference_temperature"
#variable_x<-"Minimum_air_temperature"
#variable_y<-"Mean_Precipitation"
variable_y<-"Relative_humidity"
@@ -53,7 +53,7 @@ variable_y<-"Relative_humidity"
#variable_y<-"daylength"
variable_z<-"daylength"
-#variable_x<-"Maximum_air_temperature"
+variable_x<-"Maximum_air_temperature"
### Set the size of the bins for Uniform conditional incidence and resolution in rounding the values of the bins.
### For example:
## if the size of the bin for humidity is 5 and we want the interval of he bins like ...70,75,80... we choose delta_variable<-5, round<-0
@@ -91,7 +91,7 @@ Env_laboratory_data$Difference_temperature<-Env_laboratory_data$Maximum_air_temp
###### helpful but not necessary. So this is commented
source("Function_Wavelet_analysis.R", chdir=TRUE)
-Pars_wavelet_analysis<-c("../Graphs/Power_Spectrum.tiff","../Graphs/Power_Average.tiff")
+Pars_wavelet_analysis<-c("../Graphs/Power_Spectrum.tiff","../Graphs/Power_Average.tiff","No")
Env_Pathogen_data$Cases_computed<-Env_Pathogen_data$Cases
campylobacter_data_national<-wavelet_analysis(Env_Pathogen_data, Pars_wavelet_analysis)
diff --git a/Codes/Function_Conditional_Incidence_Uniform.R b/Codes/Function_Conditional_Incidence_Uniform.R
index bc5240adb409c2b60dc8ca09194a955a754b56d3..f645ee2a1aff5fa605817268decb1014038e328c 100644
--- a/Codes/Function_Conditional_Incidence_Uniform.R
+++ b/Codes/Function_Conditional_Incidence_Uniform.R
@@ -2,7 +2,7 @@
# Analysis was done following an uniform division of the range of the environmental variables, independently of the number of observations
# This is used for the reconstruction
-Conditional_prevalence_Uniform<-function(Data_frame_1, Data_frame_2,Pars)
+Conditional_incidence_Uniform<-function(Data_frame_1, Data_frame_2,Pars)
{
with(as.list(c(Data_frame_1, Data_frame_2,Pars)), {
@@ -47,8 +47,8 @@ Env_Pathogen_data_df <- Data_frame_2[,c("PostCode","Date", variable_x,variable_y
#################### Create a data frame to fill with
-Conditional_prevalence_unif <-data.frame(character(), character(), character(), numeric(), numeric())
-colnames(Conditional_prevalence_unif) <-c(variable_x, variable_y, variable_z,"counts","residents_tot")
+Conditional_incidence_unif <-data.frame(character(), character(), character(), numeric(), numeric())
+colnames(Conditional_incidence_unif) <-c(variable_x, variable_y, variable_z,"counts","residents_tot")
@@ -114,7 +114,7 @@ for (j_z in c(j_z_min:j_z_max)){
)
colnames(data_df) <-c(variable_x, variable_y, variable_z,"counts","residents_tot") #
- Conditional_prevalence_unif <-rbind(Conditional_prevalence_unif, data_df)
+ Conditional_incidence_unif <-rbind(Conditional_incidence_unif, data_df)
}
}
}
@@ -124,8 +124,117 @@ for (j_z in c(j_z_min:j_z_max)){
}
-write.table(Conditional_prevalence_unif, paste("../Data_Base/Conditional_prevalence_",variable_x,"_",variable_y,"_",variable_z,"_",time_lag_char,info,sep=""), row.names = FALSE , sep = ",",eol = "\n")
-return(list(Conditional_prevalence_unif,breaks_x,breaks_y,breaks_z))
+write.table(Conditional_incidence_unif, paste("../Data_Base/Conditional_incidence_",variable_x,"_",variable_y,"_",variable_z,"_",time_lag_char,info,sep=""), row.names = FALSE , sep = ",",eol = "\n")
+return(list(Conditional_incidence_unif,breaks_x,breaks_y,breaks_z))
})
}
+
+
+
+
+Conditional_incidence_Uniform_two_factors<-function(Data_frame_1, Data_frame_2,Pars)
+{
+
+ with(as.list(c(Data_frame_1, Data_frame_2,Pars)), {
+
+ # Reshape the residents information to make a data frame with the postcode, year and number of residents
+
+
+ # Pars is the list of 3 variables of interest
+
+
+
+ variable_x<-Pars[1]
+ variable_y<-Pars[2]
+ info<-Pars[3]
+ delta_x<-as.numeric(Pars[4])
+ round_digit_x<-as.numeric(Pars[5])
+ delta_y<-as.numeric(Pars[6])
+ round_digit_y<-as.numeric(Pars[7])
+
+ # Select the variables of study:
+ Env_Laboratory_data_df <- Data_frame_1[,c("PostCode","Date", variable_x,variable_y,variable_z,
+ "residents")]
+ Env_Pathogen_data_df <- Data_frame_2[,c("PostCode","Date", variable_x,variable_y,variable_z,
+ "residents","Cases")]
+
+
+ ################### Divide the domain of weather variable in bins of delta size for the uniform distribution.
+
+ breaks_x <-round(seq(min(na.omit(Env_Laboratory_data_df[[variable_x]])-delta_x),
+ max(na.omit(Env_Laboratory_data_df[[variable_x]])+delta_x) , by=delta_x), round_digit_x)
+
+ breaks_y <-round(seq(min(na.omit(Env_Laboratory_data_df[[variable_y]])-delta_y),
+ max(na.omit(Env_Laboratory_data_df[[variable_y]])+delta_y) , by=delta_y), round_digit_y)
+
+
+
+
+ #################### Create a data frame to fill with
+ Conditional_incidence_unif <-data.frame(character(), character(), numeric(), numeric())
+ colnames(Conditional_incidence_unif) <-c(variable_x, variable_y,"counts","residents_tot")
+
+
+
+ Env_Pathogen_data_z <-Env_Pathogen_data_df
+
+ Env_Laboratory_data_z <-Env_Laboratory_data_df
+
+ j_y_min <- 1
+ j_y_max <-length(breaks_y)
+
+ for (j_y in c(j_y_min:j_y_max)){ # from the bins resulting above, find another variable stratification
+
+ wt_y <- match.closest (Env_Pathogen_data_z[[ variable_y]], breaks_y)
+ ww_y <-which(wt_y==j_y)
+ Env_Pathogen_data_y <-Env_Pathogen_data_z[ww_y, ]
+
+
+ wt_y_t <- match.closest (Env_Laboratory_data_z[[ variable_y]], breaks_y)
+ ww_y_t <-which(wt_y_t==j_y)
+ Env_Laboratory_data_y <-Env_Laboratory_data_z[ww_y_t, ]
+
+
+ #if (length(breaks(variable_x, delta_temp))!=0)
+ {
+ j_x_min <-1
+ j_x_max <- length( breaks_x)
+
+ for (j_x in c(j_x_min:j_x_max)){
+
+ wt_x <- match.closest (Env_Pathogen_data_y[[variable_x]], breaks_x)
+ ww_x <-which(wt_x==j_x)
+ Pathogen_stratified <-Env_Pathogen_data_y[ww_x,] #where there are cases for specific values of the 3 weather variables. residents is the number of people exposed to these situations during the entire duration of the dataset. It is the number of times the people is exposed to this caracteristics, only when a case is found. It can count a same catchment area twice. Residents total is the same irrespectively of the presence of a case
+
+
+ wt_x_t <- match.closest (Env_Laboratory_data_y[[variable_x]], breaks_x)
+ ww_x_t <-which(wt_x_t==j_x)
+ Laboratory_stratified <-Env_Laboratory_data_y[ww_x_t, ] #NA in residents at Env_Laboratory_data. we do not care if there is a case reported
+
+ Total_cases <-sum((as.numeric(na.omit(Pathogen_stratified$Cases))))
+ residents_tot <-sum((as.numeric(na.omit(Laboratory_stratified$residents)))) #number of people exposed to these weather conditions. Will be used to calculate the incidence
+
+ data_df <-data.frame (
+ breaks_x[j_x],
+ breaks_y[j_y],
+ Total_cases,
+ residents_tot
+ )
+
+ colnames(data_df) <-c(variable_x, variable_y,"counts","residents_tot") #
+ Conditional_incidence_unif <-rbind(Conditional_incidence_unif, data_df)
+ }
+ }
+ }
+
+
+ print(c(j_x, j_y,variable_x,variable_y)) # check in the output log file if x y and z cover all the conditions established in breaks()
+
+
+
+ write.table(Conditional_incidence_unif, paste("../Data_Base/Conditional_incidence_",variable_x,"_",variable_y,"_",time_lag_char,info,sep=""), row.names = FALSE , sep = ",",eol = "\n")
+ return(list(Conditional_incidence_unif,breaks_x,breaks_y))
+
+ })
+}
diff --git a/Codes/Function_Plot_Conditional_Incidence_Uniform.R b/Codes/Function_Plot_Conditional_Incidence_Uniform.R
index 22a3251f5a458ce4df960a3a7b3245db89d5a6b2..47bfcd8dbf333a0c71a251ea8b847b79e1851797 100644
--- a/Codes/Function_Plot_Conditional_Incidence_Uniform.R
+++ b/Codes/Function_Plot_Conditional_Incidence_Uniform.R
@@ -28,10 +28,8 @@ Plot_Conditional_incidence_Uniform<-function(Data_frame_1, Pars)
Conditional_incidence_unif<-subset(Conditional_incidence_unif, Conditional_incidence_unif$counts>100)
- number_days<-as.Date(last_day)-as.Date(first_day)
- number_days<-as.numeric(number_days)
- number_year<-number_days/365
- norm<-(number_year)/(1e6)
+
+ norm<-1/(1e6)
@@ -105,11 +103,200 @@ Plot_Conditional_incidence_Uniform<-function(Data_frame_1, Pars)
+Plot_Conditional_incidence_Uniform_two_factors<-function(Data_frame_1, Pars)
+{
+
+ with(as.list(c(Data_frame_1, Pars)), {
+
+ # Pars is the list of 3 variables of interest
+ # Select the variables of study:
+
+ variable_x<-Pars[1]
+ variable_y<-Pars[2]
+ time_lag_char<-Pars[3]
+
+ # Select the file with relevant information:
+ Conditional_incidence_unif <- Data_frame_1[[1]]
+
+
+ #Conditional_incidence_unif<-Conditional_incidence_unif[,-1]
+ Conditional_incidence_unif$incidence<-Conditional_incidence_unif$counts/Conditional_incidence_unif$residents_tot
+ Conditional_incidence_unif<-na.omit(Conditional_incidence_unif)
+ #For visual purposes, we are removing cases when the counts were too low resulting in exceptional large incidence
+ Conditional_incidence_unif<-subset(Conditional_incidence_unif, Conditional_incidence_unif$counts>100)
+
+
+
+ norm<-1/(1e6)
+
+
+
+ Conditional_incidence_unif$preval<-Conditional_incidence_unif$incidence/norm
+
+
+
+ ################### Choose a specific day-length ###################
+
+
+ #Conditional_incidence_unif<-na.omit(Conditional_incidence_unif)
+ ############################# General variables ###########################
+ if (variable_x=="Relative_humidity"){variable_x_char<-"Relative humidity (%)"}
+ if (variable_y=="Relative_humidity"){variable_y_char<-"Relative humidity (%)"}
+ if (variable_x=="Mean_wind_speed"){variable_x_char<-"Mean Wind Speed (knots)"}
+ if (variable_y=="Mean_wind_speed"){variable_y_char<-"Mean Wind Speed (knots)"}
+ if (variable_x=="Mean_wind_speed"){variable_x_char<-"Mean Wind Speed (knots)"}
+ if (variable_y=="Mean_wind_speed"){variable_y_char<-"Mean Wind Speed (knots)"}
+ if (variable_x=="Mean_Precipitation"){variable_x_char<-"Mean Precipitation (mm)"}
+ if (variable_y=="Mean_Precipitation"){variable_y_char<-"Mean Precipitation (mm)"}
+ if (variable_x=="daylength"){variable_x_char<-"Day-Length (hours)"}
+ if (variable_y=="daylength"){variable_y_char<-"Day-Length (hours)"}
+ if (variable_x=="Minimum_air_temperature"){variable_x_char<-expression("Minimum Air Temperature"*~degree*C)}
+ if (variable_y=="Minimum_air_temperature"){variable_y_char<-expression("Minimum Air Temperature"*~degree*C)}
+ if (variable_x=="Maximum_air_temperature"){variable_x_char<-expression("Maximum Air Temperature"*~degree*C)}
+ if (variable_y=="Maximum_air_temperature"){variable_y_char<-expression("Maximum Air Temperature"*~degree*C)}
+ if (variable_x=="Difference_temperature"){variable_x_char<-expression("Difference Maximum-Minimum Air Temperature"*~degree*C)}
+ if (variable_y=="Difference_temperature"){variable_y_char<-expression("Difference Maximum-Minimum Air Temperature"*~degree*C)}
+
+
+ Conditional_incidence_unif[,c(variable_y)]<-round(Conditional_incidence_unif[,c(variable_y)])
+
+
+ var_x_loc_plot <- ggplot(Conditional_incidence_unif,aes(x=.data[[variable_x]],y=preval,color=(.data[[variable_y]]), fill=(.data[[variable_y]])))
+
+ var_x_loc_plot <- var_x_loc_plot+ scale_fill_continuous (name =variable_y_char,type = "viridis")
+ var_x_loc_plot <- var_x_loc_plot+ scale_colour_continuous(name =variable_y_char,type = "viridis")
+
+
+ var_x_loc_plot <- var_x_loc_plot+ylab("Conditional incidence of Campylobacter Cases (per year, per million)")
+ var_x_loc_plot <- var_x_loc_plot+xlab(variable_x_char)
+
+ var_x_loc_plot <- var_x_loc_plot+ geom_point(size=3.5)
+ var_x_loc_plot <- var_x_loc_plot+ geom_line(size=1)
+
+ var_x_loc_plot <- var_x_loc_plot+ scale_x_continuous(limits=c(min(Conditional_incidence_unif[,c(variable_x)]),max(Conditional_incidence_unif[,c(variable_x)])))
+
+ var_x_loc_plot <- var_x_loc_plot+
+ theme(legend.position= c(0.125,0.75),legend.title = element_text( size = 10),
+ legend.text = element_text( size = 10),legend.background = element_blank(),
+ legend.key=element_rect(colour=NA))
+
+ var_x_loc_plot <- var_x_loc_plot+ theme(axis.title.x =element_text( colour="#990000", size=13))
+ var_x_loc_plot <- var_x_loc_plot+ theme(axis.title.y =element_text( colour="#990000", size=13))
+
+ var_x_loc_plot <- var_x_loc_plot+ theme(axis.title.x =element_text(size=13))
+ var_x_loc_plot <- var_x_loc_plot+ theme(axis.title.y =element_text(size=13))
+
+ print(var_x_loc_plot)
+
+
+
+ file_name<-paste("../Graphs/Conditional_probability_",variable_x,"_",variable_y,"_",time_lag_char,"_uniform.tiff", sep = "")
+ tiff(filename = file_name,width = 17.35, height = 17.35, units = "cm", pointsize = 9, res = 1200,compression = "lzw",antialias="default")
+
+ print(var_x_loc_plot)
+
+ dev.off()
+
+
+ })
+}
-
+Plot_Conditional_incidence_unif_constant<-function(Data_frame_1, Pars)
+{
+
+ with(as.list(c(Data_frame_1, Pars)), {
+
+ # Pars is the list of 3 variables of interest
+ # Select the variables of study:
+
+ variable<-Pars[1]
+ info<- Pars[2]
+ # Select the file with relevant information:
+ Conditional_incidence_quantiles <- Data_frame_1
+
+ norm<-1/(1e6)
+ #Conditional_incidence_quantiles<-Conditional_incidence_quantiles[,-1]
+ Conditional_incidence_quantiles$incidence<-(Conditional_incidence_quantiles$counts/Conditional_incidence_quantiles$residents_tot)/norm
+
+
+
+
+ Conditional_incidence_df<-Conditional_incidence_quantiles
+
+ ################### Choose a specific day-length ###################
+
+
+ #Conditional_incidence_quantiles<-na.omit(Conditional_incidence_quantiles)
+ ############################# General variables ###########################
+ if (variable=="Relative_humidity"){variable_x_char<-"Relative humidity (%)"}
+ if (variable=="daylength"){variable_x_char<-"Day-Length (hours)"}
+ if (variable=="Maximum_air_temperature"){variable_x_char<-expression("Maximum Air Temperature"*~degree*C)}
+
+
+
+
+ if (variable=="Maximum_air_temperature" ){
+ Conditional_incidence_plot <-
+ ggplot(Conditional_incidence_df,aes(x=Maximum_air_temperature,y=incidence))
+ Conditional_incidence_plot <- Conditional_incidence_plot+ scale_x_continuous(limits=c(10,30))
+ Conditional_incidence_plot <- Conditional_incidence_plot+ geom_point(size=3.5, color="red")+geom_line(size=1, color="red")
+
+ }
+ if (variable=="Relative_humidity" ){
+ Conditional_incidence_plot <-
+ ggplot(Conditional_incidence_df,aes(x=Relative_humidity,y=incidence))
+ Conditional_incidence_plot <- Conditional_incidence_plot+ scale_x_continuous(limits=c(60,95))
+ Conditional_incidence_plot <- Conditional_incidence_plot+ geom_point(size=3.5, color="darkblue")+geom_line(size=1, color="darkblue")
+ }
+ if (variable=="daylength"){
+ Conditional_incidence_plot <-
+ ggplot(Conditional_incidence_df,aes(x=daylength,y=incidence))
+ Conditional_incidence_plot <- Conditional_incidence_plot+ scale_x_continuous(limits=c(10,18))
+ Conditional_incidence_plot <- Conditional_incidence_plot+ geom_point(size=3.5, color="darkgreen")+geom_line(size=1, color="darkgreen")
+ }
+
+
+
+
+ Conditional_incidence_plot <- Conditional_incidence_plot+theme_bw(20)
+ Conditional_incidence_plot <- Conditional_incidence_plot+ylab("Conditional Incidence of Campylobacter Cases ")
+ Conditional_incidence_plot <- Conditional_incidence_plot+xlab(variable_x_char)
+
+ Conditional_incidence_plot <- Conditional_incidence_plot+ guides(colour="none") # guides(colour="none")
+ Conditional_incidence_plot <- Conditional_incidence_plot
+ Conditional_incidence_plot <- Conditional_incidence_plot
+
+ Conditional_incidence_plot <- Conditional_incidence_plot+
+ theme(legend.position= c(0.25,0.75),legend.title = element_text( size = 11),
+ legend.text = element_text( size = 10),legend.background = element_blank(),
+ legend.key=element_rect(colour=NA))
+
+ Conditional_incidence_plot <- Conditional_incidence_plot+ theme(axis.title.x =element_text( colour="#990000", size=15))
+ Conditional_incidence_plot <- Conditional_incidence_plot+ theme(axis.title.y =element_text( colour="#990000", size=15))
+
+
+
+ Conditional_incidence_plot <- Conditional_incidence_plot+ theme(axis.title.x =element_text(size=15))
+ Conditional_incidence_plot <- Conditional_incidence_plot+ theme(axis.title.y =element_text(size=15))
+
+ print(Conditional_incidence_plot)
+
+
+ file_name<-paste("../Graphs/Campylobacter_",info,".tiff", sep = "")
+ tiff(filename = file_name,width = 17.35, height = 17.35, units = "cm", pointsize = 12, res = 1200,compression = "lzw",antialias="default")
+
+ print(Conditional_incidence_plot)
+
+ dev.off()
+
+
+
+
+ })
+}
diff --git a/Codes/Pathogen_Linkage_fixed_time_lag.R b/Codes/Pathogen_Linkage_fixed_time_lag.R
index 0d147e8bef0325e118e3583c5add03950bb68a41..53938b2b8ecddcc18c656c7bf2ca772475dcd8dd 100644
--- a/Codes/Pathogen_Linkage_fixed_time_lag.R
+++ b/Codes/Pathogen_Linkage_fixed_time_lag.R
@@ -174,7 +174,7 @@ write.table(diagnostic_laboratory_df[-c(length(diagnostic_laboratory_df))],paste
################ Define the time lag here ###############################
-time_lag<-7
+time_lag<-90
time_lag_char<-paste(time_lag)
## This create a database for all postocode where there was at least one case