Análisis de Datos Categóricos (SOC3070)

class: center, middle, inverse, title-slide

.title[
# Análisis de Datos Categóricos (SOC3070)
]
.subtitle[
## Regresión Logística para Predicción
]
.author[
### Mauricio Bucca Profesor Asistente, Sociología UC
]

---

class: inverse, center, middle

# Evaluación de Modelos Predictivos

---

## `$\hat{\beta}$`-problems vs. `$\hat{y}$`-problems

(Mullainathan & Spiess, 2017)

.pull-left[
### `$\hat{\beta}$`-problems

- Foco en **parámetros** y su interpretación.

- Preguntas:  
  - ¿Cuál es el efecto de `$X$` sobre `$Y$`?  
  - ¿Es significativo / causal?  
  
- Objetivo: **explicación e inferencia**.

- Modelos suelen ser **transparentes** (fácil interpretar).  
]

.pull-right[
### `$\hat{y}$`-problems

- Foco en **predicciones** para nuevos casos.

- Preguntas:  
  - ¿Qué tan bien predice el modelo en datos futuros?  
  - ¿Qué modelo predice mejor?  
  
- Objetivo: **desempeño y utilidad práctica**.

- Modelos pueden ser tratados como **black-box** (no importa cómo, mientras prediga bien).  
]

---
class: inverse, center, middle

## Regresión Logística como Machine Learning

---
## Regresión Logística como Machine Learning

.middle[
.center[
![](https://www.researchgate.net/publication/376856079/figure/fig1/AS:11431281226110814@1709092685557/Famous-internet-meme-displaying-the-relation-between-statistics-Machine-Learning-and.png)
]
]
---
# Regresión Logística: el **“hello world”** de ML

- **Qué resuelve:**

- Estima la probabilidad de un evento binario (`$Y=0/1$`).  
  - Ejemplos: ¿un mail es spam? ¿alguien tendrá una aventura?

- **Cómo funciona:**

- Ajusta una combinación lineal de predictores `$X\beta$`.  
  - Pasa ese valor por la **curva sigmoide** → convierte log-odds en probabilidades (0–1).  
  - Aprende de los datos para **generalizar** a casos no observados (nuevos datos, futuro, etc.).

- **Por qué es importante:**

- Interpretable: cada predictor tiene un rol claro en la probabilidad.  
  - Eficiente: funciona bien incluso con muestras pequeñas.  
  - Fundacional: es la base sobre la cual se construyen modelos más complejos en ML.

---

## Ejemplo empírico

`$$\newcommand{\vect}[1]{\boldsymbol{#1}}$$`

Continuando con el ejemplo de clases anteriores, ajustamos el siguiente modelo:

`$$\ln \frac{p_{i}}{1-p_{i}} = \beta_{0} + \beta_{1}\text{ym}_{i} + \beta_{2}\text{male}_{i} + \beta_{3}\text{rate}_{i} + \beta_{4}\text{rate}^{2}_{i}$$`

Llamemos a este modelo, modelo `$M$`:

```
##                     Estimate Std. Error    z value   Pr(>|z|)
## (Intercept)      0.149861306 0.82300610  0.1820902 0.85551198
## ym               0.037699456 0.01813912  2.0783503 0.03767711
## factor(sex)male  0.249972343 0.19690565  1.2695032 0.20426166
## rate            -0.411820264 0.50252200 -0.8195069 0.41249725
## I(rate^2)       -0.008931688 0.07364745 -0.1212763 0.90347221
```

```
## [1] "log-likelihood: -316.105 Deviance: 632.21"
```

```
## [1] "AIC: 642.21 BIC: 664.202"
```

---
## Likelihood como función de pérdida

- En ML, entrenar un modelo = encontrar parámetros que **maximicen la verosimilitud**  
  (o equivalentemente **minimicen la pérdida**).

- En regresión logística, cada `$y_i$` es Bernoulli con probabilidad `$p_i$`:

$$
P(y_i \mid p_i) = p_i^{y_i}(1-p_i)^{1-y_i}
$$

- La probabilidad conjunta de todos los datos:

$$
\mathcal{L}(\beta) = \prod_{i=1}^{n} p_i^{y_i}(1-p_i)^{1-y_i}
$$

- Tomando logaritmos:

$$
\ell(\beta) = \sum_{i=1}^{n} \big[ y_i \ln p_i + (1-y_i) \ln(1-p_i) \big]
$$

- En ML a esto le llamamos **log-loss** o **cross-entropy loss**:

$$
\text{LogLoss} = -\ell(\beta)
$$

---
## Likelihood como función de pérdida

- Si las predicciones del modelo son correctas (`$y=1$` y `$p$` cercano a 1) → log-loss bajo.
- Si se equivocan (`$y=1$` pero `$p$` cercano a 0) → log-loss alto.
- Entrenar = **ajustar `$\beta$` para minimizar log-loss** (= maximizar likelihood).

Likelihood maximizada en nuestro ejemplo:

``` r
ll_m <- logLik(logit_affairs); print(c(Likelihood = exp(ll_m[1]), log_likelihood = ll_m[1]))
```

```
##     Likelihood log_likelihood 
##  5.217357e-138  -3.161048e+02
```

.bold[Versión ML: Log-loss promedio ]

En machine learning se suele trabajar con la **log-loss promedio** para obtener una pérdida por observación (comparable entre datasets de distinto tamaño).

$$
\text{LogLoss} = -\frac{1}{n}\,\ell(\beta)
$$

``` r
n <- length(affairsdata$everaffair_d); log_loss <- -1/n * as.numeric(ll_m[1]); log_loss
```

```
## [1] 0.5259646
```

---

## Modelos de referencia (baselines)

En ML y estadística es común usar **modelos de referencia** para poner en contexto el desempeño de un modelo más complejo.  Dos extremos:

--

- **Modelo nulo (`$M_N$`)**  
  - El más simple posible: siempre predice la probabilidad promedio global.  
  - Ventaja: parsimonioso y fácil de interpretar.  
  - Desventaja: genera las peores predicciones posibleså, ignora covariables.

``` r
# modelo nulo
logit_affairs_sex_null <- glm(everaffair_d ~ 1,family=binomial(link="logit"), data=affairsdata)
```
--

- **Modelo saturado (`$M_S$`)**  
  - El más complejo posible: ajusta un parámetro distinto para cada observación.  
  - Ventaja: fit perfecto (cero perdida).  
  - Desventaja: no generaliza → memoriza los datos.

``` r
# modelo saturado
logit_affairs_sex_sat <- glm(everaffair_d ~ factor(1:nrow(affairsdata)) ,family=binomial(link="logit"), data=affairsdata)
```

---

## Baselines vs Modelo de Interés

.center[
![dev](models.png)
]

---
## Residual Deviance y Null Deviance

En modelos logísticos, el **log-likelihood ratio** puede re-expresarse como **Deviance**, que es esencialmente una medida de **pérdida (loss)**. Dos tipos de deviance:

- **Residual Deviance:**  `$D = -2 \cdot (\log\mathcal{L}_{M} - \log \mathcal{L}_{S})$`  
  - Evalúa el ajuste de `$M$` respecto al modelo saturado (fit perfecto).

- **Null Deviance:** `$D_N = -2 \cdot (\log\mathcal{L}_{0} - \log \mathcal{L}_{S})$`  
  - Equivale al "total explicable", similar a la varianza total en OLS.

--

- **Distribución muestral:** `$D \sim \chi^{2}(df = n-k), \quad  \text{donde k es el número de parámetros en M}$`

Interpretación:

- `$D$` alto (p-value bajo) → "mal ajuste" (alto error, modelo se queda corto).  
- `$D$` bajo (p-value alto) → "buen ajuste" (más parámetros no agregan valor.

---

## Residual Deviance y Null Deviance

.center[
![dev](deviances.png)
]

---
## Residual Deviance y Log-loss promedio

```
## [1] "Null Deviance: 675.377 | Residual Deviance: 632.21"
```

``` r
# Cálculo manual de la Residual Deviance
D <- -2 * (logLik(logit_affairs)[1] - logLik(logit_affairs_sex_sat)[1])

# Relación con log-loss promedio
log_likelihood <- as.numeric(logLik(logit_affairs))
log_loss <- -1/length(affairsdata$everaffair_d) * log_likelihood
print(paste0("Log-loss promedio: ", round(log_loss,4)))
```

```
## [1] "Log-loss promedio: 0.526"
```

---
## Pseudo - `$R^2$`

--

- En modelos OLS es común medir ajuste usando el coeficiente `$R^2$`, es decir, % de varianza explicada por el modelo.

- En GLM's la varianza no es separable de la media, por tanto no se puede descomponer.

- Existe una variedad de alternativas al `$R^2$`, llamadas genéricamente pseudo - `$R^2$`. Uno de los más comunes es:

--

`$$\text{McFadden’s } R^{2} = 1 - \frac{D}{D_{0}} = 1 - \frac{(\log\mathcal{L}_{S} - \log \mathcal{L}_{M})}{ (\log\mathcal{L}_{S} - \log \mathcal{L}_{N})}$$`
 
--

.bold[Intuición:]  fracción del total del "explicable" del likelihood que es explicado por el modelo `$M$`.

- `$R^{2} \in [0,1]$`, donde 0 indica el peor fit posible y 1 indica el mejor fit posible.

---
## Pseudo - `$R^2$`

--

.bold[Residual deviance]:

``` r
R2 <- 1 - logit_affairs$deviance/logit_affairs$null.deviance; R2
```

```
## [1] 0.06391612
```

``` r
# versión automática
PseudoR2(logit_affairs, which="McFadden")
```

```
##   McFadden 
## 0.06391612
```

---
class: inverse, center, middle

# Regresión Logística como clasificador
---

## De Probabilidades a clases

- La regresión logística entrega **scores probabilísticos**:   `$\hat{p}_i = \Pr(Y_i=1\mid X_i,\hat\beta)$`

--

- Para convertir un score en decisión usamos un **umbral de corte** `$\tau$`:

$$
\hat{y}_i =
`\begin{cases}
1 & \text{si } \hat{p}_i > \tau \\
0 & \text{si } \hat{p}_i \leq \tau
\end{cases}`
$$

--

- Por defecto: `$\tau=0.5$`

- En ML, el umbral se ajusta al problema:  
  - Spam: es peor que se cuele spam → corte bajo.  
  - Subsidio: es peor dejar fuera a quien lo necesita → corte bajo.  
  - Examen de enfermedad peligrosa:
    - Detección temprana: es peor no encontrar la enfermedad → corte bajo.  
    - Confirmación: es peor asustar a alguien sano → corte alto.

.bold[👉 Un clasificador es un modelo + una regla de decisión.]

---

## Matriz de Confusión

.center[
![](https://miro.medium.com/v2/resize:fit:4800/format:webp/1*mL-nYY6MFhiG0uoR5kJaCA.jpeg)
]
---

## Matriz de Confusión

.pull-left[
![](conf_mat.png)
]

.pull-right[
- **Accuracy**: `$(TP+TN)/N$`  
  % de clasificaciones correctas

- **Misclassification Rate**: `$(FP+FN)/N$`  
  % de clasificaciones incorrectas

- **True Positive Rate (Recall / Sensitivity)**:  
  `$TP/(TP+FN)$`

- **True Negative Rate (Specificity)**:  
  `$TN/(TN+FP)$`

- **Precision (PPV)**:  
  `$TP/(TP+FP)$`

- **Prevalence**:  
  `$(TP+FN)/N$`  
]

--

📌 Estas métricas capturan distintos aspectos del desempeño de un clasificador.  
En ML es clave elegir métricas alineadas con el problema (ej: medicina ≠ marketing).

---
## Ejemplo empírico: infidelidad

Clasificamos como “Nunca infiel” a todas las personas cuya probabilidad predicha sea menor o igual a 0.5, y como “Al menos una vez” a aquellas cuya probabilidad predicha supere dicho umbral.

$$
\hat{y}_i =
`\begin{cases}
1 & \text{si } \hat{p}_i > 0.5 \quad (\text{Al menos una vez}) \\\\
0 & \text{si } \hat{p}_i \leq 0.5 \quad (\text{Nunca infiel})
\end{cases}`
$$

``` r
p_hat <- predict(logit_affairs, type = "response")
y_hat <- if_else(p_hat>0.5,1,0)

conf_mat <- confusionMatrix(
 factor(y_hat), factor(logit_affairs$model$everaffair_d),
 positive = "1", dnn = c("Predicho","Real")
)

conf_mat$table
```

```
##         Real
## Predicho   0   1
##        0 438 137
##        1  13  13
```

---
## Ejemplo empírico: infidelidad

.pull-left[
![](class_12a_files/figure-html/unnamed-chunk-12-1.png)

]

.pull-right[

Table: Métricas principales del clasificador logístico

|                  |     x|
|:-----------------|-----:|
|Sensitivity       | 0.087|
|Specificity       | 0.971|
|Precision         | 0.500|
|F1                | 0.148|
|Balanced Accuracy | 0.529|
]

---
## Umbral y predicción

- Cada umbral `$\tau$` genera un nivel de Especificidad y Sensibilidad: `$:\{\text{Specificity}(\tau), \text{Sensitivity}(\tau)\}$` para `$\tau \in [0,1]$`.

.pull-left[
![](class_12a_files/figure-html/plot_scores-1.png)
]

.pull-right[
- Umbrales bajos → más predicciones positivas
- Umbrales altos → más predicciones negativas
]

---
## Curva ROC

- La **curva ROC** conecta todos los puntos `$:\{1 - \text{Specificity}(\tau), \quad \text{Sensitivity}(\tau)\}$` para `$\tau \in [0,1]$`.

- La diagonal corresponde a un clasificador aleatorio.

- Mejor desempeño = curva más arriba a la izquierda.

.pull-left[
![](class_12a_files/figure-html/plot_scores2-1.png)
]

.pull-right[
![](class_12a_files/figure-html/plot_roc-1.png)
]

---
## AUC: Área bajo la curva

- El **AUC** resume el desempeño en un solo número: `$AUC = \int_{0}^{1}\text{Sensibilidad}(\tau)\, d(1-\text{Especificidad}(\tau))$`

- **Interpretación**: probabilidad de que el modelo asigne mayor score a un positivo que a un negativo (elegidos al azar).

.pull-left[
![](class_12a_files/figure-html/plot_auc-1.png)
]

.pull-right[
**Interpretación del AUC:**

- **AUC ≈ 0.5**: Sin poder discriminatorio (como azar)

- **AUC ≈ 0.7-0.8**: Discriminación aceptable

- **AUC > 0.9**: Discriminación excelente

- **AUC = 1**: Clasificador perfecto

👉 Nuestro modelo tiene AUC = 0.674, indicando poder discriminatorio moderado.
]

---
## Selección del umbral óptimo

- No existe un único "mejor" umbral: depende del contexto y los costos de errar.

- Criterios comunes: **Youden**, **punto más cercano a (0,1)**, etc.

Cómo extraer umbrales en R:

``` r
# Índice de Youden (maximiza TPR - FPR)
coords(roc_obj, x = "best", 
       best.method = "youden")
```

```
##   threshold specificity sensitivity
## 1  0.222457   0.5698448   0.6933333
```

``` r
# Punto más cercano a (0,1)
coords(roc_obj, x = "best", 
       best.method = "closest.topleft")
```

```
##   threshold specificity sensitivity
## 1  0.222457   0.5698448   0.6933333
```

---
## Selección del umbral óptimo

- No existe un único "mejor" umbral: depende del contexto y los costos de errar.

- Criterios comunes: **Youden**, **punto más cercano a (0,1)**, etc.

Cómo extraer umbrales en R:

``` r
# Ver todas las coordenadas
coords(roc_obj, x = "all") %>% head(n=13)
```

```
##    threshold specificity sensitivity
## 1       -Inf  0.00000000   1.0000000
## 2  0.1069478  0.00886918   1.0000000
## 3  0.1080786  0.01773836   1.0000000
## 4  0.1100684  0.03991131   0.9866667
## 5  0.1162936  0.11308204   0.9733333
## 6  0.1268852  0.14412417   0.9400000
## 7  0.1331680  0.14855876   0.9333333
## 8  0.1345351  0.18181818   0.9200000
## 9  0.1370424  0.19068736   0.9200000
## 10 0.1430321  0.23059867   0.9133333
## 11 0.1488562  0.26164080   0.8933333
## 12 0.1574665  0.31707317   0.8666667
## 13 0.1649799  0.32150776   0.8666667
```

---
class: middle

## en el próximo episodio ...

.center[
![cross](https://miro.medium.com/max/1224/1*hPDx1TYs9oMrbZRf-VYSUw.jpeg)
]

---
class: inverse, center, middle

##Hasta la próxima clase. Gracias!

Mauricio Bucca 
https://mebucca.github.io/ 
github.com/mebucca