Added more detail in fan charts.

2B_Result_Analysis.r

@@ -1,9 +1,9 @@
 library(tidyverse)
 library(gt) #For nice color coded capacity limits table.
+source("Scripts/Load_Custom_Functions/Fan_Chart_Creation_Functions.r") #Functions created to make nice graphs
 	if(!exists("SAVE_RES_LOC")){SAVE_RES_LOC <- "./Results/Primary_Simulation_Results/Main_Results/"}
 	dir.create(SAVE_RES_LOC, recursive = TRUE, showWarnings = FALSE)
 
-
 ###Process the simulations and save the main percentile results by year
 	RES <- read_csv("Results/Simulations/Kemmerer_2024_Simulation.csv")
 	RES[,"Year"] <- RES[,"Year"]
@@ -12,33 +12,6 @@ library(gt) #For nice color coded capacity limits table.
 	HIST <-	readRDS("Data/Cleaned_Data/Population_Data/RDS/Kemmerer_Diamondville_Population_Data.Rds") %>% filter(County=='Lincoln') %>% mutate(Percentile="Actual Population") %>% filter(Year>=1940)
 	######Population
 	####Fan New
-	GET_DATA <- function(RES,COL_NUM){ 
-
-	YEARS <- min(RES$Year,na.rm=TRUE):max(RES$Year,na.rm=TRUE)
-	FAN_DATA <- do.call(rbind,lapply(YEARS,function(x){quantile(as.numeric(t((RES %>% filter(Year==x))[,COL_NUM])),seq(0.01,0.99,by=0.01))}))  %>% as_tibble %>% mutate(Year=YEARS)
-	FAN_DATA <- 	rbind(FAN_DATA[1,],FAN_DATA)
-	START_VALUE <- (HIST %>% filter(Year==2024))[,COL_NUM+1] %>% as.numeric
-	FAN_DATA <- 	FAN_DATA  %>% pivot_longer(colnames(FAN_DATA %>% select(-Year)),names_to="Percentile") %>% filter(Year>2024) %>% unique
-	NUM_YEARS <- length(unique(FAN_DATA$Year) )
-	FAN_DATA$Group <- rep(c(1:49,0,rev(1:49)),NUM_YEARS)
-	FAN_DATA <- FAN_DATA %>% group_by(Year,Group)  %>% summarize(MIN=min(value),MAX=max(value))
-	TEMP <- FAN_DATA %>% filter(Year==2025) %>% mutate(Year=2024) %>% ungroup
-	TEMP[,3:4] <- START_VALUE 
-	FAN_DATA <- rbind(TEMP,FAN_DATA %>% ungroup) %>% as_tibble
-	return(FAN_DATA)
-	}
-	RES %>% pull(Sim_UUID) %>% unique %>% length()
-MAKE_GRAPH <- function(GRAPH_DATA,ALPHA=0.03,COLOR='cadetblue',LINE_WIDTH=0.75){
-	PLOT <- ggplot(data=GRAPH_DATA)
-	for(i in 1:49){
-	C_DATA <- GRAPH_DATA%>% filter(Group==i)
-	     PLOT <- 	PLOT +geom_ribbon(data=C_DATA,aes(x=Year,ymin=MIN,ymax=MAX),alpha=ALPHA,fill=COLOR)
-
-}
-	CI_90 <- rbind(GRAPH_DATA%>% filter(Group==20) %>% mutate(Interval='80%'),GRAPH_DATA%>% filter(Group==5) %>% mutate(Interval='95%'),GRAPH_DATA%>% filter(Group==0) %>% mutate(Interval='Median Prediction'))
-	PLOT <- PLOT+geom_line(aes(x=Year,y=MIN,linetype=Interval,color=Interval),linewidth=LINE_WIDTH, data=CI_90)+geom_line(aes(x=Year,y=MAX,group=Interval,linetype=Interval,color=Interval),linewidth=LINE_WIDTH ,data=CI_90)+scale_color_manual(values=c("grey50","grey80","black"))+scale_linetype_manual(values = c("solid","solid","dotdash"))
-	return(PLOT)
-}
 POP_DATA <- GET_DATA(RES,3)
 POP_PLOT <- MAKE_GRAPH(POP_DATA)
 POP_PLOT <- POP_PLOT+geom_line(data=HIST,aes(x=Year,y=Population),color='black',linewidth=0.75)+ scale_x_continuous(breaks = c(seq(1940, 2060, by = 10),2065))+ scale_y_continuous(breaks = seq(0, 35000, by = 500))+ggtitle("Kemmerer & Diamondville, Population Forecast")+ expand_limits( y = 0)+labs(color = "Prediction Interval",linetype="Prediction Interval",y="Population")+ theme_bw()+ theme(legend.position = "top",panel.grid.minor = element_blank())

2C_Upper_Bound_Result_Analysis.r

@@ -1,5 +1,6 @@
 library(tidyverse)
 library(gt) #For nice color coded capacity limits table.
+source("Scripts/Load_Custom_Functions/Fan_Chart_Creation_Functions.r") #Functions created to make nice graphs
 	if(!exists("SAVE_RES_LOC")){SAVE_RES_LOC <- "./Results/Primary_Simulation_Results/Upper_Bound_Results/"}
 	dir.create(SAVE_RES_LOC, recursive = TRUE, showWarnings = FALSE)
 
@@ -12,56 +13,29 @@ library(gt) #For nice color coded capacity limits table.
 	HIST <-	readRDS("Data/Cleaned_Data/Population_Data/RDS/Kemmerer_Diamondville_Population_Data.Rds") %>% filter(County=='Lincoln') %>% mutate(Percentile="Actual Population") %>% filter(Year>=1940)
 	######Population
 	####Fan New
-	GET_DATA <- function(RES,COL_NUM){ 
-
-	YEARS <- min(RES$Year,na.rm=TRUE):max(RES$Year,na.rm=TRUE)
-	FAN_DATA <- do.call(rbind,lapply(YEARS,function(x){quantile(as.numeric(t((RES %>% filter(Year==x))[,COL_NUM])),seq(0.01,0.99,by=0.01))}))  %>% as_tibble %>% mutate(Year=YEARS)
-	FAN_DATA <- 	rbind(FAN_DATA[1,],FAN_DATA)
-	START_VALUE <- (HIST %>% filter(Year==2024))[,COL_NUM+1] %>% as.numeric
-	FAN_DATA <- 	FAN_DATA  %>% pivot_longer(colnames(FAN_DATA %>% select(-Year)),names_to="Percentile") %>% filter(Year>2024) %>% unique
-	NUM_YEARS <- length(unique(FAN_DATA$Year) )
-	FAN_DATA$Group <- rep(c(1:49,0,rev(1:49)),NUM_YEARS)
-	FAN_DATA <- FAN_DATA %>% group_by(Year,Group)  %>% summarize(MIN=min(value),MAX=max(value))
-	TEMP <- FAN_DATA %>% filter(Year==2025) %>% mutate(Year=2024) %>% ungroup
-	TEMP[,3:4] <- START_VALUE 
-	FAN_DATA <- rbind(TEMP,FAN_DATA %>% ungroup) %>% as_tibble
-	return(FAN_DATA)
-	}
-	RES %>% pull(Sim_UUID) %>% unique %>% length()
-MAKE_GRAPH <- function(GRAPH_DATA,ALPHA=0.03,COLOR='springgreen4',LINE_WIDTH=0.75){
-	PLOT <- ggplot(data=GRAPH_DATA)
-	for(i in 1:49){
-	C_DATA <- GRAPH_DATA%>% filter(Group==i)
-	     PLOT <- 	PLOT +geom_ribbon(data=C_DATA,aes(x=Year,ymin=MIN,ymax=MAX),alpha=ALPHA,fill=COLOR)
-
-}
-	CI_90 <- rbind(GRAPH_DATA%>% filter(Group==20) %>% mutate(Interval='80%'),GRAPH_DATA%>% filter(Group==5) %>% mutate(Interval='95%'),GRAPH_DATA%>% filter(Group==0) %>% mutate(Interval='Median Prediction'))
-	PLOT <- PLOT+geom_line(aes(x=Year,y=MIN,linetype=Interval,color=Interval),linewidth=LINE_WIDTH, data=CI_90)+geom_line(aes(x=Year,y=MAX,group=Interval,linetype=Interval,color=Interval),linewidth=LINE_WIDTH ,data=CI_90)+scale_color_manual(values=c("grey50","grey80","black"))+scale_linetype_manual(values = c("solid","solid","dotdash"))
-	return(PLOT)
-}
 POP_DATA <- GET_DATA(RES,3)
-POP_PLOT <- MAKE_GRAPH(POP_DATA)
+POP_PLOT <- MAKE_GRAPH(POP_DATA,COLOR='springgreen4')
 POP_PLOT <- POP_PLOT+geom_line(data=HIST,aes(x=Year,y=Population),color='black',linewidth=0.75)+ scale_x_continuous(breaks = c(seq(1940, 2060, by = 10),2065))+ scale_y_continuous(breaks = seq(0, 35000, by = 500))+ggtitle("Kemmerer & Diamondville, Population Forecast")+ expand_limits( y = 0)+labs(color = "Prediction Interval",linetype="Prediction Interval",y="Population")+ theme_bw()+ theme(legend.position = "top",panel.grid.minor = element_blank())
 png(paste0(SAVE_RES_LOC,"Population_Fan_Chart.png"), width = 12, height = 8, units = "in", res = 600)
 POP_PLOT 
 dev.off()
 
 BIRTH_DATA <- GET_DATA(RES,4)
-BIRTH_PLOT <- MAKE_GRAPH(BIRTH_DATA)
+BIRTH_PLOT <- MAKE_GRAPH(BIRTH_DATA,COLOR='springgreen4')
 BIRTH_PLOT <- BIRTH_PLOT+geom_line(data=HIST,aes(x=Year,y=Births),color='black',linewidth=0.75)+ scale_x_continuous(breaks = c(seq(2010, 2060, by = 5),2065),limits=c(2009,2065))+ scale_y_continuous(breaks = seq(0, 35000, by = 10))+ggtitle("Kemmerer & Diamondville, Birth Forecast")+ expand_limits( y = 0)+labs(color = "Prediction Interval",linetype="Prediction Interval",y="Births")+ theme_bw()+ theme(legend.position = "top",panel.grid.minor = element_blank())
 png(paste0(SAVE_RES_LOC,"Birth_Fan_Chart.png"), width = 12, height = 8, units = "in", res = 600)
 	BIRTH_PLOT 
 dev.off()
 
 DEATH_DATA <- GET_DATA(RES,5) %>% filter(!is.na(MIN))
-DEATH_PLOT <- MAKE_GRAPH(DEATH_DATA)
+DEATH_PLOT <- MAKE_GRAPH(DEATH_DATA,COLOR='springgreen4')
 DEATH_PLOT <-	DEATH_PLOT+geom_line(data=HIST,aes(x=Year,y=Deaths),color='black',linewidth=0.75)+ scale_x_continuous(breaks = c(seq(2010, 2060, by = 5),2065),limits=c(2009,2065))+ scale_y_continuous(breaks = seq(0, 35000, by = 10))+ggtitle("Kemmerer & Diamondville, Mortality Forecast")+ expand_limits( y = 0)+labs(color = "Prediction Interval",linetype="Prediction Interval",y="Deaths")+ theme_bw()+ theme(legend.position = "top",panel.grid.minor = element_blank())
 png(paste0(SAVE_RES_LOC,"Mortality_Fan_Chart.png"), width = 12, height = 8, units = "in", res = 600)
 	DEATH_PLOT
 dev.off()
 
 MIGRATION_DATA <- GET_DATA(RES,6) %>% filter(!is.na(MIN))
-MIGRATION_PLOT <- MAKE_GRAPH(MIGRATION_DATA)
+MIGRATION_PLOT <- MAKE_GRAPH(MIGRATION_DATA,COLOR='springgreen4')
 MIGRATION_PLOT <- 	MIGRATION_PLOT+geom_line(data=HIST,aes(x=Year,y=Migration),color='black',linewidth=0.75)+ scale_x_continuous(breaks = c(seq(2010, 2060, by = 5),2065),limits=c(2009,2065))+ scale_y_continuous(breaks = seq(-1000, 1000, by = 50))+ggtitle("Kemmerer & Diamondville, Net Migration Forecast")+ expand_limits( y = 0)+labs(color = "Prediction Interval",linetype="Prediction Interval",y="Migration")+ theme_bw()+ theme(legend.position = "top",panel.grid.minor = element_blank())
 png(paste0(SAVE_RES_LOC,"Migration_Fan_Chart.png"), width = 12, height = 8, units = "in", res = 600)
 	MIGRATION_PLOT

2D_Lower_Bound_Result_Analysis.r

@@ -1,5 +1,6 @@
 library(tidyverse)
 library(gt) #For nice color coded capacity limits table.
+source("Scripts/Load_Custom_Functions/Fan_Chart_Creation_Functions.r") #Functions created to make nice graphs
 	if(!exists("SAVE_RES_LOC")){SAVE_RES_LOC <- "./Results/Primary_Simulation_Results/Lower_Bound_Results/"}
 	dir.create(SAVE_RES_LOC, recursive = TRUE, showWarnings = FALSE)
 
@@ -12,56 +13,29 @@ library(gt) #For nice color coded capacity limits table.
 	HIST <-	readRDS("Data/Cleaned_Data/Population_Data/RDS/Kemmerer_Diamondville_Population_Data.Rds") %>% filter(County=='Lincoln') %>% mutate(Percentile="Actual Population") %>% filter(Year>=1940)
 	######Population
 	####Fan New
-	GET_DATA <- function(RES,COL_NUM){ 
-
-	YEARS <- min(RES$Year,na.rm=TRUE):max(RES$Year,na.rm=TRUE)
-	FAN_DATA <- do.call(rbind,lapply(YEARS,function(x){quantile(as.numeric(t((RES %>% filter(Year==x))[,COL_NUM])),seq(0.01,0.99,by=0.01))}))  %>% as_tibble %>% mutate(Year=YEARS)
-	FAN_DATA <- 	rbind(FAN_DATA[1,],FAN_DATA)
-	START_VALUE <- (HIST %>% filter(Year==2024))[,COL_NUM+1] %>% as.numeric
-	FAN_DATA <- 	FAN_DATA  %>% pivot_longer(colnames(FAN_DATA %>% select(-Year)),names_to="Percentile") %>% filter(Year>2024) %>% unique
-	NUM_YEARS <- length(unique(FAN_DATA$Year) )
-	FAN_DATA$Group <- rep(c(1:49,0,rev(1:49)),NUM_YEARS)
-	FAN_DATA <- FAN_DATA %>% group_by(Year,Group)  %>% summarize(MIN=min(value),MAX=max(value))
-	TEMP <- FAN_DATA %>% filter(Year==2025) %>% mutate(Year=2024) %>% ungroup
-	TEMP[,3:4] <- START_VALUE 
-	FAN_DATA <- rbind(TEMP,FAN_DATA %>% ungroup) %>% as_tibble
-	return(FAN_DATA)
-	}
-	RES %>% pull(Sim_UUID) %>% unique %>% length()
-MAKE_GRAPH <- function(GRAPH_DATA,ALPHA=0.03,COLOR='firebrick2',LINE_WIDTH=0.75){
-	PLOT <- ggplot(data=GRAPH_DATA)
-	for(i in 1:49){
-	C_DATA <- GRAPH_DATA%>% filter(Group==i)
-	     PLOT <- 	PLOT +geom_ribbon(data=C_DATA,aes(x=Year,ymin=MIN,ymax=MAX),alpha=ALPHA,fill=COLOR)
-
-}
-	CI_90 <- rbind(GRAPH_DATA%>% filter(Group==20) %>% mutate(Interval='80%'),GRAPH_DATA%>% filter(Group==5) %>% mutate(Interval='95%'),GRAPH_DATA%>% filter(Group==0) %>% mutate(Interval='Median Prediction'))
-	PLOT <- PLOT+geom_line(aes(x=Year,y=MIN,linetype=Interval,color=Interval),linewidth=LINE_WIDTH, data=CI_90)+geom_line(aes(x=Year,y=MAX,group=Interval,linetype=Interval,color=Interval),linewidth=LINE_WIDTH ,data=CI_90)+scale_color_manual(values=c("grey50","grey80","black"))+scale_linetype_manual(values = c("solid","solid","dotdash"))
-	return(PLOT)
-}
 POP_DATA <- GET_DATA(RES,3)
-POP_PLOT <- MAKE_GRAPH(POP_DATA)
+POP_PLOT <- MAKE_GRAPH(POP_DATA,COLOR='firebrick2')
 POP_PLOT <- POP_PLOT+geom_line(data=HIST,aes(x=Year,y=Population),color='black',linewidth=0.75)+ scale_x_continuous(breaks = c(seq(1940, 2060, by = 10),2065))+ scale_y_continuous(breaks = seq(0, 35000, by = 500))+ggtitle("Kemmerer & Diamondville, Population Forecast")+ expand_limits( y = 0)+labs(color = "Prediction Interval",linetype="Prediction Interval",y="Population")+ theme_bw()+ theme(legend.position = "top",panel.grid.minor = element_blank())
 png(paste0(SAVE_RES_LOC,"Population_Fan_Chart.png"), width = 12, height = 8, units = "in", res = 600)
 POP_PLOT 
 dev.off()
 
 BIRTH_DATA <- GET_DATA(RES,4)
-BIRTH_PLOT <- MAKE_GRAPH(BIRTH_DATA)
+BIRTH_PLOT <- MAKE_GRAPH(BIRTH_DATA,COLOR='firebrick2')
 BIRTH_PLOT <- BIRTH_PLOT+geom_line(data=HIST,aes(x=Year,y=Births),color='black',linewidth=0.75)+ scale_x_continuous(breaks = c(seq(2010, 2060, by = 5),2065),limits=c(2009,2065))+ scale_y_continuous(breaks = seq(0, 35000, by = 10))+ggtitle("Kemmerer & Diamondville, Birth Forecast")+ expand_limits( y = 0)+labs(color = "Prediction Interval",linetype="Prediction Interval",y="Births")+ theme_bw()+ theme(legend.position = "top",panel.grid.minor = element_blank())
 png(paste0(SAVE_RES_LOC,"Birth_Fan_Chart.png"), width = 12, height = 8, units = "in", res = 600)
 	BIRTH_PLOT 
 dev.off()
 
 DEATH_DATA <- GET_DATA(RES,5) %>% filter(!is.na(MIN))
-DEATH_PLOT <- MAKE_GRAPH(DEATH_DATA)
+DEATH_PLOT <- MAKE_GRAPH(DEATH_DATA,COLOR='firebrick2')
 DEATH_PLOT <-	DEATH_PLOT+geom_line(data=HIST,aes(x=Year,y=Deaths),color='black',linewidth=0.75)+ scale_x_continuous(breaks = c(seq(2010, 2060, by = 5),2065),limits=c(2009,2065))+ scale_y_continuous(breaks = seq(0, 35000, by = 10))+ggtitle("Kemmerer & Diamondville, Mortality Forecast")+ expand_limits( y = 0)+labs(color = "Prediction Interval",linetype="Prediction Interval",y="Deaths")+ theme_bw()+ theme(legend.position = "top",panel.grid.minor = element_blank())
 png(paste0(SAVE_RES_LOC,"Mortality_Fan_Chart.png"), width = 12, height = 8, units = "in", res = 600)
 	DEATH_PLOT
 dev.off()
 
 MIGRATION_DATA <- GET_DATA(RES,6) %>% filter(!is.na(MIN))
-MIGRATION_PLOT <- MAKE_GRAPH(MIGRATION_DATA)
+MIGRATION_PLOT <- MAKE_GRAPH(MIGRATION_DATA,COLOR='firebrick2')
 MIGRATION_PLOT <- 	MIGRATION_PLOT+geom_line(data=HIST,aes(x=Year,y=Migration),color='black',linewidth=0.75)+ scale_x_continuous(breaks = c(seq(2010, 2060, by = 5),2065),limits=c(2009,2065))+ scale_y_continuous(breaks = seq(-1000, 1000, by = 50))+ggtitle("Kemmerer & Diamondville, Net Migration Forecast")+ expand_limits( y = 0)+labs(color = "Prediction Interval",linetype="Prediction Interval",y="Migration")+ theme_bw()+ theme(legend.position = "top",panel.grid.minor = element_blank())
 png(paste0(SAVE_RES_LOC,"Migration_Fan_Chart.png"), width = 12, height = 8, units = "in", res = 600)
 	MIGRATION_PLOT

Scripts/Load_Custom_Functions/Fan_Chart_Creation_Functions.r

@@ -0,0 +1,32 @@
+#Functions to create consistent fan functions from simulations
+	###Create a data set which pairs of the upper and lower confidence intervals. This allows the layers to be plotted
+GET_DATA <- function(RES,COL_NUM){ 
+	YEARS <- min(RES$Year,na.rm=TRUE):max(RES$Year,na.rm=TRUE)
+	LEVELS <- seq(0.01,0.99,by=0.005)
+	GROUPS <- floor(length(LEVELS)/2)
+	FAN_DATA <- do.call(rbind,lapply(YEARS,function(x){quantile(as.numeric(t((RES %>% filter(Year==x))[,COL_NUM])),LEVELS)}))  %>% as_tibble %>% mutate(Year=YEARS)
+	FAN_DATA <- 	rbind(FAN_DATA[1,],FAN_DATA)
+	START_VALUE <- (HIST %>% filter(Year==2024))[,COL_NUM+1] %>% as.numeric
+	FAN_DATA <- 	FAN_DATA  %>% pivot_longer(colnames(FAN_DATA %>% select(-Year)),names_to="Percentile") %>% filter(Year>2024) %>% unique
+	NUM_YEARS <- length(unique(FAN_DATA$Year) )
+	FAN_DATA$Group <- rep(c(1:GROUPS,0,rev(1:GROUPS)),NUM_YEARS)
+	FAN_DATA <- FAN_DATA %>% group_by(Year,Group)  %>% summarize(MIN=min(value),MAX=max(value))
+	TEMP <- FAN_DATA %>% filter(Year==2025) %>% mutate(Year=2024) %>% ungroup
+	TEMP[,3:4] <- START_VALUE 
+	FAN_DATA <- rbind(TEMP,FAN_DATA %>% ungroup) %>% as_tibble
+	return(FAN_DATA)
+	}
+###Create a data set which loops through many confidence bounds to create a layered fan chart. This output can be stacked on with other data later 
+MAKE_GRAPH <- function(GRAPH_DATA,ALPHA=0.03,COLOR='cadetblue',LINE_WIDTH=0.75){
+	PLOT <- ggplot(data=GRAPH_DATA)
+	LENGTH <- floor(nrow(GRAPH_DATA)/2)
+	for(i in 1:LENGTH){
+	C_DATA <- GRAPH_DATA%>% filter(Group==i)
+	     PLOT <- 	PLOT +geom_ribbon(data=C_DATA,aes(x=Year,ymin=MIN,ymax=MAX),alpha=ALPHA,fill=COLOR)
+
+}
+	CI_90 <- rbind(GRAPH_DATA%>% filter(Group==20) %>% mutate(Interval='80%'),GRAPH_DATA%>% filter(Group==5) %>% mutate(Interval='95%'),GRAPH_DATA%>% filter(Group==0) %>% mutate(Interval='Median Prediction'))
+	PLOT <- PLOT+geom_line(aes(x=Year,y=MIN,linetype=Interval,color=Interval),linewidth=LINE_WIDTH, data=CI_90)+geom_line(aes(x=Year,y=MAX,group=Interval,linetype=Interval,color=Interval),linewidth=LINE_WIDTH ,data=CI_90)+scale_color_manual(values=c("grey50","grey80","black"))+scale_linetype_manual(values = c("solid","solid","dotdash"))
+	return(PLOT)
+}
+