Segmentación y estructura latente del mercado residencial de Cali
Actividad 1 · Modelos Estadísticos para la Toma de Decisiones
Oriana Giraldo Arcia
2026-08-23
# ============================================================
# LIBRERÍAS
# ============================================================
library(paqueteMODELOS) # base de datos `vivienda`
library(dplyr) # manipulación de datos
library(tidyr) # reestructuración de datos
library(stringr) # normalización de cadenas
library(ggplot2) # gráficos
library(ggrepel) # etiquetas sin solapamiento
library(knitr) # tablas
library(kableExtra) # formato de tablas
library(e1071) # asimetría y curtosis
library(nortest) # Anderson-Darling y Lilliefors
library(naniar) # patrón de faltantes y prueba de Little
library(missMDA) # imputación por ACP regularizado
library(robustbase) # covarianza robusta MCD
library(FSA) # comparaciones múltiples de Dunn
library(scales) # formato de ejes
library(FactoMineR) # ACP, AC y ACM
library(cluster) # silueta, gap, CLARA
library(fpc) # estabilidad bootstrap de la partición
# Paquetes usados solo con `::` para evitar enmascaramientos (p. ej.
# psych::alpha sobre ggplot2::alpha). Se comprueba su disponibilidad: omitir
# esta verificación haría fallar el informe a mitad de compilación.
ns_requeridos <- c("psych", "rrcov", "leaflet", "tibble")
ns_faltan <- ns_requeridos[!vapply(ns_requeridos, requireNamespace,
logical(1), quietly = TRUE)]
if (length(ns_faltan) > 0)
stop("Faltan paquetes requeridos (véase el Anexo A1): ",
paste(ns_faltan, collapse = ", "))
select <- dplyr::select # `select` puede quedar enmascarado# ============================================================
# PALETA PASTEL Y TEMA GRÁFICO
# ============================================================
PAS_MENTA <- "#A8D5BA"
PAS_DURAZNO <- "#F7C6A3"
PAS_CIELO <- "#AEC6E8"
PAS_LILA <- "#D9C2E9"
PAS_AMBAR <- "#F6E1A0"
PAS_GRIS <- "#E0E0E0"
TINTA <- "#3A3A3A" # texto, líneas de referencia y ejes
PAL_CAT <- c(PAS_CIELO, PAS_DURAZNO, PAS_MENTA, PAS_LILA, PAS_AMBAR)
# Escalas continuas pastel: de gris muy claro a azul cielo.
escala_cont <- function(...) scale_colour_gradient(low = "#F0F0F0", high = "#7FA8D0", ...)
tema_informe <- theme_minimal(base_size = 12) +
theme(
plot.title = element_text(face = "bold", colour = TINTA, size = 13),
plot.subtitle = element_text(colour = "#6B6B6B", size = 10.5),
axis.title = element_text(colour = "#6B6B6B", size = 10.5),
panel.grid.minor = element_blank(),
panel.grid.major = element_line(colour = "#F0F0F0"),
legend.position = "bottom"
)
theme_set(tema_informe)# ============================================================
# FUNCIONES AUXILIARES
# ============================================================
# Tablas sin franjas ni encabezado coloreado: solo negrita y línea inferior.
t_kbl <- function(x, caption, digits = 3, align = NULL, scroll = FALSE, height = NULL) {
out <- kableExtra::kbl(x, caption = caption, digits = digits, format = "html",
align = align, format.args = list(big.mark = ",")) |>
kableExtra::kable_styling(bootstrap_options = c("hover", "condensed"),
full_width = FALSE, position = "center",
font_size = 13) |>
kableExtra::row_spec(0, bold = TRUE,
extra_css = "border-bottom: 1px solid #BFBFBF;")
if (scroll) out <- kableExtra::scroll_box(out, width = "100%", height = height)
out
}
# V de Cramér: tamaño del efecto asociado al estadístico chi-cuadrado.
cramer_v <- function(tabla_cont) {
if (min(dim(tabla_cont)) < 2) return(NA_real_)
pr <- suppressWarnings(chisq.test(tabla_cont))
as.numeric(sqrt(pr$statistic / (sum(tabla_cont) * (min(dim(tabla_cont)) - 1))))
}
# Índice de Rand ajustado entre dos particiones.
ari <- function(a, b) {
tc <- table(a, b)
c2 <- function(x) x * (x - 1) / 2
s_ij <- sum(c2(tc)); s_i <- sum(c2(rowSums(tc))); s_j <- sum(c2(colSums(tc)))
esperado <- s_i * s_j / c2(sum(tc))
(s_ij - esperado) / ((s_i + s_j) / 2 - esperado)
}
etiqueta_efecto <- function(x, cortes, etiquetas) as.character(cut(x, breaks = cortes, labels = etiquetas))1 Presentación del estudio
1.1 Situación y preguntas
Una empresa inmobiliaria que opera en Cali dispone de un registro de ofertas residenciales con atributos físicos, socioeconómicos y de localización. Las variables se han revisado de forma aislada, de modo que no se conoce la estructura de dependencia entre ellas ni si la oferta admite una división en grupos de producto.
El estudio se organiza en torno a tres preguntas:
| Pregunta | Herramienta | Sección |
|---|---|---|
| ¿Qué factores latentes ordenan la oferta? | Componentes principales | 4 |
| ¿La oferta admite una división en segmentos? | Conglomerados | 5 |
| ¿Cómo se asocian tipo, estrato y barrio? | Correspondencias | 6 |
Las tres herramientas dependen de la escala de medición, de los valores extremos y del tratamiento de la ausencia. La preparación de datos de la Sección 3 es por eso parte del análisis y no un paso previo.
1.2 Objetivos
Objetivo general. Describir la estructura de la oferta residencial urbana mediante técnicas multivariantes, identificando factores latentes, segmentos de producto y patrones de asociación entre atributos categóricos.
Objetivos Específicos.
- Diagnosticar la calidad del registro: duplicación, valores fuera de dominio, inconsistencias, ausencia y valores extremos.
- Identificar los ejes que resumen la variación conjunta de los atributos cuantitativos, trabajando sobre el precio unitario y no sobre el precio total.
- Construir y validar una segmentación robusta de la oferta.
- Establecer los patrones de asociación entre tipo de vivienda, estrato y barrio.
- Traducir los resultados en criterios de gestión para la empresa.
2 Los datos
2.1 Procedencia y alcance
El registro proviene del portal OLX, se obtuvo por webscraping y se distribuye
en el paquete paqueteMODELOS. Recoge ofertas de vivienda en Cali con
coordenadas geográficas.
# ============================================================
# CARGA DE LA BASE
# ============================================================
data("vivienda")
vivienda <- as.data.frame(vivienda)
dim_base <- dim(vivienda)El registro tiene 8,322 filas y 13 columnas.
Dos rasgos del dato acotan el estudio completo. El campo preciom recoge el
precio pedido por el vendedor, no el precio de cierre, de modo que describe la
oferta y no el mercado realizado. Y la captura desde un solo portal no
constituye muestreo probabilístico, así que los contrastes que siguen funcionan
como reglas de decisión reproducibles sobre el conjunto observado. Por eso cada
valor \(p\) va acompañado de un tamaño del efecto, que es la cantidad con sentido
cuando \(n\) ronda los \(10^4\).
# ============================================================
# TABLA: DICCIONARIO DE VARIABLES
# ============================================================
meta <- tibble::tribble(
~Variable, ~Escala, ~Unidad, ~Descripcion,
"id", "Nominal (identificador)", "—", "Identificador del registro de oferta",
"zona", "Nominal", "—", "Zona de la ciudad en que se ubica el inmueble",
"piso", "Ordinal", "n.º de piso", "Piso de ubicación (en principio, solo apartamentos)",
"estrato", "Ordinal", "3 a 6 (obs.)", "Estrato socioeconómico; escala oficial 1–6, observados 3 a 6",
"preciom", "Razón (continua)", "millones COP", "Precio de oferta del inmueble",
"areaconst", "Razón (continua)", "m²", "Área construida",
"parqueaderos", "Razón (discreta)", "conteo", "Número de parqueaderos",
"banios", "Razón (discreta)", "conteo", "Número de baños",
"habitaciones", "Razón (discreta)", "conteo", "Número de habitaciones",
"tipo", "Nominal", "—", "Tipo de vivienda (casa / apartamento)",
"barrio", "Nominal", "—", "Barrio de ubicación",
"longitud", "Razón (continua)", "grados dec.", "Coordenada de longitud geográfica",
"latitud", "Razón (continua)", "grados dec.", "Coordenada de latitud geográfica"
)
diccionario <- meta |>
filter(Variable %in% names(vivienda)) |>
mutate(
`Tipo en R` = vapply(Variable, function(v) class(vivienda[[v]])[1], character(1)),
Niveles = vapply(Variable, function(v) dplyr::n_distinct(vivienda[[v]], na.rm = TRUE), integer(1)),
`NA` = vapply(Variable, function(v) sum(is.na(vivienda[[v]])), integer(1)),
`% NA` = round(100 * `NA` / nrow(vivienda), 2)
) |>
select(Variable, `Tipo en R`, Escala, Unidad, Niveles, `NA`, `% NA`,
`Descripción` = Descripcion)
t_kbl(diccionario,
caption = "Diccionario de variables de la base `vivienda`: escala de medición, cardinalidad y magnitud de la ausencia.",
align = c("l","l","l","l","r","r","r","l"), scroll = TRUE)| Variable | Tipo en R | Escala | Unidad | Niveles | NA | % NA | Descripción |
|---|---|---|---|---|---|---|---|
| id | numeric | Nominal (identificador) | — | 8,319 | 3 | 0.04 | Identificador del registro de oferta |
| zona | character | Nominal | — | 5 | 3 | 0.04 | Zona de la ciudad en que se ubica el inmueble |
| piso | character | Ordinal | n.º de piso | 12 | 2,638 | 31.70 | Piso de ubicación (en principio, solo apartamentos) |
| estrato | numeric | Ordinal | 3 a 6 (obs.) | 4 | 3 | 0.04 | Estrato socioeconómico; escala oficial 1–6, observados 3 a 6 |
| preciom | numeric | Razón (continua) | millones COP | 539 | 2 | 0.02 | Precio de oferta del inmueble |
| areaconst | numeric | Razón (continua) | m² | 652 | 3 | 0.04 | Área construida |
| parqueaderos | numeric | Razón (discreta) | conteo | 10 | 1,605 | 19.29 | Número de parqueaderos |
| banios | numeric | Razón (discreta) | conteo | 11 | 3 | 0.04 | Número de baños |
| habitaciones | numeric | Razón (discreta) | conteo | 11 | 3 | 0.04 | Número de habitaciones |
| tipo | character | Nominal | — | 2 | 3 | 0.04 | Tipo de vivienda (casa / apartamento) |
| barrio | character | Nominal | — | 436 | 3 | 0.04 | Barrio de ubicación |
| longitud | numeric | Razón (continua) | grados dec. | 2,928 | 3 | 0.04 | Coordenada de longitud geográfica |
| latitud | numeric | Razón (continua) | grados dec. | 3,679 | 3 | 0.04 | Coordenada de latitud geográfica |
La escala de medición determina el papel de cada variable. Las de razón son
candidatas a activas. estrato es ordinal y entra como suplementaria, porque
tratarla como métrica supondría equidistancia entre niveles. Las coordenadas se
reservan para la representación espacial. La variable piso recibe un
tratamiento propio en la Sección 3.2.
# ============================================================
# TABLA: COBERTURA DE LA ESCALA DE ESTRATO
# ============================================================
estratos_obs <- sort(unique(na.omit(vivienda$estrato)))
t_kbl(tibble::tibble(
`Escala oficial` = paste(1:6, collapse = ", "),
`Niveles observados` = paste(estratos_obs, collapse = ", "),
`Niveles sin registros` = paste(setdiff(1:6, estratos_obs), collapse = ", "),
`Registros en estratos 1 o 2` = sum(vivienda$estrato <= 2, na.rm = TRUE)),
caption = "Cobertura efectiva de la variable `estrato`.",
align = c("c","c","c","r"))| Escala oficial | Niveles observados | Niveles sin registros | Registros en estratos 1 o 2 |
|---|---|---|---|
| 1, 2, 3, 4, 5, 6 | 3, 4, 5, 6 | 1, 2 | 0 |
El registro no contiene inmuebles de estratos 1 y 2. La consecuencia es de alcance: las conclusiones se refieren al segmento comprendido entre los estratos 3 y 6, la validación externa de los ejes factoriales se apoya en cuatro categorías, y el extremo inferior del mercado queda fuera del estudio.
2.2 Normalización de campos
# ============================================================
# NORMALIZACIÓN Y VARIABLES DERIVADAS
# ============================================================
# Distintas versiones del paquete registran las categorías en mayúscula
# sostenida o en formato título. Sin normalizar, "ZONA NORTE" y "Zona Norte"
# generarían niveles espurios en las tablas de contingencia.
viv <- vivienda |>
mutate(
zona = str_to_title(str_squish(zona)),
tipo = str_to_title(str_squish(tipo)),
barrio = str_to_title(str_squish(barrio)),
piso_num = suppressWarnings(as.integer(piso)),
estrato_num = as.numeric(estrato), # antes de factorizar
estrato = factor(estrato, levels = sort(unique(na.omit(estrato))), ordered = TRUE),
zona = factor(zona),
tipo = factor(tipo),
# `piso` no se descarta: la presencia o ausencia del campo se conserva como
# indicador del anuncio y se usará como variable suplementaria en el ACM.
piso_declarado = factor(ifelse(is.na(piso), "Piso no declarado", "Piso declarado"))
)
num_activas <- c("preciom", "areaconst", "parqueaderos", "banios", "habitaciones")# ============================================================
# TABLA: FRECUENCIAS DE LAS VARIABLES NOMINALES
# ============================================================
niveles_cat <- bind_rows(lapply(c("zona", "tipo"), function(v) {
tibble::tibble(
Variable = v,
`Categoría` = levels(viv[[v]]),
n = as.integer(table(viv[[v]])),
`%` = round(100 * as.numeric(table(viv[[v]])) / sum(!is.na(viv[[v]])), 2)
)
}))
t_kbl(niveles_cat,
caption = "Distribución de frecuencias de las variables nominales.",
align = c("l","l","r","r"))| Variable | Categoría | n | % |
|---|---|---|---|
| zona | Zona Centro | 124 | 1.49 |
| zona | Zona Norte | 1,920 | 23.08 |
| zona | Zona Oeste | 1,198 | 14.40 |
| zona | Zona Oriente | 351 | 4.22 |
| zona | Zona Sur | 4,726 | 56.81 |
| tipo | Apartamento | 5,100 | 61.31 |
| tipo | Casa | 3,219 | 38.69 |
3 Preparación y diagnóstico de los datos
3.1 Consistencia del registro
Se revisan tres familias de problemas: duplicación, valores fuera del dominio admisible e inconsistencias lógicas entre campos.
# ============================================================
# TABLA: AUDITORÍA DE INTEGRIDAD Y CONSISTENCIA
# ============================================================
n_ini <- nrow(viv)
chk <- function(cond) sum(cond, na.rm = TRUE) # cuenta ignorando NA
# Caja de coordenadas del perímetro urbano de Cali
lon_ok <- c(-76.60, -76.44)
lat_ok <- c( 3.32, 3.53)
auditoria <- tibble::tribble(
~`Verificación`, ~Casos,
"Filas completamente duplicadas", chk(duplicated(viv)),
"Identificadores `id` duplicados", chk(duplicated(viv$id)),
"`preciom` menor o igual a cero", chk(viv$preciom <= 0),
"`areaconst` menor o igual a cero", chk(viv$areaconst <= 0),
"`areaconst` inferior a 20 m² (implausible)", chk(viv$areaconst < 20),
"`habitaciones` igual a cero", chk(viv$habitaciones == 0),
"`banios` igual a cero", chk(viv$banios == 0),
"`parqueaderos` negativo", chk(viv$parqueaderos < 0),
"`estrato` fuera del rango 1–6", chk(viv$estrato_num < 1 | viv$estrato_num > 6),
"`estrato` en niveles 1 o 2 (ausentes de la base)", chk(viv$estrato_num <= 2),
"Casa con `piso` registrado", chk(viv$tipo == "Casa" & !is.na(viv$piso)),
"Apartamento sin `piso`", chk(viv$tipo == "Apartamento" & is.na(viv$piso)),
"Longitud fuera del perímetro urbano", chk(viv$longitud < lon_ok[1] | viv$longitud > lon_ok[2]),
"Latitud fuera del perímetro urbano", chk(viv$latitud < lat_ok[1] | viv$latitud > lat_ok[2])
) |>
mutate(`% del total` = round(100 * Casos / n_ini, 2))
t_kbl(auditoria,
caption = "Auditoría de integridad y consistencia lógica de la base original.",
align = c("l","r","r"))| Verificación | Casos | % del total |
|---|---|---|
| Filas completamente duplicadas | 1 | 0.01 |
Identificadores id duplicados
|
2 | 0.02 |
preciom menor o igual a cero
|
0 | 0.00 |
areaconst menor o igual a cero
|
0 | 0.00 |
areaconst inferior a 20 m² (implausible)
|
0 | 0.00 |
habitaciones igual a cero
|
66 | 0.79 |
banios igual a cero
|
45 | 0.54 |
parqueaderos negativo
|
0 | 0.00 |
estrato fuera del rango 1–6
|
0 | 0.00 |
estrato en niveles 1 o 2 (ausentes de la base)
|
0 | 0.00 |
Casa con piso registrado
|
1,965 | 23.61 |
Apartamento sin piso
|
1,381 | 16.59 |
| Longitud fuera del perímetro urbano | 0 | 0.00 |
| Latitud fuera del perímetro urbano | 0 | 0.00 |
Criterio aplicado. Los valores no positivos en preciom, areaconst,
banios y habitaciones son errores de captura y pasan a NA, para recibir el
tratamiento de la Sección 3.2. Suprimir esas filas produciría
sesgo de selección si el error dependiera de las características del inmueble.
Los duplicados exactos sí se retiran, porque duplican la masa de las categorías
en el análisis de correspondencias.
# ============================================================
# DEPURACIÓN
# ============================================================
# (i) Registros inutilizables, por dos criterios:
# a) filas sin ninguna información sustantiva;
# b) filas sin identificación (zona, tipo y barrio ausentes a la vez).
# Un registro del caso (b) no puede localizarse ni clasificarse, y degenera
# cualquier tabla de contingencia: su marginal nulo anula la frecuencia
# esperada y deja el estadístico chi-cuadrado indefinido (0/0).
vars_sustantivas <- setdiff(names(viv), c("id", "piso", "piso_num", "estrato_num",
"piso_declarado"))
sin_informacion <- rowSums(!is.na(viv[vars_sustantivas])) == 0
sin_identificacion <- is.na(viv$zona) & is.na(viv$tipo) & is.na(viv$barrio)
filas_invalidas <- sin_informacion | sin_identificacion
n_invalidas <- sum(filas_invalidas)
# (ii) Errores de captura: valores no positivos en variables cuyo dominio
# los excluye por definición.
viv <- viv[!filas_invalidas, ] |>
mutate(
preciom = ifelse(preciom <= 0, NA_real_, preciom),
areaconst = ifelse(areaconst <= 0, NA_real_, areaconst),
habitaciones = ifelse(habitaciones <= 0, NA_real_, habitaciones),
banios = ifelse(banios <= 0, NA_real_, banios),
parqueaderos = ifelse(parqueaderos < 0, NA_real_, parqueaderos)
) |>
distinct() |>
droplevels()
n_depurado <- nrow(viv)La depuración retira 3 filas inutilizables y los duplicados exactos. Quedan 8,319 registros de 8,322.
3.2 Ausencia de información
El tratamiento depende del mecanismo generador (Rubin, 1976): MCAR, MAR o MNAR. Imputar sin verificarlo contrae la varianza y puede sesgar las estimaciones.
# ============================================================
# TABLA: MAGNITUD DE LA AUSENCIA POR VARIABLE
# ============================================================
resumen_na <- tibble::tibble(
Variable = names(viv),
`Casos faltantes` = vapply(viv, function(x) sum(is.na(x)), integer(1)),
`%` = round(100 * vapply(viv, function(x) mean(is.na(x)), numeric(1)), 2)
) |>
filter(`Casos faltantes` > 0,
!Variable %in% c("piso_num", "estrato_num")) |> # derivadas
arrange(desc(`Casos faltantes`))
t_kbl(resumen_na,
caption = "Variables con valores faltantes tras la depuración, ordenadas por magnitud.",
align = c("l","r","r"))| Variable | Casos faltantes | % |
|---|---|---|
| piso | 2,635 | 31.67 |
| parqueaderos | 1,602 | 19.26 |
| habitaciones | 66 | 0.79 |
| banios | 45 | 0.54 |
# ============================================================
# FIGURA: MAGNITUD DE LA AUSENCIA
# ============================================================
naniar::gg_miss_var(viv |> select(-piso_num, -estrato_num), show_pct = TRUE) +
labs(title = "La ausencia se concentra en un número reducido de variables",
y = "% de registros con valor faltante") +
tema_informe +
theme(axis.text.y = element_text(size = 9))
Figure 3.1: Porcentaje de registros con valor faltante en cada variable de la base depurada.
# ============================================================
# TABLA: PATRONES CONJUNTOS DE AUSENCIA
# ============================================================
# Se codifica cada fila como una cadena de indicadores (1 = observado,
# 0 = faltante) y se cuentan los patrones distintos.
vars_pat <- c(num_activas, "piso")
patron_txt <- apply(!is.na(viv[vars_pat]), 1, function(r) paste0(as.integer(r), collapse = ""))
patrones <- tibble::tibble(Patrón = patron_txt) |>
count(Patrón, name = "Registros") |>
arrange(desc(Registros)) |>
mutate(`% del total` = round(100 * Registros / nrow(viv), 2),
`Variables ausentes` = vapply(Patrón, function(p)
paste(vars_pat[strsplit(p, "")[[1]] == "0"], collapse = ", "),
character(1)),
`Variables ausentes` = ifelse(`Variables ausentes` == "", "ninguna",
`Variables ausentes`))
t_kbl(head(patrones, 10),
caption = "Diez patrones conjuntos de ausencia más frecuentes (orden de las posiciones: preciom, areaconst, parqueaderos, banios, habitaciones, piso).",
align = c("l","r","r","l"), scroll = TRUE)| Patrón | Registros | % del total | Variables ausentes |
|---|---|---|---|
| 111111 | 4,787 | 57.54 | ninguna |
| 111110 | 1,901 | 22.85 | piso |
| 110111 | 863 | 10.37 | parqueaderos |
| 110110 | 692 | 8.32 | parqueaderos, piso |
| 110000 | 21 | 0.25 | parqueaderos, banios, habitaciones, piso |
| 111101 | 14 | 0.17 | habitaciones |
| 110100 | 11 | 0.13 | parqueaderos, habitaciones, piso |
| 110101 | 6 | 0.07 | parqueaderos, habitaciones |
| 111000 | 6 | 0.07 | banios, habitaciones, piso |
| 110001 | 5 | 0.06 | parqueaderos, banios, habitaciones |
3.2.1 El campo piso
# ============================================================
# TABLA: AUSENCIA DE `piso` SEGÚN TIPO DE VIVIENDA
# ============================================================
piso_tipo <- viv |>
group_by(tipo) |>
summarise(
n = n(),
`Con piso registrado` = sum(!is.na(piso)),
`Faltantes en piso` = sum(is.na(piso)),
`% faltante` = round(100 * mean(is.na(piso)), 2),
.groups = "drop"
)
t_kbl(piso_tipo,
caption = "Presencia y ausencia de `piso` según el tipo de vivienda.",
align = c("l","r","r","r","r"))| tipo | n | Con piso registrado | Faltantes en piso | % faltante |
|---|---|---|---|---|
| Apartamento | 5,100 | 3,719 | 1,381 | 27.08 |
| Casa | 3,219 | 1,965 | 1,254 | 38.96 |
Si la ausencia fuera estructural, el porcentaje de faltantes sería del 100 % en casas y cercano a cero en apartamentos. La Tabla 3.4 muestra otro comportamiento: buena parte de las casas declara piso y una fracción comparable de apartamentos no lo hace. El campo no tiene significado unívoco, así que su valor numérico no se usa.
Tratamiento elegido. En lugar de descartar la variable, se conserva el
indicador piso_declarado, que registra si el anuncio diligenció el campo. La
ausencia deja de ser un hueco y pasa a ser un atributo del anuncio, susceptible
de análisis: entra como suplementaria en el análisis de correspondencias
múltiple de la Sección 6.5 para comprobar si se asocia al tipo de
inmueble o al estrato.
3.2.2 Mecanismo en las variables cuantitativas
Prueba de Little (1988). \(H_0\): la ausencia es completamente aleatoria. El estadístico \(d^2\) sigue una \(\chi^2\) bajo \(H_0\), con \(\alpha = 0.05\).
# ============================================================
# TABLA: PRUEBA DE LITTLE PARA MCAR
# ============================================================
datos_little <- viv |> select(all_of(num_activas), estrato_num, longitud, latitud)
little <- tryCatch(naniar::mcar_test(datos_little), error = function(e) NULL)
if (!is.null(little)) {
res_little <- tibble::tibble(
`Estadístico d²` = round(little$statistic, 2),
`gl` = little$df,
`Valor p` = format.pval(little$p.value, digits = 4, eps = 1e-16),
`Patrones` = little$missing.patterns,
`Decisión` = ifelse(little$p.value < 0.05,
"Se rechaza H₀: el mecanismo no es MCAR",
"No se rechaza H₀: compatible con MCAR")
)
t_kbl(res_little, caption = "Prueba de Little para el contraste global de MCAR.",
align = c("r","r","r","r","l"))
} else {
cat("La prueba de Little no pudo calcularse por singularidad de la matriz de covarianzas.")
}| Estadístico d² | gl | Valor p | Patrones | Decisión |
|---|---|---|---|---|
| 1,458 | 44 | < 0.0000000000000001 | 8 | Se rechaza H₀: el mecanismo no es MCAR |
El contraste es global e indica que alguna variable se aparta de MCAR, sin señalar cuál. Se complementa con una prueba por variable frente al estrato, que es la variable observada sin faltantes más asociada con el resto. El tamaño del efecto es la \(V\) de Cramér, \[V=\sqrt{\frac{\chi^2}{n\,\min(r-1,\,c-1)}},\] necesaria porque con \(n\) del orden de \(10^4\) el valor \(p\) rechaza ante diferencias irrelevantes.
# ============================================================
# TABLA: MECANISMO DE AUSENCIA FRENTE AL ESTRATO
# ============================================================
# Solo se contrasta en variables con ausencia superior al 0.5%: por debajo de
# ese umbral la tabla de contingencia es degenerada y la aproximación
# chi-cuadrado deja de ser válida.
vars_na <- names(viv)[vapply(viv, function(x) mean(is.na(x)), numeric(1)) > 0.005]
vars_na <- setdiff(vars_na, c("piso", "piso_num", "estrato_num"))
mecanismo <- bind_rows(lapply(vars_na, function(v) {
tc_v <- table(is.na(viv[[v]]), viv$estrato)
if (min(dim(tc_v)) < 2 || any(rowSums(tc_v) == 0)) {
return(tibble::tibble(Variable = v,
`% faltante` = round(100 * mean(is.na(viv[[v]])), 2),
`χ²` = NA_real_, `gl` = NA_integer_,
`Valor p` = NA_character_, `V de Cramér` = NA_real_,
`Frec. esperada mínima` = NA_real_,
`Lectura` = "Tabla degenerada: prueba no aplicable"))
}
pr <- suppressWarnings(chisq.test(tc_v))
v_c <- cramer_v(tc_v)
tibble::tibble(
Variable = v,
`% faltante` = round(100 * mean(is.na(viv[[v]])), 2),
`χ²` = round(as.numeric(pr$statistic), 2),
`gl` = as.integer(pr$parameter),
`Valor p` = format.pval(pr$p.value, digits = 4, eps = 1e-16),
`V de Cramér` = round(v_c, 4),
`Frec. esperada mínima` = round(min(pr$expected), 2),
`Lectura` = ifelse(pr$p.value >= 0.05, "Compatible con MCAR",
ifelse(v_c < 0.10, "Depende del estrato, con asociación despreciable",
"Depende del estrato: mecanismo al menos MAR"))
)
}))
t_kbl(mecanismo,
caption = "Contraste del mecanismo de ausencia frente al estrato socioeconómico.",
scroll = TRUE)| Variable | % faltante | χ² | gl | Valor p | V de Cramér | Frec. esperada mínima | Lectura |
|---|---|---|---|---|---|---|---|
| parqueaderos | 19.26 | 1,518.73 | 3 | < 0.0000000000000001 | 0.427 | 279.81 | Depende del estrato: mecanismo al menos MAR |
| banios | 0.54 | 8.02 | 3 | 0.04562 | 0.031 | 7.86 | Depende del estrato, con asociación despreciable |
| habitaciones | 0.79 | 16.57 | 3 | 0.0008653 | 0.045 | 11.53 | Depende del estrato, con asociación despreciable |
La frecuencia esperada mínima es condición de validez: por debajo de 5 la aproximación \(\chi^2\) deja de sostenerse.
3.2.3 El caso de parqueaderos
# ============================================================
# TABLA: RECORRIDO OBSERVADO DE `parqueaderos`
# ============================================================
t_kbl(tibble::tibble(
`Mínimo observado` = min(viv$parqueaderos, na.rm = TRUE),
`¿Se observa el valor 0?` = ifelse(any(viv$parqueaderos == 0, na.rm = TRUE), "Sí", "No"),
`Faltantes` = sum(is.na(viv$parqueaderos)),
`% faltante` = round(100 * mean(is.na(viv$parqueaderos)), 2)),
caption = "Recorrido observado de `parqueaderos` y magnitud de su ausencia.",
align = c("r","c","r","r"))| Mínimo observado | ¿Se observa el valor 0? | Faltantes | % faltante |
|---|---|---|---|
| 1 | No | 1,602 | 19.26 |
El valor cero no aparece nunca entre los datos observados, mientras que una fracción alta de anuncios omite el campo. Cabe la hipótesis de que la omisión codifique el cero, es decir, un mecanismo MNAR.
Qué permite decidir la evidencia. Que los registros sin dato difieran en precio, área o estrato no separa MAR de MNAR: bajo MAR la probabilidad de ausencia depende justamente de las variables observadas. Los contrastes solo descartan MCAR. Distinguir MAR de MNAR exigiría conocer los valores ausentes.
# ============================================================
# TABLA: PERFIL DE LOS REGISTROS CON Y SIN DATO
# ============================================================
comparacion_pq <- viv |>
mutate(Grupo = ifelse(is.na(parqueaderos), "Sin dato de parqueaderos",
"Con dato de parqueaderos")) |>
group_by(Grupo) |>
summarise(
n = n(),
`Precio mediano` = median(preciom, na.rm = TRUE),
`Área mediana` = median(areaconst, na.rm = TRUE),
`Estrato mediano` = median(estrato_num, na.rm = TRUE),
`Baños medianos` = median(banios, na.rm = TRUE),
`% casas` = round(100 * mean(tipo == "Casa", na.rm = TRUE), 1),
.groups = "drop"
)
t_kbl(comparacion_pq,
caption = "Perfil comparativo de los registros con y sin dato de parqueaderos.",
digits = 1, scroll = TRUE)| Grupo | n | Precio mediano | Área mediana | Estrato mediano | Baños medianos | % casas |
|---|---|---|---|---|---|---|
| Con dato de parqueaderos | 6,717 | 355 | 130 | 5 | 3 | 37.0 |
| Sin dato de parqueaderos | 1,602 | 179 | 90 | 4 | 2 | 45.8 |
# ============================================================
# TABLA: CONTRASTE DE MANN-WHITNEY SOBRE EL PRECIO
# ============================================================
# Se usa Mann-Whitney por la asimetría documentada. Tamaño del efecto:
# r biserial de rangos en la convención de Kerby, r = 2U/(n1*n2) - 1, acotado
# en [-1, 1]. Con U = W de wilcox.test(sin, con), r < 0 indica que el primer
# grupo se sitúa por debajo del segundo.
grupo_sin <- viv$preciom[is.na(viv$parqueaderos)]
grupo_con <- viv$preciom[!is.na(viv$parqueaderos)]
mw <- suppressWarnings(wilcox.test(grupo_sin, grupo_con))
n1 <- sum(!is.na(grupo_sin)); n2 <- sum(!is.na(grupo_con))
r_bis <- 2 * as.numeric(mw$statistic) / (n1 * n2) - 1
t_kbl(tibble::tibble(
`Estadístico U` = round(as.numeric(mw$statistic), 0),
`Valor p` = format.pval(mw$p.value, digits = 4, eps = 1e-16),
`r biserial de rangos` = round(r_bis, 4),
`Magnitud` = etiqueta_efecto(abs(r_bis), c(-Inf, .10, .30, .50, Inf),
c("Despreciable","Pequeño","Moderado","Grande")),
`Dirección` = ifelse(r_bis < 0, "Los registros sin dato presentan precios menores",
"Los registros sin dato presentan precios mayores"),
`Lectura` = ifelse(mw$p.value < 0.05 & abs(r_bis) >= 0.10,
"La ausencia se asocia a variables observadas: se descarta MCAR",
"Sin diferencia sustantiva: compatible con MCAR")),
caption = "Contraste de Mann-Whitney del precio de oferta entre registros con y sin dato de parqueaderos.",
scroll = TRUE)| Estadístico U | Valor p | r biserial de rangos | Magnitud | Dirección | Lectura |
|---|---|---|---|---|---|
| 2,903,476 | < 0.0000000000000001 | -0.46 | Moderado | Los registros sin dato presentan precios menores | La ausencia se asocia a variables observadas: se descarta MCAR |
Tratamiento elegido: imputación, no recodificación a cero. La lectura MNAR
es plausible, pero no es contrastable, y convertir en ceros el 19 % del registro
equivale a inventar un atributo observado a partir de una conjetura. Se prefiere
el tratamiento conservador: parqueaderos se imputa junto con las demás
variables cuantitativas, de modo que el valor asignado se deduce de la
estructura de correlación observada y no de una hipótesis sobre el formulario
del portal. La recodificación a cero se conserva como escenario alternativo y
se somete a comparación explícita en la Sección 4.5: si la
solución factorial coincide bajo ambos tratamientos, la elección deja de ser
crítica.
# ============================================================
# ESCENARIO ALTERNATIVO PARA EL ANÁLISIS DE SENSIBILIDAD
# ============================================================
# Se guarda la versión con la ausencia recodificada a cero. No sustituye a la
# variable original: se usará solo para comparar soluciones factoriales.
viv$parqueaderos_cero <- ifelse(is.na(viv$parqueaderos), 0, viv$parqueaderos)
t_kbl(tibble::tibble(
`Tratamiento` = c("Imputación (principal)", "Recodificación a cero (alternativo)"),
`Valores asignados` = c(sum(is.na(viv$parqueaderos)), sum(is.na(viv$parqueaderos))),
`% del registro` = round(100 * mean(is.na(viv$parqueaderos)), 2),
`Supuesto` = c("La ausencia se explica por la estructura de correlación observada",
"Toda ausencia corresponde a la carencia del atributo")),
caption = "Los dos tratamientos considerados para la ausencia en `parqueaderos`.",
align = c("l","r","r","l"), scroll = TRUE)| Tratamiento | Valores asignados | % del registro | Supuesto |
|---|---|---|---|
| Imputación (principal) | 1,602 | 19.26 | La ausencia se explica por la estructura de correlación observada |
| Recodificación a cero (alternativo) | 1,602 | 19.26 | Toda ausencia corresponde a la carencia del atributo |
3.2.4 Imputación
Criterio. La imputación por la media es inadmisible: el ACP estima la matriz
de correlaciones, y la media contrae la varianza y atenúa las correlaciones
hacia cero. La eliminación por lista, comportamiento por defecto de prcomp,
descartaría cerca de una quinta parte del registro y sesga fuera de MCAR. Se
emplea missMDA::imputePCA() (Josse & Husson, 2016), que reconstruye los
valores ausentes con las primeras \(S\) componentes bajo regularización, con \(S\)
elegido por validación cruzada.
# ============================================================
# IMPUTACIÓN POR ACP REGULARIZADO
# ============================================================
X_num <- viv |> select(all_of(num_activas))
mascara_na <- is.na(X_num) # posiciones efectivamente imputadas
ncp_cv <- missMDA::estim_ncpPCA(X_num, ncp.min = 0, ncp.max = 4,
method.cv = "Kfold", nbsim = 20, verbose = FALSE)
imp <- missMDA::imputePCA(X_num, ncp = ncp_cv$ncp, scale = TRUE)
X_imp <- as.data.frame(imp$completeObs)
# Restricción de dominio aplicada ÚNICAMENTE a las celdas imputadas. Alterar un
# valor observado dentro del paso de imputación sería una edición silenciosa de
# los datos originales.
redondea_conteo <- function(col, mask, minimo) { col[mask] <- pmax(round(col[mask]), minimo); col }
acota_continua <- function(col, mask, minimo) { col[mask] <- pmax(col[mask], minimo); col }
X_imp$preciom <- acota_continua(X_imp$preciom, mascara_na[, "preciom"],
min(X_num$preciom, na.rm = TRUE))
X_imp$areaconst <- acota_continua(X_imp$areaconst, mascara_na[, "areaconst"],
min(X_num$areaconst, na.rm = TRUE))
X_imp$parqueaderos <- redondea_conteo(X_imp$parqueaderos, mascara_na[, "parqueaderos"], 0)
X_imp$banios <- redondea_conteo(X_imp$banios, mascara_na[, "banios"], 1)
X_imp$habitaciones <- redondea_conteo(X_imp$habitaciones, mascara_na[, "habitaciones"], 1)
ncp_elegido <- ncp_cv$ncp
# Verificación explícita de que ningún valor observado fue modificado
sin_alterar <- all(mapply(function(a, b, m) isTRUE(all.equal(a[!m], b[!m])),
X_num, X_imp, as.data.frame(mascara_na)))El número de dimensiones seleccionado es \(S = 4\). La verificación confirma que ningún valor observado se alteró: conforme.
# ============================================================
# TABLA: VERIFICACIÓN DE LA IMPUTACIÓN
# ============================================================
comparacion_imp <- bind_rows(lapply(num_activas, function(v) {
obs <- X_num[[v]][!is.na(X_num[[v]])]
imp_v <- X_imp[[v]]
tibble::tibble(
Variable = v,
`n imputados` = sum(is.na(X_num[[v]])),
`Media obs.` = mean(obs),
`Media final` = mean(imp_v),
`DE obs.` = sd(obs),
`DE final` = sd(imp_v),
`Δ% en DE` = round(100 * (sd(imp_v) - sd(obs)) / sd(obs), 2),
# Contracción que produciría la imputación por la media, en forma cerrada:
# s_media/s_obs = sqrt((n_obs - 1)/(n - 1)), porque los valores imputados
# no aportan suma de cuadrados. Es la referencia frente a la que debe
# juzgarse el poder discriminante de la columna anterior.
`Δ% bajo imputación por la media` =
round(100 * (sqrt((length(obs) - 1) / (nrow(X_num) - 1)) - 1), 3)
)
}))
t_kbl(comparacion_imp,
caption = "Verificación de la imputación: preservación de media y dispersión.",
digits = 2, scroll = TRUE)| Variable | n imputados | Media obs. | Media final | DE obs. | DE final | Δ% en DE | Δ% bajo imputación por la media |
|---|---|---|---|---|---|---|---|
| preciom | 0 | 433.90 | 433.90 | 328.67 | 328.67 | 0.00 | 0.00 |
| areaconst | 0 | 174.93 | 174.93 | 142.96 | 142.96 | 0.00 | 0.00 |
| parqueaderos | 1,602 | 1.84 | 1.75 | 1.12 | 1.07 | -4.96 | -10.14 |
| banios | 45 | 3.13 | 3.13 | 1.41 | 1.41 | 0.04 | -0.27 |
| habitaciones | 66 | 3.63 | 3.64 | 1.43 | 1.43 | -0.01 | -0.40 |
Lectura de la verificación. En banios y habitaciones la ausencia es
inferior al 1 %, así que la conservación de media y dispersión era esperable con
casi cualquier método. El caso informativo es parqueaderos, donde se imputa
cerca del 19 % de los valores: la columna de desviación estándar permite
comprobar que el procedimiento no comprimió la dispersión, que es el riesgo
específico de imputar volúmenes altos. La última columna calcula, en forma
cerrada, la contracción que habría producido la imputación por la media,
\(s_{\text{media}}/s_{\text{obs}}=\sqrt{(n_{\text{obs}}-1)/(n-1)}\), y sirve de
referencia frente a la cual juzgar el resultado.
# ============================================================
# BASE CONSOLIDADA Y PRECIO UNITARIO
# ============================================================
viv_c <- viv
viv_c[num_activas] <- X_imp
# Precio por metro cuadrado. Es la variable de valoración habitual en el
# mercado inmobiliario y separa dos efectos que `preciom` confunde: cuánto
# vale el inmueble por su tamaño y cuánto vale por su calidad y localización.
viv_c$precio_m2 <- viv_c$preciom / viv_c$areaconst3.3 Valores extremos
3.3.1 Detección univariada
Se combinan dos reglas. La de Tukey (1977) marca los valores fuera de \([Q_1 - k\,\mathrm{IQR},\; Q_3 + k\,\mathrm{IQR}]\) con \(k=1.5\) y \(k=3\). La puntuación \(z\) robusta (Iglewicz & Hoaglin, 1993) usa mediana y desviación absoluta mediana con umbral 3.5: \[z_i^{\mathrm{rob}}=\frac{0.6745\,\bigl(x_i-\mathrm{Med}(x)\bigr)}{\mathrm{MAD}(x)}.\] La regla robusta se prefiere a la clásica porque media y desviación estándar están contaminadas por los propios valores que se busca detectar.
Detalle de implementación. mad() aplica constant = 1.4826 por defecto y
ya devuelve un estimador consistente de \(\sigma\). Multiplicar además por 0.6745
duplicaría la corrección y relajaría el umbral de \(5.19\,\mathrm{MAD}\) a
\(7.69\,\mathrm{MAD}\). Se usa mad(x, constant = 1), y se reporta la MAD sin
escalar porque se anula cuando más de la mitad de los datos coincide con la
mediana, caso en que la regla no está definida.
# ============================================================
# TABLA: ATÍPICOS UNIVARIADOS
# ============================================================
detecta_atipicos <- function(x, nombre) {
q <- quantile(x, c(0.25, 0.75), na.rm = TRUE)
iqr <- q[2] - q[1]
mad_x <- mad(x, na.rm = TRUE, constant = 1) # MAD sin reescalar
z_rob <- if (mad_x > 0) 0.6745 * (x - median(x, na.rm = TRUE)) / mad_x
else rep(NA_real_, length(x))
tibble::tibble(
Variable = nombre,
Q1 = q[1],
Q3 = q[2],
IQR = iqr,
`MAD (sin escalar)` = mad_x,
`Leves (k=1.5)` = sum(x < q[1] - 1.5*iqr | x > q[2] + 1.5*iqr, na.rm = TRUE),
`Extremos (k=3)` = sum(x < q[1] - 3*iqr | x > q[2] + 3*iqr, na.rm = TRUE),
`z robusto > 3.5` = if (mad_x > 0) sum(abs(z_rob) > 3.5, na.rm = TRUE) else NA_integer_,
`% extremos` = round(100 * sum(x < q[1] - 3*iqr | x > q[2] + 3*iqr, na.rm = TRUE) / length(x), 2)
)
}
atipicos_uni <- bind_rows(lapply(num_activas, function(v) detecta_atipicos(viv_c[[v]], v)))
t_kbl(atipicos_uni,
caption = "Detección de valores atípicos univariados por los criterios de Tukey y z robusto.",
digits = 2, scroll = TRUE)| Variable | Q1 | Q3 | IQR | MAD (sin escalar) | Leves (k=1.5) | Extremos (k=3) | z robusto > 3.5 | % extremos |
|---|---|---|---|---|---|---|---|---|
| preciom | 220 | 540 | 320 | 140 | 552 | 132 | 537 | 1.59 |
| areaconst | 80 | 229 | 149 | 57 | 382 | 94 | 517 | 1.13 |
| parqueaderos | 1 | 2 | 1 | 0 | 597 | 115 | — | 1.38 |
| banios | 2 | 4 | 2 | 1 | 72 | 0 | 24 | 0.00 |
| habitaciones | 3 | 4 | 1 | 1 | 833 | 272 | 134 | 3.27 |
# ============================================================
# FIGURA: DIAGRAMAS DE CAJA
# ============================================================
viv_c |>
select(all_of(num_activas)) |>
pivot_longer(everything(), names_to = "Variable", values_to = "Valor") |>
ggplot(aes(x = "", y = Valor)) +
geom_boxplot(fill = PAS_CIELO, colour = "#7A7A7A",
outlier.colour = PAS_DURAZNO, outlier.alpha = 0.35, outlier.size = 0.9) +
facet_wrap(~ Variable, scales = "free_y", nrow = 1) +
labs(title = "La asimetría de precio y área domina la estructura univariada",
x = NULL, y = NULL) +
theme(axis.text.x = element_blank())
Figure 3.2: Distribución de las variables cuantitativas activas. En ocre, los valores atípicos por el criterio de Tukey con k = 1.5.
3.3.2 Detección multivariada
Un inmueble puede ser corriente en cada variable y atípico en la combinación,
por ejemplo un área grande con precio bajo. Como el ACP y la segmentación
operan sobre la estructura conjunta, esta es la detección que importa. La
distancia de Mahalanobis clásica sufre enmascaramiento, porque los propios
atípicos inflan el centro y la covarianza. Se usa el estimador MCD
(Rousseeuw & Van Driessen, 1999) con alpha = 0.75.
# ============================================================
# TABLA: MAHALANOBIS CLÁSICA FRENTE A ROBUSTA (MCD)
# ============================================================
X_mcd <- as.matrix(viv_c[num_activas])
mcd <- robustbase::covMcd(X_mcd, alpha = 0.75)
d2_rob <- mahalanobis(X_mcd, center = mcd$center, cov = mcd$cov)
d2_cla <- mahalanobis(X_mcd, center = colMeans(X_mcd), cov = cov(X_mcd))
p_dim <- length(num_activas)
corte <- qchisq(0.975, df = p_dim)
viv_c$atipico_mv <- d2_rob > corte
t_kbl(tibble::tibble(
Criterio = c("Mahalanobis clásica", "Mahalanobis robusta (MCD)"),
`Umbral χ²(0.975, 5)` = round(corte, 3),
`Atípicos detectados` = c(sum(d2_cla > corte), sum(d2_rob > corte)),
`% del total` = round(100 * c(sum(d2_cla > corte), sum(d2_rob > corte)) / nrow(X_mcd), 2)),
caption = "Comparación entre la detección clásica y la robusta de atípicos multivariados.",
align = c("l","r","r","r"))| Criterio | Umbral χ²(0.975, 5) | Atípicos detectados | % del total |
|---|---|---|---|
| Mahalanobis clásica | 12.83 | 751 | 9.03 |
| Mahalanobis robusta (MCD) | 12.83 | 2,591 | 31.15 |
El corte \(\chi^2_{0.975,p}\) marcaría cerca del 2.5 % de las observaciones incluso bajo normalidad multivariada exacta, que aquí se rechaza. Además MCD ajusta el modelo sobre el subconjunto más compacto, de manera que con una distribución asimétrica la cola superior queda fuera por construcción. Una proporción alta de casos señalados indica asimetría, no contaminación.
# ============================================================
# TABLA: PERFIL DE LOS ATÍPICOS MULTIVARIADOS
# ============================================================
perfil_at <- viv_c |>
group_by(Grupo = ifelse(atipico_mv, "Atípico multivariado", "Resto")) |>
summarise(
n = n(),
`Precio mediano` = median(preciom),
`Área mediana` = median(areaconst),
`Estrato modal` = names(sort(table(estrato), decreasing = TRUE))[1],
`% casas` = round(100 * mean(tipo == "Casa", na.rm = TRUE), 1),
.groups = "drop"
)
t_kbl(perfil_at,
caption = "Perfil comparativo de los registros marcados como atípicos multivariados.",
digits = 1)| Grupo | n | Precio mediano | Área mediana | Estrato modal | % casas |
|---|---|---|---|---|---|
| Atípico multivariado | 2,591 | 650 | 300 | 6 | 77.5 |
| Resto | 5,728 | 265 | 93 | 5 | 21.2 |
Tratamiento elegido: se conservan. El grupo señalado tiene precio y área muy superiores, estrato modal máximo y predominio de casas, es decir, corresponde a un segmento real de alto valor. El control de su influencia se hace por dos vías: estimación robusta del ACP en la Sección 4.5 y uso de medoides, y no de centroides, en la segmentación.
3.4 Descripción univariada
# ============================================================
# TABLA: ESTADÍSTICOS DESCRIPTIVOS
# ============================================================
resumen_estadistico <- function(x, nombre) {
tibble::tibble(
Variable = nombre,
n = sum(!is.na(x)),
Media = mean(x, na.rm = TRUE),
DE = sd(x, na.rm = TRUE),
`CV (%)` = 100 * sd(x, na.rm = TRUE) / mean(x, na.rm = TRUE),
`Mínimo` = min(x, na.rm = TRUE),
Q1 = quantile(x, 0.25, na.rm = TRUE),
Mediana = median(x, na.rm = TRUE),
Q3 = quantile(x, 0.75, na.rm = TRUE),
`Máximo` = max(x, na.rm = TRUE),
`Asimetría` = e1071::skewness(x, na.rm = TRUE, type = 2),
Curtosis = e1071::kurtosis(x, na.rm = TRUE, type = 2)
)
}
descriptivos <- bind_rows(lapply(c(num_activas, "precio_m2"), function(v) resumen_estadistico(viv_c[[v]], v)))
t_kbl(descriptivos,
caption = "Estadísticos descriptivos de las variables cuantitativas activas.",
digits = 2, scroll = TRUE)| Variable | n | Media | DE | CV (%) | Mínimo | Q1 | Mediana | Q3 | Máximo | Asimetría | Curtosis |
|---|---|---|---|---|---|---|---|---|---|---|---|
| preciom | 8,319 | 433.90 | 328.67 | 75.75 | 58.00 | 220.00 | 330.00 | 540.00 | 1,999.00 | 1.85 | 3.68 |
| areaconst | 8,319 | 174.93 | 142.96 | 81.72 | 30.00 | 80.00 | 123.00 | 229.00 | 1,745.00 | 2.69 | 12.93 |
| parqueaderos | 8,319 | 1.75 | 1.07 | 60.99 | 1.00 | 1.00 | 1.00 | 2.00 | 10.00 | 2.41 | 9.06 |
| banios | 8,319 | 3.13 | 1.41 | 45.15 | 1.00 | 2.00 | 3.00 | 4.00 | 10.00 | 0.98 | 1.14 |
| habitaciones | 8,319 | 3.64 | 1.43 | 39.27 | 1.00 | 3.00 | 3.00 | 4.00 | 10.00 | 1.81 | 4.14 |
| precio_m2 | 8,319 | 2.72 | 1.08 | 39.75 | 0.15 | 1.92 | 2.64 | 3.38 | 9.47 | 0.68 | 0.80 |
La asimetría se calcula con el estimador insesgado y la curtosis se reporta como exceso sobre la normal. Precio y área tienen asimetría positiva marcada y coeficientes de variación por encima del 100 %, así que la mediana sustituye a la media como medida de posición en el resto del informe.
3.5 Forma de la distribución
Hipótesis. \(H_0\): la variable procede de una población normal. Se aplican Anderson-Darling, con mayor sensibilidad en las colas, y Lilliefors. Se descarta Shapiro-Wilk porque su implementación en R se limita a \(n \le 5000\).
# ============================================================
# TABLA: PRUEBAS DE NORMALIDAD
# ============================================================
prueba_normalidad <- function(x, nombre) {
x <- x[!is.na(x)]
ad <- nortest::ad.test(x)
lil <- nortest::lillie.test(x)
tibble::tibble(
Variable = nombre,
`AD (p)` = format.pval(ad$p.value, digits = 3, eps = 1e-16),
`Lilliefors (p)` = format.pval(lil$p.value, digits = 3, eps = 1e-16),
`D de Lilliefors` = round(as.numeric(lil$statistic), 4),
`Asimetría` = round(e1071::skewness(x, type = 2), 3),
`Decisión (α=0.05)` = ifelse(ad$p.value < 0.05, "Se rechaza H₀", "No se rechaza H₀")
)
}
normalidad <- bind_rows(lapply(c(num_activas, "precio_m2"), function(v) prueba_normalidad(viv_c[[v]], v)))
t_kbl(normalidad,
caption = "Contraste de normalidad de las variables cuantitativas activas.",
scroll = TRUE)| Variable | AD (p) | Lilliefors (p) | D de Lilliefors | Asimetría | Decisión (α=0.05) |
|---|---|---|---|---|---|
| preciom | <0.0000000000000001 | <0.0000000000000001 | 0.168 | 1.850 | Se rechaza H₀ |
| areaconst | <0.0000000000000001 | <0.0000000000000001 | 0.183 | 2.694 | Se rechaza H₀ |
| parqueaderos | <0.0000000000000001 | <0.0000000000000001 | 0.271 | 2.410 | Se rechaza H₀ |
| banios | <0.0000000000000001 | <0.0000000000000001 | 0.203 | 0.980 | Se rechaza H₀ |
| habitaciones | <0.0000000000000001 | <0.0000000000000001 | 0.286 | 1.812 | Se rechaza H₀ |
| precio_m2 | <0.0000000000000001 | <0.0000000000000001 | 0.038 | 0.679 | Se rechaza H₀ |
Con más de ocho mil observaciones cualquier prueba de normalidad rechaza \(H_0\) ante desviaciones mínimas, así que el valor \(p\) no decide. El estadístico \(D\) de Lilliefors mide la máxima discrepancia con la normal ajustada y funciona como tamaño del efecto. La decisión se apoya en \(D\), la asimetría y la Figura 3.3.
# ============================================================
# FIGURA: GRÁFICOS CUANTIL-CUANTIL
# ============================================================
viv_c |>
select(all_of(num_activas)) |>
pivot_longer(everything(), names_to = "Variable", values_to = "Valor") |>
ggplot(aes(sample = Valor)) +
stat_qq(colour = PAS_CIELO, alpha = 0.25, size = 0.7) +
stat_qq_line(colour = PAS_DURAZNO, linewidth = 0.7) +
facet_wrap(~ Variable, scales = "free", nrow = 1) +
labs(title = "Ninguna variable activa se ajusta a la normalidad",
x = "Cuantiles teóricos", y = "Cuantiles muestrales")
Figure 3.3: Gráficos cuantil-cuantil de las variables activas frente a la distribución normal. La desviación sistemática en la cola superior es la firma de la asimetría positiva.
# ============================================================
# TABLA: COMPARACIÓN DE TRANSFORMACIONES
# ============================================================
# Log-verosimilitud perfilada de Box-Cox para el modelo de solo intercepto.
# Se calcula de forma explícita porque MASS::boxcox() reevalúa el `lm` en el
# entorno de llamada y falla dentro de una función.
perfil_bc <- function(y, lambdas = seq(-1, 1.5, by = 0.005)) {
n <- length(y); sl <- sum(log(y))
ll <- vapply(lambdas, function(l) {
z <- if (abs(l) < 1e-9) log(y) else (y^l - 1) / l
-n/2 * log(mean((z - mean(z))^2)) + (l - 1) * sl
}, numeric(1))
list(lambda = lambdas[which.max(ll)], max_ll = max(ll))
}
t_kbl(bind_rows(lapply(c("preciom", "precio_m2", "areaconst"), function(v) {
y <- viv_c[[v]]
bc <- perfil_bc(y)
z <- if (abs(bc$lambda) < 1e-9) log(y) else (y^bc$lambda - 1) / bc$lambda
tibble::tibble(
Variable = v,
`Asimetría original` = round(e1071::skewness(y, type = 2), 3),
`Asimetría con log` = round(e1071::skewness(log(y), type = 2), 3),
`Lambda óptimo` = round(bc$lambda, 3),
`Asimetría con lambda óptimo` = round(e1071::skewness(z, type = 2), 3))
})),
caption = "Efecto del logaritmo y de la transformación de Box-Cox óptima sobre la asimetría.",
digits = 3)| Variable | Asimetría original | Asimetría con log | Lambda óptimo | Asimetría con lambda óptimo |
|---|---|---|---|---|
| preciom | 1.850 | 0.245 | -0.150 | 0.015 |
| precio_m2 | 0.679 | -0.566 | 0.420 | -0.005 |
| areaconst | 2.694 | 0.514 | -0.385 | 0.088 |
Transformación elegida. Se aplica logaritmo a precio_m2 y a areaconst.
La potencia óptima reduce algo más la asimetría, pero carece de lectura
económica, mientras que en escala logarítmica las diferencias se interpretan
como variaciones porcentuales, algo necesario para leer las componentes
principales. Las variables de conteo se dejan en su escala: su asimetría es
menor y transformarlas alteraría su naturaleza discreta.
viv_c <- viv_c |>
mutate(log_precio_m2 = log(precio_m2),
log_area = log(areaconst),
log_preciom = log(preciom))3.6 Relaciones entre variables
Se usa Spearman por la asimetría, por la presencia de valores extremos que
desestabilizan a Pearson y por el carácter ordinal de estrato. Al operar sobre
rangos es invariante ante transformaciones monótonas.
# ============================================================
# FIGURA: MAPA DE CALOR DE CORRELACIONES
# ============================================================
vars_corr <- c(num_activas, "precio_m2", "estrato_num")
mat_cor <- cor(viv_c[vars_corr], method = "spearman", use = "pairwise.complete.obs")
mat_largo <- as.data.frame(as.table(mat_cor)) |>
rename(V1 = Var1, V2 = Var2, rho = Freq)
ggplot(mat_largo, aes(x = V1, y = V2, fill = rho)) +
geom_tile(colour = "white", linewidth = 0.6) +
geom_text(aes(label = sprintf("%.2f", rho)),
colour = TINTA,
size = 3.6) +
scale_fill_gradient2(low = PAS_MENTA, mid = "white", high = PAS_CIELO,
midpoint = 0, limits = c(-1, 1), name = expression(rho)) +
labs(title = "El precio se asocia sobre todo con el área construida y el estrato",
x = NULL, y = NULL) +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
panel.grid = element_blank())
Figure 3.4: Matriz de correlaciones de Spearman entre las variables cuantitativas y el estrato.
# ============================================================
# TABLA: ASOCIACIÓN DE CADA VARIABLE CON EL PRECIO
# ============================================================
cor_precio <- mat_largo |>
filter(V1 == "preciom", V2 != "preciom") |>
mutate(Magnitud = etiqueta_efecto(abs(rho), c(-Inf, .20, .40, .60, .80, Inf),
c("Muy débil","Débil","Moderada","Fuerte","Muy fuerte"))) |>
select(Variable = V2, `ρ de Spearman` = rho, Magnitud) |>
arrange(desc(abs(`ρ de Spearman`)))
t_kbl(cor_precio,
caption = "Asociación de cada variable con el precio de oferta, ordenada por magnitud.",
digits = 3, align = c("l","r","l"))| Variable | ρ de Spearman | Magnitud |
|---|---|---|
| areaconst | 0.822 | Muy fuerte |
| banios | 0.781 | Fuerte |
| parqueaderos | 0.768 | Fuerte |
| estrato_num | 0.710 | Fuerte |
| habitaciones | 0.440 | Moderada |
| precio_m2 | 0.331 | Débil |
Prueba de Kruskal-Wallis. \(H_0\): la distribución del precio es la misma en las cinco zonas. Tamaño del efecto \(\varepsilon^2 = H/(n-1)\), proporción de la variabilidad en rangos atribuible al factor (Tomczak & Tomczak, 2014).
# ============================================================
# TABLA: KRUSKAL-WALLIS DEL PRECIO SEGÚN ZONA
# ============================================================
kw_zona <- kruskal.test(preciom ~ zona, data = viv_c)
n_kw <- sum(!is.na(viv_c$preciom) & !is.na(viv_c$zona))
eps2 <- as.numeric(kw_zona$statistic) / (n_kw - 1)
t_kbl(tibble::tibble(
`Estadístico H` = round(as.numeric(kw_zona$statistic), 2),
`gl` = as.integer(kw_zona$parameter),
`Valor p` = format.pval(kw_zona$p.value, digits = 4, eps = 1e-16),
`ε²` = round(eps2, 4),
`Magnitud` = etiqueta_efecto(eps2, c(-Inf, .01, .06, .14, Inf),
c("Despreciable","Pequeño","Moderado","Grande"))),
caption = "Prueba de Kruskal-Wallis para el precio de oferta según zona de la ciudad.",
align = c("r","r","r","r","l"))| Estadístico H | gl | Valor p | ε² | Magnitud |
|---|---|---|---|---|
| 1,064 | 4 | < 0.0000000000000001 | 0.128 | Moderado |
# ============================================================
# TABLA: COMPARACIONES MÚLTIPLES DE DUNN (HOLM)
# ============================================================
dunn_zona <- FSA::dunnTest(preciom ~ zona, data = viv_c, method = "holm")
dunn_tabla <- dunn_zona$res |>
mutate(`Valor p ajustado` = format.pval(P.adj, digits = 3, eps = 1e-16),
Significativo = ifelse(P.adj < 0.05, "Sí", "No")) |>
select(`Comparación` = Comparison, Z, `Valor p ajustado`, Significativo) |>
arrange(desc(abs(Z)))
t_kbl(dunn_tabla,
caption = "Comparaciones múltiples de Dunn entre zonas, con corrección de Holm.",
digits = 3, align = c("l","r","r","c"), scroll = TRUE, height = "320px")| Comparación | Z | Valor p ajustado | Significativo |
|---|---|---|---|
| Zona Norte - Zona Oeste | -27.901 | < 0.0000000000000001 | Sí |
| Zona Oeste - Zona Oriente | 25.573 | < 0.0000000000000001 | Sí |
| Zona Oeste - Zona Sur | 24.045 | < 0.0000000000000001 | Sí |
| Zona Oriente - Zona Sur | -13.996 | < 0.0000000000000001 | Sí |
| Zona Centro - Zona Oeste | -11.458 | < 0.0000000000000001 | Sí |
| Zona Norte - Zona Sur | -9.217 | < 0.0000000000000001 | Sí |
| Zona Norte - Zona Oriente | 9.042 | < 0.0000000000000001 | Sí |
| Zona Centro - Zona Oriente | 4.511 | 0.0000194 | Sí |
| Zona Centro - Zona Sur | -3.331 | 0.00173 | Sí |
| Zona Centro - Zona Norte | -0.579 | 0.56283 | No |
La corrección de Holm controla la tasa de error por familia y domina a Bonferroni sin supuestos adicionales. El \(H_0\) de Kruskal-Wallis es la igualdad de distribuciones y no la de medianas, de modo que leer el resultado en términos de nivel de precios exige formas similares, condición que la Figura 3.5 permite verificar en escala logarítmica.
# ============================================================
# FIGURA: PRECIO POR ZONA
# ============================================================
ggplot(viv_c, aes(x = reorder(zona, log_preciom, FUN = median), y = log_preciom)) +
geom_boxplot(fill = PAS_CIELO, colour = "#7A7A7A",
outlier.colour = PAS_DURAZNO, outlier.alpha = 0.25, outlier.size = 0.8) +
stat_summary(fun = median, geom = "point", colour = "#7FA8D0", size = 2.4) +
coord_flip() +
labs(title = "La distribución del precio de oferta difiere entre zonas",
subtitle = "La comparación de niveles supone formas distribucionales semejantes",
x = NULL, y = "log(precio en millones COP)")
Figure 3.5: Distribución del precio de oferta en escala logarítmica según zona. Los puntos marcan la mediana de cada zona.
Prueba \(\chi^2\) de independencia. \(H_0\): tipo de vivienda y zona son independientes. Se acompaña de \(V\) de Cramér y de los residuos estandarizados ajustados, que bajo \(H_0\) siguen aproximadamente una \(N(0,1)\): los valores \(|r_{ij}|>2\) señalan las celdas responsables del rechazo.
# ============================================================
# TABLA: PRUEBA CHI-CUADRADO TIPO x ZONA
# ============================================================
tc <- table(viv_c$tipo, viv_c$zona)
prueba_tz <- chisq.test(tc)
t_kbl(tibble::tibble(
`Estadístico χ²` = round(as.numeric(prueba_tz$statistic), 2),
`gl` = as.integer(prueba_tz$parameter),
`Valor p` = format.pval(prueba_tz$p.value, digits = 4, eps = 1e-16),
`V de Cramér` = round(cramer_v(tc), 4),
`Frec. esperada mínima` = round(min(prueba_tz$expected), 1)),
caption = "Prueba χ² de independencia entre tipo de vivienda y zona de la ciudad.")| Estadístico χ² | gl | Valor p | V de Cramér | Frec. esperada mínima |
|---|---|---|---|---|
| 690.9 | 4 | < 0.0000000000000001 | 0.288 | 48 |
# ============================================================
# TABLA: RESIDUOS ESTANDARIZADOS AJUSTADOS
# ============================================================
t_kbl(as.data.frame.matrix(round(prueba_tz$stdres, 2)) |>
tibble::rownames_to_column("Tipo"),
caption = "Residuos estandarizados ajustados de la tabla tipo × zona.")| Tipo | Zona Centro | Zona Norte | Zona Oeste | Zona Oriente | Zona Sur |
|---|---|---|---|---|---|
| Apartamento | -9.66 | 1.12 | 18.89 | -17.15 | -5.01 |
| Casa | 9.66 | -1.12 | -18.89 | 17.15 | 5.01 |
# ============================================================
# TABLA: DECISIONES DE PREPARACIÓN
# ============================================================
t_kbl(tibble::tribble(
~Aspecto, ~Hallazgo, ~Tratamiento,
"Consistencia", "Duplicados exactos y valores no positivos en variables de dominio positivo",
"Duplicados eliminados; valores fuera de dominio convertidos en NA",
"`piso`", "La ausencia no depende del tipo de vivienda",
"El valor numérico se descarta; la declaración del campo se conserva como indicador",
"`parqueaderos`", "Se descarta MCAR; el cero nunca se observa",
"Imputación por ACP regularizado, con la recodificación a cero como escenario alternativo",
"`banios`, `habitaciones`", "Ausencia inferior al 1 %",
"Imputación conjunta con ncp por validación cruzada",
"Valores extremos", "El grupo señalado por MCD corresponde al segmento de alto valor",
"Se conservan; control con ACP robusto y con medoides",
"Distribución", "Normalidad rechazada en todas las variables",
"Logaritmo en precio unitario y área; métodos de rangos en el bivariado"),
caption = "Resumen de los hallazgos del diagnóstico y del tratamiento aplicado.",
align = c("l","l","l"), scroll = TRUE)| Aspecto | Hallazgo | Tratamiento |
|---|---|---|
| Consistencia | Duplicados exactos y valores no positivos en variables de dominio positivo | Duplicados eliminados; valores fuera de dominio convertidos en NA |
piso
|
La ausencia no depende del tipo de vivienda | El valor numérico se descarta; la declaración del campo se conserva como indicador |
parqueaderos
|
Se descarta MCAR; el cero nunca se observa | Imputación por ACP regularizado, con la recodificación a cero como escenario alternativo |
banios, habitaciones
|
Ausencia inferior al 1 % | Imputación conjunta con ncp por validación cruzada |
| Valores extremos | El grupo señalado por MCD corresponde al segmento de alto valor | Se conservan; control con ACP robusto y con medoides |
| Distribución | Normalidad rechazada en todas las variables | Logaritmo en precio unitario y área; métodos de rangos en el bivariado |
Las decisiones comparten un criterio: la ausencia de un dato es información, y
su tratamiento va después de identificar el mecanismo. La prueba de Little
rechaza MCAR para el conjunto y el examen por variable localiza el rechazo en
parqueaderos; en las demás la asociación es detectable pero de magnitud
despreciable, distinción entre significación y relevancia que el tamaño muestral
vuelve obligatoria.
4 Estructura latente: componentes principales
4.1 Planteamiento
Con \(\mathbf{X}\) la matriz \(n\times p\) centrada y reducida, el ACP busca direcciones ortonormales que maximizan sucesivamente la varianza proyectada. La solución son los vectores propios de \(\mathbf{R}=\frac{1}{n-1}\mathbf{X}^{\top}\mathbf{X}\), con valores propios \(\lambda_1\ge\dots\ge\lambda_p\ge 0\) iguales a la varianza de cada componente (Johnson & Wichern, 2007; Peña, 2002). Que \(\mathbf{R}\) sea simétrica y semidefinida positiva asegura base ortonormal y valores propios reales no negativos (Lay, 2012, cap. 7). Se trabaja sobre correlaciones porque las unidades son incomparables.
Elección del conjunto activo. La variable de valor no es el precio total
sino el precio por metro cuadrado. Incluir preciom y areaconst juntos
introduce la superficie dos veces, una directamente y otra dentro del precio, y
fuerza un primer eje de tamaño en el que todas las variables cargan con el mismo
signo. Con log(precio_m2) el tamaño y el valor unitario quedan separados, que
es la distinción que interesa a una empresa que valora inmuebles.
# ============================================================
# TABLA: ROLES DE LAS VARIABLES
# ============================================================
activas_acp <- c("log_precio_m2", "log_area", "parqueaderos", "banios", "habitaciones")
cols_acp <- c(activas_acp, "estrato", "zona", "tipo")
# Índice de filas completas. Define el conjunto sobre el que operan el ACP y la
# segmentación, y se aplica también a `viv_c` para que las etiquetas de grupo
# sean asignables fila a fila.
viv_c$.fila <- seq_len(nrow(viv_c))
idx_acp <- complete.cases(viv_c[, cols_acp])
n_excl <- sum(!idx_acp)
viv_c <- viv_c[idx_acp, ]
datos_acp <- viv_c[, cols_acp]
t_kbl(tibble::tribble(
~Variable, ~Rol, ~Motivo,
"log_precio_m2", "Activa", "Valor unitario del inmueble, en escala logarítmica",
"log_area", "Activa", "Tamaño construido, en escala logarítmica",
"parqueaderos", "Activa", "Dotación, con los faltantes imputados",
"banios", "Activa", "Dotación",
"habitaciones", "Activa", "Distribución del espacio",
"estrato", "Suplementaria", "Ordinal: como activa supondría equidistancia",
"zona", "Suplementaria", "Nominal",
"tipo", "Suplementaria", "Nominal",
"preciom", "Excluida", "Contenida en el precio unitario y en el área; su inclusión duplicaría el tamaño",
"longitud, latitud", "Excluidas", "Coordenadas, no atributos del inmueble"),
caption = "Papel de cada variable en el análisis de componentes principales.",
align = c("l","l","l"), scroll = TRUE)| Variable | Rol | Motivo |
|---|---|---|
| log_precio_m2 | Activa | Valor unitario del inmueble, en escala logarítmica |
| log_area | Activa | Tamaño construido, en escala logarítmica |
| parqueaderos | Activa | Dotación, con los faltantes imputados |
| banios | Activa | Dotación |
| habitaciones | Activa | Distribución del espacio |
| estrato | Suplementaria | Ordinal: como activa supondría equidistancia |
| zona | Suplementaria | Nominal |
| tipo | Suplementaria | Nominal |
| preciom | Excluida | Contenida en el precio unitario y en el área; su inclusión duplicaría el tamaño |
| longitud, latitud | Excluidas | Coordenadas, no atributos del inmueble |
El filtrado por casos completos excluye 0 filas. El análisis opera sobre 8,319 inmuebles. Las suplementarias se proyectan una vez construidos los ejes, de modo que no influyen en su orientación.
4.2 Condiciones previas
Esfericidad de Bartlett (1951). \(H_0\): \(\mathbf{R}=\mathbf{I}_p\), sin estructura común que resumir.
Índice KMO (Kaiser, 1974). Compara correlaciones simples y parciales; se considera aceptable por encima de 0.70.
# ============================================================
# TABLA: BARTLETT Y KMO GLOBAL
# ============================================================
X_acp <- as.matrix(datos_acp[activas_acp])
R_acp <- cor(X_acp)
bart <- psych::cortest.bartlett(R_acp, n = nrow(X_acp))
kmo <- psych::KMO(R_acp)
t_kbl(tibble::tibble(
Prueba = c("Esfericidad de Bartlett", "KMO global"),
`Estadístico` = c(round(bart$chisq, 1), round(kmo$MSA, 4)),
`gl` = c(bart$df, NA_integer_),
`Valor p` = c(format.pval(bart$p.value, digits = 4, eps = 1e-16), NA_character_),
Lectura = c(ifelse(bart$p.value < 0.05,
"Se rechaza la incorrelación: existe estructura factorial",
"No se rechaza la incorrelación: el ACP no está justificado"),
etiqueta_efecto(kmo$MSA, c(-Inf, .5, .6, .7, .8, .9, Inf),
c("Inaceptable","Pobre","Mediocre","Aceptable","Buena","Excelente")))),
caption = "Contrastes de adecuación de los datos al análisis factorial.",
scroll = TRUE)| Prueba | Estadístico | gl | Valor p | Lectura |
|---|---|---|---|---|
| Esfericidad de Bartlett | 21,135.300 | 10 | < 0.0000000000000001 | Se rechaza la incorrelación: existe estructura factorial |
| KMO global | 0.651 | — | — | Mediocre |
# ============================================================
# TABLA: KMO POR VARIABLE
# ============================================================
t_kbl(tibble::tibble(Variable = names(kmo$MSAi), `KMO individual` = round(kmo$MSAi, 4)),
caption = "Medida de adecuación muestral individual por variable.",
align = c("l","r"))| Variable | KMO individual |
|---|---|
| log_precio_m2 | 0.319 |
| log_area | 0.673 |
| parqueaderos | 0.737 |
| banios | 0.668 |
| habitaciones | 0.738 |
Con más de ocho mil filas Bartlett rechaza \(H_0\) ante correlaciones mínimas, así que funciona como condición necesaria y no como evidencia de utilidad. El índice con capacidad discriminante es el KMO, que no depende del tamaño muestral. Los valores individuales de la Tabla 4.3 no son homogéneos, y las variables por debajo de 0.70 quedan peor integradas en el factor común.
4.3 Número de componentes
# ============================================================
# TABLA: VALORES PROPIOS Y CRITERIOS DE RETENCIÓN
# ============================================================
acp <- FactoMineR::PCA(datos_acp, quali.sup = 6:8, scale.unit = TRUE,
ncp = 5, graph = FALSE)
vp <- as.data.frame(acp$eig)
names(vp) <- c("Valor propio", "% de varianza", "% acumulado")
p_var <- length(activas_acp)
# Criterio del bastón roto (Jackson, 1993): bajo reparto aleatorio de la
# varianza total entre p componentes, la proporción esperada de la componente k
# es b_k = (1/p) * sum_{i=k}^{p} 1/i. Se retienen las componentes cuya
# proporción observada supera la esperada bajo ese modelo nulo.
baston <- vapply(1:p_var, function(k) sum(1 / (k:p_var)) / p_var, numeric(1))
retencion <- tibble::tibble(
Componente = paste0("CP", 1:p_var),
`Valor propio` = vp$`Valor propio`,
`% varianza` = vp$`% de varianza`,
`% acumulado` = vp$`% acumulado`,
`Kaiser (λ > 1)` = ifelse(vp$`Valor propio` > 1, "Retener", "Descartar"),
`Bastón roto (%)` = 100 * baston,
`Bastón roto` = ifelse(vp$`% de varianza` > 100 * baston, "Retener", "Descartar")
)
k_kaiser <- sum(vp$`Valor propio` > 1)
k_baston <- sum(vp$`% de varianza` > 100 * baston)
k_80 <- which(vp$`% acumulado` >= 80)[1]
n_ret <- 2 # componentes retenidas: gobierna interpretación y segmentación
t_kbl(retencion,
caption = "Valores propios y criterios de retención de componentes.",
digits = 3, scroll = TRUE)| Componente | Valor propio | % varianza | % acumulado | Kaiser (λ > 1) | Bastón roto (%) | Bastón roto |
|---|---|---|---|---|---|---|
| CP1 | 2.834 | 56.687 | 56.69 | Retener | 45.67 | Retener |
| CP2 | 1.228 | 24.557 | 81.24 | Retener | 25.67 | Descartar |
| CP3 | 0.489 | 9.777 | 91.02 | Descartar | 15.67 | Descartar |
| CP4 | 0.289 | 5.772 | 96.79 | Descartar | 9.00 | Descartar |
| CP5 | 0.160 | 3.207 | 100.00 | Descartar | 4.00 | Descartar |
# ============================================================
# FIGURA: SEDIMENTACIÓN
# ============================================================
sedim <- tibble::tibble(Componente = 1:p_var,
Observado = vp$`% de varianza`,
`Bastón roto` = 100 * baston)
ggplot(sedim, aes(x = Componente)) +
geom_col(aes(y = Observado), fill = PAS_MENTA, width = 0.6) +
geom_line(aes(y = Observado), colour = PAS_CIELO, linewidth = 0.9) +
geom_point(aes(y = Observado), colour = PAS_CIELO, size = 2.4) +
geom_line(aes(y = `Bastón roto`), colour = PAS_DURAZNO, linewidth = 0.9) +
geom_point(aes(y = `Bastón roto`), colour = PAS_DURAZNO, size = 2.2) +
geom_hline(yintercept = 100 / p_var, linetype = "dashed", colour = "#7A7A7A") +
scale_x_continuous(breaks = 1:p_var) +
labs(title = "Retención de componentes: observado frente a modelo nulo",
subtitle = "Barras y línea azul: varianza explicada. Línea naranja: bastón roto. Discontinua: umbral de Kaiser",
x = "Componente principal", y = "% de varianza explicada")
Figure 4.1: Gráfico de sedimentación: varianza observada frente al perfil esperado bajo el modelo del bastón roto. La línea discontinua marca el umbral de Kaiser expresado en porcentaje.
Los criterios difieren: el bastón roto retiene 1, Kaiser 2 y el umbral del 80 % acumulado exige 2. Kaiser es el más permisivo y el bastón roto el más exigente, por ser el único con modelo nulo explícito.
Componentes retenidas. Se conservan dos, con estatus distinto. La primera está respaldada por los criterios formales y sostiene las conclusiones sustantivas. La segunda se retiene por su legibilidad y porque el plano factorial es el espacio en que se construye la segmentación; los hallazgos que dependan solo de ella se enuncian como hipótesis.
4.4 Contenido de los ejes
Contribución y calidad de representación: \[\mathrm{Contrib}_{jk}=\frac{v_{jk}^{2}}{\sum_{l}v_{lk}^{2}}\times 100, \qquad \cos^{2}_{jk}=\frac{c_{jk}^{2}}{\sum_{m=1}^{p}c_{jm}^{2}}.\]
# ============================================================
# TABLA: CARGAS, CONTRIBUCIONES Y COS2
# ============================================================
cargas <- as.data.frame(acp$var$coord[, 1:n_ret, drop = FALSE])
contrib <- as.data.frame(acp$var$contrib[, 1:n_ret, drop = FALSE])
coseno <- as.data.frame(acp$var$cos2[, 1:n_ret, drop = FALSE])
names(cargas) <- paste0("Coord. CP", 1:n_ret)
names(contrib) <- paste0("Contrib. CP", 1:n_ret, " (%)")
names(coseno) <- paste0("cos² CP", 1:n_ret)
tabla_var <- cbind(Variable = rownames(cargas), cargas, contrib, coseno,
`cos² acumulado` = rowSums(acp$var$cos2[, 1:n_ret, drop = FALSE]))
t_kbl(tabla_var,
caption = "Coordenadas, contribuciones y calidad de representación de las variables activas.",
digits = 3, scroll = TRUE)| Variable | Coord. CP1 | Coord. CP2 | Contrib. CP1 (%) | Contrib. CP2 (%) | cos² CP1 | cos² CP2 | cos² acumulado | |
|---|---|---|---|---|---|---|---|---|
| log_precio_m2 | log_precio_m2 | -0.305 | 0.888 | 3.272 | 64.177 | 0.093 | 0.788 | 0.881 |
| log_area | log_area | 0.921 | -0.011 | 29.899 | 0.010 | 0.847 | 0.000 | 0.848 |
| parqueaderos | parqueaderos | 0.721 | 0.495 | 18.335 | 19.945 | 0.520 | 0.245 | 0.765 |
| banios | banios | 0.875 | 0.242 | 27.041 | 4.765 | 0.766 | 0.059 | 0.825 |
| habitaciones | habitaciones | 0.780 | -0.369 | 21.454 | 11.103 | 0.608 | 0.136 | 0.744 |
# ============================================================
# FIGURA: CÍRCULO DE CORRELACIONES
# ============================================================
circulo <- data.frame(x = cos(seq(0, 2*pi, length.out = 300)),
y = sin(seq(0, 2*pi, length.out = 300)))
vec_var <- data.frame(
Variable = rownames(acp$var$coord),
CP1 = acp$var$coord[, 1],
CP2 = acp$var$coord[, 2],
cos2 = rowSums(acp$var$cos2[, 1:2])
)
ggplot() +
geom_path(data = circulo, aes(x, y), colour = "#C8C8C8", linewidth = 0.6) +
geom_hline(yintercept = 0, colour = "#C8C8C8", linewidth = 0.4) +
geom_vline(xintercept = 0, colour = "#C8C8C8", linewidth = 0.4) +
geom_segment(data = vec_var, aes(x = 0, y = 0, xend = CP1, yend = CP2, colour = cos2),
arrow = arrow(length = unit(0.18, "cm")), linewidth = 0.9) +
ggrepel::geom_text_repel(data = vec_var, aes(x = CP1, y = CP2, label = Variable),
size = 3.6, colour = TINTA) +
scale_colour_gradient(low = "#E9EEF5", high = "#7FA8D0", limits = c(0, 1), name = expression(cos^2)) +
coord_equal(xlim = c(-1.05, 1.05), ylim = c(-1.05, 1.05)) +
labs(title = "Todas las variables cargan positivamente sobre el primer eje",
x = paste0("CP1 (", round(vp$`% de varianza`[1], 1), "%)"),
y = paste0("CP2 (", round(vp$`% de varianza`[2], 1), "%)"))
Figure 4.2: Círculo de correlaciones en el primer plano factorial. La longitud del vector indica la calidad de representación; el coseno del ángulo entre dos vectores aproxima su correlación.
# ============================================================
# TABLA: COORDENADAS DE LAS CATEGORÍAS SUPLEMENTARIAS
# ============================================================
sup_coord <- as.data.frame(acp$quali.sup$coord[, 1:n_ret, drop = FALSE])
names(sup_coord) <- paste0("Coord. CP", 1:n_ret)
sup_coord$`cos²` <- rowSums(acp$quali.sup$cos2[, 1:n_ret, drop = FALSE])
sup_coord <- cbind(`Categoría` = rownames(sup_coord), sup_coord)
t_kbl(sup_coord,
caption = "Coordenadas y calidad de representación de las categorías suplementarias en el primer plano factorial.",
digits = 3, scroll = TRUE, height = "360px")| Categoría | Coord. CP1 | Coord. CP2 | cos² | |
|---|---|---|---|---|
| estrato_3 | estrato_3 | -0.323 | -1.226 | 0.948 |
| estrato_4 | estrato_4 | -0.608 | -0.320 | 0.978 |
| estrato_5 | estrato_5 | -0.008 | 0.092 | 0.585 |
| estrato_6 | estrato_6 | 0.898 | 1.111 | 0.982 |
| Zona Centro | Zona Centro | 0.527 | -1.240 | 0.969 |
| Zona Norte | Zona Norte | -0.252 | -0.302 | 0.938 |
| Zona Oeste | Zona Oeste | 0.220 | 0.873 | 0.931 |
| Zona Oriente | Zona Oriente | 0.514 | -1.762 | 0.979 |
| Zona Sur | Zona Sur | -0.006 | 0.065 | 0.642 |
| Apartamento | Apartamento | -0.841 | 0.277 | 0.997 |
| Casa | Casa | 1.332 | -0.440 | 0.997 |
# ============================================================
# FIGURA: BARICENTROS DE LAS SUPLEMENTARIAS
# ============================================================
bar_sup <- as.data.frame(acp$quali.sup$coord[, 1:2])
names(bar_sup) <- c("CP1", "CP2")
bar_sup$Categoria <- rownames(bar_sup)
bar_sup$Variable <- ifelse(grepl("^estrato", bar_sup$Categoria), "Estrato",
ifelse(grepl("^Zona", bar_sup$Categoria), "Zona", "Tipo"))
ggplot(bar_sup, aes(x = CP1, y = CP2, colour = Variable)) +
geom_hline(yintercept = 0, colour = "#C8C8C8", linewidth = 0.4) +
geom_vline(xintercept = 0, colour = "#C8C8C8", linewidth = 0.4) +
geom_point(size = 3) +
ggrepel::geom_text_repel(aes(label = Categoria), size = 3.6, show.legend = FALSE) +
scale_colour_manual(values = c(Estrato = PAS_CIELO, Zona = PAS_LILA, Tipo = PAS_DURAZNO)) +
labs(title = "Las categorías suplementarias se ordenan por gama sobre el primer eje",
x = paste0("CP1 (", round(vp$`% de varianza`[1], 1), "%)"),
y = paste0("CP2 (", round(vp$`% de varianza`[2], 1), "%)"))
Figure 4.3: Baricentros de las categorías suplementarias en el primer plano factorial. Cada punto ocupa la posición media de los individuos de esa categoría.
Las variables activas y las categorías suplementarias ocupan espacios distintos. Las primeras están en el espacio de las variables, con coordenada igual a una correlación acotada en \([-1,1]\). Las segundas están en el espacio de los individuos, con coordenada igual al baricentro de sus puntuaciones factoriales. Por eso se representan por separado.
# ============================================================
# FIGURA: NUBE DE INDIVIDUOS POR ESTRATO
# ============================================================
coord_ind <- as.data.frame(acp$ind$coord[, 1:2])
names(coord_ind) <- c("CP1", "CP2")
coord_ind$estrato <- datos_acp$estrato
ggplot(coord_ind, aes(x = CP1, y = CP2, colour = estrato)) +
geom_point(alpha = 0.30, size = 0.8) +
stat_ellipse(type = "norm", level = 0.68, linewidth = 0.9) +
scale_colour_manual(values = PAL_CAT, name = "Estrato") +
geom_hline(yintercept = 0, colour = "#C8C8C8", linewidth = 0.4) +
geom_vline(xintercept = 0, colour = "#C8C8C8", linewidth = 0.4) +
labs(title = "El estrato se ordena de forma monótona a lo largo del primer eje",
subtitle = "Los centros se desplazan sin inversiones, pero las regiones se solapan ampliamente",
x = paste0("CP1 (", round(vp$`% de varianza`[1], 1), "%)"),
y = paste0("CP2 (", round(vp$`% de varianza`[2], 1), "%)")) +
guides(colour = guide_legend(override.aes = list(alpha = 1, size = 3)))
Figure 4.4: Nube de individuos en el primer plano factorial, coloreada por estrato socioeconómico (variable suplementaria). Las elipses son regiones de concentración del 68 % bajo un modelo normal bivariado.
La Tabla 4.5 y la Figura 4.2 determinan la interpretación. Al haber sustituido el precio total por el precio unitario, el primer eje ya no puede ser un factor de tamaño puro: el valor por metro cuadrado y las variables de dimensión (área, habitaciones) responden a lógicas distintas, de modo que el eje separa calidad frente a magnitud en lugar de ordenar simplemente por escala. Los baricentros del estrato, que es suplementario, permiten verificar cuál de los dos ejes recoge el gradiente socioeconómico, y esa verificación es externa porque el estrato no intervino en la construcción del espacio.
4.5 Verificación de la solución
Se comprueban las dos fuentes de fragilidad detectadas: los inmuebles de alto
apalancamiento y el tratamiento de la ausencia en parqueaderos. El instrumento
es el coeficiente de congruencia de Tucker,
\(\varphi=\sum_i a_ib_i/\sqrt{\sum_i a_i^2\sum_i b_i^2}\), que compara vectores de
cargas con independencia de escala y signo; por encima de 0.95 indica
equivalencia (Lorenzo-Seva & ten Berge, 2006).
# ============================================================
# TABLA: SOLUCIÓN CLÁSICA FRENTE A SOLUCIÓN ROBUSTA
# ============================================================
acp_clasico <- prcomp(X_acp, scale. = TRUE)
acp_robusto <- rrcov::PcaHubert(X_acp, k = n_ret, scale = TRUE, alpha = 0.75)
congruencia <- function(a, b) abs(sum(a * b) / sqrt(sum(a^2) * sum(b^2)))
t_kbl(tibble::tibble(
Componente = paste0("CP", 1:n_ret),
`Congruencia de Tucker` = vapply(1:n_ret, function(k)
congruencia(acp_clasico$rotation[, k], rrcov::getLoadings(acp_robusto)[, k]), numeric(1))) |>
mutate(Lectura = ifelse(`Congruencia de Tucker` >= 0.95, "Equivalencia factorial",
ifelse(`Congruencia de Tucker` >= 0.85, "Similitud aceptable",
"Discrepancia relevante"))),
caption = "Congruencia entre la solución clásica y la solución robusta (ROBPCA de Hubert).",
digits = 4)| Componente | Congruencia de Tucker | Lectura |
|---|---|---|
| CP1 | 0.9313 | Similitud aceptable |
| CP2 | 0.9109 | Similitud aceptable |
# ============================================================
# TABLA: SENSIBILIDAD AL TRATAMIENTO DE `parqueaderos`
# ============================================================
# Escenario alternativo: ausencia recodificada a cero en lugar de imputada.
# `.fila` conserva la correspondencia fila a fila con `viv`.
X_alt <- viv_c[, activas_acp]
X_alt$parqueaderos <- viv$parqueaderos_cero[viv_c$.fila]
acp_alt <- prcomp(X_alt, scale. = TRUE)
t_kbl(tibble::tibble(
Componente = paste0("CP", 1:n_ret),
`% varianza (imputado)` = 100 * (acp_clasico$sdev^2 / sum(acp_clasico$sdev^2))[1:n_ret],
`% varianza (cero)` = 100 * (acp_alt$sdev^2 / sum(acp_alt$sdev^2))[1:n_ret],
`Congruencia de Tucker` = vapply(1:n_ret, function(k)
congruencia(acp_clasico$rotation[, k], acp_alt$rotation[, k]), numeric(1))) |>
mutate(Lectura = ifelse(`Congruencia de Tucker` >= 0.95,
"La solución no depende del tratamiento de la ausencia",
"La solución depende del tratamiento de la ausencia")),
caption = "Comparación de la solución factorial bajo imputación y bajo recodificación a cero.",
digits = 4, scroll = TRUE)| Componente | % varianza (imputado) | % varianza (cero) | Congruencia de Tucker | Lectura |
|---|---|---|---|---|
| CP1 | 56.69 | 53.82 | 0.9980 | La solución no depende del tratamiento de la ausencia |
| CP2 | 24.56 | 26.61 | 0.9937 | La solución no depende del tratamiento de la ausencia |
La congruencia con la solución robusta indica si la estructura está determinada por los registros extremos. La segunda tabla resuelve la decisión pendiente de la Sección 3.2.3: si ambas soluciones son congruentes, la elección entre imputar y recodificar deja de afectar a las conclusiones, y si no lo son, prevalece la solución imputada por ser la que no incorpora una conjetura sobre el formulario del portal.
# ============================================================
# TABLA: SÍNTESIS DEL ACP
# ============================================================
t_kbl(tibble::tribble(
~Elemento, ~Resultado,
"Condiciones previas", "Esfericidad rechazada y KMO global aceptable, con adecuación individual desigual",
"Conjunto activo", "Precio unitario en lugar de precio total, para no contar dos veces el tamaño",
"Dimensionalidad", "Una componente con respaldo formal; dos para la representación y la segmentación",
"Validación externa", "El estrato, suplementario, permite identificar el eje socioeconómico",
"Robustez", "Congruencia de Tucker frente a la estimación robusta y frente al escenario alternativo"),
caption = "Síntesis del análisis de componentes principales.",
align = c("l","l"), scroll = TRUE)| Elemento | Resultado |
|---|---|
| Condiciones previas | Esfericidad rechazada y KMO global aceptable, con adecuación individual desigual |
| Conjunto activo | Precio unitario en lugar de precio total, para no contar dos veces el tamaño |
| Dimensionalidad | Una componente con respaldo formal; dos para la representación y la segmentación |
| Validación externa | El estrato, suplementario, permite identificar el eje socioeconómico |
| Robustez | Congruencia de Tucker frente a la estimación robusta y frente al escenario alternativo |
5 Segmentación de la oferta
5.1 Método
Se busca una partición en \(k\) grupos que minimice la heterogeneidad interna. El criterio adoptado es el de \(k\)-medoides: en lugar del centroide, cada grupo se representa por una observación real, el medoide, que minimiza la suma de distancias al resto del grupo. La razón es directa: la nube es fuertemente asimétrica y conserva un segmento de alto valor que se decidió no eliminar (Sección 3); el centroide se desplaza hacia esa cola, el medoide no. Con \(n\) superior a ocho mil se usa CLARA (Kaufman & Rousseeuw, 2005), que aplica PAM sobre submuestras y devuelve la partición asociada al mejor conjunto de medoides.
La agregación jerárquica de Ward con consolidación queda como método de contraste, no como método principal. Todas las etapas operan sobre el mismo espacio: las 2 componentes retenidas.
5.2 Predisposición al agrupamiento
Un algoritmo de partición devuelve grupos aunque los datos no tengan estructura, de modo que este paso es el que da sentido al resto.
Estadístico de Hopkins. Con \(w_i\) la distancia de cada punto real a su vecino real más próximo y \(u_i\) la de un punto uniforme a su vecino real más próximo, \[H=\frac{\sum_{i=1}^{m}u_i}{\sum_{i=1}^{m}u_i+\sum_{i=1}^{m}w_i}.\] Sin agrupamiento \(H\approx 0.5\). Algunas implementaciones devuelven \(1-H\), así que se calcula de forma explícita.
# ============================================================
# TABLA: ESTADÍSTICO DE HOPKINS
# ============================================================
# Espacio de trabajo de toda la sección: las n_ret componentes retenidas.
puntuaciones <- as.data.frame(acp$ind$coord[, seq_len(n_ret), drop = FALSE])
names(puntuaciones) <- paste0("CP", seq_len(n_ret))
cp_vars <- paste0("CP", seq_len(n_ret))
hopkins <- function(X, m = 400) {
X <- as.matrix(X); n_x <- nrow(X); p_x <- ncol(X)
idx <- sample(n_x, m)
U <- sapply(1:p_x, function(j) runif(m, min(X[, j]), max(X[, j])))
w <- sapply(idx, function(i) { d <- sqrt(colSums((t(X) - X[i, ])^2)); min(d[-i]) })
u <- sapply(1:m, function(i) { d <- sqrt(colSums((t(X) - U[i, ])^2)); min(d) })
sum(u) / (sum(u) + sum(w))
}
# Se repite el cálculo para disponer de una distribución y no de un valor único
H_rep <- replicate(30, hopkins(puntuaciones, m = 400))
t_kbl(tibble::tibble(
Repeticiones = length(H_rep),
`H medio` = mean(H_rep),
`H mínimo` = min(H_rep),
`H máximo` = max(H_rep),
`Valor bajo H₀` = 0.5,
Lectura = ifelse(mean(H_rep) > 0.75, "Tendencia de agrupamiento marcada",
ifelse(mean(H_rep) > 0.60, "Tendencia de agrupamiento moderada",
"Sin evidencia de agrupamiento"))),
caption = "Estadístico de Hopkins calculado sobre 30 submuestras independientes.",
digits = 4, scroll = TRUE)| Repeticiones | H medio | H mínimo | H máximo | Valor bajo H₀ | Lectura |
|---|---|---|---|---|---|
| 30 | 0.9519 | 0.9427 | 0.9591 | 0.5 | Tendencia de agrupamiento marcada |
La hipótesis nula de Hopkins es la uniformidad sobre el hiperrectángulo envolvente, que es poco exigente: una distribución con núcleo denso y colas largas produce valores cercanos a 1 sin contener grupos separados. Un \(H\) alto descarta uniformidad, no establece conglomerados.
# ============================================================
# FIGURA: MATRIZ DE DISIMILARIDAD ORDENADA
# ============================================================
sub_odi <- puntuaciones[sample(nrow(puntuaciones), 400), cp_vars, drop = FALSE]
d_odi <- dist(sub_odi)
orden <- hclust(d_odi, method = "ward.D2")$order
M_odi <- as.matrix(d_odi)[orden, orden]
odi_largo <- expand.grid(i = seq_len(nrow(M_odi)), j = seq_len(ncol(M_odi)))
odi_largo$d <- as.vector(M_odi)
ggplot(odi_largo, aes(x = i, y = j, fill = d)) +
geom_raster() +
scale_fill_gradient(low = "#7FA8D0", high = "#F7F7F7", name = "Distancia") +
coord_equal() +
labs(title = "La disimilaridad ordenada no muestra bloques nítidos",
subtitle = "Gradación continua: concentración, no separación en subpoblaciones",
x = NULL, y = NULL) +
theme(axis.text = element_blank(), panel.grid = element_blank())
Figure 5.1: Matriz de disimilaridad ordenada por clasificación jerárquica sobre una submuestra aleatoria. La existencia de grupos separados se manifestaría como bloques oscuros nítidos a lo largo de la diagonal.
5.3 Elección de \(k\)
Se emplean cuatro criterios de fundamento distinto, todos evaluados sobre la misma submuestra y, para cada \(k\), sobre la misma partición de PAM.
# ============================================================
# TABLA: ÍNDICES DE CALIDAD PARA DISTINTOS k
# ============================================================
set.seed(2026)
k_max <- 8
# Los índices basados en distancias exigen n(n-1)/2 entradas. Se evalúan sobre
# una submuestra común a todos ellos.
sub_idx <- sample(nrow(puntuaciones), 1500)
sub_pun <- puntuaciones[sub_idx, cp_vars, drop = FALSE]
d_sub <- dist(sub_pun)
# Davies-Bouldin: DB = (1/k) sum_i max_{j!=i} (S_i + S_j)/M_ij, con S_i la
# dispersión media intragrupo y M_ij la distancia entre representantes.
db_index <- function(X, cl) {
X <- as.matrix(X); ks <- sort(unique(cl))
cent <- t(sapply(ks, function(g) colMeans(X[cl == g, , drop = FALSE])))
S <- sapply(ks, function(g) mean(sqrt(rowSums((X[cl == g, , drop = FALSE] -
matrix(cent[which(ks == g), ], sum(cl == g), ncol(X), byrow = TRUE))^2))))
M <- as.matrix(dist(cent))
mean(sapply(seq_along(ks), function(i) max((S[i] + S[-i]) / M[i, -i])))
}
indices_k <- bind_rows(lapply(2:k_max, function(k) {
# pamonce = 5 activa la versión optimizada del algoritmo (Reynolds et al.),
# que evita recorrer todos los pares en cada intercambio de medoide.
pm <- cluster::pam(sub_pun, k = k, metric = "euclidean", pamonce = 5)
sil <- cluster::silhouette(pm$clustering, d_sub)
# Se desactivan los índices que no se reportan: cada uno recorre la matriz
# de distancias completa y multiplica el tiempo de cómputo.
st <- fpc::cluster.stats(d_sub, pm$clustering, silhouette = FALSE,
wgap = FALSE, sepindex = FALSE, G2 = FALSE, G3 = FALSE)
tibble::tibble(
k = k,
`Distancia media al medoide` = pm$objective[2],
`Silueta media` = mean(sil[, 3]),
`Calinski-Harabasz` = st$ch,
`Davies-Bouldin` = db_index(sub_pun, pm$clustering)
)
}))
t_kbl(indices_k,
caption = "Índices de calidad de la partición por k-medoides para distintos valores de k.",
digits = 3, scroll = TRUE)| k | Distancia media al medoide | Silueta media | Calinski-Harabasz | Davies-Bouldin |
|---|---|---|---|---|
| 2 | 1.255 | 0.457 | 1,405 | 0.916 |
| 3 | 1.027 | 0.426 | 1,347 | 0.931 |
| 4 | 0.912 | 0.428 | 1,349 | 0.844 |
| 5 | 0.834 | 0.353 | 1,259 | 0.915 |
| 6 | 0.768 | 0.354 | 1,267 | 0.883 |
| 7 | 0.715 | 0.345 | 1,212 | 0.903 |
| 8 | 0.666 | 0.349 | 1,248 | 0.888 |
# ============================================================
# TABLA: ESTADÍSTICO GAP CON PAM
# ============================================================
# clusGap exige una función que devuelva una lista con el elemento `cluster`.
# Se usa CLARA y no PAM: es el mismo criterio de k-medoides que adopta el
# informe, pero evita recalcular la matriz de distancias completa en cada una
# de las K.max * (B + 1) llamadas, que es lo que vuelve inviable a PAM aquí.
claraFUN <- function(x, k) {
list(cluster = cluster::clara(x, k, samples = 10,
sampsize = min(nrow(x), 100),
pamLike = TRUE)$clustering)
}
# El gap se evalúa sobre una submuestra de la submuestra: su papel es comparar
# la inercia observada con la de una referencia sin estructura, y esa
# comparación no requiere el conjunto completo.
sub_gap <- sub_pun[sample(nrow(sub_pun), 1000), , drop = FALSE]
gap <- cluster::clusGap(sub_gap, FUN = claraFUN, K.max = k_max, B = 25)
k_gap <- cluster::maxSE(gap$Tab[, "gap"], gap$Tab[, "SE.sim"], method = "firstSEmax")
t_kbl(as.data.frame(gap$Tab) |>
tibble::rownames_to_column("k") |>
select(k, logW, E.logW, gap, SE.sim),
caption = "Estadístico gap de Tibshirani, Walther y Hastie (2001), calculado con k-medoides sobre una submuestra de 1,000 observaciones y B = 25 referencias.",
digits = 4, scroll = TRUE)| k | logW | E.logW | gap | SE.sim |
|---|---|---|---|---|
| 1 | 6.429 | 7.075 | 0.6460 | 0.0118 |
| 2 | 6.080 | 6.712 | 0.6309 | 0.0116 |
| 3 | 5.879 | 6.545 | 0.6666 | 0.0113 |
| 4 | 5.749 | 6.384 | 0.6349 | 0.0137 |
| 5 | 5.670 | 6.256 | 0.5861 | 0.0088 |
| 6 | 5.593 | 6.146 | 0.5532 | 0.0134 |
| 7 | 5.514 | 6.083 | 0.5683 | 0.0148 |
| 8 | 5.461 | 6.017 | 0.5562 | 0.0127 |
El gap compara \(\log W_k\) observado con su esperanza bajo una referencia sin estructura y retiene el menor \(k\) tal que \(\mathrm{Gap}(k)\ge\mathrm{Gap}(k+1)-s_{k+1}\), criterio que incorpora la incertidumbre de simulación en vez de tomar el máximo.
# ============================================================
# FIGURA: CRITERIOS DE SELECCIÓN DE k
# ============================================================
indices_k |>
pivot_longer(-k, names_to = "Criterio", values_to = "Valor") |>
ggplot(aes(x = k, y = Valor)) +
geom_line(colour = PAS_CIELO, linewidth = 1.1) +
geom_point(colour = "#7FA8D0", size = 2.4) +
facet_wrap(~ Criterio, scales = "free_y", ncol = 2) +
scale_x_continuous(breaks = 2:k_max) +
labs(title = "Criterios de selección del número de grupos",
subtitle = "Silueta y Calinski-Harabasz: mayor es mejor. Davies-Bouldin y distancia al medoide: menor es mejor",
x = "Número de grupos (k)", y = NULL)
Figure 5.2: Evolución conjunta de los criterios de selección del número de grupos.
# ============================================================
# TABLA: k SUGERIDO POR CADA CRITERIO
# ============================================================
resumen_k <- tibble::tibble(
Criterio = c("Silueta media (máximo)", "Calinski-Harabasz (máximo)",
"Davies-Bouldin (mínimo)", "Gap de Tibshirani"),
`k sugerido` = c(indices_k$k[which.max(indices_k$`Silueta media`)],
indices_k$k[which.max(indices_k$`Calinski-Harabasz`)],
indices_k$k[which.min(indices_k$`Davies-Bouldin`)],
k_gap)
)
t_kbl(resumen_k,
caption = "Número de grupos sugerido por cada criterio.",
align = c("l","r"))| Criterio | k sugerido |
|---|---|
| Silueta media (máximo) | 2 |
| Calinski-Harabasz (máximo) | 2 |
| Davies-Bouldin (mínimo) | 4 |
| Gap de Tibshirani | 1 |
5.4 Partición
# ============================================================
# TABLA: PARTICIÓN POR k-MEDOIDES (CLARA)
# ============================================================
k_final <- as.integer(names(sort(table(resumen_k$`k sugerido`), decreasing = TRUE))[1])
# CLARA extrae varias submuestras, aplica PAM en cada una y conserva el conjunto
# de medoides con menor distancia media sobre el total de observaciones.
clara_fit <- cluster::clara(puntuaciones[, cp_vars, drop = FALSE], k = k_final,
samples = 50, sampsize = min(nrow(puntuaciones), 200),
pamLike = TRUE, rngR = TRUE)
viv_c$cluster <- factor(clara_fit$clustering)
puntuaciones$cluster <- viv_c$cluster
medoides <- as.data.frame(clara_fit$medoids)
names(medoides) <- cp_vars
medoides$cluster <- factor(seq_len(k_final))
t_kbl(viv_c |> count(cluster, name = "n") |>
mutate(`% del total` = round(100 * n / sum(n), 2)) |>
bind_cols(round(medoides[, cp_vars], 3)),
caption = "Tamaño de cada grupo y coordenadas factoriales de su medoide.",
align = c("l","r","r","r","r"), scroll = TRUE)| cluster | n | % del total | CP1 | CP2 | |
|---|---|---|---|---|---|
| 3096 | 1 | 4,741 | 56.99 | -1.322 | -0.034 |
| 830 | 2 | 3,578 | 43.01 | 1.350 | 0.094 |
El medoide es una observación real del registro, de modo que cada grupo queda representado por un inmueble concreto y no por un promedio que puede no existir en el mercado. Esa propiedad, junto con su resistencia a los valores extremos, es la razón de la elección.
# ============================================================
# FIGURA: GRUPOS EN EL PLANO FACTORIAL
# ============================================================
# La envolvente convexa describe EXTENSIÓN, no separación: dos grupos contiguos
# de un continuo producen envolventes que se tocan.
envolventes <- puntuaciones |>
group_by(cluster) |>
slice(grDevices::chull(CP1, CP2)) |>
ungroup()
ggplot(puntuaciones, aes(x = CP1, y = CP2, colour = cluster)) +
geom_hline(yintercept = 0, colour = "#DEDEDE", linewidth = 0.4) +
geom_vline(xintercept = 0, colour = "#DEDEDE", linewidth = 0.4) +
geom_polygon(data = envolventes, aes(fill = cluster), alpha = 0.18,
colour = NA, show.legend = FALSE) +
geom_point(alpha = 0.30, size = 0.7) +
stat_ellipse(type = "norm", level = 0.68, linewidth = 0.9) +
geom_point(data = medoides, size = 4.4, shape = 21, fill = "white",
stroke = 1.5, show.legend = FALSE) +
scale_colour_manual(values = PAL_CAT, name = "Grupo") +
scale_fill_manual(values = PAL_CAT, guide = "none") +
labs(title = "Los grupos se ordenan sobre el primer plano factorial",
x = paste0("CP1 (", round(vp$`% de varianza`[1], 1), "%)"),
y = paste0("CP2 (", round(vp$`% de varianza`[2], 1), "%)")) +
guides(colour = guide_legend(override.aes = list(alpha = 1, size = 3)))
Figure 5.3: Grupos en el primer plano factorial. La envolvente delimita la extensión de cada grupo; el punto blanco marca su medoide.
5.4.1 Densidad entre grupos
La inspección visual no decide si entre los grupos media una región de baja densidad. El contraste se hace sobre la recta que une dos medoides. Con \(\mathbf{d}\) el vector unitario entre ellos, se proyecta cada individuo, \(u_i=(\mathbf{x}_i-\bar{\mathbf{x}})^{\top}\mathbf{d}\), y se examina la densidad de \(u\). El criterio no es la existencia de dos máximos, que una distribución asimétrica puede producir, sino la profundidad de la depresión entre ellos, \[\nu=1-\frac{f(\text{mín})}{f(\text{modo menor})}\in[0,1],\] junto con la densidad en el punto de corte.
# ============================================================
# FIGURA: DENSIDAD EN LA DIRECCIÓN DE SEPARACIÓN
# ============================================================
M <- as.matrix(puntuaciones[, cp_vars, drop = FALSE])
C <- as.matrix(medoides[, cp_vars, drop = FALSE])
d_vec <- C[2, ] - C[1, ]
d_vec <- d_vec / sqrt(sum(d_vec^2))
u <- as.numeric(sweep(M, 2, colMeans(M)) %*% d_vec)
u_cent <- as.numeric(sweep(C, 2, colMeans(M)) %*% d_vec)
u_frontera <- mean(u_cent[1:2])
ggplot(tibble::tibble(u = u, Grupo = puntuaciones$cluster), aes(x = u)) +
geom_density(aes(fill = Grupo, colour = Grupo), alpha = 0.35, linewidth = 0.6) +
geom_density(colour = TINTA, linewidth = 1.0) +
geom_vline(xintercept = u_frontera, colour = "#C98F63", linetype = "dashed", linewidth = 0.9) +
geom_vline(xintercept = u_cent, colour = "#9A9A9A", linetype = "dotted", linewidth = 0.6) +
scale_fill_manual(values = PAL_CAT) +
scale_colour_manual(values = PAL_CAT) +
labs(title = "Densidad a lo largo de la dirección que separa los grupos",
subtitle = "Curva oscura: densidad de toda la oferta. Discontinua: punto medio entre medoides",
x = "Proyección sobre la recta que une los medoides", y = "Densidad")
Figure 5.4: Densidad de los individuos proyectados sobre la recta que une los dos medoides principales. La línea discontinua marca el punto medio entre ambos.
# ============================================================
# TABLA: DIAGNÓSTICO DE MULTIMODALIDAD
# ============================================================
dens_u <- density(u, n = 1024)
# Un valle es un mínimo local interior. Se identifican los cambios de signo de
# la primera diferencia para contar modos y valles de forma explícita.
dy <- diff(dens_u$y)
maximos <- which(dy[-length(dy)] > 0 & dy[-1] <= 0) + 1
minimos <- which(dy[-length(dy)] < 0 & dy[-1] >= 0) + 1
if (length(maximos) >= 2) {
ord <- order(dens_u$y[maximos], decreasing = TRUE)[1:2]
m1 <- sort(maximos[ord])
entre <- minimos[minimos > m1[1] & minimos < m1[2]]
valle <- if (length(entre) > 0) 1 - min(dens_u$y[entre]) / min(dens_u$y[m1]) else 0
} else {
valle <- 0
}
f_frontera <- approx(dens_u$x, dens_u$y, u_frontera)$y
t_kbl(tibble::tibble(
`Modos detectados` = length(maximos),
`Profundidad del valle ν` = round(valle, 4),
`Densidad en la frontera` = round(f_frontera, 4),
`Densidad máxima` = round(max(dens_u$y), 4),
`% de la densidad máxima` = round(100 * f_frontera / max(dens_u$y), 1),
Lectura = if (length(maximos) < 2) {
"Unimodal: no existe región de baja densidad entre los grupos"
} else if (valle < 0.10) {
"Dos máximos sin valle apreciable (ν < 0.10): compatible con un continuo asimétrico"
} else if (valle < 0.33) {
"Bimodalidad débil: existe una depresión poco marcada"
} else {
"Bimodalidad marcada: existe una región de baja densidad entre los grupos"
}),
caption = "Diagnóstico de multimodalidad en la dirección de separación.",
digits = 4, scroll = TRUE)| Modos detectados | Profundidad del valle ν | Densidad en la frontera | Densidad máxima | % de la densidad máxima | Lectura |
|---|---|---|---|---|---|
| 4 | 0 | 0.1869 | 0.3467 | 53.9 | Dos máximos sin valle apreciable (ν < 0.10): compatible con un continuo asimétrico |
5.5 Validación
5.5.1 Cohesión y separación
Con \(a(i)\) la disimilaridad media dentro del grupo y \(b(i)\) la mínima frente a otro grupo, el ancho de silueta es \(s(i)=\bigl(b(i)-a(i)\bigr)/\max\{a(i),b(i)\}\), acotado en \([-1,1]\). Los valores negativos señalan observaciones mejor situadas en otro grupo (Rousseeuw, 1987).
# ============================================================
# TABLA: ANCHO DE SILUETA POR CONGLOMERADO
# ============================================================
# La silueta se evalúa sobre las mismas coordenadas en que se trazó la
# frontera: calcularla en otro espacio no mediría la cohesión de esta partición.
sil_idx <- sample(nrow(puntuaciones), 3000)
sil_fin <- cluster::silhouette(as.integer(puntuaciones$cluster[sil_idx]),
dist(puntuaciones[sil_idx, cp_vars, drop = FALSE]))
t_kbl(as.data.frame(sil_fin[, 1:3]) |>
group_by(Grupo = cluster) |>
summarise(`n en la submuestra` = n(),
`Silueta media` = mean(sil_width),
`% con silueta < 0` = round(100 * mean(sil_width < 0), 2),
.groups = "drop") |>
mutate(`Interpretación` = etiqueta_efecto(`Silueta media`,
c(-Inf, .25, .50, .70, Inf),
c("Estructura débil o ausente","Estructura razonable",
"Estructura sólida","Estructura muy marcada"))),
caption = "Ancho de silueta por grupo (submuestra de 3,000 observaciones).",
digits = 4, scroll = TRUE)| Grupo | n en la submuestra | Silueta media | % con silueta < 0 | Interpretación |
|---|---|---|---|---|
| 1 | 1,777 | 0.6059 | 0.00 | Estructura sólida |
| 2 | 1,223 | 0.2083 | 20.36 | Estructura débil o ausente |
5.5.2 Estabilidad por remuestreo
# ============================================================
# TABLA: ESTABILIDAD BOOTSTRAP (JACCARD)
# ============================================================
# `claraCBI` de fpc reproduce el criterio adoptado: en cada réplica vuelve a
# aplicar k-medoides sobre los datos remuestreados. No hace falta una interfaz
# propia, a diferencia de lo que ocurriría con una partición jerárquica
# consolidada, que clusterboot no puede reconstruir con su método por defecto.
boot_idx <- sample(nrow(puntuaciones), 2000)
# `usepam = FALSE` hace que claraCBI use CLARA y no PAM en cada réplica, que es
# el algoritmo efectivamente empleado en la partición reportada y el único
# viable con B = 50 réplicas.
cb <- fpc::clusterboot(puntuaciones[boot_idx, cp_vars, drop = FALSE],
B = 50, clustermethod = fpc::claraCBI,
k = k_final, usepam = FALSE, samples = 20,
seed = 2026, count = FALSE)
t_kbl(tibble::tibble(
Grupo = seq_along(cb$bootmean),
`Jaccard medio` = cb$bootmean,
Disoluciones = cb$bootbrd,
Recuperaciones = cb$bootrecover,
Lectura = etiqueta_efecto(cb$bootmean, c(-Inf, .5, .6, .75, .85, Inf),
c("Inestable","Dudoso","Con patrón","Estable","Muy estable"))),
caption = "Estabilidad de los grupos por remuestreo bootstrap (B = 50).",
digits = 4, scroll = TRUE)| Grupo | Jaccard medio | Disoluciones | Recuperaciones | Lectura |
|---|---|---|---|---|
| 1 | 0.9582 | 0 | 50 | Muy estable |
| 2 | 0.9462 | 0 | 50 | Muy estable |
El coeficiente de Jaccard mide la coincidencia entre cada grupo original y su homólogo en cada réplica. Como referencia (Hennig, 2007): por debajo de 0.60 el grupo no es fiable, entre 0.60 y 0.75 hay patrón sin delimitación nítida, y por encima de 0.85 el grupo es estable. Conviene recordar que el remuestreo mide reproducibilidad y no separación: el corte de un continuo denso se reproduce sin que exista discontinuidad.
5.5.3 Contraste con otros criterios de agregación
# ============================================================
# TABLA: CONCORDANCIA CON OTROS MÉTODOS
# ============================================================
# El ACP de reporte se estimó con ncp = 5. Agrupar sobre las cinco componentes
# equivaldría, por invariancia de la distancia euclídea ante rotaciones
# ortogonales, a agrupar sobre las variables estandarizadas. Se reestima el ACP
# truncado en las componentes retenidas para que el contraste sea comparable.
acp_clust <- FactoMineR::PCA(datos_acp, quali.sup = 6:8, scale.unit = TRUE,
ncp = n_ret, graph = FALSE)
stopifnot(ncol(acp_clust$ind$coord) == n_ret,
nrow(acp_clust$ind$coord) == nrow(viv_c))
hcpc <- FactoMineR::HCPC(acp_clust, nb.clust = k_final, kk = 100, consol = TRUE,
graph = FALSE, min = 2, max = k_max)
km_ref <- kmeans(puntuaciones[, cp_vars, drop = FALSE], centers = k_final, nstart = 25)
comp <- tibble::tibble(
`Comparación` = c("k-medoides frente a Ward consolidado",
"k-medoides frente a k-medias directo"),
ARI = c(ari(viv_c$cluster, hcpc$data.clust$clust),
ari(viv_c$cluster, km_ref$cluster)))
t_kbl(comp |>
mutate(Lectura = etiqueta_efecto(ARI, c(-Inf, .20, .40, .60, .80, Inf),
c("Concordancia nula","Concordancia leve","Concordancia moderada",
"Concordancia sustancial: mismos grupos con frontera desplazada",
"Concordancia casi perfecta"))),
caption = "Concordancia entre la partición adoptada y dos criterios de agregación alternativos.",
digits = 4, scroll = TRUE)| Comparación | ARI | Lectura |
|---|---|---|
| k-medoides frente a Ward consolidado | 0.4444 | Concordancia moderada |
| k-medoides frente a k-medias directo | 0.6779 | Concordancia sustancial: mismos grupos con frontera desplazada |
Los tres métodos operan sobre las mismas coordenadas y con el mismo \(k\), de modo que cualquier discrepancia es atribuible al criterio de agregación. Una concordancia alta indicaría grupos bien delimitados; una concordancia parcial indica que los métodos reconocen la misma estructura pero sitúan la frontera en lugares algo distintos, comportamiento propio de un continuo denso.
5.5.4 Concordancia con variables externas
# ============================================================
# TABLA: ARI FRENTE A VARIABLES CATEGÓRICAS
# ============================================================
t_kbl(tibble::tibble(
`Comparación` = c("Grupo frente a estrato", "Grupo frente a zona",
"Grupo frente a tipo de vivienda"),
ARI = c(ari(viv_c$cluster, viv_c$estrato),
ari(viv_c$cluster, viv_c$zona),
ari(viv_c$cluster, viv_c$tipo))) |>
mutate(Lectura = etiqueta_efecto(ARI, c(-Inf, .05, .20, .40, Inf),
c("Concordancia nula","Concordancia leve",
"Concordancia moderada","Concordancia alta"))),
caption = "Índice de Rand ajustado entre la partición obtenida y las variables categóricas.",
digits = 4, scroll = TRUE)| Comparación | ARI | Lectura |
|---|---|---|
| Grupo frente a estrato | 0.0602 | Concordancia leve |
| Grupo frente a zona | 0.0139 | Concordancia nula |
| Grupo frente a tipo de vivienda | 0.3078 | Concordancia moderada |
5.6 Perfil de los grupos
# ============================================================
# TABLA: PERFIL DE LOS CONGLOMERADOS
# ============================================================
t_kbl(viv_c |>
group_by(Grupo = cluster) |>
summarise(n = n(),
`Precio mediano` = median(preciom),
`Precio por m2` = median(precio_m2),
`Área mediana` = median(areaconst),
Habitaciones = median(habitaciones),
`Baños` = median(banios),
Parqueaderos = median(parqueaderos),
`Estrato modal` = names(sort(table(estrato), decreasing = TRUE))[1],
`% casas` = round(100 * mean(tipo == "Casa"), 1),
`Zona predominante` = names(sort(table(zona), decreasing = TRUE))[1],
.groups = "drop"),
caption = "Perfil de los grupos según los atributos originales (medianas).",
digits = 1, scroll = TRUE)| Grupo | n | Precio mediano | Precio por m2 | Área mediana | Habitaciones | Baños | Parqueaderos | Estrato modal | % casas | Zona predominante |
|---|---|---|---|---|---|---|---|---|---|---|
| 1 | 4,741 | 240 | 2.8 | 85 | 3 | 2 | 1 | 5 | 15.7 | Zona Sur |
| 2 | 3,578 | 560 | 2.3 | 250 | 4 | 4 | 2 | 6 | 69.1 | Zona Sur |
Los grupos se construyeron para maximizar la separación en estas variables, así que el rechazo de \(H_0\) en la tabla siguiente está garantizado por construcción y no prueba que los grupos existan, cosa que deciden la silueta y el bootstrap. La pregunta que sí responde es cuáles atributos separan más.
# ============================================================
# TABLA: CAPACIDAD DISCRIMINANTE DE CADA VARIABLE
# ============================================================
t_kbl(bind_rows(lapply(c(num_activas, "precio_m2"), function(v) {
kw <- kruskal.test(viv_c[[v]] ~ viv_c$cluster)
e2 <- as.numeric(kw$statistic) / (nrow(viv_c) - 1)
tibble::tibble(
Variable = v,
`Estadístico H` = round(as.numeric(kw$statistic), 1),
gl = as.integer(kw$parameter),
`Valor p` = format.pval(kw$p.value, digits = 3, eps = 1e-16),
`ε²` = round(e2, 4),
Magnitud = etiqueta_efecto(e2, c(-Inf, .01, .06, .14, Inf),
c("Despreciable","Pequeño","Moderado","Grande"))
)
})) |> arrange(desc(`ε²`)),
caption = "Capacidad discriminante de cada variable entre grupos, ordenada por tamaño del efecto.",
scroll = TRUE)| Variable | Estadístico H | gl | Valor p | ε² | Magnitud |
|---|---|---|---|---|---|
| areaconst | 5,547.5 | 1 | <0.0000000000000001 | 0.667 | Grande |
| banios | 5,081.3 | 1 | <0.0000000000000001 | 0.611 | Grande |
| preciom | 4,003.2 | 1 | <0.0000000000000001 | 0.481 | Grande |
| habitaciones | 3,816.5 | 1 | <0.0000000000000001 | 0.459 | Grande |
| parqueaderos | 3,360.7 | 1 | <0.0000000000000001 | 0.404 | Grande |
| precio_m2 | 259.6 | 1 | <0.0000000000000001 | 0.031 | Pequeño |
# ============================================================
# FIGURA: PERFIL ESTANDARIZADO
# ============================================================
viv_c |>
select(cluster, all_of(c(num_activas, "precio_m2"))) |>
mutate(across(all_of(c(num_activas, "precio_m2")), ~ as.numeric(scale(.)))) |>
pivot_longer(-cluster, names_to = "Variable", values_to = "z") |>
group_by(cluster, Variable) |>
summarise(z_medio = mean(z), .groups = "drop") |>
ggplot(aes(x = Variable, y = z_medio, fill = cluster)) +
geom_col(position = position_dodge(width = 0.8), width = 0.72) +
geom_hline(yintercept = 0, colour = "#7A7A7A", linewidth = 0.5) +
scale_fill_manual(values = PAL_CAT, name = "Grupo") +
labs(title = "Perfil estandarizado de cada grupo",
subtitle = "El valor cero corresponde a la media general de la oferta",
x = NULL, y = "Puntuación estandarizada media")
Figure 5.5: Perfil de cada conglomerado en las variables activas, expresadas en puntuaciones estandarizadas para permitir la comparación entre escalas.
5.7 Distribución en el territorio
# ============================================================
# FIGURA: MAPA ESTÁTICO
# ============================================================
mapa_df <- viv_c |> select(longitud, latitud, cluster, zona) |> tidyr::drop_na()
ggplot(mapa_df, aes(x = longitud, y = latitud, colour = cluster)) +
geom_point(alpha = 0.35, size = 0.7) +
scale_colour_manual(values = PAL_CAT, name = "Grupo") +
coord_quickmap() +
labs(title = "Ambos segmentos coexisten en toda la ciudad",
subtitle = "La segmentación por gama no reproduce una división geográfica",
x = "Longitud", y = "Latitud") +
guides(colour = guide_legend(override.aes = list(alpha = 1, size = 3)))
Figure 5.6: Distribución espacial de los grupos sobre las coordenadas geográficas de los inmuebles.
# ============================================================
# FIGURA: MAPA INTERACTIVO
# ============================================================
set.seed(2026)
mapa_sub <- mapa_df[sample(nrow(mapa_df), min(2000, nrow(mapa_df))), ]
paleta <- leaflet::colorFactor(PAL_CAT[seq_len(k_final)], domain = mapa_sub$cluster)
leaflet::leaflet(mapa_sub) |>
leaflet::addProviderTiles(leaflet::providers$CartoDB.Positron) |>
leaflet::addCircleMarkers(
lng = ~longitud, lat = ~latitud, radius = 3, stroke = FALSE,
fillOpacity = 0.6, color = ~paleta(cluster),
popup = ~paste0("Grupo: ", cluster, "<br>Zona: ", zona)) |>
leaflet::addLegend("bottomright", pal = paleta, values = ~cluster,
title = "Grupo", opacity = 0.8)Figure 5.7: Mapa interactivo de la oferta segmentada sobre una muestra aleatoria, para mantener acotado el tamaño del documento.
# ============================================================
# TABLA: COMPOSICIÓN ZONAL DE CADA CONGLOMERADO
# ============================================================
tabla_cz <- table(viv_c$cluster, viv_c$zona)
prueba_cz <- chisq.test(tabla_cz)
t_kbl(as.data.frame.matrix(round(100 * prop.table(tabla_cz, 1), 1)) |>
tibble::rownames_to_column("Grupo"),
caption = "Distribución porcentual de cada grupo entre las zonas de la ciudad.",
scroll = TRUE)| Grupo | Zona Centro | Zona Norte | Zona Oeste | Zona Oriente | Zona Sur |
|---|---|---|---|---|---|
| 1 | 1.1 | 26.3 | 10.6 | 3.5 | 58.4 |
| 2 | 2.0 | 18.8 | 19.5 | 5.1 | 54.7 |
# ============================================================
# TABLA: ASOCIACIÓN CONGLOMERADO x ZONA
# ============================================================
t_kbl(tibble::tibble(
`Estadístico χ²` = round(as.numeric(prueba_cz$statistic), 1),
gl = as.integer(prueba_cz$parameter),
`Valor p` = format.pval(prueba_cz$p.value, digits = 3, eps = 1e-16),
`V de Cramér` = round(cramer_v(tabla_cz), 4)),
caption = "Asociación entre la segmentación obtenida y la zona de la ciudad.")| Estadístico χ² | gl | Valor p | V de Cramér |
|---|---|---|---|
| 188 | 4 | <0.0000000000000001 | 0.15 |
5.8 Lectura de la segmentación
La combinación de resultados permite decidir qué es la partición obtenida. La matriz de disimilaridad ordenada, la densidad a lo largo de la dirección que une los medoides, el ancho de silueta y la concordancia con otros criterios de agregación apuntan en la misma dirección: se trata de una división operativa de un mercado continuo, no del descubrimiento de subpoblaciones separadas. Eso es coherente con el resultado del ACP, donde la variación se resume en muy pocas direcciones.
La consecuencia práctica no es negativa. Una división reproducible y de perfiles contrastados permite organizar líneas de producto, políticas de precio y equipos comerciales, siempre que se acepte que los inmuebles situados cerca de la frontera admiten cualquiera de las dos clasificaciones, y que su asignación es una decisión comercial y no un dato del mercado.
Al haberse construido el espacio sobre el precio unitario, los grupos no se ordenan solo por precio total. Las tablas de perfil permiten comprobar si la separación responde al valor por metro cuadrado, al tamaño o a la dotación, que es la distinción con implicaciones para la valoración. La concordancia con la zona y con el estrato indica hasta qué punto la segmentación coincide con la geografía del mercado o la corta transversalmente.
6 Asociación entre atributos categóricos
6.1 Planteamiento
Para una tabla de contingencia con frecuencias relativas \(p_{ij}\), el análisis descompone la inercia total \[\Phi^{2}=\frac{\chi^{2}}{n}=\sum_{i,j}\frac{(p_{ij}-p_{i\cdot}p_{\cdot j})^{2}}{p_{i\cdot}p_{\cdot j}}\] en ejes ortogonales, con la distancia \(\chi^2\) entre perfiles: \[d^{2}_{\chi^2}(i,i')=\sum_{j}\frac{1}{p_{\cdot j}}\left(\frac{p_{ij}}{p_{i\cdot}}-\frac{p_{i'j}}{p_{i'\cdot}}\right)^{2}.\] La métrica pondera cada columna por el inverso de su masa, de modo que las categorías poco frecuentes pesan igual que las frecuentes. De ahí la equivalencia distribucional: fusionar categorías de perfil idéntico no altera la solución (Greenacre, 2017; Lebart et al., 2006). El número de dimensiones no triviales es \(\min(r-1, c-1)\).
Este apartado se organiza en torno al barrio, que es la unidad de mercado con la que trabaja el equipo comercial. El cruce con el estrato sustituye al cruce con la zona, porque la zona agrupa barrios muy heterogéneos y ya quedó representada en el ACP.
6.2 Estado de la variable barrio
# ============================================================
# TABLA: AUDITORÍA DE `barrio`
# ============================================================
fb <- sort(table(viv_c$barrio), decreasing = TRUE)
# Un barrio debería pertenecer a una sola zona: la relación es jerárquica.
barrio_zona <- viv_c |> count(barrio, zona) |> count(barrio, name = "zonas_distintas")
t_kbl(tibble::tribble(
~Indicador, ~Valor, ~`Consecuencia`,
"Categorías totales", as.character(length(fb)),
"Cardinalidad muy alta frente al tamaño muestral",
"Categorías con 20 registros o menos", as.character(sum(fb <= 20)),
"Frecuencias esperadas insuficientes para la aproximación χ²",
"Categorías con 5 registros o menos", as.character(sum(fb <= 5)),
"Masa mínima con contribución desproporcionada a la inercia",
"Categorías para cubrir el 50 % de la oferta", as.character(which(cumsum(fb)/sum(fb) >= 0.5)[1]),
"Concentración fuerte en pocas categorías",
"Barrios asignados a más de una zona", as.character(sum(barrio_zona$zonas_distintas > 1)),
"La relación barrio–zona no es funcional",
"Barrios asignados a tres zonas", as.character(sum(barrio_zona$zonas_distintas == 3)),
"Inconsistencia grave de georreferenciación"),
caption = "Auditoría de calidad de la variable `barrio`.",
align = c("l","r","l"), scroll = TRUE)| Indicador | Valor | Consecuencia |
|---|---|---|
| Categorías totales | 407 | Cardinalidad muy alta frente al tamaño muestral |
| Categorías con 20 registros o menos | 324 | Frecuencias esperadas insuficientes para la aproximación χ² |
| Categorías con 5 registros o menos | 240 | Masa mínima con contribución desproporcionada a la inercia |
| Categorías para cubrir el 50 % de la oferta | 16 | Concentración fuerte en pocas categorías |
| Barrios asignados a más de una zona | 93 | La relación barrio–zona no es funcional |
| Barrios asignados a tres zonas | 14 | Inconsistencia grave de georreferenciación |
La revisión detecta cuatro problemas simultáneos: variantes de escritura de una misma unidad, denominaciones truncadas, nombres de zona registrados como barrio e inconsistencia jerárquica, con barrios asignados a más de una zona.
Tratamiento elegido. El umbral no se fija por frecuencia mínima sino por cobertura: se retienen los barrios más frecuentes hasta cubrir el 60 % de la oferta, criterio que se justifica solo, porque fija de antemano qué proporción del mercado queda representada de forma nominal. El resto se agrupa en una categoría residual que reúne perfiles heterogéneos: su baricentro tiende al origen por composición y no se interpreta. La equivalencia distribucional autoriza fusionar categorías de perfil idéntico, condición que estas no cumplen, de modo que la agrupación es un recurso pragmático y así se declara.
# ============================================================
# TABLA: AGRUPACIÓN DE `barrio` POR COBERTURA
# ============================================================
cobertura_objetivo <- 0.60
orden_acum <- cumsum(fb) / sum(fb)
n_retenidos <- which(orden_acum >= cobertura_objetivo)[1]
barrios_frec <- names(fb)[seq_len(n_retenidos)]
# Se excluyen las denominaciones que en realidad son nombres de zona
barrios_frec <- setdiff(barrios_frec, levels(viv_c$zona))
viv_c <- viv_c |>
mutate(barrio_agr = factor(ifelse(barrio %in% barrios_frec, barrio, "Resto de barrios")))
t_kbl(tibble::tibble(
`Cobertura objetivo` = paste0(100 * cobertura_objetivo, " %"),
`Barrios retenidos` = length(barrios_frec),
`Frecuencia del menor retenido` = as.integer(min(fb[barrios_frec])),
`Cobertura alcanzada` = paste0(round(100 * mean(viv_c$barrio %in% barrios_frec), 1), " %"),
`Registros en la residual` = sum(viv_c$barrio_agr == "Resto de barrios")),
caption = "Resultado de la agrupación de `barrio` por criterio de cobertura.",
align = rep("r", 5), scroll = TRUE)| Cobertura objetivo | Barrios retenidos | Frecuencia del menor retenido | Cobertura alcanzada | Registros en la residual |
|---|---|---|---|---|
| 60 % | 26 | 60 | 59.4 % | 3,380 |
6.3 Tipo de vivienda y estrato
# ============================================================
# TABLA: PERFILES FILA DE TIPO x ESTRATO
# ============================================================
tc_te <- table(viv_c$tipo, viv_c$estrato)
t_kbl(as.data.frame.matrix(round(100 * prop.table(tc_te, 1), 1)) |>
tibble::rownames_to_column("Tipo"),
caption = "Distribución porcentual de cada tipo de vivienda entre los estratos observados.")| Tipo | 3 | 4 | 5 | 6 |
|---|---|---|---|---|
| Apartamento | 12.5 | 27.5 | 34.6 | 25.3 |
| Casa | 25.3 | 22.5 | 30.6 | 21.6 |
# ============================================================
# TABLA: DESCOMPOSICIÓN DE LA ASOCIACIÓN TIPO x ESTRATO
# ============================================================
prueba_te <- chisq.test(tc_te)
inercia_te <- as.numeric(prueba_te$statistic) / sum(tc_te)
t_kbl(tibble::tibble(
`Chi cuadrado` = round(as.numeric(prueba_te$statistic), 2),
gl = as.integer(prueba_te$parameter),
`Valor p` = format.pval(prueba_te$p.value, digits = 4, eps = 1e-16),
`Inercia total` = round(inercia_te, 5),
`V de Cramér` = round(sqrt(inercia_te / (min(dim(tc_te)) - 1)), 4),
`Dimensiones no triviales` = min(dim(tc_te)) - 1,
`Frecuencia esperada mínima` = round(min(prueba_te$expected), 1)),
caption = "Asociación entre tipo de vivienda y estrato socioeconómico.",
scroll = TRUE)| Chi cuadrado | gl | Valor p | Inercia total | V de Cramér | Dimensiones no triviales | Frecuencia esperada mínima |
|---|---|---|---|---|---|---|
| 224.3 | 3 | < 0.0000000000000001 | 0.027 | 0.164 | 1 | 562.2 |
Como el tipo de vivienda tiene dos categorías, \(\min(r-1,c-1)=1\): la tabla posee una sola dimensión no trivial, que recoge el 100 % de la inercia. No procede mapa en el plano. El análisis se reduce a ordenar los estratos sobre una recta que va del predominio de apartamentos al de casas.
# ============================================================
# FIGURA: EJE FACTORIAL TIPO x ESTRATO
# ============================================================
ca_te <- FactoMineR::CA(as.data.frame.matrix(tc_te), ncp = 1, graph = FALSE)
coord_est <- data.frame(
Estrato = paste("Estrato", rownames(ca_te$col$coord)),
Coord = ca_te$col$coord[, 1],
Masa = ca_te$call$marge.col)
ggplot(coord_est, aes(x = Coord, y = 0)) +
geom_hline(yintercept = 0, colour = "#CFCFCF", linewidth = 0.5) +
geom_point(aes(size = Masa), colour = PAS_CIELO, alpha = 0.95) +
ggrepel::geom_text_repel(aes(label = Estrato), size = 3.6, nudge_y = 0.02, colour = TINTA) +
scale_size_continuous(range = c(4, 11), name = "Masa") +
scale_y_continuous(limits = c(-0.05, 0.06), breaks = NULL) +
labs(title = "Composición tipológica de la oferta según el estrato",
subtitle = "Negativo: predominio de apartamentos. Positivo: predominio de casas",
x = "Coordenada factorial", y = NULL)
Figure 6.1: Coordenadas de los estratos sobre el único eje factorial de la tabla tipo × estrato. El tamaño del punto es proporcional a la masa del estrato.
# ============================================================
# TABLA: CONTRIBUCIONES POR ESTRATO
# ============================================================
t_kbl(data.frame(
Estrato = paste("Estrato", rownames(ca_te$col$coord)),
Masa = round(ca_te$call$marge.col, 4),
Coordenada = round(ca_te$col$coord[, 1], 4),
`Contribución (%)` = round(ca_te$col$contrib[, 1], 2),
`Coseno cuadrado` = round(ca_te$col$cos2[, 1], 4),
check.names = FALSE) |>
arrange(desc(`Contribución (%)`)),
caption = "Masa, coordenada y contribución de cada estrato al eje factorial.")| Estrato | Masa | Coordenada | Contribución (%) | Coseno cuadrado | |
|---|---|---|---|---|---|
| 3 | Estrato 3 | 0.175 | 0.356 | 81.98 | 1 |
| 4 | Estrato 4 | 0.256 | -0.095 | 8.62 | 1 |
| 6 | Estrato 6 | 0.239 | -0.075 | 5.02 | 1 |
| 5 | Estrato 5 | 0.331 | -0.060 | 4.38 | 1 |
6.4 Barrio y estrato
El cruce anterior se agota en un eje. Para trabajar en el plano se analiza la tabla barrio agrupado por estrato, que admite \(\min(r-1,3)\) dimensiones y es la que interesa al equipo comercial.
# ============================================================
# TABLA: DESCOMPOSICIÓN DE LA INERCIA BARRIO x ESTRATO
# ============================================================
tc_be <- table(viv_c$barrio_agr, viv_c$estrato)
prueba_be <- suppressWarnings(chisq.test(tc_be))
inercia_be <- as.numeric(prueba_be$statistic) / sum(tc_be)
ca_be <- FactoMineR::CA(as.data.frame.matrix(tc_be), graph = FALSE)
eig_be <- as.data.frame(ca_be$eig)
names(eig_be) <- c("Valor propio", "% de inercia", "% acumulado")
eig_be <- tibble::rownames_to_column(eig_be, "Eje")
t_kbl(eig_be,
caption = "Descomposición de la inercia en la tabla barrio × estrato.",
digits = 4, scroll = TRUE)| Eje | Valor propio | % de inercia | % acumulado |
|---|---|---|---|
| dim 1 | 0.5844 | 58.15 | 58.15 |
| dim 2 | 0.2568 | 25.55 | 83.70 |
| dim 3 | 0.1639 | 16.30 | 100.00 |
# ============================================================
# TABLA: CONTRASTE DE INDEPENDENCIA BARRIO x ESTRATO
# ============================================================
t_kbl(tibble::tibble(
`Chi cuadrado` = round(as.numeric(prueba_be$statistic), 1),
gl = as.integer(prueba_be$parameter),
`Valor p` = format.pval(prueba_be$p.value, digits = 4, eps = 1e-16),
`Inercia total` = round(inercia_be, 4),
`V de Cramér` = round(sqrt(inercia_be / (min(dim(tc_be)) - 1)), 4),
`Frecuencia esperada mínima` = round(min(prueba_be$expected), 2)),
caption = "Contraste de independencia entre barrio agrupado y estrato.",
scroll = TRUE)| Chi cuadrado | gl | Valor p | Inercia total | V de Cramér | Frecuencia esperada mínima |
|---|---|---|---|---|---|
| 8,361 | 78 | < 0.0000000000000001 | 1.005 | 0.579 | 10.48 |
La frecuencia esperada mínima condiciona la validez de la aproximación. Si algún valor queda por debajo de 5, el estadístico \(\chi^2\) pierde fiabilidad, aunque la descomposición factorial de la inercia sigue siendo un descriptor legítimo de la tabla observada.
# ============================================================
# FIGURA: MAPA FACTORIAL BARRIO x ESTRATO
# ============================================================
puntos_be <- rbind(
data.frame(Etiqueta = rownames(ca_be$row$coord), ca_be$row$coord[, 1:2],
Elemento = "Barrio", Masa = ca_be$call$marge.row),
data.frame(Etiqueta = paste("Estrato", rownames(ca_be$col$coord)),
ca_be$col$coord[, 1:2], Elemento = "Estrato", Masa = ca_be$call$marge.col))
names(puntos_be)[2:3] <- c("Dim1", "Dim2")
ggplot(puntos_be, aes(x = Dim1, y = Dim2, colour = Elemento, size = Masa)) +
geom_hline(yintercept = 0, colour = "#DEDEDE", linewidth = 0.4) +
geom_vline(xintercept = 0, colour = "#DEDEDE", linewidth = 0.4) +
geom_point(alpha = 0.9) +
ggrepel::geom_text_repel(aes(label = Etiqueta), size = 3.2, max.overlaps = 25,
colour = TINTA, show.legend = FALSE) +
scale_colour_manual(values = c(Barrio = PAS_CIELO, Estrato = PAS_DURAZNO)) +
scale_size_continuous(range = c(2, 8), guide = "none") +
labs(title = "Los barrios se ordenan por composición socioeconómica",
subtitle = "El tamaño del punto es proporcional a la masa de la categoría",
x = paste0("Eje 1 (", round(eig_be$`% de inercia`[1], 1), "%)"),
y = paste0("Eje 2 (", round(eig_be$`% de inercia`[2], 1), "%)"))
Figure 6.2: Mapa factorial simétrico de la tabla barrio × estrato. La proximidad entre un barrio y un estrato indica sobrerrepresentación relativa.
# ============================================================
# TABLA: CONTRIBUCIONES EN EL PLANO BARRIO x ESTRATO
# ============================================================
t_kbl(rbind(
data.frame(Categoria = rownames(ca_be$row$coord), Elemento = "Barrio",
`Contribución eje 1` = ca_be$row$contrib[, 1],
`Contribución eje 2` = ca_be$row$contrib[, 2],
`Coseno en el plano` = rowSums(ca_be$row$cos2[, 1:2]), check.names = FALSE),
data.frame(Categoria = paste("Estrato", rownames(ca_be$col$coord)), Elemento = "Estrato",
`Contribución eje 1` = ca_be$col$contrib[, 1],
`Contribución eje 2` = ca_be$col$contrib[, 2],
`Coseno en el plano` = rowSums(ca_be$col$cos2[, 1:2]), check.names = FALSE)) |>
arrange(desc(`Contribución eje 1`)),
caption = "Contribuciones y calidad de representación en el plano barrio × estrato.",
digits = 3, scroll = TRUE, height = "340px")| Categoria | Elemento | Contribución eje 1 | Contribución eje 2 | Coseno en el plano | |
|---|---|---|---|---|---|
| 6 | Estrato 6 | Estrato | 74.781 | 0.005 | 0.995 |
| Pance | Pance | Barrio | 23.081 | 0.061 | 0.988 |
| Ciudad Jardín | Ciudad Jardín | Barrio | 22.738 | 0.317 | 0.997 |
| 4 | Estrato 4 | Estrato | 15.040 | 2.344 | 0.501 |
| Santa Teresita | Santa Teresita | Barrio | 12.253 | 0.133 | 0.994 |
| Valle Del Lili | Valle Del Lili | Barrio | 7.773 | 12.706 | 0.795 |
| Normandía | Normandía | Barrio | 6.573 | 0.103 | 0.999 |
| 5 | Estrato 5 | Estrato | 6.126 | 23.530 | 0.612 |
| Resto de barrios | Resto de barrios | Barrio | 5.386 | 37.117 | 0.975 |
| 3 | Estrato 3 | Estrato | 4.053 | 74.121 | 0.968 |
| Los Cristales | Los Cristales | Barrio | 3.774 | 0.476 | 0.990 |
| Parcelaciones Pance | Parcelaciones Pance | Barrio | 3.770 | 0.002 | 0.985 |
| El Caney | El Caney | Barrio | 2.104 | 0.817 | 0.368 |
| Cristales | Cristales | Barrio | 1.977 | 0.278 | 0.984 |
| La Flora | La Flora | Barrio | 1.367 | 10.313 | 0.543 |
| Aguacatal | Aguacatal | Barrio | 1.153 | 0.378 | 0.844 |
| El Refugio | El Refugio | Barrio | 1.129 | 0.613 | 0.461 |
| Ciudad 2000 | Ciudad 2000 | Barrio | 1.052 | 0.021 | 0.261 |
| Prados Del Norte | Prados Del Norte | Barrio | 1.029 | 1.257 | 0.758 |
| Caney | Caney | Barrio | 0.812 | 0.586 | 0.389 |
| La Hacienda | La Hacienda | Barrio | 0.775 | 3.908 | 0.669 |
| El Limonar | El Limonar | Barrio | 0.666 | 1.968 | 0.830 |
| El Peñon | El Peñon | Barrio | 0.535 | 0.474 | 0.848 |
| Nueva Tequendama | Nueva Tequendama | Barrio | 0.418 | 1.436 | 0.873 |
| Brisas De Los | Brisas De Los | Barrio | 0.403 | 15.390 | 0.954 |
| Versalles | Versalles | Barrio | 0.388 | 1.659 | 0.817 |
| Urbanización La Flora | Urbanización La Flora | Barrio | 0.334 | 2.639 | 0.482 |
| Quintas De Don | Quintas De Don | Barrio | 0.281 | 2.081 | 0.555 |
| Santa Isabel | Santa Isabel | Barrio | 0.126 | 0.873 | 0.828 |
| Acopi | Acopi | Barrio | 0.094 | 0.148 | 0.933 |
| El Ingenio | El Ingenio | Barrio | 0.008 | 4.247 | 0.469 |
6.5 Análisis de correspondencias múltiple
El ACM extiende el análisis a \(Q\) variables mediante la tabla disyuntiva completa. Su inercia total es \((J-Q)/Q\), cantidad que depende solo del diseño de la codificación y no de la asociación entre variables, de modo que los porcentajes brutos subestiman la estructura recogida. Se reportan la corrección de Benzécri, que retiene los ejes con \(\lambda>1/Q\) y los reescala, y el ajuste de Greenacre, menos optimista, que además descuenta la inercia de los bloques diagonales de la matriz de Burt (Greenacre, 2017, cap. 19).
Las variables activas son tipo, estrato y barrio agrupado. Entran como
suplementarias la zona, el grupo obtenido en la segmentación y el indicador
piso_declarado, lo que permite comprobar si la calidad del anuncio se asocia a
algún perfil de mercado.
# ============================================================
# TABLA: CANTIDADES DE LA CORRECCIÓN DE INERCIA
# ============================================================
datos_acm <- viv_c |>
mutate(estrato_lab = factor(paste("Estrato", estrato))) |>
select(tipo, estrato_lab, barrio_agr, zona, cluster, piso_declarado) |>
tidyr::drop_na() |>
droplevels()
Q_acm <- 3
J_acm <- sum(sapply(datos_acm[, 1:3], nlevels))
# El número de dimensiones no triviales es exactamente J - Q. FactoMineR trunca
# `$eig` al valor de `ncp`, así que se solicita el número completo para que las
# correcciones no se calculen sobre una suma parcial.
ncp_acm <- J_acm - Q_acm
acm <- FactoMineR::MCA(datos_acm, quali.sup = 4:6, graph = FALSE, ncp = ncp_acm)
lambda <- acm$eig[, 1]
t_kbl(tibble::tibble(
Concepto = c("Variables activas (Q)", "Categorías activas (J)",
"Dimensiones esperadas (J - Q)", "Dimensiones devueltas",
"Suma de valores propios", "Inercia teórica (J - Q)/Q",
"Suma de cuadrados de los valores propios",
"Cota inferior por Cauchy-Schwarz"),
Valor = c(Q_acm, J_acm, J_acm - Q_acm, length(lambda),
round(sum(lambda), 4), round((J_acm - Q_acm)/Q_acm, 4),
round(sum(lambda^2), 4), round((J_acm - Q_acm)/Q_acm^2, 4))),
caption = "Verificación de las cantidades que intervienen en la corrección de la inercia.",
align = c("l","r"))| Concepto | Valor |
|---|---|
| Variables activas (Q) | 3.000 |
| Categorías activas (J) | 33.000 |
| Dimensiones esperadas (J - Q) | 30.000 |
| Dimensiones devueltas | 30.000 |
| Suma de valores propios | 10.000 |
| Inercia teórica (J - Q)/Q | 10.000 |
| Suma de cuadrados de los valores propios | 3.589 |
| Cota inferior por Cauchy-Schwarz | 3.333 |
# ============================================================
# TABLA: INERCIA BRUTA Y CORREGIDA
# ============================================================
lam_b <- ifelse(lambda > 1/Q_acm, ((Q_acm/(Q_acm - 1)) * (lambda - 1/Q_acm))^2, 0)
pct_benzecri <- 100 * lam_b / sum(lam_b)
I_aj <- (Q_acm/(Q_acm - 1)) * (sum(lambda^2) - (J_acm - Q_acm)/Q_acm^2)
# Por Cauchy-Schwarz la inercia ajustada es no negativa por construcción; un
# valor negativo indicaría un vector de valores propios incompleto.
greenacre_valido <- I_aj > 0
pct_greenacre <- if (greenacre_valido) 100 * lam_b / I_aj else NA_real_
t_kbl(tibble::tibble(
Eje = paste("Dim", seq_along(lambda)),
`Valor propio` = lambda,
`% bruto` = acm$eig[, 2],
`% Benzécri` = pct_benzecri,
`% Greenacre` = pct_greenacre) |>
filter(`Valor propio` > 1/Q_acm | dplyr::row_number() <= 4),
caption = "Inercia del análisis múltiple: porcentajes brutos y corregidos.",
digits = 4, scroll = TRUE)| Eje | Valor propio | % bruto | % Benzécri | % Greenacre |
|---|---|---|---|---|
| Dim 1 | 0.5888 | 5.888 | 49.611 | 38.230 |
| Dim 2 | 0.5392 | 5.392 | 32.203 | 24.816 |
| Dim 3 | 0.4691 | 4.691 | 14.012 | 10.797 |
| Dim 4 | 0.4074 | 4.074 | 4.174 | 3.217 |
| Dim 5 | 0.3333 | 3.333 | 0.000 | 0.000 |
| Dim 6 | 0.3333 | 3.333 | 0.000 | 0.000 |
| Dim 7 | 0.3333 | 3.333 | 0.000 | 0.000 |
| Dim 8 | 0.3333 | 3.333 | 0.000 | 0.000 |
| Dim 9 | 0.3333 | 3.333 | 0.000 | 0.000 |
| Dim 10 | 0.3333 | 3.333 | 0.000 | 0.000 |
| Dim 11 | 0.3333 | 3.333 | 0.000 | 0.000 |
| Dim 12 | 0.3333 | 3.333 | 0.000 | 0.000 |
# ============================================================
# FIGURA: MAPA FACTORIAL DEL ACM
# ============================================================
cat_act <- as.data.frame(acm$var$coord[, 1:2])
cat_act$Etiqueta <- rownames(cat_act)
cat_act$Contrib <- acm$var$contrib[, 1] + acm$var$contrib[, 2]
cat_act$Elemento <- "Activa"
cat_sup <- as.data.frame(acm$quali.sup$coord[, 1:2])
cat_sup$Etiqueta <- rownames(cat_sup)
cat_sup$Contrib <- NA_real_
cat_sup$Elemento <- "Suplementaria"
cat_sup$Etiqueta[cat_sup$Etiqueta %in% levels(datos_acm$cluster)] <-
paste("Grupo", cat_sup$Etiqueta[cat_sup$Etiqueta %in% levels(datos_acm$cluster)])
puntos_acm <- rbind(cat_act, cat_sup)
names(puntos_acm)[1:2] <- c("Dim1", "Dim2")
# Se etiquetan las categorías de tipo y estrato, las suplementarias y los ocho
# barrios de mayor contribución, para no saturar el mapa.
claves <- c(levels(datos_acm$tipo), levels(datos_acm$estrato_lab))
top_barrios <- cat_act$Etiqueta[order(cat_act$Contrib, decreasing = TRUE)]
top_barrios <- setdiff(top_barrios, claves)[1:8]
etiquetadas <- puntos_acm |>
filter(Etiqueta %in% c(claves, top_barrios) | Elemento == "Suplementaria")
subtitulo <- if (greenacre_valido) {
paste0("Inercia ajustada de Greenacre: ", round(pct_greenacre[1], 1),
" % y ", round(pct_greenacre[2], 1), " % en los dos primeros ejes")
} else {
paste0("Inercia corregida de Benzécri: ", round(pct_benzecri[1], 1),
" % y ", round(pct_benzecri[2], 1), " % en los dos primeros ejes")
}
ggplot(puntos_acm, aes(x = Dim1, y = Dim2, colour = Elemento)) +
geom_hline(yintercept = 0, colour = "#DEDEDE", linewidth = 0.4) +
geom_vline(xintercept = 0, colour = "#DEDEDE", linewidth = 0.4) +
geom_point(size = 2.6, alpha = 0.9) +
ggrepel::geom_text_repel(data = etiquetadas, aes(label = Etiqueta), size = 3.2,
max.overlaps = 30, colour = TINTA, show.legend = FALSE) +
scale_colour_manual(values = c(Activa = PAS_CIELO, Suplementaria = PAS_DURAZNO)) +
labs(title = "Estructura conjunta de tipo, estrato y barrio",
subtitle = subtitulo, x = "Dimensión 1", y = "Dimensión 2")
Figure 6.3: Mapa factorial del análisis de correspondencias múltiple. Se etiquetan tipo, estrato, los barrios de mayor contribución y las categorías suplementarias.
# ============================================================
# TABLA: COORDENADAS Y CONTRIBUCIONES EN EL ACM
# ============================================================
t_kbl(data.frame(
`Categoría` = rownames(acm$var$coord),
`Coord. Dim1` = acm$var$coord[, 1],
`Coord. Dim2` = acm$var$coord[, 2],
`Contrib. Dim1 (%)` = acm$var$contrib[, 1],
`Contrib. Dim2 (%)` = acm$var$contrib[, 2],
`cos² plano` = rowSums(acm$var$cos2[, 1:2]),
check.names = FALSE) |>
arrange(desc(`Contrib. Dim1 (%)`)),
caption = "Coordenadas y contribuciones de las categorías activas en el análisis múltiple.",
digits = 3, scroll = TRUE, height = "360px")| Categoría | Coord. Dim1 | Coord. Dim2 | Contrib. Dim1 (%) | Contrib. Dim2 (%) | cos² plano | |
|---|---|---|---|---|---|---|
| Estrato 6 | Estrato 6 | 1.660 | -0.094 | 37.266 | 0.130 | 0.868 |
| Pance | Pance | 1.996 | -0.297 | 11.173 | 0.270 | 0.212 |
| Ciudad Jardín | Ciudad Jardín | 1.750 | -0.348 | 10.800 | 0.466 | 0.211 |
| Estrato 4 | Estrato 4 | -0.689 | 0.420 | 6.872 | 2.791 | 0.224 |
| Santa Teresita | Santa Teresita | 1.928 | 0.590 | 6.653 | 0.680 | 0.133 |
| Resto de barrios | Resto de barrios | -0.400 | -0.786 | 3.686 | 15.527 | 0.533 |
| Normandía | Normandía | 1.854 | 0.642 | 3.603 | 0.471 | 0.073 |
| Estrato 3 | Estrato 3 | -0.566 | -1.617 | 3.168 | 28.236 | 0.621 |
| Valle Del Lili | Valle Del Lili | -0.655 | 1.180 | 2.947 | 10.437 | 0.251 |
| Estrato 5 | Estrato 5 | -0.367 | 0.597 | 2.524 | 7.282 | 0.243 |
| Los Cristales | Los Cristales | 1.416 | 0.682 | 2.100 | 0.532 | 0.047 |
| Parcelaciones Pance | Parcelaciones Pance | 2.048 | -0.706 | 1.741 | 0.226 | 0.035 |
| Cristales | Cristales | 1.392 | 0.653 | 1.095 | 0.263 | 0.024 |
| El Caney | El Caney | -0.827 | 0.558 | 0.973 | 0.483 | 0.026 |
| Aguacatal | Aguacatal | 0.928 | 0.135 | 0.639 | 0.015 | 0.012 |
| Ciudad 2000 | Ciudad 2000 | -0.971 | -0.405 | 0.616 | 0.117 | 0.013 |
| El Refugio | El Refugio | -0.787 | 0.661 | 0.506 | 0.390 | 0.015 |
| La Flora | La Flora | -0.447 | 1.120 | 0.501 | 3.433 | 0.067 |
| Prados Del Norte | Prados Del Norte | -0.696 | 0.972 | 0.419 | 0.892 | 0.022 |
| Caney | Caney | -0.770 | 0.771 | 0.355 | 0.389 | 0.013 |
| El Limonar | El Limonar | -0.605 | 0.398 | 0.336 | 0.159 | 0.009 |
| El Peñon | El Peñon | 0.905 | 0.974 | 0.334 | 0.423 | 0.013 |
| La Hacienda | La Hacienda | -0.529 | 0.945 | 0.317 | 1.102 | 0.024 |
| Brisas De Los | Brisas De Los | -0.715 | -2.325 | 0.285 | 3.294 | 0.059 |
| Casa | Casa | -0.110 | -0.681 | 0.267 | 11.102 | 0.301 |
| Nueva Tequendama | Nueva Tequendama | -0.628 | 0.644 | 0.196 | 0.225 | 0.007 |
| Apartamento | Apartamento | 0.070 | 0.430 | 0.168 | 7.007 | 0.301 |
| Versalles | Versalles | -0.548 | 1.185 | 0.145 | 0.741 | 0.015 |
| Urbanización La Flora | Urbanización La Flora | -0.469 | 1.154 | 0.124 | 0.821 | 0.016 |
| Quintas De Don | Quintas De Don | -0.441 | 1.250 | 0.096 | 0.848 | 0.016 |
| Acopi | Acopi | -0.210 | 0.106 | 0.047 | 0.013 | 0.001 |
| Santa Isabel | Santa Isabel | -0.327 | 0.812 | 0.046 | 0.314 | 0.006 |
| El Ingenio | El Ingenio | -0.012 | 0.782 | 0.000 | 0.922 | 0.015 |
# ============================================================
# TABLA: VALORES-TEST
# ============================================================
# El valor-test contrasta si la coordenada de una categoría difiere de cero:
# bajo independencia se distribuye aproximadamente N(0,1), de modo que
# |valor-test| > 1.96 señala categorías significativamente desplazadas.
t_kbl(data.frame(
`Categoría` = rownames(acm$var$v.test),
`v.test Dim1` = acm$var$v.test[, 1],
`v.test Dim2` = acm$var$v.test[, 2],
check.names = FALSE) |>
mutate(`Significativa en Dim1` = ifelse(abs(`v.test Dim1`) > 1.96, "Sí", "No")) |>
arrange(desc(abs(`v.test Dim1`))),
caption = "Valores-test de las categorías activas sobre los dos primeros ejes.",
digits = 2, scroll = TRUE, height = "340px")| Categoría | v.test Dim1 | v.test Dim2 | Significativa en Dim1 | |
|---|---|---|---|---|
| Estrato 6 | Estrato 6 | 84.82 | -4.79 | Sí |
| Pance | Pance | 41.56 | -6.18 | Sí |
| Ciudad Jardín | Ciudad Jardín | 41.14 | -8.17 | Sí |
| Estrato 4 | Estrato 4 | -36.84 | 22.46 | Sí |
| Santa Teresita | Santa Teresita | 31.77 | 9.72 | Sí |
| Resto de barrios | Resto de barrios | -30.20 | -59.32 | Sí |
| Estrato 3 | Estrato 3 | -23.75 | -67.84 | Sí |
| Estrato 5 | Estrato 5 | -23.54 | 38.26 | Sí |
| Normandía | Normandía | 23.22 | 8.04 | Sí |
| Valle Del Lili | Valle Del Lili | -22.20 | 39.98 | Sí |
| Los Cristales | Los Cristales | 17.73 | 8.54 | Sí |
| Parcelaciones Pance | Parcelaciones Pance | 16.05 | -5.53 | Sí |
| Cristales | Cristales | 12.75 | 5.98 | Sí |
| El Caney | El Caney | -12.11 | 8.17 | Sí |
| Aguacatal | Aguacatal | 9.76 | 1.41 | Sí |
| Ciudad 2000 | Ciudad 2000 | -9.57 | -3.99 | Sí |
| La Flora | La Flora | -8.78 | 21.98 | Sí |
| El Refugio | El Refugio | -8.69 | 7.29 | Sí |
| Casa | Casa | -8.00 | -49.36 | Sí |
| Apartamento | Apartamento | 8.00 | 49.36 | Sí |
| Prados Del Norte | Prados Del Norte | -7.91 | 11.04 | Sí |
| Caney | Caney | -7.26 | 7.27 | Sí |
| El Limonar | El Limonar | -7.08 | 4.66 | Sí |
| El Peñon | El Peñon | 7.03 | 7.57 | Sí |
| La Hacienda | La Hacienda | -6.89 | 12.30 | Sí |
| Brisas De Los | Brisas De Los | -6.50 | -21.16 | Sí |
| Nueva Tequendama | Nueva Tequendama | -5.39 | 5.53 | Sí |
| Versalles | Versalles | -4.63 | 10.03 | Sí |
| Urbanización La Flora | Urbanización La Flora | -4.29 | 10.57 | Sí |
| Quintas De Don | Quintas De Don | -3.78 | 10.73 | Sí |
| Acopi | Acopi | -2.67 | 1.34 | Sí |
| Santa Isabel | Santa Isabel | -2.62 | 6.52 | Sí |
| El Ingenio | El Ingenio | -0.17 | 11.27 | No |
6.6 Lectura del apartado
El primer cruce ordena los estratos según la composición tipológica de su oferta. La tabla de contribuciones indica cuáles estratos sostienen el eje: una contribución baja significa que el perfil de ese estrato coincide con el marginal, es decir, que su mezcla de casas y apartamentos reproduce la del conjunto. La asociación puede ser estadísticamente clara y proceder, aun así, de unas pocas categorías.
El cruce barrio por estrato es el que aporta contenido comercial. Los dos primeros ejes concentran la mayor parte de la inercia, de modo que el mapa aproxima bien la tabla, y los barrios se ordenan según su composición socioeconómica. Los barrios situados junto al origen tienen perfil próximo al promedio de la ciudad y no admiten interpretación por posición; los alejados concentran un estrato determinado y son los que definen la estructura.
En el análisis múltiple interesan tres proyecciones suplementarias. La zona
permite comprobar si la estructura por barrios reproduce la división
administrativa o la corta. Los grupos de la segmentación indican, por su
distancia al origen, en qué medida el producto y el territorio son dimensiones
distintas del mercado: coordenadas próximas al origen significan que ambos
segmentos tienen perfiles categóricos parecidos y que la segmentación no es
geográfica. Y el indicador piso_declarado responde a una pregunta operativa:
si se aleja del origen, la calidad de diligenciamiento del anuncio se concentra
en cierto tipo de inmueble o de barrio, lo que convierte un defecto de captura
en un patrón identificable y corregible.
# ============================================================
# TABLA: SÍNTESIS DEL APARTADO
# ============================================================
t_kbl(tibble::tribble(
~Elemento, ~Resultado,
"`barrio`", "Cuatro defectos concurrentes; se agrupa por cobertura del 60 % y se usa con esa reserva",
"Tipo por estrato", "Una sola dimensión no trivial; ordena los estratos por composición tipológica",
"Barrio por estrato", "Dos ejes concentran la mayor parte de la inercia y ordenan los barrios socioeconómicamente",
"Inercia del ACM", "Los porcentajes brutos subestiman la estructura; se reportan Benzécri y Greenacre",
"Suplementarias", "Zona, grupo de segmentación e indicador de piso declarado proyectados sobre el mapa"),
caption = "Síntesis del análisis de correspondencias.",
align = c("l","l"), scroll = TRUE)| Elemento | Resultado |
|---|---|
barrio
|
Cuatro defectos concurrentes; se agrupa por cobertura del 60 % y se usa con esa reserva |
| Tipo por estrato | Una sola dimensión no trivial; ordena los estratos por composición tipológica |
| Barrio por estrato | Dos ejes concentran la mayor parte de la inercia y ordenan los barrios socioeconómicamente |
| Inercia del ACM | Los porcentajes brutos subestiman la estructura; se reportan Benzécri y Greenacre |
| Suplementarias | Zona, grupo de segmentación e indicador de piso declarado proyectados sobre el mapa |
7 Discusión y cierre
7.1 Integración de los tres análisis
# ============================================================
# TABLA: INTEGRACIÓN DE RESULTADOS
# ============================================================
t_kbl(tibble::tribble(
~Pregunta, ~Herramienta, ~Respuesta,
"¿Qué factores ordenan la oferta?", "Componentes principales",
"Un eje principal construido sobre valor unitario y dimensión, validado externamente por el estrato",
"¿La oferta admite una división en segmentos?", "k-medoides",
"Una división operativa de un mercado continuo, reproducible por remuestreo",
"¿La división coincide con el territorio?", "Segmentación y correspondencias",
"Se responde con el índice de Rand ajustado y con la posición de los grupos en el mapa factorial",
"¿Cómo se asocian tipo, estrato y barrio?", "Correspondencias simple y múltiple",
"Una ordenación socioeconómica de los barrios, resumida en dos ejes"),
caption = "Integración de los resultados de las tres técnicas.",
align = c("l","l","l"), scroll = TRUE)| Pregunta | Herramienta | Respuesta |
|---|---|---|
| ¿Qué factores ordenan la oferta? | Componentes principales | Un eje principal construido sobre valor unitario y dimensión, validado externamente por el estrato |
| ¿La oferta admite una división en segmentos? | k-medoides | Una división operativa de un mercado continuo, reproducible por remuestreo |
| ¿La división coincide con el territorio? | Segmentación y correspondencias | Se responde con el índice de Rand ajustado y con la posición de los grupos en el mapa factorial |
| ¿Cómo se asocian tipo, estrato y barrio? | Correspondencias simple y múltiple | Una ordenación socioeconómica de los barrios, resumida en dos ejes |
7.2 Alcance de los resultados
Conviene fijar hasta dónde llegan las afirmaciones anteriores. El registro describe ofertas publicadas y no transacciones cerradas, de modo que el precio analizado es el pedido por el vendedor e incorpora expectativas y margen de negociación; todo lo dicho se refiere a la estructura de la oferta. La captura desde un solo portal no constituye muestreo probabilístico, así que los contrastes funcionan como reglas de decisión sobre el conjunto observado y no como base de generalización. La ausencia de los estratos 1 y 2 restringe el alcance al segmento comprendido entre los estratos 3 y 6. Y la división en grupos es un corte reproducible de un continuo, no el descubrimiento de subpoblaciones: los tamaños de efecto que la describen se calculan sobre las mismas variables que la generaron, de modo que ordenan la contribución relativa de cada atributo sin probar nada sobre la estructura del mercado.
7.3 Conclusiones
El análisis de componentes principales muestra que la variación de la oferta se resume en muy pocas direcciones. Al construir el espacio sobre el precio por metro cuadrado en lugar del precio total, el primer eje ya no puede ser un factor de tamaño puro: la Tabla 4.5 muestra cómo se reparten el valor unitario y las variables de dimensión, que es la distinción relevante para valorar. El estrato socioeconómico, que no participó en la construcción del espacio, reproduce el ordenamiento de uno de los ejes: esa coincidencia es una validación externa, porque muestra que la estructura recuperada con atributos físicos y precio coincide con una clasificación administrativa independiente.
La segmentación por \(k\)-medoides produce una división reproducible. Los cuatro diagnósticos aplicados deciden si además está naturalmente separada: los bloques en la matriz de disimilaridad ordenada, el valle de densidad entre grupos, el ancho de silueta y la concordancia entre criterios de agregación. La conclusión no es que la segmentación carezca de valor, sino que debe leerse como lo que es: una estratificación operativa de un mercado continuo. Los inmuebles próximos a la frontera admiten cualquiera de las dos clasificaciones, y decidir en qué segmento se ofertan es una decisión comercial con consecuencias de precio y de público, no una lectura del dato. La elección del medoide como representante refuerza esta interpretación, porque cada grupo queda descrito por un inmueble que existe en el registro y no por un promedio que puede no corresponder a ningún producto real.
El análisis de correspondencias aporta la dimensión territorial que las dos técnicas anteriores no capturan. Los barrios se ordenan según su composición socioeconómica, con unos pocos definiendo los ejes y una mayoría situada cerca del origen por tener perfil próximo al promedio de la ciudad. La proyección de los grupos de segmentación como categorías suplementarias permite responder directamente a la pregunta de si producto y territorio son la misma cosa: cuanto más cerca del origen queden, más independiente resulta la segmentación por gama respecto de la geografía, con la implicación operativa de que una organización comercial por zonas no captura la estructura del producto.
Tres resultados del diagnóstico de datos tienen valor propio para la empresa. La ausencia en el campo de parqueaderos no es aleatoria y se concentra en los inmuebles de menor gama, de modo que el formulario del portal está produciendo información sesgada; registrar explícitamente el valor cero resolvería el problema sin coste. El campo de piso carece de significado unívoco, pero su declaración o su omisión constituyen un atributo del anuncio que sí admite análisis. Y la variable de barrio presenta cuatro defectos concurrentes que impiden usarla como variable activa, cuando es precisamente la unidad con la que opera el equipo comercial: adoptar un catálogo cerrado de barrios con validación de pertenencia a zona es la corrección con mayor retorno analítico de todas las posibles.
Por último, las dos decisiones metodológicas discutibles se sometieron a verificación explícita. La solución factorial se comparó con su versión robusta, para comprobar que no está determinada por el segmento de alto valor que se optó por conservar, y con el escenario alternativo de tratamiento de la ausencia. El coeficiente de congruencia de Tucker permite decidir en ambos casos si las conclusiones dependen de la decisión adoptada, que es la condición para presentarlas como resultados y no como consecuencias de una elección del analista.
8 Referencias
Bartlett, M. S. (1951). The effect of standardization on a \(\chi^2\) approximation in factor analysis. Biometrika, 38(3/4), 337–344.
Box, G. E. P., & Cox, D. R. (1964). An analysis of transformations. Journal of the Royal Statistical Society: Series B, 26(2), 211–252.
Greenacre, M. (2017). Correspondence Analysis in Practice (3.ª ed.). CRC Press.
Hennig, C. (2007). Cluster-wise assessment of cluster stability. Computational Statistics & Data Analysis, 52(1), 258–271.
Hubert, M., Rousseeuw, P. J., & Vanden Branden, K. (2005). ROBPCA: A new approach to robust principal component analysis. Technometrics, 47(1), 64–79.
Husson, F., Lê, S., & Pagès, J. (2017). Exploratory Multivariate Analysis by Example Using R (2.ª ed.). CRC Press.
Iglewicz, B., & Hoaglin, D. C. (1993). How to Detect and Handle Outliers. ASQC Quality Press.
Jackson, D. A. (1993). Stopping rules in principal components analysis. Ecology, 74(8), 2204–2214.
Johnson, R. A., & Wichern, D. W. (2007). Applied Multivariate Statistical Analysis (6.ª ed.). Pearson.
Josse, J., & Husson, F. (2016). missMDA: A package for handling missing values in multivariate data analysis. Journal of Statistical Software, 70(1), 1–31.
Kaiser, H. F. (1974). An index of factorial simplicity. Psychometrika, 39(1), 31–36.
Kaufman, L., & Rousseeuw, P. J. (2005). Finding Groups in Data: An Introduction to Cluster Analysis. Wiley. [Cap. 2: PAM; cap. 3: CLARA].
Lay, D. C. (2012). Álgebra lineal y sus aplicaciones (4.ª ed.). Pearson. [Cap. 7: §7.1 diagonalización de matrices simétricas; §7.4 descomposición en valores singulares; §7.5 componentes principales].
Lebart, L., Morineau, A., & Piron, M. (2006). Statistique exploratoire multidimensionnelle (4.ª ed.). Dunod.
Little, R. J. A. (1988). A test of missing completely at random for multivariate data with missing values. Journal of the American Statistical Association, 83(404), 1198–1202.
Lorenzo-Seva, U., & ten Berge, J. M. F. (2006). Tucker’s congruence coefficient as a meaningful index of factor similarity. Methodology, 2(2), 57–64.
Peña, D. (2002). Análisis de datos multivariantes. McGraw-Hill.
Rousseeuw, P. J. (1987). Silhouettes: A graphical aid to the interpretation and validation of cluster analysis. Journal of Computational and Applied Mathematics, 20, 53–65.
Rousseeuw, P. J., & Van Driessen, K. (1999). A fast algorithm for the minimum covariance determinant estimator. Technometrics, 41(3), 212–223.
Rubin, D. B. (1976). Inference and missing data. Biometrika, 63(3), 581–592.
Tibshirani, R., Walther, G., & Hastie, T. (2001). Estimating the number of clusters in a data set via the gap statistic. Journal of the Royal Statistical Society: Series B, 63(2), 411–423.
Tomczak, M., & Tomczak, E. (2014). The need to report effect size estimates revisited. Trends in Sport Sciences, 1(21), 19–25.
Tukey, J. W. (1977). Exploratory Data Analysis. Addison-Wesley.
9 Anexos
9.1 A1. Preparación del entorno
# ── Paquete de datos del curso (una sola vez) ───────────────────────────────
# install.packages("devtools")
# devtools::install_github("centromagis/paqueteMODELOS", force = TRUE)
paquetes <- c("dplyr", "tidyr", "stringr", "ggplot2", "ggrepel", "knitr",
"kableExtra", "e1071", "nortest", "naniar", "missMDA",
"robustbase", "FSA", "scales", "FactoMineR", "cluster", "fpc",
"tibble", "psych", "rrcov", "leaflet")
faltantes <- setdiff(paquetes, rownames(installed.packages()))
if (length(faltantes) > 0) install.packages(faltantes)9.2 A2. Secuencia de tratamiento de los datos
# ── 1. Carga, normalización de cadenas y construcción de `piso_declarado` ────
# ── 2. Depuración: duplicados, filas sin identificación y valores fuera de
# dominio convertidos en NA ──────────────────────────────────────────────
# ── 3. Mecanismo de ausencia: prueba de Little y chi cuadrado por variable
# frente al estrato, con V de Cramér como tamaño del efecto ─────────────
# ── 4. Imputación conjunta por ACP regularizado, incluida `parqueaderos`;
# escenario alternativo con la ausencia recodificada a cero ─────────────
# ── 5. Precio por metro cuadrado y logaritmos ───────────────────────────────
# ── 6. Valores extremos: Tukey, z robusto y MCD. Se conservan ───────────────
# ── 7. ACP sobre precio unitario, área y dotación ───────────────────────────
# ── 8. Segmentación por k-medoides y validación ─────────────────────────────
# ── 9. Correspondencias sobre tipo, estrato y barrio ────────────────────────9.3 A3. Registro de decisiones
# ============================================================
# TABLA: DECISIONES Y ALTERNATIVAS CONSIDERADAS
# ============================================================
t_kbl(tibble::tribble(
~`Decisión`, ~`Alternativa considerada`, ~Motivo,
"Precio por metro cuadrado como variable activa", "Precio total junto al área",
"Incluir ambos cuenta el tamaño dos veces y fuerza un primer eje de tamaño",
"Imputación de `parqueaderos`", "Recodificación de la ausencia a cero",
"La lectura MNAR no es contrastable; se prefiere no inventar ceros y se comprueba la sensibilidad",
"`piso_declarado` como indicador", "Exclusión de la variable",
"La omisión del campo es un atributo del anuncio y admite análisis",
"k-medoides con CLARA", "Ward consolidado o k-medias",
"El medoide resiste la cola de alto valor que se decidió conservar; el centroide no",
"Agrupación de `barrio` por cobertura", "Umbral de frecuencia mínima",
"Fija de antemano qué proporción del mercado queda representada nominalmente",
"Correcciones de Benzécri y Greenacre", "Porcentajes brutos de inercia",
"La inercia total del ACM depende del diseño de codificación, no de la asociación"),
caption = "Decisiones metodológicas y alternativas consideradas.",
align = c("l","l","l"), scroll = TRUE)| Decisión | Alternativa considerada | Motivo |
|---|---|---|
| Precio por metro cuadrado como variable activa | Precio total junto al área | Incluir ambos cuenta el tamaño dos veces y fuerza un primer eje de tamaño |
Imputación de parqueaderos
|
Recodificación de la ausencia a cero | La lectura MNAR no es contrastable; se prefiere no inventar ceros y se comprueba la sensibilidad |
piso_declarado como indicador
|
Exclusión de la variable | La omisión del campo es un atributo del anuncio y admite análisis |
| k-medoides con CLARA | Ward consolidado o k-medias | El medoide resiste la cola de alto valor que se decidió conservar; el centroide no |
Agrupación de barrio por cobertura
|
Umbral de frecuencia mínima | Fija de antemano qué proporción del mercado queda representada nominalmente |
| Correcciones de Benzécri y Greenacre | Porcentajes brutos de inercia | La inercia total del ACM depende del diseño de codificación, no de la asociación |
9.4 A4. Entorno de cómputo
# ============================================================
# TABLA: ENTORNO DE CÓMPUTO
# ============================================================
sesion <- sessionInfo()
t_kbl(tibble::tibble(
Elemento = c("Versión de R", "Plataforma", "Sistema operativo", "Semilla global"),
Valor = c(sesion$R.version$version.string, sesion$platform,
sesion$running, "2026")),
caption = "Entorno de cómputo empleado en la compilación del informe.",
align = c("l","l"))| Elemento | Valor |
|---|---|
| Versión de R | R version 4.6.1 (2026-06-24 ucrt) |
| Plataforma | x86_64-w64-mingw32/x64 |
| Sistema operativo | Windows 11 x64 (build 26200) |
| Semilla global | 2026 |
Oriana Giraldo Arcia