act2mod1opt8
Mariano Bautista A01656904
2025-02-18
Ejemplo en clase : Poblacion
Instalar paquetes y llamar librerias
#install.packages("forecast")
library(forecast)## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo#install.packages("tidyverse")
library(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.1 ✔ tibble 3.2.1
## ✔ lubridate 1.9.4 ✔ tidyr 1.3.1
## ✔ purrr 1.0.4## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorsInstalar base de datos
#file.choose()
# Instalar y cargar el paquete necesario (si no lo has hecho previamente)
library(readxl)
# Seleccionar el archivo usando un cuadro de diálogo
ruta_archivo <- file.choose()
# Leer el archivo Excel seleccionado
poblacion <- read.csv(ruta_archivo)
# Si es necesario, convertir a data frame
poblacion <- as.data.frame(poblacion)Entender la base de datos
summary(poblacion)## state year population
## Length:6020 Min. :1900 Min. : 43000
## Class :character 1st Qu.:1930 1st Qu.: 901483
## Mode :character Median :1960 Median : 2359000
## Mean :1960 Mean : 3726003
## 3rd Qu.:1990 3rd Qu.: 4541883
## Max. :2019 Max. :39512223str(poblacion)## '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 ...head(poblacion)## state year population
## 1 AK 1950 135000
## 2 AK 1951 158000
## 3 AK 1952 189000
## 4 AK 1953 205000
## 5 AK 1954 215000
## 6 AK 1955 222000Serie de Tiempo en Texas
poblacion_texas <- poblacion %>% filter(state=="TX")
ggplot(poblacion_texas, aes(x=year, y=population)) +
geom_line() +
labs(title = "Poblacion de Texas", x ="Año",
y = "Poblacion")
ts_texas <- ts(poblacion_texas$population, start=1900, frequency = 1) # Serie de Tiempo Anual
# ts_texas <- ts(poblacion_texas$population, start=1900, frequency = 4) # Serie de Tiempo Trimestral
# ts_texas <- ts(poblacion_texas$population, start=1900, frequency = 12) # Serie de Tiempo mensual
arima_texas <- auto.arima(ts_texas)
summary(arima_texas)## Series: ts_texas
## ARIMA(0,2,2)
##
## Coefficients:
## ma1 ma2
## -0.5950 -0.1798
## s.e. 0.0913 0.0951
##
## sigma^2 = 1.031e+10: log likelihood = -1527.14
## AIC=3060.28 AICc=3060.5 BIC=3068.6
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 12147.62 99818.31 59257.39 0.1046163 0.5686743 0.2672197
## ACF1
## Training set -0.02136734pronostico_texas <- forecast(arima_texas, level = 95, h = 10)
pronostico_texas## Point Forecast Lo 95 Hi 95
## 2020 29398472 29199487 29597457
## 2021 29806827 29463665 30149990
## 2022 30215183 29742956 30687410
## 2023 30623538 30024100 31222977
## 2024 31031894 30303359 31760429
## 2025 31440249 30579246 32301253
## 2026 31848605 30851090 32846119
## 2027 32256960 31118581 33395339
## 2028 32665316 31381587 33949044
## 2029 33073671 31640070 34507272plot(pronostico_texas, main = "Poblacion en Texas")
Ejercicio en Clase Lunes 17: Mapa
Instalar paquetes y llamar librerias
#install.packages("forecast")
library(forecast)
#install.packages("tidyverse")
library(tidyverse)
#install.packages("maps")
library(maps)##
## Adjuntando el paquete: 'maps'## The following object is masked from 'package:purrr':
##
## mapInstalar Base de Datos
#file.choose()
# Instalar y cargar el paquete necesario (si no lo has hecho previamente)
library(readxl)
# Seleccionar el archivo usando un cuadro de diálogo
ruta_archivo <- file.choose()
# Leer el archivo Excel seleccionado
poblacion <- read.csv(ruta_archivo)
# Si es necesario, convertir a data frame
poblacion <- as.data.frame(poblacion)Entender la Base de Datos
summary(poblacion)## state year population
## Length:6020 Min. :1900 Min. : 43000
## Class :character 1st Qu.:1930 1st Qu.: 901483
## Mode :character Median :1960 Median : 2359000
## Mean :1960 Mean : 3726003
## 3rd Qu.:1990 3rd Qu.: 4541883
## Max. :2019 Max. :39512223str(poblacion)## '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 ...head(poblacion)## state year population
## 1 AK 1950 135000
## 2 AK 1951 158000
## 3 AK 1952 189000
## 4 AK 1953 205000
## 5 AK 1954 215000
## 6 AK 1955 222000Serie de Tiempo en Texas
poblacion_texas <- poblacion %>% filter(state=="TX")
ggplot(poblacion_texas, aes(x=year, y=population)) +
geom_line() +
labs(title = "Poblacion de Texas", x ="Año",
y = "Poblacion")
ts_texas <- ts(poblacion_texas$population, start=1900, frequency = 1) # Serie de Tiempo Anual
# ts_texas <- ts(poblacion_texas$population, start=1900, frequency = 4) # Serie de Tiempo Trimestral
# ts_texas <- ts(poblacion_texas$population, start=1900, frequency = 12) # Serie de Tiempo mensual
arima_texas <- auto.arima(ts_texas)
summary(arima_texas)## Series: ts_texas
## ARIMA(0,2,2)
##
## Coefficients:
## ma1 ma2
## -0.5950 -0.1798
## s.e. 0.0913 0.0951
##
## sigma^2 = 1.031e+10: log likelihood = -1527.14
## AIC=3060.28 AICc=3060.5 BIC=3068.6
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 12147.62 99818.31 59257.39 0.1046163 0.5686743 0.2672197
## ACF1
## Training set -0.02136734pronostico_texas <- forecast(arima_texas, level = 95, h = 10)
pronostico_texas## Point Forecast Lo 95 Hi 95
## 2020 29398472 29199487 29597457
## 2021 29806827 29463665 30149990
## 2022 30215183 29742956 30687410
## 2023 30623538 30024100 31222977
## 2024 31031894 30303359 31760429
## 2025 31440249 30579246 32301253
## 2026 31848605 30851090 32846119
## 2027 32256960 31118581 33395339
## 2028 32665316 31381587 33949044
## 2029 33073671 31640070 34507272plot(pronostico_texas, main = "Poblacion en Texas")
Crear un Mapa
# Crear un mapa de EUA por decada, con un gradiente verde-rojo de la poblacion por estado, desde 1950 hasta 2050
map(database = "state")
map(database = "state", regions = "Texas", col = "red", fill = TRUE, add = TRUE)
map(database = "state", regions = "New York", col = "green", fill = TRUE, add = TRUE)
Generar Pronostico por Cada Estado
warning=FALSESuponiendo que en ‘poblacion’ tienes columnas:
state (abreviatura, ej. “TX”), year (año), population (valor)
poblacion_extended <- poblacion %>%
group_by(state) %>%
arrange(year) %>%
# Para cada estado, creamos un df con valores observados + pronóstico
do({
df <- .
# Años mínimos y máximos en tus datos
min_year <- min(df$year)
max_year <- max(df$year)
# Creamos la serie de tiempo anual
ts_pop <- ts(df$population, start = min_year, frequency = 1)
# Determinamos cuántos años faltan para llegar a 2050
horizon <- 2050 - max_year
# Si el dataset llega hasta antes de 2050, hacemos forecast
if(horizon > 0){
fit <- auto.arima(ts_pop)
fc <- forecast(fit, h = horizon)
# Data frame con los datos pronosticados
years_forecast <- (max_year + 1):2050
df_forecast <- data.frame(
state = unique(df$state),
year = years_forecast,
population = as.numeric(fc$mean)
)
# Unimos histórico + forecast
df_all <- bind_rows(
# Histórico (columnas relevantes)
df %>% select(state, year, population),
# Futuro
df_forecast
)
} else {
# Si ya tenemos datos hasta 2050 o más, no pronosticamos
df_all <- df
}
df_all
}) %>%
ungroup()Convertir abreviaturas de estado (p.ej., “TX”) a nombres completos
df_state_names <- data.frame(
state_abb = state.abb,
state_name = tolower(state.name),
stringsAsFactors = FALSE
)
poblacion_full <- poblacion_extended %>%
left_join(df_state_names, by = c("state" = "state_abb"))Unir la geometría de los estados con los datos de población
4.1) Obtenemos las coordenadas de polígonos con ‘map_data(“state”)’
states_map <- map_data("state")4.2) Hacemos left_join para heredar la columna ‘population’ (por estado y año)
map_data_joined <- states_map %>%
left_join(poblacion_full, by = c("region" = "state_name"))## Warning in left_join(., poblacion_full, by = c(region = "state_name")): Detected an unexpected many-to-many relationship between `x` and `y`.
## ℹ Row 1 of `x` matches multiple rows in `y`.
## ℹ Row 102 of `y` matches multiple rows in `x`.
## ℹ If a many-to-many relationship is expected, set `relationship =
## "many-to-many"` to silence this warning.Filtrar décadas entre 1950 y 2050
map_data_decadas <- map_data_joined %>%
mutate(year = as.numeric(year)) %>%
filter(year >= 1950, year <= 2050, year %% 10 == 0)Graficar con ggplot2 (rojo = población baja, verde = población alta)Filtrar décadas entre 1950 y 2050
ggplot(map_data_decadas, aes(x = long,
y = lat,
group = group,
fill = population)) +
geom_polygon(color = "black", size = 0.1) +
# Gradiente de rojo a verde
scale_fill_gradient(low = "red", high = "green",
na.value = "grey90") +
facet_wrap(~ year) +
coord_fixed(1.3) +
labs(title = "Población de EUA por Estado (1950 - 2050, por década)",
fill = "Población estimada") +
theme_void() +
theme(legend.position = "right",
strip.text = element_text(face = "bold"))## 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.
1. Actividad 2. leche saborizada hersey
library(forecast)
library(tidyverse)
library(forecast)
library(ggplot2)
library(readxl)imagen_path <- file.choose()
#Cargar la imagen con knitr::include_graphics
knitr::include_graphics(imagen_path)
# Seleccionar el archivo usando un cuadro de diálogo
ruta_archivo <- file.choose()
# Leer el archivo Excel seleccionado
ventas <- read_excel(ruta_archivo)
# Si es necesario, convertir a data frame
ventas <- as.data.frame(ventas)
# Ver las primeras filas del archivo
head(ventas)## Ventas
## 1 25520.51
## 2 23740.11
## 3 26253.58
## 4 25868.43
## 5 27072.87
## 6 27150.503. Aritma
ts_ventas <- ts(ventas$Ventas, start = c(2017,1), frequency = 12)
autoplot(ts_ventas) + labs(title= "Ventas de leche saborizada Hershey`s", x = "tiempo", y = "Miles de dolares")
arima_ventas <- auto.arima(ts_ventas)
summary(arima_ventas)## Series: ts_ventas
## ARIMA(1,0,0)(1,1,0)[12] with drift
##
## Coefficients:
## ar1 sar1 drift
## 0.6383 -0.5517 288.8979
## s.e. 0.1551 0.2047 14.5026
##
## sigma^2 = 202701: log likelihood = -181.5
## AIC=371 AICc=373.11 BIC=375.72
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 25.22158 343.864 227.17 0.08059932 0.7069542 0.06491044 0.2081026pronostico_ventas <- forecast(arima_ventas, level=95, h=12)
pronostico_ventas## Point Forecast Lo 95 Hi 95
## Jan 2020 35498.90 34616.48 36381.32
## Feb 2020 34202.17 33155.28 35249.05
## Mar 2020 36703.01 35596.10 37809.92
## Apr 2020 36271.90 35141.44 37402.36
## May 2020 37121.98 35982.07 38261.90
## Jun 2020 37102.65 35958.90 38246.40
## Jul 2020 37151.04 36005.73 38296.34
## Aug 2020 38564.64 37418.70 39710.58
## Sep 2020 38755.22 37609.03 39901.42
## Oct 2020 39779.02 38632.72 40925.32
## Nov 2020 38741.63 37595.28 39887.97
## Dec 2020 38645.86 37499.50 39792.22autoplot(pronostico_ventas) + labs(title= "pronostico de ventas 2020 de leche saborizada Hershey`s", x= "tiempo", y="miles de dolares")
ventas$mes <- 1:36
regresion_ventas <- lm(Ventas ~ mes, data = ventas)
summary(regresion_ventas)##
## Call:
## lm(formula = Ventas ~ mes, data = ventas)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2075.79 -326.41 33.74 458.40 1537.04
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 24894.67 275.03 90.52 <2e-16 ***
## mes 298.37 12.96 23.02 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 808 on 34 degrees of freedom
## Multiple R-squared: 0.9397, Adjusted R-squared: 0.9379
## F-statistic: 529.8 on 1 and 34 DF, p-value: < 2.2e-16siguiente_año <- data.frame(mes=37:38)
prediccion_regresion <- predict(regresion_ventas, siguiente_año)
prediccion_regresion## 1 2
## 35934.49 36232.86plot(ventas$mes, ventas$Ventas, main = "Pronostico de ventas 2020 de leche saborizada Hershey`s", xlab = "Tiempo", ylab = "MIles de dolares")
abline(regresion_ventas, col="blue")
points(siguiente_año$mes, prediccion_regresion, col = "red")
predicciones_reales <- predict(regresion_ventas, ventas)
predicciones_reales## 1 2 3 4 5 6 7 8
## 25193.04 25491.42 25789.79 26088.16 26386.54 26684.91 26983.28 27281.66
## 9 10 11 12 13 14 15 16
## 27580.03 27878.40 28176.78 28475.15 28773.52 29071.90 29370.27 29668.64
## 17 18 19 20 21 22 23 24
## 29967.02 30265.39 30563.76 30862.14 31160.51 31458.88 31757.26 32055.63
## 25 26 27 28 29 30 31 32
## 32354.00 32652.38 32950.75 33249.12 33547.50 33845.87 34144.25 34442.62
## 33 34 35 36
## 34740.99 35039.37 35337.74 35636.11MAPE <- mean(abs((ventas$Ventas - predicciones_reales)/ ventas$Ventas))*100
MAPE## [1] 2.0112973. concusiones
El mejor modelo que se adapta a la serie es el
SSARIMA con un MAPE
de 0.71%, comparado con la Regresion lineal que su MAPE es de
2.01%.
Para el siguiente año, la proyeccion de ventas es la siguiente:
| Mes y Año | Escenario Esperado | Escenario Pesimista | Escenario Optimista |
|---|---|---|---|
| Jan 2020 | 35498.90 | 34616.48 | 36381.32 |
| Feb 2020 | 34202.17 | 33155.28 | 35249.05 |
| Mar 2020 | 36703.01 | 35596.10 | 37809.92 |
| Apr 2020 | 36271.90 | 35141.44 | 37402.36 |
| May 2020 | 37121.98 | 35982.07 | 38261.90 |
| Jun 2020 | 37102.65 | 35958.90 | 38246.40 |
| Jul 2020 | 37151.04 | 36005.73 | 38296.34 |
| Aug 2020 | 38564.64 | 37418.70 | 39710.58 |
| Sep 2020 | 38755.22 | 37609.03 | 39901.42 |
| Oct 2020 | 39779.02 | 38632.72 | 40925.32 |
| Nov 2020 | 38741.63 | 37595.28 | 39887.97 |
| Dec 2020 | 38645.86 | 37499.50 | 39792.22 |
# Seleccionar el archivo usando un cuadro de diálogo
ruta_archivo <- file.choose()
# Leer el archivo Excel seleccionado
ventas_por_anio <- read_csv(ruta_archivo)## Rows: 48 Columns: 3
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): mes
## dbl (2): anio, ventas
##
## ℹ 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.# Si es necesario, convertir a data frame
ventas_por_anio <- as.data.frame(ventas_por_anio)
# Ver las primeras filas del archivo
head(ventas_por_anio)## mes anio ventas
## 1 Enero 2017 25520.51
## 2 Febrero 2017 23740.11
## 3 Marzo 2017 26253.58
## 4 Abril 2017 25868.43
## 5 Mayo 2017 27072.87
## 6 Junio 2017 27150.50ggplot(ventas_por_anio, aes(x=mes, y = ventas,
col=as.factor(anio), group = anio)) +
geom_line() + labs(title = "ventas de leche saborizada Hershey`s por año", x= "mes", y= "Miles de dolares")
#uestra recomendacion seria realizar campañas punlicitarias par aumentar el consumo a principios de año