library(tidyverse) library(forecast) library(lmtest) source("Scripts/Functions.r") #source("Scripts/Load_Wyoming_Web_Data.r") DF <- FRED_GET('WYLINC3POP','LN_POP') %>% select(-YEAR) TS <- 1000*ts(DF,start=c(1970),end=c(2024),frequency=1) BC <- BoxCox.lambda(TS) MODEL <- auto.arima(TS, lambda = BC) forecast(MODEL,h=20) ARMA_POP <- forecast(MODEL,h=35,fan=TRUE,bootstrap = TRUE,npaths=100000,biasadj=FALSE) ARMA_POP plot(ARMA_POP,main="Lincoln County Population Forecast",xlab="Year",ylab="Population") help("forecast") ####Employment to pop ratio EMP <- FRED_GET('LAUCN560230000000005','EMP') %>% inner_join(FRED_GET('WYLINC3POP','LN_POP')) %>% mutate(LN_POP=1000*LN_POP) 2900/3300 EMP <- EMP %>% mutate(RATIO=LN_POP/EMP) ggplot(aes(x=YEAR,y=RATIO),data=EMP)+geom_line() AVG_POP_RATIO <- mean(EMP$RATIO) SD_POP_RATIO <- sd(EMP$RATIO) #####City level #See data http://eadiv.state.wy.us/pop/ KEM <- c(3273,3523,3688,3667,3626,3637,3611,3388,3156,3040,3020,3029,2989,2959,2976,2963,2910,2807,2729,2690,2657,2608,2575,2561,2574,2579,2603,2640,2679,2692,2642,2597,2551,2575,2578,2554,2544,2499,2457,2435,2413,2445,2445,2404,2378) DIAMOND <- c(1000,1070,1114,1101,1082,1078,1063,991,916,876,864,863,847,835,835,827,808,774,748,732,705,695,690,689,695,700,710,723,738,745,731,704,677,667,652,629,613,586,559,540,523,526,527,521,517) AREA_POP <- KEM+DIAMOND LN <- c(12177,13254,14031,14110,14111,14319,14384,13658,12875,12552,12625,12975,13124,13329,13759,14073,14206,14099,14114,14338,14621,14697,14858,15117,15539,15917,16429,17013,17629,18082,18083,17946,17822,18148,18346,18473,18766,18899,19042,19379,19658,20174,20690,20909,21000) NO_CITY <- c(10095,10392,10747,10944,11043) YEAR <- 1980:2024 DATA <- cbind(YEAR,LN,AREA_POP) %>% as_tibble %>% rename("Lincoln County"=LN,"Kemmerer & Diamondvile"=AREA_POP) #Old daata:Period Ends in 1970 #See in part http://eadiv.state.wy.us/demog_data/cntycity_hist.htm LN_OLD <- c(12487,10894,10286,9023,9018,8640) LN_OLD <- cbind(seq(1920,1970,by=10),LN_OLD) %>% as_tibble colnames(LN_OLD) <- c("YEAR","LN_POP") KEM_OLD <- c(843,1517,1884,2026,1667,2028,2292) KEM_OLD <- cbind(seq(1910,1970,by=10),KEM_OLD) %>% as_tibble colnames(KEM_OLD) <- c("YEAR","KEM_POP") DIAMOND_OLD <- c(696,726,812,586,415,398,485) DIAMOND_OLD <- cbind(seq(1910,1970,by=10),DIAMOND_OLD) %>% as_tibble colnames(DIAMOND_OLD) <- c("YEAR","DIA_POP") OLD_DATA <- inner_join(KEM_OLD,DIAMOND_OLD) %>% full_join(LN_OLD) GRAPH_DATA <- ts(DATA %>% select(-YEAR),start=c(1980),end=c(2024),frequency=1) png("Population.png") plot(GRAPH_DATA ,main="Regional Population Trends",type="b",lwd=4,col="blue") dev.off() lines(GRAPH_DATA,col="blue") ?plot ############## TS <- EMP %>% select(EMP) %>% ts(start=c(1990),end=c(2024),frequency=1) BC <- BoxCox.lambda(TS) MODEL2 <- auto.arima(TS, lambda = BC) MODEL2 ARMA_EMP <- forecast(MODEL2,h=35,fan=TRUE,bootstrap = TRUE,npaths=100000,biasadj=FALSE) plot(ARMA_EMP,main="Lincoln County Employment",xlab="Year",ylab="Employed") ########## DATA <- FRED_GET('BPPRIV056023','PRIV_HOUSING') %>% inner_join(FRED_GET('WYLINC3POP','LN_POP')) %>% inner_join(FRED_GET('ATNHPIUS56023A','HOUSE_PRICE_INDEX')) %>% left_join(FRED_GET('LAUCN560230000000005','EMP')) %>%left_join(FRED_GET('DCOILWTICO','WTI')) %>%left_join(FRED_GET('PCOALAUUSDM','COAL')) %>% mutate(LN_POP=1000*LN_POP) %>% select(-YEAR) %>% ts() %>% log() %>% diff() grangertest(PRIV_HOUSING~HOUSE_PRICE_INDEX,data=DATA,order=1) grangertest(HOUSE_PRICE_INDEX~PRIV_HOUSING,data=DATA,order=2) grangertest(LN_POP~PRIV_HOUSING,data=DATA,order=1) grangertest(HOUSE_PRICE_INDEX~COAL,data=DATA,order=2) grangertest(LN_POP~WTI,data=DATA,order=1) grangertest(COAL~WTI,data=DATA,order=2) grangertest(WTI~COAL,data=DATA,order=2) grangertest(LN_POP~COAL,data=DATA,order=1) grangertest(EMP~LN_POP,data=DATA,order=2) grangertest(PRIV_HOUSING~LN_POP,data=DATA,order=2) ## Signficant grangertest(EMP~LN_POP,data=DATA,order=2) ## Signficant grangertest(EMP~PRIV_HOUSING,data=DATA,order=2) ##Signficant grangertest(EMP~LN_POP,data=DATA,order=2) ##Signficant grangertest(PRIV_HOUSING~COAL,data=DATA,order=1) ##Signficant grangertest(PRIV_HOUSING~COAL,data=DATA,order=2) ##Signficant grangertest(PRIV_HOUSING~COAL,data=DATA,order=3) ##Signficant grangertest(EMP~COAL,data=DATA,order=1) ##Signficant grangertest(LN_POP~EMP,data=DATA,order=2) grangertest(PRIV_HOUSING~EMP,data=DATA,order=2) plot(DATA) library(BVAR) MOD <- bvar(DATA,lags=1, n_draw=40000) opt_irf <- bv_irf(horizon = 25, identification = TRUE) plot(irf(MOD,opt_irf,conf_bands = c(0.05, 0.1,0.15)),area=TRUE) head(DATA) DATA # ATNHPIUS56023A Housing price index # Median income MHIWY56023A052NCEN # Employed Persons in Lincoln LAUCN560230000000005 #####Plan and ideas #1) Review IMPLAN for industry multipliers #2) Review IMPLAN for employment to population multipliers (imparted) #3) Find a list of all planned new projects #4) Use the IMPLAN multipliers for each sector to estimate total change #5) Develop survey to estimate likelihood of new projects #6) Compare to the ARMA percentile #7) Adjust the ARMA up assuming some of these outputs are known. #8) Occupancy rate from IMPLAN as a housing cap when projecting #9) Model housing construciton rate (Maybe) #10) Employment rate by age in IMPLAN ####Other ideas, develop larger plan? Maybe look at decline in other industries as a proportion of employment ###Seperate out Kemmer and Diamondville? http://eadiv.state.wy.us/pop/wyc&sc40.htm