From e1fc4206b4ae7bf9432460f53a3a51a899e6f1ff Mon Sep 17 00:00:00 2001 From: Alex Gebben Work Date: Tue, 9 Dec 2025 13:10:49 -0700 Subject: [PATCH] Creatd histogram and capacity tables --- 2B_Result_Analysis.r | 31 ++++++++++++++++++++++++++++--- Prelim_Process.sh | 1 + 2 files changed, 29 insertions(+), 3 deletions(-) diff --git a/2B_Result_Analysis.r b/2B_Result_Analysis.r index e502462..cb5735f 100644 --- a/2B_Result_Analysis.r +++ b/2B_Result_Analysis.r @@ -1,4 +1,5 @@ library(tidyverse) +library(gt) #For nice color coded capacity limits table. ###Process the simulations and save the main percentile results by year RES <- read_csv("Results/Simulations/Kemmerer_2024_Simulation.csv") RES[,"Year"] <- RES[,"Year"] @@ -41,16 +42,40 @@ POP_PLOT BIRTH_DATA <- GET_DATA(RES,4) BIRTH_PLOT <- MAKE_GRAPH(BIRTH_DATA) -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 = 10),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()) +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()) BIRTH_PLOT DEATH_DATA <- GET_DATA(RES,5) %>% filter(!is.na(MIN)) DEATH_PLOT <- MAKE_GRAPH(DEATH_DATA) -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 = 10),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()) +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()) DEATH_PLOT MIGRATION_DATA <- GET_DATA(RES,6) %>% filter(!is.na(MIN)) MIGRATION_PLOT <- MAKE_GRAPH(MIGRATION_DATA) -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 = 10),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()) +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()) MIGRATION_PLOT +#####Key year table + AVG_VALUES <- KEY %>% group_by(Year) %>% summarize(MED=median(Population),MEAN=mean(Population)) +AVG_VALUES <- rbind(AVG_VALUES[,1:2]%>% rename(Value=MED) %>% mutate('Summary Stat.'="Median"),AVG_VALUES[,c(1,3)] %>% rename(Value=MEAN) %>% mutate('Summary Stat.'="Mean")) + +HISTOGRAM <- ggplot(KEY, aes(x = Population,group=-Year,Color=Year,fill=Year)) + geom_histogram(alpha=0.3,bins=800)+geom_vline(data = AVG_VALUES, aes(xintercept = Value,group=`Summary Stat.`,color = `Summary Stat.`), size = 0.75)+scale_color_manual(values=c("red","black","black"))+ facet_grid(rows=vars(Year))+ scale_x_continuous(breaks = c(seq(0, 10000, by = 500)))+ theme_bw()+ theme(legend.position = "top",panel.grid.minor = element_blank())+ylab("Number of Simulation")+guides(fill= guide_legend(nrow = 1)) +HISTOGRAM +#rm(KEY_YEARS) +POP_LEVELS <- seq(2000,6000,100) +YEARS <- c(2030,2035,2045,2055,2065) +for(i in YEARS ){ + KEY <- RES %>% filter(Year==i ) %>% pull(Population) + ECDF <- ecdf(KEY) + ECDF_VALUES <- ECDF(POP_LEVELS) + if(!exists("KEY_YEARS")){KEY_YEARS<- ECDF_VALUES} else{KEY_YEARS<- cbind(KEY_YEARS,ECDF_VALUES)} + +} +colnames(KEY_YEARS) <- YEARS +rownames(KEY_YEARS) <- POP_LEVELS +PLOT_GREEN <- "forestgreen" +PLOT_YELLOW <- "yellow" +PLOT_RED <- "red" +KEY_YEARS <- KEY_YEARS%>% as.data.frame +Capacity_Risk <- KEY_YEARS%>% gt(rownames_to_stub = TRUE,caption="Year") %>% data_color( fn = scales::col_numeric( palette = c(PLOT_RED, PLOT_YELLOW, PLOT_GREEN), domain = c(0, 1) ) ) %>% fmt_percent( decimals = 1, drop_trailing_zeros = FALSE) %>% tab_stubhead(label =c("Capacity")) +gtsave( data = Capacity_Risk , filename = "my_table.html") diff --git a/Prelim_Process.sh b/Prelim_Process.sh index 05cf832..e52952c 100644 --- a/Prelim_Process.sh +++ b/Prelim_Process.sh @@ -18,3 +18,4 @@ Rscript "./Scripts/2G_Single_Age_Sex_ARIMA_Models.r" Rscript "./Scripts/2E_Move_Current_Demographic_Data_to_Current_Year.r" #Produce final results for either the simulation, or information for the report, but not anything used in later stages of the simulation Rscript "./Scripts/3A_Population_Pyramid.r" +#wkhtmltopdf --disable-smart-shrinking --no-stop-slow-scripts --enable-local-file-access --page-width 85mm --page-height 328mm my_table.html output.pdf