R/popandvac.R

In arranhamlet/popvac_package: Generate population level vaccination data from a csv template.

#' Generates population and vaccination arrays
#' Takes the initial information to generate population and vaccination arrays for the desired countries. Only inputs the EPI coverage for vaccination and sets it up for the addition of vaccination campaigns.
#' 
#' @param epi This is a csv containing the routine infant immunization coverage for all the countries you want to generate population vaccination coverage for.
#' @param popdf The csv file generated by the function csv_to_adm2pop \code{\link{csv_to_adm2pop}}
#' @param countryiso This is a vector of the ISO codes for all countries you want to generate population vaccination coverage for.
#' #@param shp0 This is a country level shapefule that contains all the countries you want to generate population vaccination coverage for.
#' #@param shp1 This is a first administrative division shapefule that contains all the countries you want to generate population vaccination coverage for.
#' @param shp2 This is a second administrative division shapefule that contains all the countries you want to generate population vaccination coverage for.
#' @param scenario This changes the underlying assumptions of the population size. Base does nothing, "hiCov" reduces the population by 25 percent, "loCov" increases the population by 25 percent.
#' @param GAVI.switch This implements different vaccination campaigns based on the column "scenario" in the campaigncsv parameter. "-all.epi" ignores the contribution of epi, "5year.mass" implements five year mass vaccination campaigns,"epi90" does...something?,"+fut.epi" includes future epi coverage.
#' @param skew This affects how vaccination campaigns are applied. If the skew is set to -1 then campaigns are applied to the whole population, regardless of prior vaccination status. If the skew is set to 0 then campaigns are applied only to the unvaccinated population.
#' @examples
#' popandvac(epi,countryiso,shp1,shp2)
#' @keywords internal

# popdf<-read.csv("//fi--didenas1.dide.local/yf/Arran/PhD/Code/age_distribution_generator/data/WPP2015_INT_F3_Population_By_Age_Annual_Single_Medium.csv")

popandvac<-function(epi,popdf,countryiso,shplist,scenario="base",GAVI.switch="base",skew=0){
  
  #Release shapefiles
  shp2_v1.0<-shplist$shp2_v1.0
  shp2_v2.0<-shplist$shp2_v2.0
  shp2_v2.8<-shplist$shp2_v2.8

  #Subset population to country and remove years 2051:2100
  countrypop<-popdf[popdf$ISO %in% countryiso,]
  # no2051to2100<-which(countrypop$Date %in% 2051:2100)
  countryearpop<-if(range(as.numeric(countrypop$Date))[2]>2050) countrypop[-which(countrypop$Date>2050),] else countrypop

  #Find the full names of the countries
  countryiso_names<-as.character(as.character(shp2_v2.8$NAME_0[shp2_v2.8$ISO %in% countryiso]))
  
  #Subset shapefiles to included ISOs
  shp2_v2.8sub<-shp2_v2.8[shp2_v2.8$ISO %in% countryiso,]
  
  # dat2<-data.frame(shp2_v2.8sub[order(shp2_v2.8sub$ID_2),])
  dat2<-data.frame(shp2_v2.8sub)
  
  #Set up ID's for locations consistent with previous formatting 
  dat2$adm0_adm1<-paste(dat2$ISO, dat2$ID_1, sep="")
  dat2$adm0_adm2<-paste(dat2$ISO, dat2$ID_2, sep="")
  
  #Set different lengths
  l1<-length(table(dat2$adm0_adm2))
  l2<-length(table(dat2$year))
  l3<-ncol(dat2)-9-8 

  #Setting the ISO identification and names to those present in 1950 (should be all of them)
  dn1.0<-dat2$ISO
  dn1.1<-paste(dat2$ISO,dat2$ID_1,sep="")
  dn1.2<-if(length(dat2$SP_ID)==0) dat2$adm0_adm2 else dat2$SP_ID
  dn1.0.name<-dat2$NAME_0
  dn1.1.name<-dat2$NAME_1
  dn1.2.name<-dat2$NAME_2
  # dn2<-1940:2050
  dn2<-1950:2050
  dn3<-paste("a",0:100,sep="")
  
  #Format gadm1_adm1
  gadm1_adm1$adm0_adm1<-paste(gadm1_adm1$adm0.new,gadm1_adm1$adm1.new, sep="")
  mm1<-match(dn1.1, gadm1_adm1$adm0_adm1)
  dn1.1_gadm1<-gadm1_adm1$adm1.old[mm1]
  #Format gadm1_adm2
  gadm1_adm2$adm0_adm2<-paste(gadm1_adm2$adm0.new,gadm1_adm2$adm2.new, sep="")
  mm2<-match(dn1.2, gadm1_adm2$adm0_adm2)
  dn1.2_gadm1<-gadm1_adm2$adm2.old[mm2]

  dndf<-data.frame("adm0"=dn1.0,"adm1"=dn1.1,"adm2"=dn1.2,"dn1.0.name"=dn1.0.name,"dn1.1.name"=dn1.1.name,"dn1.2.name"=dn1.2.name,
                   "dn1.1_gadm1"=dn1.1_gadm1,"dn1.2_gadm1"=dn1.2_gadm1)
  

  pop1.adm2<-as.matrix(countryearpop[,which(names(countryearpop) %in% paste("X",0:100,sep=""))])
  
  
  dim(pop1.adm2)<-c(length(dn1.0),length(dn2),length(dn3))
  
  # changing the underlying assumptions of the population size:
  pop1.adm2.base<-pop1.adm2

  if(scenario=="base") {
    tag<-""
    pop1.adm2<-pop1.adm2.base
  } else if(scenario == "hiCov") {
    tag<-"_highCoverage"
    pop1.adm2<-0.75*pop1.adm2.base  # for the low population, i.e. high coverage scenario.
  } else if(scenario == "loCov") {
    pop.adm2<-1.25*pop1.adm2.base  # for the high population, i.e. low coverage scenario.
    tag<-"_lowCoverage"
  } else if(scenario=="M.base") {
    tag<-"_Moreau"
    pop1.adm2<-pop1.adm2.base
  } else if(scenario=="M.hiCov") {
    tag<-"_Moreau_highCoverage"
    pop1.adm2<-0.75*pop1.adm2.base
  } else if(scenario=="M.loCov") {
    tag<-"_Moreau_lowCoverage"
    pop1.adm2<-1.25*pop1.adm2.base
  }
  
  vacc.cov<-rep(0,length(pop1.adm2))
  dim(vacc.cov)<-dim(pop1.adm2)
  
  ## implement EPI coverage:
  if(GAVI.switch != "-all.epi") {
    epi<-epi[epi$ISO_code %in% countryiso,]#1:40]
    epi[is.na(epi)] = 0
    y.end<-2049
    if(GAVI.switch %in% c("epi90","5year.mass")) y.end = 2013
    for(y in 1980:y.end) {
      cc<-which(names(epi)==paste("X",y,sep=""))
      for(rr in 1:nrow(epi)) {
        cnt<-epi$ISO_code[rr]
        cover<-epi[rr,cc]
        if(cover>0) vacc.cov<-add.vacc.campaign(dndf,pop1.adm2,vacc.cov,year=y,country=cnt, coverage = cover/100, agemin = 0, agemax = 0, skew=skew)    
      } 
    } 
  }

  exportvals<-list(dndf,vacc.cov,dn2,dn3,shp2_v2.8sub,pop1.adm2)
  names(exportvals)<-c("dndf","vacc.cov","dn2","dn3","shp2_v2.8sub","pop1.adm2")
  exportvals

}