Actividad 2

Ejercicio Clase: Población

Instalar paquetes y librerías

#install.packages("forecast")
#install.packages("tidyverse")
#install.packages("dplyr")
#install.packages("ggplot2")
#install.packages("sf")

library(forecast)
library(tidyverse)
library(dplyr)
library(ggplot2)
library(sf)

Importar base de datos

poblacion <- read.csv("population.csv")

Entender la base de datos

# Analisis Descriptivo
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.   :39512223

# Tipos de Variables
str(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 ...

Series de Tiempo de la Población en Texas

poblacion_TX <- poblacion %>%
  filter(state == "TX")

ggplot(poblacion_TX, aes(x=year, y=population)) +
  geom_line() +
  labs(title =  "Población Texas", x = "Año", y = "Población")

ts_TX <- ts(poblacion_TX$population, start = 1900, frequency = 1) # ST Anual

#ARIMA
arima_TX <- auto.arima(ts_TX)

summary(arima_TX)

## Series: ts_TX 
## 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.02136734

#Pronostico del ARIMA
pronostico_ArimaTX <- forecast(arima_TX, level = 95, h=10)
print(pronostico_ArimaTX)

##      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 34507272

#Grafico del Pronostico
plot(pronostico_ArimaTX, main = "Poblacion en Texas 1990 - 2029 (Pronostico)")

Ejercicio Clase: Mapas

Importar Librerias y Paquetes

#install.packages("maps")

library(maps)

Analisis Descriptivo y Crear un mapa

# Analisis Descriptivo
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.   :39512223

# Tipos de Variables
str(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 ...

#Crear un mapa de EUA
map(database = "state")
map(database = "state", regions = "Texas", col = "red", fill = TRUE, add = TRUE)
map(database = "state", regions = "New York", col = "blue", fill = TRUE, add = TRUE)
title(main = "Pronostico de Poblacion de EUA")

Mapas Poblacionales de EUA por década 1950 - 2050

# Dataset a modificar
poblacion_extendida <- poblacion

# Identificamos los estados únicos
estados <- unique(poblacion$state)

# Para cada estado, entrenamos un modelo ARIMA y hacemos forecast hasta 2050
for(st in estados){
  
  # Filtrar datos de ese estado y ordenarlos por año
  datos_st <- poblacion %>%
    filter(state == st) %>%
    arrange(year)
  
  # Determinar el último año disponible
  ultimo_anio <- max(datos_st$year)
  
  # Crear serie de tiempo
  ts_st <- ts(datos_st$population,
              start = min(datos_st$year),
              end   = ultimo_anio,
              frequency = 1)  # Anual
  
  # Entrenar el modelo ARIMA
  modelo_st <- auto.arima(ts_st)
  
  # Calcular cuántos años hay que pronosticar
  # (solo hacemos forecast si el ultimo_anio es < 2050)
  h_years <- 2050 - ultimo_anio
  
  if(h_years > 0){
    # Hacemos el pronóstico
    pronostico <- forecast(modelo_st, h = h_years)
    
    # Creamos un data frame con los resultados
    anios_pronostico <- (ultimo_anio + 1):2050
    poblacion_pronosticada <- as.numeric(pronostico$mean)
    
    df_forecast <- data.frame(
      state = st,
      year  = anios_pronostico,
      population = poblacion_pronosticada
    )
    
    # Agregamos filas con la población pronosticada
    poblacion_extendida <- rbind(poblacion_extendida, df_forecast)
  }
}

# 3) Función para graficar el mapa para un año dado ------------
plot_map <- function(year) {
  
  # Filtramos los datos para el año solicitado
  data_year <- poblacion_extendida %>%
    filter(year == !!year)
  
  # Cargar datos geográficos de EE.UU.
  states_map <- map_data("state")
  
  # Necesitamos relacionar la sigla (state) con el nombre en minúsculas
  # R trae dos vectores auxiliares: state.abb (siglas) y state.name (nombres completos)
  data_year <- data_year %>%
    mutate(region = tolower(state.name[match(state, state.abb)])) %>%
    right_join(states_map, by = "region")
  
  # Graficar
  ggplot(data_year, aes(x = long, y = lat, group = group, fill = population)) +
    geom_polygon(color = "black") +
    # Gradiente verde -> rojo
    scale_fill_gradient(
      low = "green",   # color mínimo
      high = "red",    # color máximo
      name = "Población"
    ) +
    labs(
      title = paste("Población por Estado en", year)
    ) +
    theme_void() +
    theme(
      legend.position = "right",
      plot.title = element_text(size = 16, face = "bold")
    )
}

# 4) Graficar el mapa cada década entre 1900 y 2050 -------------
for(year in seq(1900, 2050, by = 10)) {
  print(plot_map(year))
}

Actividad 2: Leche Saborizada Hershey’s

Importar librerias y paquetes

#install.packages("readxl")

library(readxl)

Importar Bases de Datos

ventas <- read_excel("Ventas_Históricas_Lechitas.xlsx")

1. Modelo ARIMA

ts_ventas <- ts(ventas$Ventas, start = c(2017,1), frequency = 12) # ST Mensual

autoplot(ts_ventas) +
  labs(title = "Ventas de Leche Saborizada Hershey's", x = "Tiempo", y = "Miles de USD")

#ARIMA
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.2081026

# Generar pronostico de los proximos 12 meses (Año 2020)
pronostico_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.22

# Graficar pronostico ARIMA
autoplot(pronostico_ventas) +
  labs(title = "Pronostico de Ventas 2020 de Leche Saborizada Hershey's",
       x = "Tiempo", y = "Miles de USD")

2. Modelo Regresión Lineal

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-16

siguiente_anio <- data.frame(mes = 37:48)
prediccion_regresion <- predict(regresion_ventas, siguiente_anio)
prediccion_regresion

##        1        2        3        4        5        6        7        8 
## 35934.49 36232.86 36531.23 36829.61 37127.98 37426.35 37724.73 38023.10 
##        9       10       11       12 
## 38321.47 38619.85 38918.22 39216.59

plot(ventas$mes, ventas$Ventas, 
     main = "Pronostico de Ventas 2020 de Leche Saborizada Hershey's",
     xlim = c(1,50),
     ylim = range(c(ventas$Ventas,prediccion_regresion)),
     xlab = "Tiempo",
     ylab = "Miles USD") +
abline(regresion_ventas, col = "blue") +
points(siguiente_anio$mes, prediccion_regresion, col = "red")

## integer(0)

predicciones_reales <- predict(regresion_ventas, ventas)
MAPE_reg <- mean(abs((ventas$Ventas - predicciones_reales)/ ventas$Ventas))*100

MAPE_reg

## [1] 2.011297

3. Elección del Modelo y Conclusión

# MAPE ARIMA
MAPE_arima <- 0.7069542
MAPE_arima

## [1] 0.7069542

# MAPE REG
MAPE_reg

## [1] 2.011297

Conclusion:
La predicción del SARIMA es más acertada que la proyección de la regresión lineal. Al contar con un MAPE menor (0.71) con respecto a la Regresion Lineal (2.01), se elige el SARIMA como el mejor modelo.

Por lo que se puede esperar que las ventas de Hersheys incrementen en el 2020 de la siguiente manera.

Escenarios:

Mes y Año	Realista	Pesimista	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

ventas_anio <- read.csv("ventas_por_anio.csv")

ggplot(ventas_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 USD") +
  scale_x_continuous(breaks = 1:12)

Recomendacion:
- Aumentar el consumo de leche saborizada Hershey’s en el primer semestre del año mediante campañas publicitarias.