options(repos = c(CRAN = "https://cloud.r-project.org"))
#install.packages("readxl")
library(readxl)
#install.packages("dplyr")
library(dplyr)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union#install.packages("tidyr")
library(tidyr)
#install.packages("plm")
library(plm)##
## Attaching package: 'plm'## The following objects are masked from 'package:dplyr':
##
## between, lag, leadEste archivo es un panel de datos empresa-año: cada fila es una combinación de una empresa (cusip) y un año (2012–2021). Es un tipo de dataset muy usado en econometría de innovación (gasto en I+D vs. patentes), popularizado por estudios como los de Hall, Griliches y Hausman sobre patentes y R&D. Diccionario de variables: cusip: identificador de la empresa (no es una variable predictiva, es un ID). merger: 1 si la empresa tuvo una fusión importante ese año, 0 si no. employ: empleados, en miles. return: retorno de la acción, en %. patents: patentes solicitadas en el año. patentsg: patentes concedidas en el año — esta es nuestra variable objetivo (Y). stckpr: precio de la acción. rnd: gasto en I+D, en millones de dólares corrientes. rndeflt: gasto en I+D, en millones de dólares deflactados (precios constantes de 1972). rndstck: “stock” acumulado de I+D (inversión acumulada, no solo del año). sales: ventas, en millones de dólares corrientes. sic: código de industria a 4 dígitos. year: año de la observación.
# Cargamos la base de datos
df1 <- read_xls("C:\\Users\\me\\Desktop\\reto\\bds\\PATENT 3.xls")# Relación con variables explicativas
vars_num <- df1 %>% select(employ, return, patents, stckpr, rnd, rndeflt, rndstck, sales, patentsg)
round(cor(vars_num, use = "complete.obs"), 2)## employ return patents stckpr rnd rndeflt rndstck sales patentsg
## employ 1.00 0.04 0.71 0.38 0.86 0.90 0.85 0.81 0.74
## return 0.04 1.00 0.08 0.31 0.08 0.09 0.07 0.05 0.09
## patents 0.71 0.08 1.00 0.45 0.55 0.61 0.55 0.50 0.94
## stckpr 0.38 0.31 0.45 1.00 0.40 0.45 0.39 0.33 0.47
## rnd 0.86 0.08 0.55 0.40 1.00 0.98 0.99 0.79 0.59
## rndeflt 0.90 0.09 0.61 0.45 0.98 1.00 0.96 0.78 0.64
## rndstck 0.85 0.07 0.55 0.39 0.99 0.96 1.00 0.78 0.59
## sales 0.81 0.05 0.50 0.33 0.79 0.78 0.78 1.00 0.55
## patentsg 0.74 0.09 0.94 0.47 0.59 0.64 0.59 0.55 1.00# Convertimos nuestra bdd a datos panel
df_modelo <- df1 %>%
select(cusip, year, patentsg, rnd, employ, sales, patents, merger)
pdf <- pdata.frame(df_modelo, index = c("cusip", "year"))Variable dependiente (Y): patentsg, las patentes
otorgadas/concedidas a la empresa ese año (no las
solicitadas). Variables explicativas elegidas: rnd (gasto
en I+D), employ (empleados), sales (ventas),
patents (patentes solicitadas ese año) y
merger (si hubo fusión). La idea es explicar cuántas
patentes le conceden a una empresa (patentsg) a
partir de cuánto invierte, qué tan grande es, y cuánto solicitó. A
continuación se presenta un modelo de regresión multiple para tener algo
con qué comparar nuestros siguientes modelos.
modelo_referencia <- lm(patentsg ~ rnd + employ + sales + patents + merger, data = pdf)
summary(modelo_referencia)##
## Call:
## lm(formula = patentsg ~ rnd + employ + sales + patents + merger,
## data = pdf)
##
## Residuals:
## Min 1Q Median 3Q Max
## -266.57 -3.55 -1.49 0.79 644.42
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.3832054 0.5906000 2.342 0.01927 *
## rnd 0.0265024 0.0090625 2.924 0.00349 **
## employ 0.1359658 0.0282759 4.809 1.62e-06 ***
## sales 0.0003127 0.0002654 1.178 0.23884
## patents 0.9653663 0.0112721 85.642 < 2e-16 ***
## merger -2.4621552 4.0822979 -0.603 0.54648
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 25.5 on 2230 degrees of freedom
## (24 observations deleted due to missingness)
## Multiple R-squared: 0.8976, Adjusted R-squared: 0.8973
## F-statistic: 3908 on 5 and 2230 DF, p-value: < 2.2e-16Este modelo muestra una R cuadrada de .89, lo que significa que las variables escogidas explican el 89% del comportamiento de la variable dependiente. Pasamos a los modelos Pooling, Within y Random
form <- patentsg ~ rnd + employ + sales + patents + merger
pooling <- plm(form, data = pdf, model = "pooling")
within <- plm(form, data = pdf, model = "within")
random <- plm(form, data = pdf, model = "random")Resumen de Pooling
summary(pooling)## Pooling Model
##
## Call:
## plm(formula = form, data = pdf, model = "pooling")
##
## Unbalanced Panel: n = 225, T = 8-10, N = 2236
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -266.56790 -3.54922 -1.49325 0.78669 644.42008
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## (Intercept) 1.38320543 0.59059998 2.3420 0.019267 *
## rnd 0.02650238 0.00906254 2.9244 0.003486 **
## employ 0.13596581 0.02827591 4.8085 1.622e-06 ***
## sales 0.00031272 0.00026542 1.1782 0.238838
## patents 0.96536633 0.01127212 85.6420 < 2.2e-16 ***
## merger -2.46215521 4.08229788 -0.6031 0.546484
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 14152000
## Residual Sum of Squares: 1449500
## R-Squared: 0.89758
## Adj. R-Squared: 0.89735
## F-statistic: 3908.43 on 5 and 2230 DF, p-value: < 2.22e-16Resumen de Within / Efectos Fijos
summary(within)## Oneway (individual) effect Within Model
##
## Call:
## plm(formula = form, data = pdf, model = "within")
##
## Unbalanced Panel: n = 225, T = 8-10, N = 2236
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -223.29533 -1.88554 -0.31819 1.51373 262.80141
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## rnd -0.1309985 0.0125503 -10.4379 < 2.2e-16 ***
## employ -0.0663792 0.0607919 -1.0919 0.2750038
## sales -0.0013453 0.0003577 -3.7609 0.0001742 ***
## patents 0.0670273 0.0188908 3.5481 0.0003969 ***
## merger 1.8352766 3.5413607 0.5182 0.6043476
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 715600
## Residual Sum of Squares: 586560
## R-Squared: 0.18033
## Adj. R-Squared: 0.086759
## F-statistic: 88.2654 on 5 and 2006 DF, p-value: < 2.22e-16Resumen de Random / Efectos Aleatorios
summary(random)## Oneway (individual) effect Random Effect Model
## (Swamy-Arora's transformation)
##
## Call:
## plm(formula = form, data = pdf, model = "random")
##
## Unbalanced Panel: n = 225, T = 8-10, N = 2236
##
## Effects:
## var std.dev share
## idiosyncratic 292.4 17.1 1
## individual 0.0 0.0 0
## theta:
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0 0 0 0 0 0
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -266.56790 -3.54922 -1.49325 0.78669 644.42008
##
## Coefficients:
## Estimate Std. Error z-value Pr(>|z|)
## (Intercept) 1.38320543 0.59059998 2.3420 0.019179 *
## rnd 0.02650238 0.00906254 2.9244 0.003451 **
## employ 0.13596581 0.02827591 4.8085 1.52e-06 ***
## sales 0.00031272 0.00026542 1.1782 0.238712
## patents 0.96536633 0.01127212 85.6420 < 2.2e-16 ***
## merger -2.46215521 4.08229788 -0.6031 0.546422
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 14152000
## Residual Sum of Squares: 1449500
## R-Squared: 0.89758
## Adj. R-Squared: 0.89735
## Chisq: 19542.1 on 5 DF, p-value: < 2.22e-16Pruebas de especificación entre los diferentes modelos
pFtest(within, pooling)##
## F test for individual effects
##
## data: form
## F = 13.175, df1 = 224, df2 = 2006, p-value < 2.2e-16
## alternative hypothesis: significant effectsphtest(random, within)##
## Hausman Test
##
## data: form
## chisq = 110.59, df = 5, p-value < 2.2e-16
## alternative hypothesis: one model is inconsistentefectos_fijos <- fixef(within)
head(efectos_fijos)## 800 4626 4671 7500 7603 20753
## 41.2082685 5.9277947 1.4333872 0.8330017 0.5298820 2.0559662Pasamos a la prueba de pronóstico
escenario_2022 <- df_modelo %>%
filter(year == 2021, complete.cases(across(c(rnd, employ, sales, patents, merger))))
fe_df <- data.frame(cusip = names(efectos_fijos), efecto_fijo = as.numeric(efectos_fijos))
# pdata.frame guarda cusip como texto; igualamos el tipo para poder unir
fe_df$cusip <- as.numeric(as.character(fe_df$cusip))
coefs <- coef(within)
pronostico_2022 <- escenario_2022 %>%
inner_join(fe_df, by = "cusip") %>%
mutate(
patentsg_pronostico_2022 =
efecto_fijo +
coefs["rnd"] * rnd +
coefs["employ"] * employ +
coefs["sales"] * sales +
coefs["patents"] * patents +
coefs["merger"] * merger
) %>%
select(cusip, patentsg_2021 = patentsg, patentsg_pronostico_2022)
cat("Empresas con pronóstico calculado:", nrow(pronostico_2022), "de", n_distinct(df_modelo$cusip), "\n")## Empresas con pronóstico calculado: 223 de 226knitr::kable(head(pronostico_2022, 15), digits = 1)| cusip | patentsg_2021 | patentsg_pronostico_2022 |
|---|---|---|
| 800 | 70 | 38.4 |
| 4626 | 7 | 4.5 |
| 4671 | 1 | 1.2 |
| 7500 | 1 | 0.7 |
| 7603 | 0 | 0.2 |
| 20753 | 0 | 1.6 |
| 21367 | 0 | 1.0 |
| 23519 | 6 | 14.5 |
| 29069 | 0 | 0.8 |
| 38213 | 2 | 1.9 |
| 54303 | 2 | 0.9 |
| 67131 | 0 | 2.7 |
| 67383 | 9 | 3.5 |
| 74077 | 18 | 19.5 |
| 77491 | 0 | 0.3 |
summary(pronostico_2022$patentsg_pronostico_2022)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -13.6016 0.8866 3.2913 22.7981 19.5489 772.9755cat("Pronósticos negativos (sin sentido para un conteo):", sum(pronostico_2022$patentsg_pronostico_2022 < 0), "\n")## Pronósticos negativos (sin sentido para un conteo): 8Conclusiones finales para el ejercicio 1:
Las pruebas de especificación (pFtest,
phtest) coinciden en que efectos fijos es el modelo
correcto para estos datos, ya que sí hay características propias de cada
empresa que importan y que están correlacionadas con las variables
explicativas, así que un modelo pooling (que ignora esa identidad de
empresa) queda descartado y un modelo de efectos aleatorios tampoco es
consistente frente al de efectos fijos. El pronóstico de patentes
otorgadas para 2022 se hizo bajo el supuesto de que cada empresa
mantiene en 2022 los mismos valores de I+D, empleo, ventas, patentes
solicitadas y fusión que tuvo en 2021, usando el efecto fijo propio de
cada empresa más los coeficientes del modelo de efectos fijos. Como
patentsg es un conteo y este es un modelo lineal, algunos
pronósticos salen negativos, lo cual no tiene sentido para un número de
patentes pero es una limitación esperada de usar un modelo lineal sobre
una variable de conteo, no un error del código.
Si tuvieras que invertir en alguna sub-categoría, ¿en cuál lo harías? Justifica ampliamente tu respuesta.
df2 <- read.csv("C:\\Users\\me\\Desktop\\reto\\bds\\market_sizes_raw.csv")raw <- read_excel("C:\\Users\\me\\Desktop\\reto\\bds\\Market_sizes.xlsx", sheet = "Statistics Data", skip = 5, col_types = "text")
names(raw)[1] <- "Geography"
raw <- raw %>% filter(!is.na(Geography), Geography == "Mexico")
categorias_producto <- c("Bath and Shower", "Deodorants", "Depilatories", "Fragrances",
"Hair Care", "Men's Grooming", "Skin Care", "Sun Care")
raw <- raw %>% filter(Category %in% categorias_producto)
datos <- raw %>%
select(Category, `2011`:`2025`) %>%
pivot_longer(-Category, names_to = "year", values_to = "value") %>%
mutate(year = as.integer(year), value = suppressWarnings(as.numeric(value)))
knitr::kable(datos %>% filter(year %in% c(2011, 2015, 2019, 2020, 2021, 2025)) %>%
pivot_wider(names_from = year, values_from = value),
digits = 0, caption = "Valores de mercado (MXN millones), años seleccionados")| Category | 2011 | 2015 | 2019 | 2020 | 2021 | 2025 |
|---|---|---|---|---|---|---|
| Bath and Shower | 8410 | 10815 | 12960 | 14634 | 15469 | 21342 |
| Deodorants | 9151 | 12379 | 16745 | 13835 | 15091 | 21824 |
| Depilatories | 802 | 1106 | 1516 | 1491 | 1564 | 1873 |
| Fragrances | 18728 | 22568 | 28483 | 25516 | 31801 | 56742 |
| Hair Care | 24854 | 30416 | 37247 | 36951 | 39555 | 56290 |
| Men’s Grooming | 18672 | 24774 | 32693 | 29669 | 32939 | 49508 |
| Skin Care | 25650 | 30964 | 40638 | 42196 | 47981 | 71830 |
| Sun Care | 1145 | 1693 | 2444 | 2043 | 2294 | 4326 |
Aquí cada subcategoría (Category) hace el papel de
“empresa” y cada año el papel de periodo de tiempo, igual que en el
ejemplo de clase (empresa, anio). Guardamos
una copia numérica del año (anio_num) porque
pdata.frame convierte la columna que se usa como índice de
tiempo, y para usar el año como variable explicativa continua (y poder
pronosticar años futuros como 2026) necesitamos esa copia numérica
aparte.
datos$anio_num <- datos$year
pdf_piel <- pdata.frame(datos, index = c("Category", "year"))Variable dependiente (Y): value (tamaño de mercado, en
MXN millones). Variable explicativa (X): anio_num (el año),
igual que en el ejemplo de clase donde ventas se explicaba
con publicidad.
modelo_referencia_piel <- lm(value ~ anio_num, data = pdf_piel)
summary(modelo_referencia_piel)##
## Call:
## lm(formula = value ~ anio_num, data = pdf_piel)
##
## Residuals:
## Min 1Q Median 3Q Max
## -30062 -10733 -673 13521 39895
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2937049.2 670768.7 -4.379 2.60e-05 ***
## anio_num 1466.2 332.4 4.411 2.29e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 15730 on 118 degrees of freedom
## Multiple R-squared: 0.1415, Adjusted R-squared: 0.1343
## F-statistic: 19.46 on 1 and 118 DF, p-value: 2.287e-05form_piel <- value ~ anio_num
pooling_piel <- plm(form_piel, data = pdf_piel, model = "pooling")
within_piel <- plm(form_piel, data = pdf_piel, model = "within")
random_piel <- plm(form_piel, data = pdf_piel, model = "random")Resumen de Pooling
summary(pooling_piel)## Pooling Model
##
## Call:
## plm(formula = form_piel, data = pdf_piel, model = "pooling")
##
## Balanced Panel: n = 8, T = 15, N = 120
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -30062.20 -10733.24 -672.88 13521.05 39895.00
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## (Intercept) -2937049.22 670768.71 -4.3786 2.598e-05 ***
## anio_num 1466.17 332.39 4.4110 2.287e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 3.4018e+10
## Residual Sum of Squares: 2.9203e+10
## R-Squared: 0.14155
## Adj. R-Squared: 0.13427
## F-statistic: 19.4565 on 1 and 118 DF, p-value: 2.2866e-05Resumen de Within / Efectos Fijos
summary(within_piel)## Oneway (individual) effect Within Model
##
## Call:
## plm(formula = form_piel, data = pdf_piel, model = "within")
##
## Balanced Panel: n = 8, T = 15, N = 120
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -9759.9 -3014.6 -1523.3 2555.2 19311.3
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## anio_num 1466.17 116.12 12.626 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 8.168e+09
## Residual Sum of Squares: 3352800000
## R-Squared: 0.58952
## Adj. R-Squared: 0.55993
## F-statistic: 159.413 on 1 and 111 DF, p-value: < 2.22e-16Resumen de Random / Efectos Aleatorios
summary(random_piel)## Oneway (individual) effect Random Effect Model
## (Swamy-Arora's transformation)
##
## Call:
## plm(formula = form_piel, data = pdf_piel, model = "random")
##
## Balanced Panel: n = 8, T = 15, N = 120
##
## Effects:
## var std.dev share
## idiosyncratic 30205839 5496 0.11
## individual 244180653 15626 0.89
## theta: 0.9096
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -11596.07 -3143.95 -709.43 1841.40 21172.87
##
## Coefficients:
## Estimate Std. Error z-value Pr(>|z|)
## (Intercept) -2937049.22 234403.59 -12.530 < 2.2e-16 ***
## anio_num 1466.17 116.12 12.626 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 8379500000
## Residual Sum of Squares: 3564300000
## R-Squared: 0.57464
## Adj. R-Squared: 0.57104
## Chisq: 159.413 on 1 DF, p-value: < 2.22e-16pFtest(within_piel, pooling_piel)##
## F test for individual effects
##
## data: form_piel
## F = 122.26, df1 = 7, df2 = 111, p-value < 2.2e-16
## alternative hypothesis: significant effectsphtest(random_piel, within_piel)##
## Hausman Test
##
## data: form_piel
## chisq = 4.9585e-12, df = 1, p-value = 1
## alternative hypothesis: one model is inconsistentLa prueba F rechaza que pooling sea válido: sí hay diferencias propias de cada subcategoría que pooling ignora. La prueba de Hausman, en cambio, no rechaza que random effects sea consistente (el estadístico sale prácticamente en cero), así que aquí el modelo elegido es random effects — a diferencia de patentes, donde ganaba efectos fijos.
random effects solo nos da un intercepto y una pendiente
comunes para todas las subcategorías (igual que en el ejemplo de clase,
donde el pronóstico de ventas usaba un único intercepto y
una única pendiente para las 4 empresas). Aquí eso no sirve para
diferenciar entre subcategorías, porque el año 2026 es el mismo para
todas. Como la prueba de Hausman mostró que efectos fijos y random
effects prácticamente coinciden en este caso (la pendiente del año es
idéntica en ambos), usamos fixef() — la misma herramienta
ya usada en el ejercicio de patentes — para obtener el nivel propio de
cada subcategoría y así sí poder pronosticar un valor distinto para cada
una.
efectos_fijos_piel <- fixef(within_piel)
efectos_fijos_piel## Bath and Shower Deodorants Depilatories Fragrances Hair Care
## -2945068 -2943789 -2957351 -2928228 -2921814
## Men's Grooming Skin Care Sun Care
## -2927314 -2916465 -2956364pendiente_piel <- coef(within_piel)["anio_num"]
futuro_piel <- expand.grid(Category = names(efectos_fijos_piel), anio_num = 2026:2028)
futuro_piel$efecto_fijo <- efectos_fijos_piel[as.character(futuro_piel$Category)]
futuro_piel$pronostico <- futuro_piel$efecto_fijo + pendiente_piel * futuro_piel$anio_num
knitr::kable(futuro_piel %>% select(Category, anio_num, pronostico) %>%
pivot_wider(names_from = anio_num, values_from = pronostico),
digits = 0, caption = "Pronóstico de valor de mercado 2026-2028 (MXN millones)")| Category | 2026 | 2027 | 2028 |
|---|---|---|---|
| Bath and Shower | 25383 | 26849 | 28315 |
| Deodorants | 26662 | 28128 | 29594 |
| Depilatories | 13099 | 14565 | 16031 |
| Fragrances | 42222 | 43688 | 45155 |
| Hair Care | 48637 | 50103 | 51569 |
| Men’s Grooming | 43136 | 44602 | 46068 |
| Skin Care | 53985 | 55451 | 56917 |
| Sun Care | 14086 | 15552 | 17019 |
Skin Care. Este modelo asume que todas las subcategorías crecen la misma cantidad de pesos cada año (una sola pendiente común de aproximadamente 1,466 millones de MXN anuales) y solo las distingue por su nivel propio (el efecto fijo de cada una), así que la comparación más honesta que este modelo permite hacer es de nivel, no de tasa de crecimiento. Bajo esa comparación, Skin Care ya es la subcategoría con el valor de mercado más alto en 2025 (cerca de MXN 71,830 millones, muy por encima de Fragrances y Hair Care, las siguientes más grandes) y también la que tiene el efecto fijo más alto del modelo, por lo que el pronóstico 2026-2028 la vuelve a dejar como la más grande. Vale la pena decirlo con honestidad: como el modelo obliga a todas las subcategorías a crecer al mismo ritmo en pesos, el pronóstico de Skin Care para 2028 (alrededor de MXN 56,900 millones) queda incluso por debajo de su valor real de 2025, lo que sugiere que su crecimiento real ha sido más rápido que el promedio que el modelo le está asignando, y que subcategorías pequeñas como Depilatories o Sun Care terminan con pronósticos desproporcionadamente altos frente a su tamaño real — una limitación clara de asumir una sola pendiente para las 8 subcategorías. Aun con esa limitación, la conclusión de invertir en Skin Care se sostiene porque no depende del pronóstico lineal en sí, sino de que ya es, por un margen amplio, la subcategoría más grande y de mayor nivel del mercado de cuidado personal en México, y todo indica que su ritmo de crecimiento real es al menos tan alto como el de las demás subcategorías grandes.
url_gdp <- "https://raw.githubusercontent.com/datasets/gdp/master/data/gdp.csv"
url_infl <- "https://raw.githubusercontent.com/datasets/inflation/master/data/inflation-consumer.csv"
url_pop <- "https://raw.githubusercontent.com/datasets/population/master/data/population.csv"
gdp <- read.csv(url_gdp)
infl <- read.csv(url_infl)
pop <- read.csv(url_pop)
paises <- c("CHN", "JPN", "KOR")
nombres_paises <- c(CHN = "China", JPN = "Japón", KOR = "Corea del Sur")
gdp <- gdp %>% filter(Country.Code %in% paises) %>% select(Country.Code, Year, gdp_usd = Value)
pop <- pop %>% filter(Country.Code %in% paises) %>% select(Country.Code, Year, poblacion = Value)
infl <- infl %>% filter(Country.Code %in% paises) %>% select(Country.Code, Year, inflacion_pct = Inflation)
datos_bm <- gdp %>%
inner_join(pop, by = c("Country.Code", "Year")) %>%
inner_join(infl, by = c("Country.Code", "Year")) %>%
filter(Year >= 1990, Year <= 2023) %>%
mutate(
pais = nombres_paises[Country.Code],
gdp_percapita_usd = gdp_usd / poblacion
) %>%
rename(anio = Year) %>%
select(pais, Country.Code, anio, gdp_usd, poblacion, gdp_percapita_usd, inflacion_pct) %>%
arrange(pais, anio)
cat("Filas:", nrow(datos_bm), " | Países:", n_distinct(datos_bm$pais), " | Años:", min(datos_bm$anio), "-", max(datos_bm$anio), "\n")## Filas: 102 | Países: 3 | Años: 1990 - 2023knitr::kable(datos_bm %>% filter(anio %in% c(1990, 2000, 2010, 2020, 2023)), digits = 1,
caption = "Muestra de los datos importados")| pais | Country.Code | anio | gdp_usd | poblacion | gdp_percapita_usd | inflacion_pct |
|---|---|---|---|---|---|---|
| China | CHN | 1990 | 3.608579e+11 | 1135185000 | 317.9 | 5.7 |
| China | CHN | 2000 | 1.211332e+12 | 1262645000 | 959.4 | 2.1 |
| China | CHN | 2010 | 6.087192e+12 | 1337705000 | 4550.5 | 6.9 |
| China | CHN | 2020 | 1.468774e+13 | 1411100000 | 10408.7 | 0.5 |
| China | CHN | 2023 | 1.779478e+13 | 1410710000 | 12614.1 | -0.6 |
| Corea del Sur | KOR | 1990 | 2.833658e+11 | 42869283 | 6610.0 | 10.1 |
| Corea del Sur | KOR | 2000 | 5.761794e+11 | 47008111 | 12257.0 | 1.0 |
| Corea del Sur | KOR | 2010 | 1.143672e+12 | 49554112 | 23079.3 | 2.7 |
| Corea del Sur | KOR | 2020 | 1.644313e+12 | 51836239 | 31721.3 | 1.6 |
| Corea del Sur | KOR | 2023 | 1.712793e+12 | 51712619 | 33121.4 | 2.1 |
| Japón | JPN | 1990 | 3.185905e+12 | 123478000 | 25801.4 | 2.6 |
| Japón | JPN | 2000 | 4.968359e+12 | 126843000 | 39169.4 | -1.3 |
| Japón | JPN | 2010 | 5.759072e+12 | 128070000 | 44968.2 | -1.9 |
| Japón | JPN | 2020 | 5.055587e+12 | 126261000 | 40040.8 | 0.9 |
| Japón | JPN | 2023 | 4.212945e+12 | 124516650 | 33834.4 | 3.8 |
Panel balanceado: 3 países × 34 años = 102 observaciones, sin datos faltantes.
pdf_bm <- pdata.frame(datos_bm, index = c("pais", "anio"))
modelo_referencia_bm <- lm(gdp_percapita_usd ~ inflacion_pct, data = pdf_bm)
summary(modelo_referencia_bm)##
## Call:
## lm(formula = gdp_percapita_usd ~ inflacion_pct, data = pdf_bm)
##
## Residuals:
## Min 1Q Median 3Q Max
## -29089 -8392 2748 9522 23319
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 26913.7 1455.2 18.495 < 2e-16 ***
## inflacion_pct -2413.5 322.9 -7.474 3.02e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 12190 on 100 degrees of freedom
## Multiple R-squared: 0.3584, Adjusted R-squared: 0.352
## F-statistic: 55.86 on 1 and 100 DF, p-value: 3.02e-11form_bm <- gdp_percapita_usd ~ inflacion_pct
pooling_bm <- plm(form_bm, data = pdf_bm, model = "pooling")
within_bm <- plm(form_bm, data = pdf_bm, model = "within")
random_bm <- plm(form_bm, data = pdf_bm, model = "random")Resumen de Pooling
summary(pooling_bm)## Pooling Model
##
## Call:
## plm(formula = form_bm, data = pdf_bm, model = "pooling")
##
## Balanced Panel: n = 3, T = 34, N = 102
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -29088.7 -8391.5 2748.3 9521.7 23318.5
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## (Intercept) 26913.65 1455.17 18.4952 < 2.2e-16 ***
## inflacion_pct -2413.48 322.93 -7.4738 3.02e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 2.3176e+10
## Residual Sum of Squares: 1.487e+10
## R-Squared: 0.35839
## Adj. R-Squared: 0.35197
## F-statistic: 55.8571 on 1 and 100 DF, p-value: 3.0202e-11Resumen de Within / Efectos Fijos
summary(within_bm)## Oneway (individual) effect Within Model
##
## Call:
## plm(formula = form_bm, data = pdf_bm, model = "within")
##
## Balanced Panel: n = 3, T = 34, N = 102
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -13310.366 -4410.235 -26.126 3237.947 14535.611
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## inflacion_pct -849.76 179.68 -4.7294 7.559e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 4181900000
## Residual Sum of Squares: 3404800000
## R-Squared: 0.18583
## Adj. R-Squared: 0.1609
## F-statistic: 22.3675 on 1 and 98 DF, p-value: 7.559e-06Resumen de Random / Efectos Aleatorios
summary(random_bm)## Oneway (individual) effect Random Effect Model
## (Swamy-Arora's transformation)
##
## Call:
## plm(formula = form_bm, data = pdf_bm, model = "random")
##
## Balanced Panel: n = 3, T = 34, N = 102
##
## Effects:
## var std.dev share
## idiosyncratic 34742880 5894 0.627
## individual 20630819 4542 0.373
## theta: 0.7828
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -13724.627 -5046.723 -19.408 4771.156 14505.859
##
## Coefficients:
## Estimate Std. Error z-value Pr(>|z|)
## (Intercept) 23222.90 2958.95 7.8484 4.215e-15 ***
## inflacion_pct -946.06 193.71 -4.8839 1.040e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 5078300000
## Residual Sum of Squares: 4100300000
## R-Squared: 0.19258
## Adj. R-Squared: 0.18451
## Chisq: 23.852 on 1 DF, p-value: 1.0403e-06En los tres modelos el coeficiente de inflación es negativo y significativo: más inflación se asocia con menor PIB per cápita, consistente con la idea económica de que la inestabilidad de precios va de la mano de menor desarrollo relativo en esta muestra.
pFtest(within_bm, pooling_bm)##
## F test for individual effects
##
## data: form_bm
## F = 165, df1 = 2, df2 = 98, p-value < 2.2e-16
## alternative hypothesis: significant effectsEl p-value es prácticamente 0: se rechaza pooling. Hay efectos propios de cada país que pooling está ignorando, así que efectos fijos le gana a pooling.
phtest(random_bm, within_bm)##
## Hausman Test
##
## data: form_bm
## chisq = 1.7694, df = 1, p-value = 0.1835
## alternative hypothesis: one model is inconsistentAquí el resultado es distinto al del ejercicio de patentes: el p-value es 0.18, no significativo, así que no se rechaza que random effects sea consistente. A diferencia del ejercicio de patentes, aquí sí conviene usar random effects, porque es más eficiente que efectos fijos y la prueba de Hausman no encontró evidencia de que esté sesgado. El modelo elegido es efectos aleatorios (random).
Random effects tiene un solo intercepto y una sola pendiente para los tres países, así que pronosticar es tan simple como en el ejemplo original de clase: intercepto + pendiente * inflación. El único supuesto que necesitamos es un valor de inflación futura por país; usamos el último dato observado (2023) como escenario “si la inflación se mantiene igual”.
intercepto_bm <- coef(random_bm)["(Intercept)"]
pendiente_bm <- coef(random_bm)["inflacion_pct"]
inflacion_2023 <- datos_bm %>% filter(anio == 2023) %>% select(pais, inflacion_pct)
pronostico_bm <- expand.grid(pais = unique(datos_bm$pais), anio = 2024:2026) %>%
left_join(inflacion_2023, by = "pais") %>%
mutate(gdp_percapita_pronostico = intercepto_bm + pendiente_bm * inflacion_pct)
knitr::kable(pronostico_bm %>% pivot_wider(names_from = anio, values_from = gdp_percapita_pronostico,
id_cols = pais),
digits = 0, caption = "PIB per cápita pronosticado 2024-2026 (US$, supone inflación = la de 2023)")| pais | 2024 | 2025 | 2026 |
|---|---|---|---|
| China | 23774 | 23774 | 23774 |
| Corea del Sur | 21270 | 21270 | 21270 |
| Japón | 19634 | 19634 | 19634 |
Como el modelo de random effects solo usa la inflación para explicar el PIB per cápita y le pega el mismo intercepto/pendiente a los tres países, el pronóstico 2024-2026 no mueve mucho el nivel de cada país frente a su valor de 2023, ya que el modelo no tiene forma de “saber” que China ha venido creciendo de forma sostenida si esa tendencia no pasa por la inflación; esto es una limitación real del modelo y vale la pena decirla, no esconderla.
El PIB per cápita de China, Japón y Corea del Sur tiene diferencias de nivel tan grandes que un modelo pooled, que ignora la identidad del país, queda claramente rechazado por la prueba F (p < 0.001). A diferencia del ejercicio de patentes, aquí la prueba de Hausman (p = 0.18) no rechaza efectos aleatorios, así que el modelo elegido es random effects, más eficiente que efectos fijos porque aprovecha también la variación entre países y no solo dentro de cada país en el tiempo. La relación estimada entre inflación y PIB per cápita es negativa y significativa en los tres modelos: cada punto porcentual adicional de inflación se asocia con, aproximadamente, US$946 menos de PIB per cápita según el modelo de random effects. El ajuste del modelo es modesto (R² ≈ 0.19), lo que tiene sentido porque el desarrollo económico depende de muchos factores más (capital humano, instituciones, comercio, tecnología) que no están en este modelo de una sola variable, y el pronóstico 2024-2026 hereda esa misma limitación: como el “shock” de inflación no cambia mucho de un año a otro, el pronóstico tampoco se mueve mucho, lo cual no es un error del código sino una consecuencia esperable de haber usado solo una variable explicativa. Por último, esta correlación negativa entre inflación y PIB per cápita no debe leerse como que “bajar la inflación garantiza” un mayor PIB per cápita: con 3 países y una sola variable, es evidencia sugestiva, no una prueba causal.