Este documento presenta el procedimiento completo para adaptar el modelo hedónico de precios a la canasta efectiva de autos usados del IPC. La idea central es dejar de trabajar con el universo general de avisos y pasar a una muestra restringida a los modelos que efectivamente pertenecen a la canasta. Esto mejora la representatividad del ejercicio y permite defender mejor la metodología.
El documento está organizado por bloques. En cada uno se explica:
El flujo del trabajo sigue la siguiente lógica:
En este bloque se cargan los paquetes necesarios para trabajar con la base, limpiar texto, estimar modelos, calcular errores robustos y revisar diagnósticos. También se incluye una función sencilla para instalar paquetes en caso de que alguno no esté disponible.
Este bloque debe ejecutarse primero porque todo el resto del procedimiento depende de estas librerías.
Conviene comprobar que todas las librerías carguen sin errores. Si alguna falla, el código posterior también fallará.
rm(list = ls())
###############################################################
# BLOQUE 0: LIBRERÍAS
###############################################################
paquetes <- c(
"tidyverse", "readxl", "stringi", "janitor", "broom",
"lmtest", "sandwich", "car", "tseries", "MASS"
)
instalar_si_falta <- function(pkgs){
for(p in pkgs){
if(!requireNamespace(p, quietly = TRUE)) install.packages(p)
}
}
instalar_si_falta(paquetes)
library(tidyverse)
library(readxl)
library(stringi)
library(janitor)
library(broom)
library(lmtest)
library(sandwich)
library(car)
library(tseries)
library(MASS)
options(scipen = 999)Aquí se leen dos archivos:
Además, se limpian los nombres de columnas para evitar problemas por espacios en blanco. Este paso es importante porque varios errores ocurren simplemente porque el nombre de una columna no coincide exactamente con el que se está usando en el código.
Hay que verificar:
MODELO.###############################################################
# BLOQUE 1: CARGA DE DATOS
###############################################################
archivo_base <- "C:/Users/mnorellanam/Documents/Base_A_2022_02_MARZO.xlsx"
archivo_canasta <- "C:/Users/mnorellanam/Documents/Base_B.xlsx"
df <- read_excel(archivo_base, sheet = "Datos MercadoLibre")
canasta <- read_excel(archivo_canasta, sheet = "Hoja1")
names(df) <- trimws(names(df))
names(canasta) <- trimws(names(canasta))
cat("Filas base principal:", nrow(df), "\n")
cat("Filas canasta:", nrow(canasta), "\n")
cat("\nColumnas base principal:\n")
print(names(df))
cat("\nColumnas canasta:\n")
print(names(canasta))Este bloque reúne todas las funciones de apoyo del procedimiento. La lógica es dejar en un solo lugar las funciones que se reutilizan más adelante.
Las funciones incluidas cumplen distintos propósitos:
normalizar_texto(): limpia texto general, como marcas,
títulos o comunas.normalizar_modelo_catalogo(): hace una limpieza más
agresiva sobre nombres de modelos, eliminando expresiones como
NEW, HB, SEDAN, etc.extraer_modelo_titulo(): intenta obtener el modelo
desde el título del aviso.calcular_vif(): revisa multicolinealidad.imprimir_resumen_modelo(): muestra coeficientes
clásicos y robustos.diagnosticos_modelo(): calcula pruebas de normalidad,
heterocedasticidad, especificación y medidas de ajuste.graficos_residuos(): genera gráficos diagnósticos.eliminar_outliers_iterativo(): depura observaciones
influyentes de manera iterativa.Lo importante aquí es que todas las funciones queden definidas correctamente. Si alguna de ellas presenta error, el flujo posterior no va a ejecutarse.
###############################################################
# BLOQUE 2: FUNCIONES AUXILIARES
###############################################################
normalizar_texto <- function(texto){
texto <- as.character(texto)
texto <- trimws(tolower(texto))
texto <- stringi::stri_trans_general(texto, "Latin-ASCII")
texto <- gsub("[[:punct:]]", " ", texto)
texto <- gsub("\\s+", " ", texto)
texto <- trimws(texto)
texto[texto == ""] <- NA_character_
texto
}
normalizar_modelo_catalogo <- function(texto){
texto <- as.character(texto)
texto <- trimws(toupper(texto))
texto <- stringi::stri_trans_general(texto, "Latin-ASCII")
texto <- gsub("-", " ", texto)
texto <- gsub("/", " ", texto)
texto <- gsub("\\.", " ", texto)
texto <- gsub("\\bALL NEW\\b", "", texto)
texto <- gsub("\\bNEW\\b", "", texto)
texto <- gsub("\\bNUEVA\\b", "", texto)
texto <- gsub("\\bNUEVO\\b", "", texto)
texto <- gsub("\\bHB\\b", "", texto)
texto <- gsub("\\bH B\\b", "", texto)
texto <- gsub("\\bSEDAN\\b", "", texto)
texto <- gsub("\\bSPORT\\b", "", texto)
texto <- gsub("\\bNB\\b", "", texto)
texto <- gsub("\\bMT\\b", "", texto)
texto <- gsub("\\bAT\\b", "", texto)
texto <- gsub("\\bMEC\\b", "", texto)
texto <- gsub("\\bAUT\\b", "", texto)
texto <- gsub("\\bLX\\b", "", texto)
texto <- gsub("\\bEX\\b", "", texto)
texto <- gsub("\\bGL\\b", "", texto)
texto <- gsub("\\bGLS\\b", "", texto)
texto <- gsub("I10", "I 10", texto)
texto <- gsub("I20", "I 20", texto)
texto <- gsub("CX([0-9])", "CX \\1", texto)
texto <- gsub("CX\\s*([0-9]{2})", "CX \\1", texto)
texto <- gsub("\\s+", " ", texto)
texto <- trimws(texto)
texto[texto == ""] <- NA_character_
texto
}
extraer_modelo_titulo <- function(titulo, marca_std){
if(is.na(titulo) || is.na(marca_std)) return(NA_character_)
txt <- normalizar_texto(titulo)
marca_txt <- normalizar_texto(marca_std)
txt <- gsub(paste0("^", marca_txt, "\\s+"), "", txt)
palabras <- unlist(strsplit(txt, " "))
if(length(palabras) == 0) return(NA_character_)
if(length(palabras) == 1) return(toupper(palabras[1]))
modelo <- paste(palabras[1:min(2, length(palabras))], collapse = " ")
modelo <- toupper(modelo)
modelo
}
calcular_vif <- function(modelo_ols){
mm <- model.matrix(modelo_ols)
mm <- as.data.frame(mm)
if("(Intercept)" %in% names(mm)) mm[["(Intercept)"]] <- NULL
if(ncol(mm) < 2){
return(data.frame(Variable = names(mm), VIF = NA_real_))
}
vif_vals <- car::vif(modelo_ols)
if(is.matrix(vif_vals)) vif_vals <- vif_vals[, 1]
data.frame(Variable = names(vif_vals), VIF = as.numeric(vif_vals))
}
imprimir_resumen_modelo <- function(modelo_ols, robust_type = "HC1"){
cat("\n--- Resumen OLS clásico ---\n")
print(summary(modelo_ols))
cat(paste0("\n--- Coeficientes con errores robustos ", robust_type, " ---\n"))
tabla_robusta <- lmtest::coeftest(
modelo_ols,
vcov = sandwich::vcovHC(modelo_ols, type = robust_type)
)
print(tabla_robusta)
invisible(tabla_robusta)
}
diagnosticos_modelo <- function(modelo_ols, nombre_modelo = "Modelo"){
residuos <- resid(modelo_ols)
pred <- fitted(modelo_ols)
jb <- tryCatch(tseries::jarque.bera.test(residuos), error = function(e) NULL)
muestra_shapiro <- residuos
if(length(residuos) > 5000){
set.seed(123)
muestra_shapiro <- sample(residuos, 5000)
}
sh <- tryCatch(shapiro.test(muestra_shapiro), error = function(e) NULL)
bp_p <- tryCatch(lmtest::bptest(modelo_ols)$p.value, error = function(e) NA_real_)
white_p <- tryCatch(lmtest::bptest(modelo_ols, ~ fitted(modelo_ols) + I(fitted(modelo_ols)^2))$p.value,
error = function(e) NA_real_)
reset_p <- tryCatch(lmtest::resettest(modelo_ols, power = 2:3, type = "fitted")$p.value,
error = function(e) NA_real_)
vif_tabla <- tryCatch(calcular_vif(modelo_ols), error = function(e) data.frame(Variable = NA, VIF = NA))
max_vif <- suppressWarnings(max(vif_tabla$VIF, na.rm = TRUE))
if(is.infinite(max_vif)) max_vif <- NA_real_
rmse <- sqrt(mean((model.response(model.frame(modelo_ols)) - pred)^2))
mae <- mean(abs(model.response(model.frame(modelo_ols)) - pred))
cat("\n", strrep("=", 70), "\n", sep = "")
cat("DIAGNÓSTICOS -", nombre_modelo, "\n")
cat(strrep("=", 70), "\n", sep = "")
cat("\nNormalidad\n")
cat("Jarque-Bera p-value:", ifelse(is.null(jb), NA, round(jb$p.value, 6)), "\n")
cat("Shapiro-Wilk p-value:", ifelse(is.null(sh), NA, round(sh$p.value, 6)), "\n")
cat("\nHeterocedasticidad\n")
cat("Breusch-Pagan p-value:", bp_p, "\n")
cat("White p-value:", white_p, "\n")
cat("\nEspecificación\n")
cat("RESET Ramsey p-value:", reset_p, "\n")
cat("\nMulticolinealidad (top 15 VIF)\n")
print(head(vif_tabla[order(-vif_tabla$VIF), ], 15))
cat("\nMétricas de ajuste\n")
cat("R2:", round(summary(modelo_ols)$r.squared, 4), "\n")
cat("R2 ajustado:", round(summary(modelo_ols)$adj.r.squared, 4), "\n")
cat("AIC:", round(AIC(modelo_ols), 2), "\n")
cat("BIC:", round(BIC(modelo_ols), 2), "\n")
cat("RMSE log:", round(rmse, 4), "\n")
cat("MAE log:", round(mae, 4), "\n")
list(
jb_p = ifelse(is.null(jb), NA, jb$p.value),
shapiro_p = ifelse(is.null(sh), NA, sh$p.value),
bp_p = bp_p,
white_p = white_p,
reset_p = reset_p,
max_vif = max_vif,
rmse = rmse,
mae = mae,
vif_tabla = vif_tabla
)
}
graficos_residuos <- function(modelo_ols, titulo = "Diagnóstico de residuos"){
residuos <- resid(modelo_ols)
pred <- fitted(modelo_ols)
oldpar <- par(no.readonly = TRUE)
on.exit(par(oldpar))
par(mfrow = c(1, 3))
hist(residuos, breaks = 30, main = paste("Histograma -", titulo), xlab = "Residuos")
qqnorm(residuos, main = paste("QQ-plot -", titulo))
qqline(residuos, col = 2)
plot(pred, residuos,
main = paste("Residuos vs Ajustados -", titulo),
xlab = "Ajustados", ylab = "Residuos")
abline(h = 0, col = 2, lty = 2)
}
eliminar_outliers_iterativo <- function(df_base, formula_modelo,
max_iter = 20,
max_pct_eliminar = 0.15,
min_flags = 2,
verbose = TRUE){
df_work <- df_base
n_inicial <- nrow(df_base)
max_eliminar <- floor(n_inicial * max_pct_eliminar)
historial <- data.frame()
for(i in 1:max_iter){
modelo <- lm(formula_modelo, data = df_work)
cooks <- cooks.distance(modelo)
rstud <- rstudent(modelo)
hats <- hatvalues(modelo)
umbral_cook <- 4 / nrow(df_work)
umbral_hat <- 2 * length(coef(modelo)) / nrow(df_work)
flags <- (abs(rstud) > 3) + (cooks > umbral_cook) + (hats > umbral_hat)
idx_eliminar <- which(flags >= min_flags)
if(verbose){
cat("\nIteración", i,
"- obs actuales:", nrow(df_work),
"- candidatas a eliminar:", length(idx_eliminar), "\n")
}
historial <- bind_rows(
historial,
data.frame(
iteracion = i,
n_actual = nrow(df_work),
n_eliminar = length(idx_eliminar),
cook_umbral = umbral_cook,
hat_umbral = umbral_hat
)
)
if(length(idx_eliminar) == 0) break
if((n_inicial - (nrow(df_work) - length(idx_eliminar))) > max_eliminar){
if(verbose) cat("Se alcanzó el máximo % de eliminación permitido.\n")
break
}
df_work <- df_work[-idx_eliminar, , drop = FALSE]
}
modelo_final <- lm(formula_modelo, data = df_work)
list(
df_final = df_work,
modelo_final_ols = modelo_final,
modelo_final_rob = lmtest::coeftest(modelo_final, vcov = sandwich::vcovHC(modelo_final, type = "HC1")),
historial = historial
)
}En este punto se limpian las variables cuantitativas principales del modelo: precio y kilometraje. Se convierten a formato numérico, se eliminan filas con datos faltantes en variables esenciales y se aplican filtros para evitar valores extremos o poco plausibles.
Luego se construyen LogPrecio y LogKm, ya
que la especificación hedónica se trabaja en forma logarítmica.
Conviene comprobar:
Precio_num y Km_num no tengan
demasiados NA,###############################################################
# BLOQUE 3: LIMPIEZA BASE PRINCIPAL
###############################################################
df <- df %>%
mutate(
Precio_num = as.numeric(gsub("\\.", "", as.character(Precio))),
Km_num = as.numeric(gsub("[^0-9]", "", as.character(Kilometraje)))
)
df <- df %>%
filter(!is.na(Precio_num), !is.na(Km_num), !is.na(Título), !is.na(Marca), !is.na(Comuna))
df <- df %>%
filter(Precio_num >= 1000000, Km_num >= 1500)
p_km <- quantile(df$Km_num, probs = c(0.01, 0.99), na.rm = TRUE)
p_pr <- quantile(df$Precio_num, probs = c(0.01, 0.99), na.rm = TRUE)
df <- df %>%
filter(
between(Km_num, p_km[1], p_km[2]),
between(Precio_num, p_pr[1], p_pr[2])
) %>%
mutate(
LogPrecio = log(Precio_num),
LogKm = log(Km_num)
)Aquí se preparan las variables de texto que luego se usarán para el
cruce con la canasta. Primero se normalizan marca, título y comuna.
Después se intenta extraer el modelo a partir del título del aviso.
Finalmente, se construye una llave compuesta
Marca_Modelo.
Esta etapa es fundamental porque el cruce con la canasta depende de que la base principal tenga una representación razonable del modelo del vehículo.
Conviene revisar una muestra de Marca_Modelo para
confirmar que la extracción desde el título sea coherente.
###############################################################
# BLOQUE 4: VARIABLES DE MARCA / MODELO
###############################################################
df <- df %>%
mutate(
Marca_std = toupper(normalizar_texto(Marca)),
Titulo_norm = normalizar_texto(Título),
Comuna_norm = normalizar_texto(Comuna)
)
df$Modelo_std <- mapply(extraer_modelo_titulo, df$Título, df$Marca_std)
df$Modelo_std <- toupper(df$Modelo_std)
df <- df %>%
mutate(
Marca_Modelo = paste(Marca_std, Modelo_std),
Marca_Modelo = gsub("\\s+", " ", Marca_Modelo),
Marca_Modelo = trimws(Marca_Modelo)
)Este bloque es uno de los más importantes del procedimiento. Aquí se
limpia la columna MODELO de la canasta y se normalizan los
nombres tanto en la canasta como en la base principal. Luego se aplica
una tabla de homologación para corregir diferencias de escritura entre
ambas fuentes.
Este paso es indispensable porque, en la práctica, los nombres rara
vez coinciden de forma exacta. Por ejemplo, puede aparecer
HYUNDAI GRAND I-10 en una fuente y
HYUNDAI GRAND I10 en otra.
Si el filtro por canasta posterior reconoce muy pocos modelos, casi siempre el problema está aquí. En ese caso, hay que ampliar o corregir la tabla de homologación.
###############################################################
# BLOQUE 4.1: CANASTA REAL + HOMOLOGACIÓN
###############################################################
canasta <- canasta %>%
mutate(
MODELO_base = as.character(MODELO),
MODELO_std = normalizar_modelo_catalogo(MODELO_base)
) %>%
filter(!is.na(MODELO_std), MODELO_std != "") %>%
distinct(MODELO_std, .keep_all = TRUE)
df <- df %>%
mutate(
Marca_Modelo_std = normalizar_modelo_catalogo(Marca_Modelo)
)
homologacion <- c(
"CHERY TIGGO 2" = "CHERY TIGGO 2",
"CHERY TIGGO 3" = "CHERY TIGGO 3",
"CHERY TIGGO 8" = "CHERY TIGGO 8",
"CHEVROLET CAPTIVA" = "CHEVROLET CAPTIVA",
"CHEVROLET GROOVE" = "CHEVROLET GROOVE",
"CHEVROLET SAIL NB" = "CHEVROLET SAIL",
"CHEVROLET SAIL" = "CHEVROLET SAIL",
"CHEVROLET TRACKER" = "CHEVROLET TRACKER",
"FORD EXPLORER" = "FORD EXPLORER",
"FORD ECOSPORT" = "FORD ECOSPORT",
"FORD NUEVA ECOSPORT" = "FORD ECOSPORT",
"FORD TERRITORY" = "FORD TERRITORY",
"HYUNDAI ACCENT" = "HYUNDAI ACCENT",
"HYUNDAI GRAND I10" = "HYUNDAI GRAND I 10",
"HYUNDAI GRAND I 10" = "HYUNDAI GRAND I 10",
"KIA MORNING" = "KIA MORNING",
"KIA RIO 4" = "KIA RIO 4",
"KIA SOLUTO" = "KIA SOLUTO",
"MAZDA CX30" = "MAZDA CX 30",
"MAZDA CX 30" = "MAZDA CX 30",
"MAZDA NEW MAZDA2" = "MAZDA 2",
"MAZDA 2" = "MAZDA 2",
"MG ZX" = "MG ZX",
"MG ZS" = "MG ZS",
"NISSAN KICKS" = "NISSAN KICKS",
"NISSAN MARCH" = "NISSAN MARCH",
"NISSAN QASHQAI" = "NISSAN QASHQAI",
"NISSAN VERSA" = "NISSAN VERSA",
"PEUGEOT 301" = "PEUGEOT 301",
"SUZUKI BALENO" = "SUZUKI BALENO",
"SUZUKI SWIFT" = "SUZUKI SWIFT",
"TOYOTA COROLLA CROSS" = "TOYOTA COROLLA CROSS",
"TOYOTA RAIZE" = "TOYOTA RAIZE",
"TOYOTA RAV4" = "TOYOTA RAV4",
"TOYOTA YARIS" = "TOYOTA YARIS",
"VOLKSWAGEN GOL" = "VOLKSWAGEN GOL",
"VOLKSWAGEN NUEVO GOL" = "VOLKSWAGEN GOL",
"VOLKSWAGEN POLO" = "VOLKSWAGEN POLO",
"VOLKSWAGEN T CROSS" = "VOLKSWAGEN T CROSS",
"VOLKSWAGEN T-CROSS" = "VOLKSWAGEN T CROSS"
)
df <- df %>%
mutate(
Marca_Modelo_homologado = ifelse(
Marca_Modelo_std %in% names(homologacion),
homologacion[Marca_Modelo_std],
Marca_Modelo_std
)
)
canasta <- canasta %>%
mutate(
MODELO_homologado = ifelse(
MODELO_std %in% names(homologacion),
homologacion[MODELO_std],
MODELO_std
)
)Una vez homologados los nombres, se filtra la base principal y se conservan únicamente las observaciones cuyo modelo pertenezca a la canasta.
Este paso define el universo efectivo de estimación. A partir de aquí, el modelo ya no se estima sobre todos los autos usados disponibles, sino sobre los autos relevantes para el IPC.
Es importante mirar:
###############################################################
# BLOQUE 5: FILTRO POR CANASTA
###############################################################
modelos_canasta <- unique(canasta$MODELO_homologado)
df_canasta <- df %>%
filter(Marca_Modelo_homologado %in% modelos_canasta)
cat("\n", strrep("=", 70), "\n", sep = "")
cat("FILTRO POR CANASTA\n")
cat(strrep("=", 70), "\n", sep = "")
cat("Observaciones base total :", nrow(df), "\n")
cat("Observaciones base canasta :", nrow(df_canasta), "\n")
cat("Modelos encontrados :", length(unique(df_canasta$Marca_Modelo_homologado)), "\n")
cat("\nFrecuencia por modelo encontrado:\n")
print(sort(table(df_canasta$Marca_Modelo_homologado), decreasing = TRUE))
modelos_no_encontrados <- setdiff(modelos_canasta, unique(df_canasta$Marca_Modelo_homologado))
cat("\nModelos de canasta no encontrados:\n")
print(modelos_no_encontrados)En este bloque se asigna una zona geográfica a cada observación según
su comuna. La variable Zona se utilizará como control
espacial en las especificaciones hedónicas.
La agrupación por zonas permite capturar diferencias sistemáticas de precios asociadas a localización sin sobrecargar el modelo con demasiadas categorías geográficas.
Hay que verificar si quedan comunas sin asignación. Si ocurre, el mapa de comunas debe ampliarse.
###############################################################
# BLOQUE 6: ZONAS
###############################################################
mapa_comuna_zona <- c(
"las condes" = "Oriente",
"vitacura" = "Oriente",
"providencia" = "Oriente",
"lo barnechea" = "Oriente",
"la reina" = "Oriente",
"nunoa" = "Oriente",
"penalolen" = "Oriente",
"macul" = "Oriente",
"la florida" = "Oriente",
"puente alto" = "Oriente",
"san jose de maipo" = "Oriente",
"colina" = "Norte",
"lampa" = "Norte",
"renca" = "Norte",
"quilicura" = "Norte",
"conchali" = "Norte",
"huechuraba" = "Norte",
"independencia" = "Norte",
"recoleta" = "Norte",
"san bernardo" = "Sur",
"el bosque" = "Sur",
"la cisterna" = "Sur",
"la granja" = "Sur",
"san ramon" = "Sur",
"la pintana" = "Sur",
"calera de tango" = "Sur",
"buin" = "Sur",
"paine" = "Sur",
"padre hurtado" = "Sur",
"penaflor" = "Sur",
"talagante" = "Sur",
"curacavi" = "Sur",
"melipilla" = "Sur",
"cerro navia" = "Poniente",
"estacion central" = "Poniente",
"maipu" = "Poniente",
"pudahuel" = "Poniente",
"quinta normal" = "Poniente",
"cerrillos" = "Poniente",
"santiago" = "Poniente",
"san joaquin" = "Poniente",
"san miguel" = "Poniente",
"pedro aguirre cerda" = "Poniente"
)
df_canasta$Zona <- unname(mapa_comuna_zona[df_canasta$Comuna_norm])
df_canasta <- df_canasta %>% filter(!is.na(Zona))La especificación con dummies por modelo es mucho más exigente que las demás. Si algunos modelos tienen muy pocas observaciones, el modelo puede volverse inestable y producir problemas de leverage o errores robustos poco confiables.
Por eso se define un umbral mínimo de observaciones por modelo antes de estimar la especificación C.
Si este umbral es demasiado bajo, la estimación puede presentar problemas numéricos. Si es demasiado alto, la muestra puede reducirse demasiado. Hay que buscar un equilibrio.
###############################################################
# BLOQUE 7: FILTRO DE REPRESENTATIVIDAD PARA MODELO C
###############################################################
MIN_OBS_MODELO <- 12
freq_modelo <- sort(table(df_canasta$Marca_Modelo_homologado), decreasing = TRUE)
modelos_validos <- names(freq_modelo[freq_modelo >= MIN_OBS_MODELO])
df_modelo <- df_canasta %>%
filter(Marca_Modelo_homologado %in% modelos_validos)
cat("\nModelos válidos canasta (>= ", MIN_OBS_MODELO, " obs): ",
length(modelos_validos), "\n", sep = "")
cat("Observaciones para especificación C:", nrow(df_modelo), "\n")Aquí se estiman tres modelos hedónicos alternativos:
Las especificaciones A y B se estiman sobre la base filtrada por canasta. La especificación C se estima sobre una versión más restringida, que además exige representatividad mínima por modelo.
Conviene comparar:
###############################################################
# BLOQUE 8: ESPECIFICACIONES
###############################################################
formula_A <- LogPrecio ~ LogKm + factor(Zona)
formula_B <- LogPrecio ~ LogKm + factor(Zona) + factor(Marca_std)
formula_C <- LogPrecio ~ LogKm + factor(Zona) + factor(Marca_Modelo_homologado)
resultados <- list()
for(esp in c("A", "B", "C")){
cat("\n", strrep("#", 80), "\n", sep = "")
cat("ESTIMANDO MODELO HEDÓNICO — ESPECIFICACIÓN", esp, "\n")
cat(strrep("#", 80), "\n", sep = "")
if(esp == "A"){
base_uso <- df_canasta
formula_modelo <- formula_A
} else if(esp == "B"){
base_uso <- df_canasta
formula_modelo <- formula_B
} else {
base_uso <- df_modelo
formula_modelo <- formula_C
}
modelo_ols <- lm(formula_modelo, data = base_uso)
tabla_robusta <- imprimir_resumen_modelo(modelo_ols, robust_type = "HC1")
diag <- diagnosticos_modelo(modelo_ols, nombre_modelo = paste("Especificación", esp))
graficos_residuos(modelo_ols, titulo = paste("Especificación", esp))
resultados[[esp]] <- list(
df = base_uso,
ols = modelo_ols,
robusta = tabla_robusta,
diag = diag
)
}Este bloque resume las métricas principales de las tres especificaciones para facilitar la comparación. Se reportan tamaño muestral, ajuste, criterios de información y pruebas diagnósticas.
No conviene elegir el modelo solo por R². También hay que mirar RESET, heterocedasticidad, VIF y estabilidad general.
###############################################################
# BLOQUE 9: RESUMEN COMPARATIVO
###############################################################
cat("\n", strrep("=", 70), "\n", sep = "")
cat("CUADRO COMPARATIVO — ESPECIFICACIONES A, B, C\n")
cat(strrep("=", 70), "\n", sep = "")
cat(sprintf("%-25s %10s %10s %10s\n", "", "A", "B", "C"))
cat(strrep("-", 60), "\n", sep = "")
imprimir_fila <- function(nombre, f){
vals <- sapply(c("A", "B", "C"), f)
cat(sprintf("%-25s %10s %10s %10s\n", nombre, vals[1], vals[2], vals[3]))
}
imprimir_fila("N observaciones", function(m) as.character(nrow(resultados[[m]]$df)))
imprimir_fila("R² ajustado", function(m) sprintf("%.4f", summary(resultados[[m]]$ols)$adj.r.squared))
imprimir_fila("AIC", function(m) sprintf("%.1f", AIC(resultados[[m]]$ols)))
imprimir_fila("BIC", function(m) sprintf("%.1f", BIC(resultados[[m]]$ols)))
imprimir_fila("RMSE (log)", function(m) sprintf("%.4f", resultados[[m]]$diag$rmse))
imprimir_fila("RESET p", function(m) sprintf("%.4f", resultados[[m]]$diag$reset_p))
imprimir_fila("Breusch-Pagan p", function(m) sprintf("%.4f", resultados[[m]]$diag$bp_p))
imprimir_fila("White p", function(m) sprintf("%.4f", resultados[[m]]$diag$white_p))
imprimir_fila("Max VIF", function(m) sprintf("%.2f", resultados[[m]]$diag$max_vif))Este bloque elimina observaciones problemáticas de manera iterativa usando criterios de influencia. Luego vuelve a revisar la representatividad de los modelos después de la depuración, ya que algunos pueden quedar con muy pocas observaciones.
Hay que observar:
###############################################################
# BLOQUE 10: DEPURACIÓN DE OUTLIERS
###############################################################
USAR_DEPURACION <- TRUE
if(USAR_DEPURACION){
dep <- eliminar_outliers_iterativo(
df_base = df_modelo,
formula_modelo = formula_C,
max_iter = 25,
max_pct_eliminar = 0.15,
min_flags = 2,
verbose = TRUE
)
df_final_dep <- dep$df_final
modelo_final_ols <- dep$modelo_final_ols
tabla_final_rob <- dep$modelo_final_rob
historial_dep <- dep$historial
freq_post_dep <- sort(table(df_final_dep$Marca_Modelo_homologado), decreasing = FALSE)
print(freq_post_dep)
MIN_OBS_FINAL <- 10
modelos_validos_final <- names(freq_post_dep[freq_post_dep >= MIN_OBS_FINAL])
df_final_dep <- df_final_dep %>%
filter(Marca_Modelo_homologado %in% modelos_validos_final)
modelo_final_ols <- lm(formula_C, data = df_final_dep)
tabla_final_rob <- lmtest::coeftest(
modelo_final_ols,
vcov = sandwich::vcovHC(modelo_final_ols, type = "HC1")
)
} else {
df_final_dep <- df_modelo
modelo_final_ols <- lm(formula_C, data = df_final_dep)
tabla_final_rob <- lmtest::coeftest(
modelo_final_ols,
vcov = sandwich::vcovHC(modelo_final_ols, type = "HC1")
)
}Aquí se presenta la estimación final elegida, junto con sus errores robustos, diagnósticos y distribución final por modelo.
Además, se listan los modelos de la canasta que no fueron encontrados en la base, lo que permite documentar la cobertura efectiva del ejercicio.
Hay que verificar que el modelo final mantenga:
###############################################################
# BLOQUE FINAL: MODELO ELEGIDO
###############################################################
cat("\n", strrep("=", 80), "\n", sep = "")
cat("MODELO FINAL ELEGIDO — ESPECIFICACIÓN C DEPURADA\n")
cat(strrep("=", 80), "\n", sep = "")
print(summary(modelo_final_ols))
print(tabla_final_rob)
diag_final <- diagnosticos_modelo(modelo_final_ols, nombre_modelo = "Modelo final elegido C")
graficos_residuos(modelo_final_ols, titulo = "Modelo final elegido C")
cat("\nDistribución final por modelo:\n")
print(head(sort(table(df_final_dep$Marca_Modelo_homologado), decreasing = TRUE), 30))
cat("\nModelos de canasta no encontrados:\n")
print(setdiff(unique(canasta$MODELO_homologado), unique(df$Marca_Modelo_homologado)))La adaptación del modelo hedónico a la canasta de autos usados requiere una etapa explícita de homologación entre la nomenclatura de la base de avisos y la nomenclatura de la canasta oficial. Ese es el corazón del procedimiento. Una vez resuelto ese punto, el resto del trabajo consiste en asegurar limpieza, representatividad y estabilidad econométrica.
Este enfoque permite estimar un modelo más defendible para fines de ajuste de calidad en el IPC, ya que se centra en los bienes efectivamente relevantes para la medición y no en el universo general de autos usados.