Librerias y paquetes

library(maps)
library(forecast)
library(tidyverse)
library(maps)
library(ggplot2)
library(shiny)

Obtención de base de datos

Historico de poblaciones

pop <- read.csv("historical_state_population_by_year.csv", header = FALSE)
pop <- rename(pop, "State" = "V1", "Year" = "V2", "Population" = "V3")
str(pop)
## 'data.frame':    6020 obs. of  3 variables:
##  $ State     : chr  "AK" "AK" "AK" "AK" ...
##  $ Year      : int  1950 1951 1952 1953 1954 1955 1956 1957 1958 1959 ...
##  $ Population: int  135000 158000 189000 205000 215000 222000 224000 231000 224000 224000 ...

Relacion nombre del estado con codico de 2 letras

state_name <- read_csv("us_state_two_letter_abbreviation.csv")
## Rows: 59 Columns: 2
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (2): U.S. States and Territories, Two-Letter Abbreviation
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
state_name$`U.S. States and Territories`= str_to_lower(state_name$`U.S. States and Territories`)
state_name <- rename(state_name, "states" = "U.S. States and Territories")
state_name <- rename(state_name, "State" = "Two-Letter Abbreviation")

Coordenadas

state_data <- ggplot2::map_data('state')
state_data <- fortify(state_data)
head(state_data)
##        long      lat group order  region subregion
## 1 -87.46201 30.38968     1     1 alabama      <NA>
## 2 -87.48493 30.37249     1     2 alabama      <NA>
## 3 -87.52503 30.37249     1     3 alabama      <NA>
## 4 -87.53076 30.33239     1     4 alabama      <NA>
## 5 -87.57087 30.32665     1     5 alabama      <NA>
## 6 -87.58806 30.32665     1     6 alabama      <NA>

Columna de años

years <- seq(2020, 2070, by = 1)

Nuevo México

pop_nm <- filter(pop, pop$State == "NM")
ts_nm <- ts(data = pop_nm$Population, start = c(1900,1), frequency = 1)
ts_nm
## Time Series:
## Start = 1900 
## End = 2019 
## Frequency = 1 
##   [1]  196000  206000  218000  229000  242000  255000  269000  283000  298000
##  [10]  314000  329000  334000  338000  344000  348000  351000  354000  361000
##  [19]  382000  362000  363000  370000  375000  382000  390000  396000  403000
##  [28]  410000  416000  420000  427000  436000  441000  449000  461000  475000
##  [37]  489000  503000  513000  523000  531000  506000  502000  534000  527000
##  [46]  537000  561000  582000  604000  644000  689000  717000  735000  756000
##  [55]  763000  785000  806000  847000  886000  919000  954000  965000  979000
##  [64]  989000 1006000 1012000 1007000 1000000  994000 1011000 1017055 1053737
##  [73] 1078697 1105529 1131309 1159944 1189295 1215720 1238034 1284722 1309400
##  [82] 1332748 1363823 1394361 1416717 1438361 1462729 1478520 1490337 1503901
##  [91] 1519933 1547115 1580750 1614937 1653329 1682417 1706151 1722939 1733535
## [100] 1739844 1821204 1831690 1855309 1877574 1903808 1932274 1962137 1990070
## [109] 2010662 2036802 2064552 2080450 2087309 2092273 2089568 2089291 2091630
## [118] 2091784 2092741 2096829

Modelo ARIMA

arima_nm <- auto.arima(ts_nm)
arima_nm
## Series: ts_nm 
## ARIMA(0,2,1) 
## 
## Coefficients:
##           ma1
##       -0.6390
## s.e.   0.0852
## 
## sigma^2 = 152075089:  log likelihood = -1278.75
## AIC=2561.5   AICc=2561.6   BIC=2567.04
summary(arima_nm)
## Series: ts_nm 
## ARIMA(0,2,1) 
## 
## Coefficients:
##           ma1
##       -0.6390
## s.e.   0.0852
## 
## sigma^2 = 152075089:  log likelihood = -1278.75
## AIC=2561.5   AICc=2561.6   BIC=2567.04
## 
## Training set error measures:
##                     ME     RMSE      MAE       MPE      MAPE     MASE
## Training set -183.5503 12176.75 7780.673 0.0198346 0.9789513 0.450605
##                    ACF1
## Training set 0.01632203

Generar pronóstico

pronostico_nm <- forecast(arima_nm, level = c(95), h = 51)
pronostico_nm
##      Point Forecast     Lo 95   Hi 95
## 2020        2099505 2075334.9 2123675
## 2021        2102181 2061360.1 2143002
## 2022        2104857 2046558.4 2163155
## 2023        2107533 2030502.9 2184562
## 2024        2110209 2013135.1 2207282
## 2025        2112884 1994478.4 2231291
## 2026        2115560 1974577.8 2256543
## 2027        2118236 1953482.5 2282990
## 2028        2120912 1931240.2 2310584
## 2029        2123588 1907895.5 2339281
## 2030        2126264 1883489.1 2369039
## 2031        2128940 1858058.5 2399821
## 2032        2131616 1831637.7 2431594
## 2033        2134292 1804258.0 2464326
## 2034        2136968 1775947.8 2497988
## 2035        2139644 1746733.5 2532554
## 2036        2142320 1716639.6 2568000
## 2037        2144995 1685688.4 2604303
## 2038        2147671 1653900.9 2641442
## 2039        2150347 1621296.6 2679398
## 2040        2153023 1587893.9 2718153
## 2041        2155699 1553709.7 2757689
## 2042        2158375 1518760.1 2797990
## 2043        2161051 1483060.2 2839042
## 2044        2163727 1446624.3 2880829
## 2045        2166403 1409466.0 2923340
## 2046        2169079 1371598.0 2966559
## 2047        2171755 1333032.4 3010477
## 2048        2174431 1293780.7 3055080
## 2049        2177106 1253854.1 3100359
## 2050        2179782 1213262.9 3146302
## 2051        2182458 1172017.2 3192899
## 2052        2185134 1130126.5 3240142
## 2053        2187810 1087600.1 3288020
## 2054        2190486 1044446.6 3336525
## 2055        2193162 1000674.6 3385649
## 2056        2195838  956292.2 3435384
## 2057        2198514  911307.1 3485720
## 2058        2201190  865726.9 3536653
## 2059        2203866  819558.8 3588172
## 2060        2206542  772809.7 3640273
## 2061        2209217  725486.5 3692948
## 2062        2211893  677595.6 3746191
## 2063        2214569  629143.3 3799995
## 2064        2217245  580135.7 3854355
## 2065        2219921  530578.8 3909263
## 2066        2222597  480478.1 3964716
## 2067        2225273  429839.3 4020707
## 2068        2227949  378667.6 4077230
## 2069        2230625  326968.5 4134281
## 2070        2233301  274746.8 4191855
plot(pronostico_nm)

pronostico_nm_dec <- data.frame(pronostico_nm)
pronostico_nm_dec <- pronostico_nm_dec %>% 
                        add_column(Decade = years) %>% 
                        filter(row_number() %in% c(11, 21, 31, 41, 51)) %>% 
                        add_column(State = "NM") %>%
                        select(-Lo.95, -Hi.95)
pronostico_nm_dec
##      Point.Forecast Decade State
## 2030        2126264   2030    NM
## 2040        2153023   2040    NM
## 2050        2179782   2050    NM
## 2060        2206542   2060    NM
## 2070        2233301   2070    NM

Arizona

pop_az <- filter(pop, pop$State == "AZ")
ts_az <- ts(data = pop_az$Population, start = c(1900,1), frequency = 1)
ts_az
## Time Series:
## Start = 1900 
## End = 2019 
## Frequency = 1 
##   [1]  124000  131000  138000  144000  151000  158000  167000  176000  186000
##  [10]  196000  206000  212000  217000  236000  253000  263000  282000  311000
##  [19]  320000  329000  340000  351000  360000  371000  382000  393000  403000
##  [28]  414000  422000  430000  434000  429000  426000  426000  428000  434000
##  [37]  443000  453000  466000  484000  499000  490000  524000  692000  610000
##  [46]  594000  616000  653000  690000  714000  756000  785000  842000  894000
##  [55]  933000  987000 1053000 1125000 1193000 1261000 1321000 1407000 1471000
##  [64] 1521000 1556000 1584000 1614000 1646000 1682000 1737000 1775399 1895814
##  [73] 2008291 2124438 2223196 2284847 2346157 2425197 2515316 2635571 2737774
##  [82] 2810107 2889861 2968925 3067135 3183538 3308262 3437103 3535183 3622185
##  [91] 3679056 3762394 3867333 3993390 4147561 4306908 4432308 4552207 4667277
## [100] 4778332 5160586 5273477 5396255 5510364 5652404 5839077 6029141 6167681
## [109] 6280362 6343154 6407172 6472643 6554978 6632764 6730413 6829676 6941072
## [118] 7044008 7158024 7278717

Modelo ARIMA

arima_az <- auto.arima(ts_az)
arima_az
## Series: ts_az 
## ARIMA(0,2,1) 
## 
## Coefficients:
##           ma1
##       -0.7815
## s.e.   0.0620
## 
## sigma^2 = 1.504e+09:  log likelihood = -1414.15
## AIC=2832.3   AICc=2832.4   BIC=2837.84
summary(arima_az)
## Series: ts_az 
## ARIMA(0,2,1) 
## 
## Coefficients:
##           ma1
##       -0.7815
## s.e.   0.0620
## 
## sigma^2 = 1.504e+09:  log likelihood = -1414.15
## AIC=2832.3   AICc=2832.4   BIC=2837.84
## 
## Training set error measures:
##                    ME     RMSE      MAE       MPE     MAPE      MASE       ACF1
## Training set 3812.388 38291.02 20396.69 0.3291156 1.521771 0.3286796 0.04665136

Generar pronóstico

pronostico_az <- forecast(arima_az, level = c(95), h = 51)
pronostico_az
##      Point Forecast   Lo 95    Hi 95
## 2020        7385734 7309728  7461739
## 2021        7492750 7372942  7612559
## 2022        7599767 7437644  7761889
## 2023        7706783 7501558  7912008
## 2024        7813800 7563987  8063613
## 2025        7920817 7624669  8216964
## 2026        8027833 7683509  8372157
## 2027        8134850 7740482  8529217
## 2028        8241866 7795599  8688134
## 2029        8348883 7848885  8848881
## 2030        8455900 7900375  9011424
## 2031        8562916 7950108  9175724
## 2032        8669933 7998122  9341744
## 2033        8776949 8044455  9509443
## 2034        8883966 8089147  9678785
## 2035        8990982 8132232  9849733
## 2036        9097999 8173747 10022251
## 2037        9205016 8213723 10196308
## 2038        9312032 8252193 10371871
## 2039        9419049 8289187 10548911
## 2040        9526065 8324733 10727398
## 2041        9633082 8358858 10907306
## 2042        9740099 8391588 11088609
## 2043        9847115 8422948 11271282
## 2044        9954132 8452961 11455302
## 2045       10061148 8481650 11640647
## 2046       10168165 8509035 11827295
## 2047       10275182 8535138 12015226
## 2048       10382198 8559977 12204420
## 2049       10489215 8583571 12394858
## 2050       10596231 8605939 12586523
## 2051       10703248 8627098 12779398
## 2052       10810265 8647063 12973466
## 2053       10917281 8665851 13168711
## 2054       11024298 8683477 13365118
## 2055       11131314 8699956 13562673
## 2056       11238331 8715302 13761360
## 2057       11345348 8729528 13961167
## 2058       11452364 8742649 14162079
## 2059       11559381 8754676 14364086
## 2060       11666397 8765621 14567173
## 2061       11773414 8775498 14771330
## 2062       11880431 8784317 14976544
## 2063       11987447 8792089 15182805
## 2064       12094464 8798826 15390101
## 2065       12201480 8804538 15598423
## 2066       12308497 8809234 15807759
## 2067       12415513 8812926 16018101
## 2068       12522530 8815622 16229438
## 2069       12629547 8817333 16441761
## 2070       12736563 8818066 16655061
plot(pronostico_az)

pronostico_az_dec <- data.frame(pronostico_az)
pronostico_az_dec <- pronostico_az_dec %>% 
                        add_column(Decade = years) %>% 
                        filter(row_number() %in% c(11, 21, 31, 41, 51)) %>% 
                        add_column(State = "AZ") %>%
                        select(-Lo.95, -Hi.95)
pronostico_az_dec
##      Point.Forecast Decade State
## 2030        8455900   2030    AZ
## 2040        9526065   2040    AZ
## 2050       10596231   2050    AZ
## 2060       11666397   2060    AZ
## 2070       12736563   2070    AZ

California

pop_ca <- filter(pop, pop$State == "CA")
ts_ca <- ts(data = pop_ca$Population, start = c(1900,1), frequency = 1)
ts_ca
## Time Series:
## Start = 1900 
## End = 2019 
## Frequency = 1 
##   [1]  1490000  1550000  1623000  1702000  1792000  1893000  1976000  2054000
##   [9]  2161000  2282000  2406000  2534000  2668000  2811000  2934000  3008000
##  [17]  3071000  3171000  3262000  3339000  3554000  3795000  3991000  4270000
##  [25]  4541000  4730000  4929000  5147000  5344000  5531000  5711000  5824000
##  [33]  5894000  5963000  6060000  6175000  6341000  6528000  6656000  6785000
##  [41]  6950000  7237000  7735000  8506000  8945000  9344000  9559000  9832000
##  [49] 10064000 10337000 10677000 11134000 11635000 12251000 12746000 13133000
##  [57] 13713000 14264000 14880000 15467000 15870000 16497000 17072000 17668000
##  [65] 18151000 18585000 18858000 19176000 19394000 19711000 19971069 20345939
##  [73] 20585469 20868728 21173865 21537849 21935909 22352396 22835958 23256880
##  [81] 23800800 24285933 24820009 25360026 25844393 26441109 27102237 27777158
##  [89] 28464249 29218164 29950111 30414114 30875920 31147208 31317179 31493525
##  [97] 31780829 32217708 32682794 33145121 33987977 34479458 34871843 35253159
## [105] 35574576 35827943 36021202 36250311 36604337 36961229 37319502 37638369
## [113] 37948800 38260787 38596972 38918045 39167117 39358497 39461588 39512223

Modelo ARIMA

arima_ca <- auto.arima(ts_ca)
arima_ca
## Series: ts_ca 
## ARIMA(0,2,0) 
## 
## sigma^2 = 9.971e+09:  log likelihood = -1525.79
## AIC=3053.57   AICc=3053.61   BIC=3056.34
summary(arima_ca)
## Series: ts_ca 
## ARIMA(0,2,0) 
## 
## sigma^2 = 9.971e+09:  log likelihood = -1525.79
## AIC=3053.57   AICc=3053.61   BIC=3056.34
## 
## Training set error measures:
##                     ME     RMSE      MAE        MPE     MAPE      MASE
## Training set -88.02939 99017.45 65793.72 0.06388803 0.560492 0.2059178
##                    ACF1
## Training set -0.1174913

Generar pronóstico

pronostico_ca <- forecast(arima_ca, level = c(95), h = 51)
pronostico_ca
##      Point Forecast      Lo 95    Hi 95
## 2020       39562858 39367149.6 39758566
## 2021       39613493 39175875.8 40051110
## 2022       39664128 38931854.3 40396402
## 2023       39714763 38642824.0 40786702
## 2024       39765398 38313985.8 41216810
## 2025       39816033 37949094.0 41682972
## 2026       39866668 37551015.2 42182321
## 2027       39917303 37122028.2 42712578
## 2028       39967938 36664000.2 43271876
## 2029       40018573 36178497.2 43858649
## 2030       40069208 35666856.8 44471559
## 2031       40119843 35130238.7 45109447
## 2032       40170478 34569661.1 45771295
## 2033       40221113 33986027.0 46456199
## 2034       40271748 33380144.7 47163351
## 2035       40322383 32752743.0 47892023
## 2036       40373018 32104483.6 48641552
## 2037       40423653 31435970.6 49411335
## 2038       40474288 30747758.6 50200817
## 2039       40524923 30040359.0 51009487
## 2040       40575558 29314245.3 51836871
## 2041       40626193 28569857.7 52682528
## 2042       40676828 27807606.7 53546049
## 2043       40727463 27027876.2 54427050
## 2044       40778098 26231026.6 55325169
## 2045       40828733 25417396.8 56240069
## 2046       40879368 24587306.3 57171430
## 2047       40930003 23741057.1 58118949
## 2048       40980638 22878935.0 59082341
## 2049       41031273 22001211.4 60061335
## 2050       41081908 21108144.0 61055672
## 2051       41132543 20199978.4 62065108
## 2052       41183178 19276948.5 63089407
## 2053       41233813 18339277.9 64128348
## 2054       41284448 17387180.2 65181716
## 2055       41335083 16420860.1 66249306
## 2056       41385718 15440513.7 67330922
## 2057       41436353 14446329.2 68426377
## 2058       41486988 13438487.6 69535488
## 2059       41537623 12417162.7 70658083
## 2060       41588258 11382522.0 71793994
## 2061       41638893 10334727.0 72943059
## 2062       41689528  9273933.2 74105123
## 2063       41740163  8200291.0 75280035
## 2064       41790798  7113945.5 76467651
## 2065       41841433  6015037.0 77667829
## 2066       41892068  4903701.4 78880435
## 2067       41942703  3780070.2 80105336
## 2068       41993338  2644270.7 81342405
## 2069       42043973  1496426.6 82591519
## 2070       42094608   336657.5 83852558
plot(pronostico_ca)

pronostico_ca_dec <- data.frame(pronostico_ca)
pronostico_ca_dec <- pronostico_ca_dec %>% 
                        add_column(Decade = years) %>% 
                        filter(row_number() %in% c(11, 21, 31, 41, 51)) %>% 
                        add_column(State = "CA") %>%
                        select(-Lo.95, -Hi.95)

Alabama

pop_al <- filter(pop, pop$State == "AL")
ts_al <- ts(data = pop_al$Population, start = c(1900,1), frequency = 1)
ts_al
## Time Series:
## Start = 1900 
## End = 2019 
## Frequency = 1 
##   [1] 1830000 1907000 1935000 1957000 1978000 2012000 2045000 2058000 2070000
##  [10] 2108000 2150000 2185000 2217000 2279000 2336000 2345000 2359000 2361000
##  [19] 2343000 2337000 2359000 2402000 2434000 2464000 2501000 2534000 2567000
##  [28] 2609000 2640000 2644000 2647000 2649000 2653000 2661000 2685000 2719000
##  [37] 2743000 2762000 2787000 2814000 2845000 2902000 2941000 2902000 2802000
##  [46] 2775000 2911000 2942000 2969000 3000000 3058000 3059000 3068000 3053000
##  [55] 3014000 3050000 3071000 3109000 3163000 3204000 3274000 3316000 3323000
##  [64] 3358000 3395000 3443000 3464000 3458000 3446000 3440000 3444354 3497076
##  [73] 3539400 3579780 3626499 3678814 3735139 3780403 3831836 3866248 3900368
##  [82] 3918531 3925266 3934102 3951820 3972523 3991569 4015264 4023844 4030222
##  [91] 4048508 4091025 4139269 4193114 4232965 4262731 4290403 4320281 4351037
## [100] 4369862 4452173 4467634 4480089 4503491 4530729 4569805 4628981 4672840
## [109] 4718206 4757938 4785437 4799069 4815588 4830081 4841799 4852347 4863525
## [118] 4874486 4887681 4903185

Modelo ARIMA

arima_al <- auto.arima(ts_al)
arima_al
## Series: ts_al 
## ARIMA(0,1,1) with drift 
## 
## Coefficients:
##          ma1     drift
##       0.4763  25997.59
## s.e.  0.0763   3156.13
## 
## sigma^2 = 556157014:  log likelihood = -1366.1
## AIC=2738.19   AICc=2738.4   BIC=2746.53
summary(arima_al)
## Series: ts_al 
## ARIMA(0,1,1) with drift 
## 
## Coefficients:
##          ma1     drift
##       0.4763  25997.59
## s.e.  0.0763   3156.13
## 
## sigma^2 = 556157014:  log likelihood = -1366.1
## AIC=2738.19   AICc=2738.4   BIC=2746.53
## 
## Training set error measures:
##                     ME     RMSE      MAE          MPE      MAPE      MASE
## Training set -133.9332 23286.33 15525.23 -0.006537242 0.5080397 0.5118891
##                     ACF1
## Training set 0.009996564

Generar pronóstico

pronostico_al <- forecast(arima_al, level = c(95), h = 51)
pronostico_al
##      Point Forecast   Lo 95   Hi 95
## 2020        4925961 4879739 4972183
## 2021        4951959 4869540 5034377
## 2022        4977956 4870955 5084957
## 2023        5003954 4877046 5130862
## 2024        5029951 4885861 5174042
## 2025        5055949 4896517 5215381
## 2026        5081947 4908525 5255368
## 2027        5107944 4921581 5294307
## 2028        5133942 4935478 5332405
## 2029        5159939 4950073 5369806
## 2030        5185937 4965255 5406619
## 2031        5211935 4980943 5442926
## 2032        5237932 4997073 5478792
## 2033        5263930 5013591 5514269
## 2034        5289927 5030455 5549400
## 2035        5315925 5047629 5584220
## 2036        5341923 5065085 5618760
## 2037        5367920 5082797 5653043
## 2038        5393918 5100742 5687093
## 2039        5419915 5118903 5720927
## 2040        5445913 5137263 5754562
## 2041        5471910 5155808 5788013
## 2042        5497908 5174524 5821292
## 2043        5523906 5193400 5854411
## 2044        5549903 5212427 5887379
## 2045        5575901 5231595 5920207
## 2046        5601898 5250896 5952901
## 2047        5627896 5270322 5985470
## 2048        5653894 5289867 6017921
## 2049        5679891 5309524 6050259
## 2050        5705889 5329288 6082490
## 2051        5731886 5349153 6114620
## 2052        5757884 5369115 6146653
## 2053        5783882 5389169 6178594
## 2054        5809879 5409312 6210447
## 2055        5835877 5429539 6242215
## 2056        5861874 5449847 6273902
## 2057        5887872 5470232 6305512
## 2058        5913870 5490691 6337048
## 2059        5939867 5511223 6368512
## 2060        5965865 5531823 6399907
## 2061        5991862 5552489 6431236
## 2062        6017860 5573219 6462501
## 2063        6043858 5594011 6493704
## 2064        6069855 5614863 6524848
## 2065        6095853 5635772 6555934
## 2066        6121850 5656736 6586964
## 2067        6147848 5677755 6617941
## 2068        6173846 5698826 6648865
## 2069        6199843 5719947 6679739
## 2070        6225841 5741117 6710564
plot(pronostico_al)

pronostico_al_dec <- data.frame(pronostico_al)
pronostico_al_dec <- pronostico_al_dec %>% 
                        add_column(Decade = years) %>% 
                        filter(row_number() %in% c(11, 21, 31, 41, 51)) %>% 
                        add_column(State = "AL") %>%
                        select(-Lo.95, -Hi.95)
pronostico_al_dec
##      Point.Forecast Decade State
## 2030        5185937   2030    AL
## 2040        5445913   2040    AL
## 2050        5705889   2050    AL
## 2060        5965865   2060    AL
## 2070        6225841   2070    AL

Nueva York

pop_ny <- filter(pop, pop$State == "NY")
ts_ny <- ts(data = pop_ny$Population, start = c(1900,1), frequency = 1)
ts_ny
## Time Series:
## Start = 1900 
## End = 2019 
## Frequency = 1 
##   [1]  7283000  7449000  7612000  7771000  7927000  8084000  8289000  8499000
##   [9]  8714000  8935000  9137000  9249000  9361000  9473000  9585000  9700000
##  [17]  9848000  9993000  9936000 10252000 10282000 10416000 10589000 10752000
##  [25] 10953000 11186000 11257000 11174000 11599000 12171000 12647000 12848000
##  [33] 13001000 13126000 13253000 13375000 13481000 13511000 13512000 13523000
##  [41] 13456000 13267000 13002000 12807000 12628000 12495000 13398000 13982000
##  [49] 14497000 14892000 14865000 14890000 15192000 15527000 15814000 15966000
##  [57] 16112000 16374000 16601000 16685000 16838000 17061000 17301000 17461000
##  [65] 17589000 17734000 17843000 17935000 18051000 18105000 18241391 18357982
##  [73] 18339400 18177063 18049775 18003485 17940541 17812602 17680589 17583838
##  [81] 17566754 17567734 17589738 17686905 17745684 17791672 17833419 17868848
##  [89] 17941309 17983086 18002855 18029532 18082032 18140894 18156652 18150928
##  [97] 18143805 18143184 18159175 18196601 19001780 19082838 19137800 19175939
## [105] 19171567 19132610 19104631 19132335 19212436 19307066 19399878 19499241
## [113] 19572932 19624447 19651049 19654666 19633428 19589572 19530351 19453561

Modelo ARIMA

arima_ny <- auto.arima(ts_ny)
arima_ny
## Series: ts_ny 
## ARIMA(1,2,1) 
## 
## Coefficients:
##          ar1      ma1
##       0.5065  -0.9697
## s.e.  0.0855   0.0271
## 
## sigma^2 = 2.326e+10:  log likelihood = -1575.62
## AIC=3157.25   AICc=3157.46   BIC=3165.56
summary(arima_ny)
## Series: ts_ny 
## ARIMA(1,2,1) 
## 
## Coefficients:
##          ar1      ma1
##       0.5065  -0.9697
## s.e.  0.0855   0.0271
## 
## sigma^2 = 2.326e+10:  log likelihood = -1575.62
## AIC=3157.25   AICc=3157.46   BIC=3165.56
## 
## Training set error measures:
##                     ME     RMSE      MAE         MPE      MAPE      MASE
## Training set -14563.74 149933.4 82036.53 -0.08319142 0.5752184 0.5840548
##                      ACF1
## Training set -0.009336835

Generar pronóstico

pronostico_ny <- forecast(arima_ny, level = c(95), h = 51)
pronostico_ny
##      Point Forecast    Lo 95    Hi 95
## 2020       19441372 19142479 19740265
## 2021       19461901 18913884 20009918
## 2022       19499001 18722834 20275168
## 2023       19544494 18560730 20528257
## 2024       19594237 18419928 20768546
## 2025       19646134 18294652 20997616
## 2026       19699121 18180758 21217483
## 2027       19752660 18075313 21430006
## 2028       19806478 17976225 21636732
## 2029       19860439 17881980 21838897
## 2030       19914471 17791466 22037476
## 2031       19968539 17703845 22233234
## 2032       20022626 17618478 22426774
## 2033       20076722 17534869 22618576
## 2034       20130823 17452622 22809025
## 2035       20184927 17371422 22998431
## 2036       20239031 17291010 23187052
## 2037       20293136 17211176 23375097
## 2038       20347242 17131741 23562742
## 2039       20401347 17052559 23750136
## 2040       20455453 16973503 23937403
## 2041       20509559 16894467 24124650
## 2042       20563664 16815359 24311970
## 2043       20617770 16736100 24499440
## 2044       20671876 16656621 24687131
## 2045       20725982 16576862 24875101
## 2046       20780087 16496771 25063404
## 2047       20834193 16416303 25252084
## 2048       20888299 16335416 25441182
## 2049       20942405 16254076 25630733
## 2050       20996510 16172251 25820770
## 2051       21050616 16089912 26011320
## 2052       21104722 16007036 26202407
## 2053       21158828 15923601 26394055
## 2054       21212933 15839586 26586281
## 2055       21267039 15754974 26779104
## 2056       21321145 15669750 26972540
## 2057       21375251 15583900 27166601
## 2058       21429356 15497412 27361301
## 2059       21483462 15410274 27556650
## 2060       21537568 15322477 27752659
## 2061       21591674 15234012 27949335
## 2062       21645780 15144871 28146688
## 2063       21699885 15055048 28344722
## 2064       21753991 14964536 28543446
## 2065       21808097 14873330 28742864
## 2066       21862203 14781425 28942980
## 2067       21916308 14688817 29143799
## 2068       21970414 14595503 29345325
## 2069       22024520 14501479 29547561
## 2070       22078626 14406742 29750509
plot(pronostico_ny)

pronostico_ny_dec <- data.frame(pronostico_ny)
pronostico_ny_dec <- pronostico_ny_dec %>% 
                        add_column(Decade = years) %>% 
                        filter(row_number() %in% c(11, 21, 31, 41, 51)) %>% 
                        add_column(State = "NY") %>%
                        select(-Lo.95, -Hi.95)
pronostico_ny_dec
##      Point.Forecast Decade State
## 2030       19914471   2030    NY
## 2040       20455453   2040    NY
## 2050       20996510   2050    NY
## 2060       21537568   2060    NY
## 2070       22078626   2070    NY

Mapa

Función Mapa

map_pop <- function(df, state_data, decade = 2030){
  
  df_decade <- df %>% 
    filter(df$Decade == decade)
  
  # Function for setting the aesthetics of the plot
  my_theme <- function () { 
    theme_bw() + theme(axis.text = element_text(size = 14),
                       axis.title = element_text(size = 14),
                       strip.text = element_text(size = 14),
                       panel.grid.major = element_blank(), 
                       panel.grid.minor = element_blank(),
                       panel.background = element_blank(), 
                       legend.position = "bottom",
                       panel.border = element_blank(), 
                       strip.background = element_rect(fill = 'white', colour = 'white'))
  }
  
  
  # Add the data the user wants to see to the geographical world data
  state_data['pop'] <- df_decade$Point.Forecast[match(state_data$region, df_decade$states)]

  
  # Specify the plot for the world map
  library(RColorBrewer)
  g <- ggplot() + 
    geom_polygon(data = state_data, size = 1,
                                    aes(x = long, y = lat, fill = pop, group = group)) + 
    scale_fill_gradientn(colours = brewer.pal(5, "RdBu"), na.value = 'white',limits = c(580479, 49816250)) + 
    my_theme()
  
  return(g)
}

Dataset Mapas

df_map <- rbind(pronostico_al_dec,pronostico_az_dec)
df_map <- rbind(df_map,pronostico_ca_dec)
df_map <- rbind(df_map,pronostico_nm_dec)
df_map <- rbind(df_map,pronostico_ny_dec)

df_map
##       Point.Forecast Decade State
## 2030         5185937   2030    AL
## 2040         5445913   2040    AL
## 2050         5705889   2050    AL
## 2060         5965865   2060    AL
## 2070         6225841   2070    AL
## 20301        8455900   2030    AZ
## 20401        9526065   2040    AZ
## 20501       10596231   2050    AZ
## 20601       11666397   2060    AZ
## 20701       12736563   2070    AZ
## 20302       40069208   2030    CA
## 20402       40575558   2040    CA
## 20502       41081908   2050    CA
## 20602       41588258   2060    CA
## 20702       42094608   2070    CA
## 20303        2126264   2030    NM
## 20403        2153023   2040    NM
## 20503        2179782   2050    NM
## 20603        2206542   2060    NM
## 20703        2233301   2070    NM
## 20304       19914471   2030    NY
## 20404       20455453   2040    NY
## 20504       20996510   2050    NY
## 20604       21537568   2060    NY
## 20704       22078626   2070    NY
df_map <- left_join(df_map,state_name,by = "State")

Mapa 2030

map_pop(df_map, state_data)
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

map_pop(df_map, state_data, decade = 2040)

map_pop(df_map, state_data, decade = 2050)

map_pop(df_map, state_data, decade = 2060)

map_pop(df_map, state_data, decade = 2070)

Aplicación en Shiny

Puede encontrar el app de shiny para generar el forecast cualquiera de los 51 estados, utilizando su código de 2 letras. Igualmente se despliega el mapa de la densidad de la población en USA de a cuerdo a los pronosticos realizados para las siguientes 5 décadas.
App de shiny

Shiny applications not supported in static R Markdown documents
LS0tCnRpdGxlOiAiQWN0aXZpZGFkIEludGVncmFkb3JhIgphdXRob3I6ICJVemllbCBDYXJkZWxhcyAtIEEwMTc0NjA1MCIKZGF0ZTogIjIwMjQtMDItMjgiCm91dHB1dDogCiAgaHRtbF9kb2N1bWVudDoKICAgIHRoZW1lOiBkYXJrbHkKICAgIGhpZ2hsaWdodDogYnJlZXplZGFyawogICAgbnVtYmVyX3NlY3Rpb25zOiBGQUxTRQogICAgdG9jOiBUUlVFCiAgICB0b2NfZGVwdGg6IDQKICAgIHRvY19mbG9hdDoKICAgICAgY29sbGFwc2VkOiBUUlVFCiAgICAgIHNtb290aF9zY3JvbGw6IFRSVUUKICAgIGNvZGVfZG93bmxvYWQ6IFRSVUUKLS0tCgpgYGB7ciBzZXR1cCwgaW5jbHVkZT1GQUxTRX0Ka25pdHI6Om9wdHNfY2h1bmskc2V0KGVjaG8gPSBUUlVFKQpgYGAKCiFbXShodHRwczovL21lZGlhLmxpY2RuLmNvbS9kbXMvaW1hZ2UvRDU2MTJBUUhJSDlXY2dBemhjQS9hcnRpY2xlLWNvdmVyX2ltYWdlLXNocmlua183MjBfMTI4MC8wLzE2ODEzOTg4OTExNjY/ZT0yMTQ3NDgzNjQ3JnY9YmV0YSZ0PTA5aG9uTVZ3TkYwaUpCdWRxWnpPdldJeUg1SndvU3NKUE8wMU5qUUNTYjApCgoKIyBMaWJyZXJpYXMgeSBwYXF1ZXRlcwpgYGB7ciBtZXNzYWdlPUZBTFNFLCB3YXJuaW5nPUZBTFNFfQpsaWJyYXJ5KG1hcHMpCmxpYnJhcnkoZm9yZWNhc3QpCmxpYnJhcnkodGlkeXZlcnNlKQpsaWJyYXJ5KG1hcHMpCmxpYnJhcnkoZ2dwbG90MikKbGlicmFyeShzaGlueSkKYGBgCgoKIyBPYnRlbmNpw7NuIGRlIGJhc2UgZGUgZGF0b3MKCiMjIEhpc3RvcmljbyBkZSBwb2JsYWNpb25lcwpgYGB7ciBtZXNzYWdlPUZBTFNFLCB3YXJuaW5nPUZBTFNFfQpwb3AgPC0gcmVhZC5jc3YoImhpc3RvcmljYWxfc3RhdGVfcG9wdWxhdGlvbl9ieV95ZWFyLmNzdiIsIGhlYWRlciA9IEZBTFNFKQpgYGAKCmBgYHtyfQpwb3AgPC0gcmVuYW1lKHBvcCwgIlN0YXRlIiA9ICJWMSIsICJZZWFyIiA9ICJWMiIsICJQb3B1bGF0aW9uIiA9ICJWMyIpCmBgYAoKYGBge3J9CnN0cihwb3ApCmBgYAoKIyMgUmVsYWNpb24gbm9tYnJlIGRlbCBlc3RhZG8gY29uIGNvZGljbyBkZSAyIGxldHJhcwpgYGB7cn0Kc3RhdGVfbmFtZSA8LSByZWFkX2NzdigidXNfc3RhdGVfdHdvX2xldHRlcl9hYmJyZXZpYXRpb24uY3N2IikKc3RhdGVfbmFtZSRgVS5TLiBTdGF0ZXMgYW5kIFRlcnJpdG9yaWVzYD0gc3RyX3RvX2xvd2VyKHN0YXRlX25hbWUkYFUuUy4gU3RhdGVzIGFuZCBUZXJyaXRvcmllc2ApCnN0YXRlX25hbWUgPC0gcmVuYW1lKHN0YXRlX25hbWUsICJzdGF0ZXMiID0gIlUuUy4gU3RhdGVzIGFuZCBUZXJyaXRvcmllcyIpCnN0YXRlX25hbWUgPC0gcmVuYW1lKHN0YXRlX25hbWUsICJTdGF0ZSIgPSAiVHdvLUxldHRlciBBYmJyZXZpYXRpb24iKQpgYGAKCiMjIENvb3JkZW5hZGFzCmBgYHtyfQpzdGF0ZV9kYXRhIDwtIGdncGxvdDI6Om1hcF9kYXRhKCdzdGF0ZScpCnN0YXRlX2RhdGEgPC0gZm9ydGlmeShzdGF0ZV9kYXRhKQpoZWFkKHN0YXRlX2RhdGEpCmBgYAoKIyMgQ29sdW1uYSBkZSBhw7FvcwpgYGB7cn0KeWVhcnMgPC0gc2VxKDIwMjAsIDIwNzAsIGJ5ID0gMSkKYGBgCgojIE51ZXZvIE3DqXhpY28KCmBgYHtyfQpwb3Bfbm0gPC0gZmlsdGVyKHBvcCwgcG9wJFN0YXRlID09ICJOTSIpCmBgYAoKYGBge3J9CnRzX25tIDwtIHRzKGRhdGEgPSBwb3Bfbm0kUG9wdWxhdGlvbiwgc3RhcnQgPSBjKDE5MDAsMSksIGZyZXF1ZW5jeSA9IDEpCnRzX25tCmBgYAoKIyMgTW9kZWxvIEFSSU1BCgpgYGB7cn0KYXJpbWFfbm0gPC0gYXV0by5hcmltYSh0c19ubSkKYXJpbWFfbm0Kc3VtbWFyeShhcmltYV9ubSkKYGBgCgojIyBHZW5lcmFyIHByb27Ds3N0aWNvCgpgYGB7cn0KcHJvbm9zdGljb19ubSA8LSBmb3JlY2FzdChhcmltYV9ubSwgbGV2ZWwgPSBjKDk1KSwgaCA9IDUxKQpwcm9ub3N0aWNvX25tCnBsb3QocHJvbm9zdGljb19ubSkKYGBgCgpgYGB7cn0KcHJvbm9zdGljb19ubV9kZWMgPC0gZGF0YS5mcmFtZShwcm9ub3N0aWNvX25tKQpwcm9ub3N0aWNvX25tX2RlYyA8LSBwcm9ub3N0aWNvX25tX2RlYyAlPiUgCiAgICAgICAgICAgICAgICAgICAgICAgIGFkZF9jb2x1bW4oRGVjYWRlID0geWVhcnMpICU+JSAKICAgICAgICAgICAgICAgICAgICAgICAgZmlsdGVyKHJvd19udW1iZXIoKSAlaW4lIGMoMTEsIDIxLCAzMSwgNDEsIDUxKSkgJT4lIAogICAgICAgICAgICAgICAgICAgICAgICBhZGRfY29sdW1uKFN0YXRlID0gIk5NIikgJT4lCiAgICAgICAgICAgICAgICAgICAgICAgIHNlbGVjdCgtTG8uOTUsIC1IaS45NSkKcHJvbm9zdGljb19ubV9kZWMKYGBgCgojIEFyaXpvbmEKCmBgYHtyfQpwb3BfYXogPC0gZmlsdGVyKHBvcCwgcG9wJFN0YXRlID09ICJBWiIpCmBgYAoKYGBge3J9CnRzX2F6IDwtIHRzKGRhdGEgPSBwb3BfYXokUG9wdWxhdGlvbiwgc3RhcnQgPSBjKDE5MDAsMSksIGZyZXF1ZW5jeSA9IDEpCnRzX2F6CmBgYAoKIyMgTW9kZWxvIEFSSU1BCgpgYGB7cn0KYXJpbWFfYXogPC0gYXV0by5hcmltYSh0c19heikKYXJpbWFfYXoKc3VtbWFyeShhcmltYV9heikKYGBgCgojIyBHZW5lcmFyIHByb27Ds3N0aWNvCgpgYGB7cn0KcHJvbm9zdGljb19heiA8LSBmb3JlY2FzdChhcmltYV9heiwgbGV2ZWwgPSBjKDk1KSwgaCA9IDUxKQpwcm9ub3N0aWNvX2F6CnBsb3QocHJvbm9zdGljb19heikKYGBgCgpgYGB7cn0KcHJvbm9zdGljb19hel9kZWMgPC0gZGF0YS5mcmFtZShwcm9ub3N0aWNvX2F6KQpwcm9ub3N0aWNvX2F6X2RlYyA8LSBwcm9ub3N0aWNvX2F6X2RlYyAlPiUgCiAgICAgICAgICAgICAgICAgICAgICAgIGFkZF9jb2x1bW4oRGVjYWRlID0geWVhcnMpICU+JSAKICAgICAgICAgICAgICAgICAgICAgICAgZmlsdGVyKHJvd19udW1iZXIoKSAlaW4lIGMoMTEsIDIxLCAzMSwgNDEsIDUxKSkgJT4lIAogICAgICAgICAgICAgICAgICAgICAgICBhZGRfY29sdW1uKFN0YXRlID0gIkFaIikgJT4lCiAgICAgICAgICAgICAgICAgICAgICAgIHNlbGVjdCgtTG8uOTUsIC1IaS45NSkKcHJvbm9zdGljb19hel9kZWMKYGBgCgojIENhbGlmb3JuaWEKCmBgYHtyfQpwb3BfY2EgPC0gZmlsdGVyKHBvcCwgcG9wJFN0YXRlID09ICJDQSIpCmBgYAoKYGBge3J9CnRzX2NhIDwtIHRzKGRhdGEgPSBwb3BfY2EkUG9wdWxhdGlvbiwgc3RhcnQgPSBjKDE5MDAsMSksIGZyZXF1ZW5jeSA9IDEpCnRzX2NhCmBgYAoKIyMgTW9kZWxvIEFSSU1BCgpgYGB7cn0KYXJpbWFfY2EgPC0gYXV0by5hcmltYSh0c19jYSkKYXJpbWFfY2EKc3VtbWFyeShhcmltYV9jYSkKYGBgCgojIyBHZW5lcmFyIHByb27Ds3N0aWNvCgpgYGB7cn0KcHJvbm9zdGljb19jYSA8LSBmb3JlY2FzdChhcmltYV9jYSwgbGV2ZWwgPSBjKDk1KSwgaCA9IDUxKQpwcm9ub3N0aWNvX2NhCnBsb3QocHJvbm9zdGljb19jYSkKYGBgCgpgYGB7cn0KcHJvbm9zdGljb19jYV9kZWMgPC0gZGF0YS5mcmFtZShwcm9ub3N0aWNvX2NhKQpwcm9ub3N0aWNvX2NhX2RlYyA8LSBwcm9ub3N0aWNvX2NhX2RlYyAlPiUgCiAgICAgICAgICAgICAgICAgICAgICAgIGFkZF9jb2x1bW4oRGVjYWRlID0geWVhcnMpICU+JSAKICAgICAgICAgICAgICAgICAgICAgICAgZmlsdGVyKHJvd19udW1iZXIoKSAlaW4lIGMoMTEsIDIxLCAzMSwgNDEsIDUxKSkgJT4lIAogICAgICAgICAgICAgICAgICAgICAgICBhZGRfY29sdW1uKFN0YXRlID0gIkNBIikgJT4lCiAgICAgICAgICAgICAgICAgICAgICAgIHNlbGVjdCgtTG8uOTUsIC1IaS45NSkKYGBgCgojIEFsYWJhbWEKCmBgYHtyfQpwb3BfYWwgPC0gZmlsdGVyKHBvcCwgcG9wJFN0YXRlID09ICJBTCIpCmBgYAoKYGBge3J9CnRzX2FsIDwtIHRzKGRhdGEgPSBwb3BfYWwkUG9wdWxhdGlvbiwgc3RhcnQgPSBjKDE5MDAsMSksIGZyZXF1ZW5jeSA9IDEpCnRzX2FsCmBgYAoKIyMgTW9kZWxvIEFSSU1BCgpgYGB7cn0KYXJpbWFfYWwgPC0gYXV0by5hcmltYSh0c19hbCkKYXJpbWFfYWwKc3VtbWFyeShhcmltYV9hbCkKYGBgCgojIyBHZW5lcmFyIHByb27Ds3N0aWNvCgpgYGB7cn0KcHJvbm9zdGljb19hbCA8LSBmb3JlY2FzdChhcmltYV9hbCwgbGV2ZWwgPSBjKDk1KSwgaCA9IDUxKQpwcm9ub3N0aWNvX2FsCnBsb3QocHJvbm9zdGljb19hbCkKYGBgCgpgYGB7cn0KcHJvbm9zdGljb19hbF9kZWMgPC0gZGF0YS5mcmFtZShwcm9ub3N0aWNvX2FsKQpwcm9ub3N0aWNvX2FsX2RlYyA8LSBwcm9ub3N0aWNvX2FsX2RlYyAlPiUgCiAgICAgICAgICAgICAgICAgICAgICAgIGFkZF9jb2x1bW4oRGVjYWRlID0geWVhcnMpICU+JSAKICAgICAgICAgICAgICAgICAgICAgICAgZmlsdGVyKHJvd19udW1iZXIoKSAlaW4lIGMoMTEsIDIxLCAzMSwgNDEsIDUxKSkgJT4lIAogICAgICAgICAgICAgICAgICAgICAgICBhZGRfY29sdW1uKFN0YXRlID0gIkFMIikgJT4lCiAgICAgICAgICAgICAgICAgICAgICAgIHNlbGVjdCgtTG8uOTUsIC1IaS45NSkKcHJvbm9zdGljb19hbF9kZWMKYGBgCgojIE51ZXZhIFlvcmsKCmBgYHtyfQpwb3BfbnkgPC0gZmlsdGVyKHBvcCwgcG9wJFN0YXRlID09ICJOWSIpCmBgYAoKYGBge3J9CnRzX255IDwtIHRzKGRhdGEgPSBwb3BfbnkkUG9wdWxhdGlvbiwgc3RhcnQgPSBjKDE5MDAsMSksIGZyZXF1ZW5jeSA9IDEpCnRzX255CmBgYAoKIyMgTW9kZWxvIEFSSU1BCgpgYGB7cn0KYXJpbWFfbnkgPC0gYXV0by5hcmltYSh0c19ueSkKYXJpbWFfbnkKc3VtbWFyeShhcmltYV9ueSkKYGBgCgojIyBHZW5lcmFyIHByb27Ds3N0aWNvCgpgYGB7cn0KcHJvbm9zdGljb19ueSA8LSBmb3JlY2FzdChhcmltYV9ueSwgbGV2ZWwgPSBjKDk1KSwgaCA9IDUxKQpwcm9ub3N0aWNvX255CnBsb3QocHJvbm9zdGljb19ueSkKYGBgCgpgYGB7cn0KcHJvbm9zdGljb19ueV9kZWMgPC0gZGF0YS5mcmFtZShwcm9ub3N0aWNvX255KQpwcm9ub3N0aWNvX255X2RlYyA8LSBwcm9ub3N0aWNvX255X2RlYyAlPiUgCiAgICAgICAgICAgICAgICAgICAgICAgIGFkZF9jb2x1bW4oRGVjYWRlID0geWVhcnMpICU+JSAKICAgICAgICAgICAgICAgICAgICAgICAgZmlsdGVyKHJvd19udW1iZXIoKSAlaW4lIGMoMTEsIDIxLCAzMSwgNDEsIDUxKSkgJT4lIAogICAgICAgICAgICAgICAgICAgICAgICBhZGRfY29sdW1uKFN0YXRlID0gIk5ZIikgJT4lCiAgICAgICAgICAgICAgICAgICAgICAgIHNlbGVjdCgtTG8uOTUsIC1IaS45NSkKcHJvbm9zdGljb19ueV9kZWMKYGBgCgojIE1hcGEKCiMjIEZ1bmNpw7NuIE1hcGEKYGBge3IsIGZ1bmNpb24gbWFwYX0KbWFwX3BvcCA8LSBmdW5jdGlvbihkZiwgc3RhdGVfZGF0YSwgZGVjYWRlID0gMjAzMCl7CiAgCiAgZGZfZGVjYWRlIDwtIGRmICU+JSAKICAgIGZpbHRlcihkZiREZWNhZGUgPT0gZGVjYWRlKQogIAogICMgRnVuY3Rpb24gZm9yIHNldHRpbmcgdGhlIGFlc3RoZXRpY3Mgb2YgdGhlIHBsb3QKICBteV90aGVtZSA8LSBmdW5jdGlvbiAoKSB7IAogICAgdGhlbWVfYncoKSArIHRoZW1lKGF4aXMudGV4dCA9IGVsZW1lbnRfdGV4dChzaXplID0gMTQpLAogICAgICAgICAgICAgICAgICAgICAgIGF4aXMudGl0bGUgPSBlbGVtZW50X3RleHQoc2l6ZSA9IDE0KSwKICAgICAgICAgICAgICAgICAgICAgICBzdHJpcC50ZXh0ID0gZWxlbWVudF90ZXh0KHNpemUgPSAxNCksCiAgICAgICAgICAgICAgICAgICAgICAgcGFuZWwuZ3JpZC5tYWpvciA9IGVsZW1lbnRfYmxhbmsoKSwgCiAgICAgICAgICAgICAgICAgICAgICAgcGFuZWwuZ3JpZC5taW5vciA9IGVsZW1lbnRfYmxhbmsoKSwKICAgICAgICAgICAgICAgICAgICAgICBwYW5lbC5iYWNrZ3JvdW5kID0gZWxlbWVudF9ibGFuaygpLCAKICAgICAgICAgICAgICAgICAgICAgICBsZWdlbmQucG9zaXRpb24gPSAiYm90dG9tIiwKICAgICAgICAgICAgICAgICAgICAgICBwYW5lbC5ib3JkZXIgPSBlbGVtZW50X2JsYW5rKCksIAogICAgICAgICAgICAgICAgICAgICAgIHN0cmlwLmJhY2tncm91bmQgPSBlbGVtZW50X3JlY3QoZmlsbCA9ICd3aGl0ZScsIGNvbG91ciA9ICd3aGl0ZScpKQogIH0KICAKICAKICAjIEFkZCB0aGUgZGF0YSB0aGUgdXNlciB3YW50cyB0byBzZWUgdG8gdGhlIGdlb2dyYXBoaWNhbCB3b3JsZCBkYXRhCiAgc3RhdGVfZGF0YVsncG9wJ10gPC0gZGZfZGVjYWRlJFBvaW50LkZvcmVjYXN0W21hdGNoKHN0YXRlX2RhdGEkcmVnaW9uLCBkZl9kZWNhZGUkc3RhdGVzKV0KCiAgCiAgIyBTcGVjaWZ5IHRoZSBwbG90IGZvciB0aGUgd29ybGQgbWFwCiAgbGlicmFyeShSQ29sb3JCcmV3ZXIpCiAgZyA8LSBnZ3Bsb3QoKSArIAogICAgZ2VvbV9wb2x5Z29uKGRhdGEgPSBzdGF0ZV9kYXRhLCBzaXplID0gMSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgYWVzKHggPSBsb25nLCB5ID0gbGF0LCBmaWxsID0gcG9wLCBncm91cCA9IGdyb3VwKSkgKyAKICAgIHNjYWxlX2ZpbGxfZ3JhZGllbnRuKGNvbG91cnMgPSBicmV3ZXIucGFsKDUsICJSZEJ1IiksIG5hLnZhbHVlID0gJ3doaXRlJyxsaW1pdHMgPSBjKDU4MDQ3OSwgNDk4MTYyNTApKSArIAogICAgbXlfdGhlbWUoKQogIAogIHJldHVybihnKQp9CmBgYAoKCiMjIERhdGFzZXQgTWFwYXMKYGBge3J9CmRmX21hcCA8LSByYmluZChwcm9ub3N0aWNvX2FsX2RlYyxwcm9ub3N0aWNvX2F6X2RlYykKZGZfbWFwIDwtIHJiaW5kKGRmX21hcCxwcm9ub3N0aWNvX2NhX2RlYykKZGZfbWFwIDwtIHJiaW5kKGRmX21hcCxwcm9ub3N0aWNvX25tX2RlYykKZGZfbWFwIDwtIHJiaW5kKGRmX21hcCxwcm9ub3N0aWNvX255X2RlYykKCmRmX21hcApgYGAKCmBgYHtyfQpkZl9tYXAgPC0gbGVmdF9qb2luKGRmX21hcCxzdGF0ZV9uYW1lLGJ5ID0gIlN0YXRlIikKYGBgCgojIyBNYXBhIDIwMzAKYGBge3J9Cm1hcF9wb3AoZGZfbWFwLCBzdGF0ZV9kYXRhKQpgYGAKCmBgYHtyfQptYXBfcG9wKGRmX21hcCwgc3RhdGVfZGF0YSwgZGVjYWRlID0gMjA0MCkKYGBgCgpgYGB7cn0KbWFwX3BvcChkZl9tYXAsIHN0YXRlX2RhdGEsIGRlY2FkZSA9IDIwNTApCmBgYAoKYGBge3J9Cm1hcF9wb3AoZGZfbWFwLCBzdGF0ZV9kYXRhLCBkZWNhZGUgPSAyMDYwKQpgYGAKCmBgYHtyfQptYXBfcG9wKGRmX21hcCwgc3RhdGVfZGF0YSwgZGVjYWRlID0gMjA3MCkKYGBgCgojIEFwbGljYWNpw7NuIGVuIFNoaW55ClB1ZWRlIGVuY29udHJhciBlbCBhcHAgZGUgc2hpbnkgcGFyYSBnZW5lcmFyIGVsICoqZm9yZWNhc3QqKiBjdWFscXVpZXJhIGRlIGxvcyA1MSBlc3RhZG9zLCB1dGlsaXphbmRvIHN1IGPDs2RpZ28gZGUgMiBsZXRyYXMuIElndWFsbWVudGUgc2UgZGVzcGxpZWdhIGVsIG1hcGEgZGUgbGEgZGVuc2lkYWQgZGUgbGEgcG9ibGFjacOzbiBlbiBVU0EgZGUgYSBjdWVyZG8gYSBsb3MgcHJvbm9zdGljb3MgcmVhbGl6YWRvcyBwYXJhIGxhcyBzaWd1aWVudGVzIDUgZMOpY2FkYXMuICAKW0FwcCBkZSBzaGlueV0oaHR0cHM6Ly85bWR2YzktamVzMHMwdXppZWwtY2FyZGVsYXMwcDByZXouc2hpbnlhcHBzLmlvL1VTQV9kZW1vZ3JhcGhpY3MvKQoKYGBge3IgZWNobz1GQUxTRX0KdWkgPC0gZmx1aWRQYWdlKG5hdmJhclBhZ2UoIlByb27Ds3N0aWNvcyBQcm95ZWN0YWRvcyIsCiAgICAgICAgICAgICAgICAgICAgICAgICAgIHRhYlBhbmVsKCJQcm9ub3N0aWNvcyBwb3IgZXN0YWRvIiwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgc2lkZWJhclBhbmVsKHRhZ3MkaDQoIkluZ3Jlc2EgbG9zIHNpZ3VpZW50ZXMgZGF0b3M6IiksCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICB0ZXh0SW5wdXQoImVzdGFkbyIsICJFc2NyaWJlIGVsIGPDs2RpZ28gZGUgMiBsZXRyYXMgZGVsIGVzdGFkbyBxdWUgcXVpZXJlczoiLCJISSIpLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgc2VsZWN0SW5wdXQoImNvbmZpYW56YSIsICJOaXZlbCBkZSBDb25maWFiaWxpZGFkIiwgYyg5MCw5NSw5OSkpLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgc2xpZGVySW5wdXQoInBlcmlvZG9zIiwgIsK/Q3XDoW50b3MgcGVyaW9kb3MgcXVpZXJlcyBwcm9ub3N0aWNhcj8iLCAxLDEwMSw1KQogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICApLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBtYWluUGFuZWwoaDMoIkVsIHByb25vc3RpY28gcGFyYSBlbCBlc3RhZG8gcXVlIHNlbGVjY2lvbsOzIGVzOiIpLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgcGxvdE91dHB1dCgicHJvbm9zdGljbyIpCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICkKICAgICAgICAgICAgICAgICAgICAgICAgICAgKSwgCiAgICAgICAgICAgICAgICAgICAgICAgICAgIHRhYlBhbmVsKCJNYXBhcyIsIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBzaWRlYmFyUGFuZWwodGFncyRoNCgiSW5ncmVzYSBsb3Mgc2lndWllbnRlcyBkYXRvczoiKSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHNlbGVjdElucHV0KCJkZWNhZGEiLCAiU2VsZWNjaW9uYSB1bmEgRGVjYWRhOiIsIGMoMjAzMCwyMDQwLDIwNTAsMjA2MCwyMDcwKSkKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgKSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgbWFpblBhbmVsKGgzKCJFbCBtYXBhIGRlIGxhIGRlbW9ncmFmw61hIGRlIFVTQSBwYXJhIGxhIGRlZGFjYSBzZWxlY2Npb25hZGEgZXM6IiksCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBwbG90T3V0cHV0KCJtYXBhIikKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgKQogICAgICAgICAgICAgICAgICAgICAgICAgICApKSkKCnNlcnZlciA8LSBmdW5jdGlvbihpbnB1dCwgb3V0cHV0KXsKICAKICBvdXRwdXQkcHJvbm9zdGljbyA8LSByZW5kZXJQbG90KHsKICAgIHBvcCA8LSByZWFkLmNzdigiaGlzdG9yaWNhbF9zdGF0ZV9wb3B1bGF0aW9uX2J5X3llYXIuY3N2IiwgaGVhZGVyID0gRkFMU0UpCiAgICBwb3AgPC0gcmVuYW1lKHBvcCwgIlN0YXRlIiA9ICJWMSIsICJZZWFyIiA9ICJWMiIsICJQb3B1bGF0aW9uIiA9ICJWMyIpCiAgICBzdGF0ZV9wb3AgPC0gZmlsdGVyKHBvcCwgcG9wJFN0YXRlID09IHRvdXBwZXIoaW5wdXQkZXN0YWRvKSkKICAgIAogICAgaWYgKHRvdXBwZXIoaW5wdXQkZXN0YWRvKSA9PSAiSEkiIHwgdG91cHBlcihpbnB1dCRlc3RhZG8pID09ICJBSyIpIHsKICAgICAgdHMgPC0gdHMoZGF0YSA9IHN0YXRlX3BvcCRQb3B1bGF0aW9uLCBzdGFydCA9IGMoMTk1MCwxKSwgZnJlcXVlbmN5ID0gMSkKICAgIH0gZWxzZSB7CiAgICAgIHRzIDwtIHRzKGRhdGEgPSBzdGF0ZV9wb3AkUG9wdWxhdGlvbiwgc3RhcnQgPSBjKDE5MDAsMSksIGZyZXF1ZW5jeSA9IDEpCiAgICB9CiAgICAKICAgIGFyaW1hIDwtIGF1dG8uYXJpbWEodHMpCiAgICBwcm9ub3N0aWNvIDwtIGZvcmVjYXN0KGFyaW1hLCBsZXZlbCA9IGMoYXMubnVtZXJpYyhpbnB1dCRjb25maWFuemEpKSwgaCA9IGlucHV0JHBlcmlvZG9zKQogICAgcGxvdChwcm9ub3N0aWNvKQogIH0pCiAgCiAgb3V0cHV0JG1hcGEgPC0gcmVuZGVyUGxvdCh7CiAgICBwb3AgPC0gcmVhZC5jc3YoImhpc3RvcmljYWxfc3RhdGVfcG9wdWxhdGlvbl9ieV95ZWFyLmNzdiIsIGhlYWRlciA9IEZBTFNFKQogICAgcG9wIDwtIHJlbmFtZShwb3AsICJTdGF0ZSIgPSAiVjEiLCAiWWVhciIgPSAiVjIiLCAiUG9wdWxhdGlvbiIgPSAiVjMiKQogICAgc3RhdGVzIDwtIHVuaXF1ZShwb3AkU3RhdGUpCiAgICBwb3AgPC0gcG9wICU+JSAKICAgICAgcGl2b3Rfd2lkZXIoCiAgICAgICAgbmFtZXNfZnJvbSA9IFN0YXRlLAogICAgICAgIHZhbHVlc19mcm9tID0gUG9wdWxhdGlvbikgJT4lIAogICAgICBhcnJhbmdlKFllYXIpCiAgICAKICAgIAogICAgcG9wMSA8LSBzdWJzZXQocG9wLCBzZWxlY3QgPSAtWWVhcikKICAgIHBvcDIgPC0gcG9wICU+JSBmaWx0ZXIoWWVhciA+PSAxOTUwKQogICAgcG9wMiA8LSBzdWJzZXQocG9wMiwgc2VsZWN0ID0gLVllYXIpCiAgICB5ZWFycyA8LSBzZXEoMjAyMCwgMjA3MCwgYnkgPSAxKQogICAgcmVzdWx0YWRvcyA8LSBkYXRhLmZyYW1lKHllYXIgPSB5ZWFycykKICAgIAogICAgZm9yIChpIGluIHN0YXRlcykgewogICAgICBpZiAoaSA9PSAiQUsiIHwgaSA9PSAiSEkiKSB7CiAgICAgICAgdHMgPC0gdHMoZGF0YT1wb3AyW2ldLCBzdGFydCA9IGMoMTk1MCksIGZyZXF1ZW5jeSA9IDEpCiAgICAgIH0gZWxzZSB7CiAgICAgICAgdHMgPC0gdHMoZGF0YT1wb3AxW2ldLCBzdGFydCA9IGMoMTkwMCksIGZyZXF1ZW5jeSA9IDEpCiAgICAgIH0KICAgICAgYXJpbWEgPC0gYXV0by5hcmltYSh0cykKICAgICAgcHJvbm9zdGljbyA8LSBmb3JlY2FzdChhcmltYSwgbGV2ZWwgPSBjKDk1KSwgaCA9IDUxKQogICAgICBkIDwtIGFzLmRhdGEuZnJhbWUocHJvbm9zdGljbykKICAgICAgcmVzdWx0YWRvc1tpXSA8LSBkJGBQb2ludCBGb3JlY2FzdGAKICAgIH0KICAgIAogICAgcmVzdWx0YWRvcyA8LSByZXN1bHRhZG9zICU+JSAKICAgICAgcGl2b3RfbG9uZ2VyKAogICAgICAgIGNvbHMgPSBBSzpXWSwgCiAgICAgICAgbmFtZXNfdG8gPSAiU3RhdGUiLCAKICAgICAgICB2YWx1ZXNfdG8gPSAicG9wIgogICAgICApCiAgICAKICAgIHN0YXRlX2RhdGEgPC0gZ2dwbG90Mjo6bWFwX2RhdGEoJ3N0YXRlJykKICAgIHN0YXRlX2RhdGEgPC0gZm9ydGlmeShzdGF0ZV9kYXRhKQogICAgCiAgICBzdGF0ZV9uYW1lIDwtIHJlYWRfY3N2KCJ1c19zdGF0ZV90d29fbGV0dGVyX2FiYnJldmlhdGlvbi5jc3YiKQogICAgc3RhdGVfbmFtZSRgVS5TLiBTdGF0ZXMgYW5kIFRlcnJpdG9yaWVzYD0gc3RyX3RvX2xvd2VyKHN0YXRlX25hbWUkYFUuUy4gU3RhdGVzIGFuZCBUZXJyaXRvcmllc2ApCiAgICBzdGF0ZV9uYW1lIDwtIHJlbmFtZShzdGF0ZV9uYW1lLCAic3RhdGVzIiA9ICJVLlMuIFN0YXRlcyBhbmQgVGVycml0b3JpZXMiKQogICAgc3RhdGVfbmFtZSA8LSByZW5hbWUoc3RhdGVfbmFtZSwgIlN0YXRlIiA9ICJUd28tTGV0dGVyIEFiYnJldmlhdGlvbiIpCiAgICAKICAgIGRmX21hcCA8LSBsZWZ0X2pvaW4ocmVzdWx0YWRvcyxzdGF0ZV9uYW1lLGJ5ID0gIlN0YXRlIikKICAgIAogICAgbWFwX3BvcCA8LSBmdW5jdGlvbihkZiwgc3RhdGVfZGF0YSwgZGVjYWRlID0gMjAzMCl7CiAgICAgIGRmX2RlY2FkZSA8LSBkZiAlPiUgCiAgICAgICAgZmlsdGVyKGRmJHllYXIgPT0gYXMubnVtZXJpYyhkZWNhZGUpKQogICAgICBteV90aGVtZSA8LSBmdW5jdGlvbiAoKSB7IAogICAgICAgIHRoZW1lX2J3KCkgKyB0aGVtZShheGlzLnRleHQgPSBlbGVtZW50X3RleHQoc2l6ZSA9IDE0KSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgYXhpcy50aXRsZSA9IGVsZW1lbnRfdGV4dChzaXplID0gMTQpLAogICAgICAgICAgICAgICAgICAgICAgICAgICBzdHJpcC50ZXh0ID0gZWxlbWVudF90ZXh0KHNpemUgPSAxNCksCiAgICAgICAgICAgICAgICAgICAgICAgICAgIHBhbmVsLmdyaWQubWFqb3IgPSBlbGVtZW50X2JsYW5rKCksIAogICAgICAgICAgICAgICAgICAgICAgICAgICBwYW5lbC5ncmlkLm1pbm9yID0gZWxlbWVudF9ibGFuaygpLAogICAgICAgICAgICAgICAgICAgICAgICAgICBwYW5lbC5iYWNrZ3JvdW5kID0gZWxlbWVudF9ibGFuaygpLCAKICAgICAgICAgICAgICAgICAgICAgICAgICAgbGVnZW5kLnBvc2l0aW9uID0gImJvdHRvbSIsCiAgICAgICAgICAgICAgICAgICAgICAgICAgIHBhbmVsLmJvcmRlciA9IGVsZW1lbnRfYmxhbmsoKSwgCiAgICAgICAgICAgICAgICAgICAgICAgICAgIHN0cmlwLmJhY2tncm91bmQgPSBlbGVtZW50X3JlY3QoZmlsbCA9ICd3aGl0ZScsIGNvbG91ciA9ICd3aGl0ZScpKQogICAgICB9CiAgICAgIHN0YXRlX2RhdGFbJ3BvcCddIDwtIGRmX2RlY2FkZSRwb3BbbWF0Y2goc3RhdGVfZGF0YSRyZWdpb24sIGRmX2RlY2FkZSRzdGF0ZXMpXQogICAgICBnIDwtIGdncGxvdCgpICsgCiAgICAgICAgZ2VvbV9wb2x5Z29uKGRhdGEgPSBzdGF0ZV9kYXRhLCBzaXplID0gMSwKICAgICAgICAgICAgICAgICAgICAgYWVzKHggPSBsb25nLCB5ID0gbGF0LCBmaWxsID0gcG9wLCBncm91cCA9IGdyb3VwKSkgKyAKICAgICAgICBzY2FsZV9maWxsX2dyYWRpZW50bihjb2xvdXJzID0gYnJld2VyLnBhbCg1LCAiUmRCdSIpLCBuYS52YWx1ZSA9ICd3aGl0ZScsbGltaXRzID0gYyg1ODA0NzksIDQ5ODE2MjUwKSkgKyAKICAgICAgICBteV90aGVtZSgpCiAgICAgIHJldHVybihnKQogICAgfQogICAgCiAgICBtYXBfcG9wKGRmX21hcCwgc3RhdGVfZGF0YSwgZGVjYWRlID0gaW5wdXQkZGVjYWRhKQogICAgCiAgfSkKICAKICAKICAKfQoKc2hpbnlBcHAodWk9dWksIHNlcnZlcj1zZXJ2ZXIpCmBgYAoK