Regresión Lineal Simple¶

La regresión lineal no solo ajusta una línea a los datos: también permite hacer inferencia sobre los coeficientes, evaluar si la relación es estadísticamente significativa, y cuantificar la incertidumbre de las predicciones. Este notebook conecta los conceptos de pruebas de hipótesis e intervalos de confianza con el modelo de regresión.

Librerías¶

In [ ]:

Copied!





import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import scipy.stats as stats
import statsmodels.formula.api as smf
import statsmodels.api as sm

plt.rcParams['figure.figsize'] = (9, 4)
plt.rcParams['axes.spines.top'] = False
plt.rcParams['axes.spines.right'] = False

np.random.seed(42)
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import scipy.stats as stats
import statsmodels.formula.api as smf
import statsmodels.api as sm

plt.rcParams['figure.figsize'] = (9, 4)
plt.rcParams['axes.spines.top'] = False
plt.rcParams['axes.spines.right'] = False

np.random.seed(42)

1. El Modelo de Regresión Lineal Simple¶

El modelo poblacional es:

$$Y_i = \beta_0 + \beta_1 X_i + \varepsilon_i$$

donde:

$\beta_0$: intercepto (valor esperado de $Y$ cuando $X=0$)
$\beta_1$: pendiente (cambio esperado en $Y$ por unidad de aumento en $X$)
$\varepsilon_i \sim \mathcal{N}(0, \sigma^2)$: error aleatorio, independiente e idénticamente distribuido

Los supuestos del modelo (necesarios para que la inferencia sea válida):

Supuesto	Descripción
Linealidad	La relación entre $X$ e $Y$ es lineal
Independencia	Los errores $\varepsilon_i$ son independientes entre sí
Homocedasticidad	La varianza de $\varepsilon_i$ es constante para todo $X$
Normalidad	Los errores se distribuyen Normal

💡 El método OLS (Ordinary Least Squares) estima $\hat{\beta}_0$ y $\hat{\beta}_1$ minimizando la suma de cuadrados de los residuos: $\sum_{i=1}^{n}(y_i - \hat{y}_i)^2$.

In [ ]:

Copied!





# Dataset simulado: relación entre años de experiencia y salario
np.random.seed(42)
n = 80
experiencia = np.random.uniform(1, 20, n)
salario = 800_000 + 120_000 * experiencia + np.random.normal(0, 250_000, n)

df = pd.DataFrame({'experiencia': experiencia, 'salario': salario})

print(df.describe().round(0))
# Dataset simulado: relación entre años de experiencia y salario
np.random.seed(42)
n = 80
experiencia = np.random.uniform(1, 20, n)
salario = 800_000 + 120_000 * experiencia + np.random.normal(0, 250_000, n)

df = pd.DataFrame({'experiencia': experiencia, 'salario': salario})

print(df.describe().round(0))

In [ ]:

Copied!





# Visualización inicial
fig, ax = plt.subplots()
ax.scatter(df['experiencia'], df['salario']/1e6, alpha=0.6, color='steelblue', edgecolor='white', s=60)
ax.set_title('Salario vs Años de Experiencia')
ax.set_xlabel('Años de experiencia')
ax.set_ylabel('Salario (millones CLP)')
plt.tight_layout()
plt.show()
# Visualización inicial
fig, ax = plt.subplots()
ax.scatter(df['experiencia'], df['salario']/1e6, alpha=0.6, color='steelblue', edgecolor='white', s=60)
ax.set_title('Salario vs Años de Experiencia')
ax.set_xlabel('Años de experiencia')
ax.set_ylabel('Salario (millones CLP)')
plt.tight_layout()
plt.show()

2. Ajuste con statsmodels y Lectura del `summary()`¶

statsmodels es la librería central para inferencia en regresión. A diferencia de scikit-learn, entrega directamente todos los estadísticos de inferencia.

In [ ]:

Copied!

# Ajuste del modelo con fórmula
modelo = smf.ols(formula='salario ~ experiencia', data=df).fit()
print(modelo.summary())
# Ajuste del modelo con fórmula
modelo = smf.ols(formula='salario ~ experiencia', data=df).fit()
print(modelo.summary())

Cómo leer el `summary()`¶

El output se divide en tres bloques:

Bloque 1 — Información general del modelo:

Campo	Qué indica
`R-squared`	Proporción de varianza de $Y$ explicada por el modelo (0 a 1)
`Adj. R-squared`	$R^2$ penalizado por número de predictores
`F-statistic`	Prueba de significancia global del modelo
`Prob (F-statistic)`	Valor-p de la prueba F
`AIC / BIC`	Criterios de información para comparar modelos

Bloque 2 — Tabla de coeficientes:

Columna	Qué indica
`coef`	Estimación de $\hat{\beta}_j$
`std err`	Error estándar del coeficiente
`t`	Estadístico t de la prueba $H_0: \beta_j = 0$
`P>\\|t\\|`	Valor-p bilateral
`[0.025 0.975]`	IC al 95% para el coeficiente

Bloque 3 — Diagnósticos:

Campo	Qué indica
`Omnibus / Prob(Omnibus)`	Prueba de normalidad de residuos
`Durbin-Watson`	Detección de autocorrelación en residuos (valor ~2 es ideal)
`Jarque-Bera`	Prueba alternativa de normalidad
`Cond. No.`	Número de condición (multicolinealidad, relevante en regresión múltiple)

In [ ]:

Copied!





# Acceso programático a los resultados clave
print("=== Coeficientes ===")
print(f"β₀ (Intercepto):  {modelo.params['Intercept']:,.0f} CLP")
print(f"β₁ (Experiencia): {modelo.params['experiencia']:,.0f} CLP por año adicional")

print("\n=== Bondad de ajuste ===")
print(f"R²:               {modelo.rsquared:.4f}")
print(f"R² ajustado:      {modelo.rsquared_adj:.4f}")
print(f"F-statistic:      {modelo.fvalue:.4f}")
print(f"Prob (F):         {modelo.f_pvalue:.6f}")
# Acceso programático a los resultados clave
print("=== Coeficientes ===")
print(f"β₀ (Intercepto):  {modelo.params['Intercept']:,.0f} CLP")
print(f"β₁ (Experiencia): {modelo.params['experiencia']:,.0f} CLP por año adicional")

print("\n=== Bondad de ajuste ===")
print(f"R²:               {modelo.rsquared:.4f}")
print(f"R² ajustado:      {modelo.rsquared_adj:.4f}")
print(f"F-statistic:      {modelo.fvalue:.4f}")
print(f"Prob (F):         {modelo.f_pvalue:.6f}")

3. Prueba t sobre los Coeficientes¶

Para cada coeficiente $\beta_j$, statsmodels realiza automáticamente la prueba:

$$H_0: \beta_j = 0 \qquad H_1: \beta_j \neq 0$$

El estadístico es:

$$t = \frac{\hat{\beta}_j}{SE(\hat{\beta}_j)} \sim t_{n-2}$$

$\beta_1 = 0$ significaría que $X$ no tiene efecto lineal sobre $Y$ → el modelo no aporta información
$\beta_0 = 0$ rara vez tiene interpretación sustantiva (depende del contexto)

💡 La prueba t de $\beta_1$ y la prueba F global son equivalentes en regresión simple. Divergen en regresión múltiple.

In [ ]:

Copied!





# Tabla de inferencia de los coeficientes
tabla_coef = pd.DataFrame({
    'Coeficiente': modelo.params,
    'Std. Error':  modelo.bse,
    't-stat':      modelo.tvalues,
    'p-value':     modelo.pvalues,
    'Significativo (α=0.05)': modelo.pvalues < 0.05
})
print(tabla_coef.round(4))
# Tabla de inferencia de los coeficientes
tabla_coef = pd.DataFrame({
    'Coeficiente': modelo.params,
    'Std. Error':  modelo.bse,
    't-stat':      modelo.tvalues,
    'p-value':     modelo.pvalues,
    'Significativo (α=0.05)': modelo.pvalues < 0.05
})
print(tabla_coef.round(4))

In [ ]:

Copied!





# Visualización de la distribución t y el estadístico observado para β₁
df_residual = modelo.df_resid
t_obs = modelo.tvalues['experiencia']
p_obs = modelo.pvalues['experiencia']
t_crit = stats.t.ppf(0.975, df=df_residual)

x = np.linspace(-6, 6, 500)
fig, ax = plt.subplots()
ax.plot(x, stats.t.pdf(x, df_residual), 'steelblue', lw=2)

# Región de rechazo
for signo in [1, -1]:
    x_fill = x[x >= t_crit] if signo == 1 else x[x <= -t_crit]
    ax.fill_between(x_fill, stats.t.pdf(x_fill, df_residual), color='salmon', alpha=0.5)

ax.axvline(t_obs, color='darkred', lw=2.5, linestyle='-',
           label=f't observado = {t_obs:.2f}  (p = {p_obs:.4f})')
ax.axvline(-t_crit, color='black', lw=1.5, linestyle='--', label=f't crítico = ±{t_crit:.2f}')
ax.axvline( t_crit, color='black', lw=1.5, linestyle='--')

ax.set_title('Prueba t para β₁ (experiencia)  —  H₀: β₁ = 0')
ax.set_xlabel('t')
ax.set_ylabel('Densidad')
ax.legend(fontsize=9)
plt.tight_layout()
plt.show()
# Visualización de la distribución t y el estadístico observado para β₁
df_residual = modelo.df_resid
t_obs = modelo.tvalues['experiencia']
p_obs = modelo.pvalues['experiencia']
t_crit = stats.t.ppf(0.975, df=df_residual)

x = np.linspace(-6, 6, 500)
fig, ax = plt.subplots()
ax.plot(x, stats.t.pdf(x, df_residual), 'steelblue', lw=2)

# Región de rechazo
for signo in [1, -1]:
    x_fill = x[x >= t_crit] if signo == 1 else x[x <= -t_crit]
    ax.fill_between(x_fill, stats.t.pdf(x_fill, df_residual), color='salmon', alpha=0.5)

ax.axvline(t_obs, color='darkred', lw=2.5, linestyle='-',
           label=f't observado = {t_obs:.2f}  (p = {p_obs:.4f})')
ax.axvline(-t_crit, color='black', lw=1.5, linestyle='--', label=f't crítico = ±{t_crit:.2f}')
ax.axvline( t_crit, color='black', lw=1.5, linestyle='--')

ax.set_title('Prueba t para β₁ (experiencia)  —  H₀: β₁ = 0')
ax.set_xlabel('t')
ax.set_ylabel('Densidad')
ax.legend(fontsize=9)
plt.tight_layout()
plt.show()

4. Intervalos de Confianza para los Coeficientes¶

El IC al $(1-\alpha)\%$ para $\beta_j$ es:

$$IC(\beta_j) = \hat{\beta}_j \pm t_{\alpha/2,\ n-2} \cdot SE(\hat{\beta}_j)$$

Nos indica el rango plausible del efecto real de $X$ sobre $Y$, con una determinada confianza.

In [ ]:

Copied!





# IC para los coeficientes
ic_coefs = modelo.conf_int(alpha=0.05)
ic_coefs.columns = ['IC 2.5%', 'IC 97.5%']
ic_coefs['Estimación'] = modelo.params
print("=== IC 95% para los coeficientes ===")
print(ic_coefs[['Estimación', 'IC 2.5%', 'IC 97.5%']].round(0))
# IC para los coeficientes
ic_coefs = modelo.conf_int(alpha=0.05)
ic_coefs.columns = ['IC 2.5%', 'IC 97.5%']
ic_coefs['Estimación'] = modelo.params
print("=== IC 95% para los coeficientes ===")
print(ic_coefs[['Estimación', 'IC 2.5%', 'IC 97.5%']].round(0))

In [ ]:

Copied!





# Visualización de IC para β₁ en distintos niveles de confianza
beta1 = modelo.params['experiencia']
se_beta1 = modelo.bse['experiencia']
niveles = [0.80, 0.90, 0.95, 0.99]
colores = ['#2ecc71', '#3498db', '#9b59b6', '#e74c3c']

fig, ax = plt.subplots(figsize=(9, 3))
for i, (nivel, color) in enumerate(zip(niveles, colores)):
    alpha_i = 1 - nivel
    t_i = stats.t.ppf(1 - alpha_i/2, df=df_residual)
    li_i = beta1 - t_i * se_beta1
    ls_i = beta1 + t_i * se_beta1
    ax.plot([li_i, ls_i], [i, i], color=color, lw=5, alpha=0.7,
            label=f'IC {int(nivel*100)}%: [{li_i:,.0f} , {ls_i:,.0f}]')
    ax.plot(beta1, i, 'o', color='black', ms=6, zorder=5)

ax.axvline(0, color='gray', lw=1.5, linestyle='--', label='β₁ = 0')
ax.set_yticks(range(len(niveles)))
ax.set_yticklabels([f'{int(n*100)}%' for n in niveles])
ax.set_title('IC para β₁ (efecto de experiencia sobre salario) — distintos niveles de confianza')
ax.set_xlabel('CLP por año adicional de experiencia')
ax.legend(fontsize=8, loc='upper left')
plt.tight_layout()
plt.show()
# Visualización de IC para β₁ en distintos niveles de confianza
beta1 = modelo.params['experiencia']
se_beta1 = modelo.bse['experiencia']
niveles = [0.80, 0.90, 0.95, 0.99]
colores = ['#2ecc71', '#3498db', '#9b59b6', '#e74c3c']

fig, ax = plt.subplots(figsize=(9, 3))
for i, (nivel, color) in enumerate(zip(niveles, colores)):
    alpha_i = 1 - nivel
    t_i = stats.t.ppf(1 - alpha_i/2, df=df_residual)
    li_i = beta1 - t_i * se_beta1
    ls_i = beta1 + t_i * se_beta1
    ax.plot([li_i, ls_i], [i, i], color=color, lw=5, alpha=0.7,
            label=f'IC {int(nivel*100)}%: [{li_i:,.0f} , {ls_i:,.0f}]')
    ax.plot(beta1, i, 'o', color='black', ms=6, zorder=5)

ax.axvline(0, color='gray', lw=1.5, linestyle='--', label='β₁ = 0')
ax.set_yticks(range(len(niveles)))
ax.set_yticklabels([f'{int(n*100)}%' for n in niveles])
ax.set_title('IC para β₁ (efecto de experiencia sobre salario) — distintos niveles de confianza')
ax.set_xlabel('CLP por año adicional de experiencia')
ax.legend(fontsize=8, loc='upper left')
plt.tight_layout()
plt.show()

5. IC para Predicción vs IC para la Media¶

Dado un valor $X = x^*$, podemos construir dos tipos de intervalos muy distintos:

Intervalo	¿Para qué?	Fórmula
IC para la media (`mean_ci`)	Rango plausible del valor promedio de $Y$ para $X = x^*$	Más estrecho
IC para predicción (`obs_ci`)	Rango plausible de una observación individual para $X = x^*$	Más ancho (incluye variabilidad individual)

El IC para predicción siempre es más ancho porque agrega la variabilidad del error individual $\varepsilon$ al IC de la media:

$$IC_{pred} = \hat{y} \pm t_{\alpha/2} \cdot \hat{\sigma} \sqrt{1 + \frac{1}{n} + \frac{(x^* - \bar{x})^2}{\sum(x_i - \bar{x})^2}}$$

$$IC_{media} = \hat{y} \pm t_{\alpha/2} \cdot \hat{\sigma} \sqrt{\frac{1}{n} + \frac{(x^* - \bar{x})^2}{\sum(x_i - \bar{x})^2}}$$

In [ ]:

Copied!





# Predicciones con intervalos para un rango de valores de X
x_nuevo = pd.DataFrame({'experiencia': np.linspace(1, 20, 200)})
predicciones = modelo.get_prediction(x_nuevo)
resumen_pred = predicciones.summary_frame(alpha=0.05)

print("Columnas disponibles en summary_frame:")
print(resumen_pred.columns.tolist())
print()
print(resumen_pred.head().round(0))
# Predicciones con intervalos para un rango de valores de X
x_nuevo = pd.DataFrame({'experiencia': np.linspace(1, 20, 200)})
predicciones = modelo.get_prediction(x_nuevo)
resumen_pred = predicciones.summary_frame(alpha=0.05)

print("Columnas disponibles en summary_frame:")
print(resumen_pred.columns.tolist())
print()
print(resumen_pred.head().round(0))

In [ ]:

Copied!





# Visualización de ambos intervalos
fig, ax = plt.subplots(figsize=(10, 5))

x_plot = x_nuevo['experiencia']

# Datos observados
ax.scatter(df['experiencia'], df['salario']/1e6, alpha=0.5, color='steelblue',
           edgecolor='white', s=50, label='Datos observados', zorder=5)

# Línea de regresión
ax.plot(x_plot, resumen_pred['mean']/1e6, color='black', lw=2, label='Línea ajustada')

# IC para la media
ax.fill_between(x_plot,
                resumen_pred['mean_ci_lower']/1e6,
                resumen_pred['mean_ci_upper']/1e6,
                alpha=0.4, color='steelblue', label='IC 95% para la media')

# IC para predicción individual
ax.fill_between(x_plot,
                resumen_pred['obs_ci_lower']/1e6,
                resumen_pred['obs_ci_upper']/1e6,
                alpha=0.15, color='orange', label='IC 95% para predicción individual')

ax.set_title('Regresión lineal: IC para la media vs IC para predicción individual')
ax.set_xlabel('Años de experiencia')
ax.set_ylabel('Salario (millones CLP)')
ax.legend(fontsize=9)
plt.tight_layout()
plt.show()
# Visualización de ambos intervalos
fig, ax = plt.subplots(figsize=(10, 5))

x_plot = x_nuevo['experiencia']

# Datos observados
ax.scatter(df['experiencia'], df['salario']/1e6, alpha=0.5, color='steelblue',
           edgecolor='white', s=50, label='Datos observados', zorder=5)

# Línea de regresión
ax.plot(x_plot, resumen_pred['mean']/1e6, color='black', lw=2, label='Línea ajustada')

# IC para la media
ax.fill_between(x_plot,
                resumen_pred['mean_ci_lower']/1e6,
                resumen_pred['mean_ci_upper']/1e6,
                alpha=0.4, color='steelblue', label='IC 95% para la media')

# IC para predicción individual
ax.fill_between(x_plot,
                resumen_pred['obs_ci_lower']/1e6,
                resumen_pred['obs_ci_upper']/1e6,
                alpha=0.15, color='orange', label='IC 95% para predicción individual')

ax.set_title('Regresión lineal: IC para la media vs IC para predicción individual')
ax.set_xlabel('Años de experiencia')
ax.set_ylabel('Salario (millones CLP)')
ax.legend(fontsize=9)
plt.tight_layout()
plt.show()

In [ ]:

Copied!





# Predicción puntual para un caso específico
x_caso = pd.DataFrame({'experiencia': [5, 10, 15]})
pred_caso = modelo.get_prediction(x_caso).summary_frame(alpha=0.05)

resultado = pd.DataFrame({
    'Experiencia': [5, 10, 15],
    'Predicción':   pred_caso['mean'].values,
    'IC Media [inf]':  pred_caso['mean_ci_lower'].values,
    'IC Media [sup]':  pred_caso['mean_ci_upper'].values,
    'IC Pred [inf]':   pred_caso['obs_ci_lower'].values,
    'IC Pred [sup]':   pred_caso['obs_ci_upper'].values,
})
print(resultado.round(0).to_string(index=False))
# Predicción puntual para un caso específico
x_caso = pd.DataFrame({'experiencia': [5, 10, 15]})
pred_caso = modelo.get_prediction(x_caso).summary_frame(alpha=0.05)

resultado = pd.DataFrame({
    'Experiencia': [5, 10, 15],
    'Predicción':   pred_caso['mean'].values,
    'IC Media [inf]':  pred_caso['mean_ci_lower'].values,
    'IC Media [sup]':  pred_caso['mean_ci_upper'].values,
    'IC Pred [inf]':   pred_caso['obs_ci_lower'].values,
    'IC Pred [sup]':   pred_caso['obs_ci_upper'].values,
})
print(resultado.round(0).to_string(index=False))

6. Visualización de Residuos¶

Antes de confiar en la inferencia, es fundamental verificar los supuestos del modelo a través de los residuos $e_i = y_i - \hat{y}_i$.

Cuatro gráficos diagnósticos estándar:

In [ ]:

Copied!





residuos = modelo.resid
y_hat    = modelo.fittedvalues

fig, axes = plt.subplots(2, 2, figsize=(13, 9))

# 1. Residuos vs valores ajustados (homocedasticidad)
ax = axes[0, 0]
ax.scatter(y_hat/1e6, residuos/1e6, alpha=0.6, color='steelblue', edgecolor='white', s=50)
ax.axhline(0, color='red', lw=1.5, linestyle='--')
ax.set_title('Residuos vs Valores ajustados')
ax.set_xlabel('Valores ajustados (millones CLP)')
ax.set_ylabel('Residuos (millones CLP)')

# 2. Q-Q plot (normalidad de residuos)
ax = axes[0, 1]
sm.qqplot(residuos, line='s', ax=ax, alpha=0.6, color='steelblue')
ax.set_title('Q-Q Plot de residuos')

# 3. Histograma de residuos
ax = axes[1, 0]
ax.hist(residuos/1e6, bins=20, color='steelblue', edgecolor='white', density=True, alpha=0.8)
x_range = np.linspace(residuos.min(), residuos.max(), 200)
ax.plot(x_range/1e6, stats.norm.pdf(x_range, residuos.mean(), residuos.std())*1e6,
        'red', lw=2, label='Normal teórica')
ax.set_title('Distribución de residuos')
ax.set_xlabel('Residuos (millones CLP)')
ax.set_ylabel('Densidad')
ax.legend()

# 4. Residuos vs X (linealidad)
ax = axes[1, 1]
ax.scatter(df['experiencia'], residuos/1e6, alpha=0.6, color='seagreen', edgecolor='white', s=50)
ax.axhline(0, color='red', lw=1.5, linestyle='--')
ax.set_title('Residuos vs X (experiencia)')
ax.set_xlabel('Años de experiencia')
ax.set_ylabel('Residuos (millones CLP)')

plt.suptitle('Diagnóstico de residuos', fontsize=13, y=1.01)
plt.tight_layout()
plt.show()
residuos = modelo.resid
y_hat    = modelo.fittedvalues

fig, axes = plt.subplots(2, 2, figsize=(13, 9))

# 1. Residuos vs valores ajustados (homocedasticidad)
ax = axes[0, 0]
ax.scatter(y_hat/1e6, residuos/1e6, alpha=0.6, color='steelblue', edgecolor='white', s=50)
ax.axhline(0, color='red', lw=1.5, linestyle='--')
ax.set_title('Residuos vs Valores ajustados')
ax.set_xlabel('Valores ajustados (millones CLP)')
ax.set_ylabel('Residuos (millones CLP)')

# 2. Q-Q plot (normalidad de residuos)
ax = axes[0, 1]
sm.qqplot(residuos, line='s', ax=ax, alpha=0.6, color='steelblue')
ax.set_title('Q-Q Plot de residuos')

# 3. Histograma de residuos
ax = axes[1, 0]
ax.hist(residuos/1e6, bins=20, color='steelblue', edgecolor='white', density=True, alpha=0.8)
x_range = np.linspace(residuos.min(), residuos.max(), 200)
ax.plot(x_range/1e6, stats.norm.pdf(x_range, residuos.mean(), residuos.std())*1e6,
        'red', lw=2, label='Normal teórica')
ax.set_title('Distribución de residuos')
ax.set_xlabel('Residuos (millones CLP)')
ax.set_ylabel('Densidad')
ax.legend()

# 4. Residuos vs X (linealidad)
ax = axes[1, 1]
ax.scatter(df['experiencia'], residuos/1e6, alpha=0.6, color='seagreen', edgecolor='white', s=50)
ax.axhline(0, color='red', lw=1.5, linestyle='--')
ax.set_title('Residuos vs X (experiencia)')
ax.set_xlabel('Años de experiencia')
ax.set_ylabel('Residuos (millones CLP)')

plt.suptitle('Diagnóstico de residuos', fontsize=13, y=1.01)
plt.tight_layout()
plt.show()

¿Qué buscar en cada gráfico?

Gráfico	Supuesto que verifica	Patrón problemático
Residuos vs ajustados	Homocedasticidad	Patrón de abanico (varianza crece con $\hat{y}$)
Q-Q plot	Normalidad	Puntos alejados de la diagonal
Histograma	Normalidad	Distribución muy asimétrica o multimodal
Residuos vs X	Linealidad	Patrón curvilíneo

7. Ejemplo Integrador¶

Contexto: Una empresa de logística quiere modelar el costo de envío en función de la distancia recorrida. Se dispone de 100 registros de envíos recientes.

In [ ]:

Copied!





np.random.seed(55)
n = 100
distancia = np.random.uniform(10, 500, n)  # km
costo = 5000 + 80 * distancia + np.random.normal(0, 8000, n)  # CLP

df_log = pd.DataFrame({'distancia': distancia, 'costo': costo})

# Ajuste
modelo_log = smf.ols('costo ~ distancia', data=df_log).fit()
print(modelo_log.summary())
np.random.seed(55)
n = 100
distancia = np.random.uniform(10, 500, n)  # km
costo = 5000 + 80 * distancia + np.random.normal(0, 8000, n)  # CLP

df_log = pd.DataFrame({'distancia': distancia, 'costo': costo})

# Ajuste
modelo_log = smf.ols('costo ~ distancia', data=df_log).fit()
print(modelo_log.summary())

In [ ]:

Copied!





# Interpretación de resultados
b0 = modelo_log.params['Intercept']
b1 = modelo_log.params['distancia']
ic_b1 = modelo_log.conf_int().loc['distancia']
p_b1  = modelo_log.pvalues['distancia']
r2    = modelo_log.rsquared

print("=== Interpretación del modelo ===")
print(f"β₀ = {b0:,.0f} CLP  →  Costo base (distancia = 0)")
print(f"β₁ = {b1:,.0f} CLP  →  Cada km adicional incrementa el costo en {b1:,.0f} CLP")
print(f"IC 95% β₁: [{ic_b1[0]:,.0f} , {ic_b1[1]:,.0f}] CLP/km")
print(f"p-value β₁: {p_b1:.6f}  →  {'Significativo' if p_b1 < 0.05 else 'No significativo'} al 5%")
print(f"R² = {r2:.4f}  →  El modelo explica el {r2*100:.1f}% de la variabilidad del costo")
# Interpretación de resultados
b0 = modelo_log.params['Intercept']
b1 = modelo_log.params['distancia']
ic_b1 = modelo_log.conf_int().loc['distancia']
p_b1  = modelo_log.pvalues['distancia']
r2    = modelo_log.rsquared

print("=== Interpretación del modelo ===")
print(f"β₀ = {b0:,.0f} CLP  →  Costo base (distancia = 0)")
print(f"β₁ = {b1:,.0f} CLP  →  Cada km adicional incrementa el costo en {b1:,.0f} CLP")
print(f"IC 95% β₁: [{ic_b1[0]:,.0f} , {ic_b1[1]:,.0f}] CLP/km")
print(f"p-value β₁: {p_b1:.6f}  →  {'Significativo' if p_b1 < 0.05 else 'No significativo'} al 5%")
print(f"R² = {r2:.4f}  →  El modelo explica el {r2*100:.1f}% de la variabilidad del costo")

In [ ]:

Copied!





# Predicción para nuevos envíos
nuevos_envios = pd.DataFrame({'distancia': [50, 150, 300, 450]})
pred_log = modelo_log.get_prediction(nuevos_envios).summary_frame(alpha=0.05)

resultado_log = pd.DataFrame({
    'Distancia (km)':  [50, 150, 300, 450],
    'Costo estimado':   pred_log['mean'].round(0),
    'IC Media inf':     pred_log['mean_ci_lower'].round(0),
    'IC Media sup':     pred_log['mean_ci_upper'].round(0),
    'IC Pred inf':      pred_log['obs_ci_lower'].round(0),
    'IC Pred sup':      pred_log['obs_ci_upper'].round(0),
})
print(resultado_log.to_string(index=False))
# Predicción para nuevos envíos
nuevos_envios = pd.DataFrame({'distancia': [50, 150, 300, 450]})
pred_log = modelo_log.get_prediction(nuevos_envios).summary_frame(alpha=0.05)

resultado_log = pd.DataFrame({
    'Distancia (km)':  [50, 150, 300, 450],
    'Costo estimado':   pred_log['mean'].round(0),
    'IC Media inf':     pred_log['mean_ci_lower'].round(0),
    'IC Media sup':     pred_log['mean_ci_upper'].round(0),
    'IC Pred inf':      pred_log['obs_ci_lower'].round(0),
    'IC Pred sup':      pred_log['obs_ci_upper'].round(0),
})
print(resultado_log.to_string(index=False))

In [ ]:

Copied!





# Visualización final
x_plot = pd.DataFrame({'distancia': np.linspace(10, 500, 300)})
pred_plot = modelo_log.get_prediction(x_plot).summary_frame(alpha=0.05)

fig, axes = plt.subplots(1, 2, figsize=(14, 5))

# Panel izquierdo: regresión con bandas
ax = axes[0]
ax.scatter(df_log['distancia'], df_log['costo']/1000, alpha=0.5, color='steelblue',
           edgecolor='white', s=50, label='Datos')
ax.plot(x_plot['distancia'], pred_plot['mean']/1000, 'black', lw=2, label='Ajuste')
ax.fill_between(x_plot['distancia'],
                pred_plot['mean_ci_lower']/1000, pred_plot['mean_ci_upper']/1000,
                alpha=0.4, color='steelblue', label='IC 95% media')
ax.fill_between(x_plot['distancia'],
                pred_plot['obs_ci_lower']/1000, pred_plot['obs_ci_upper']/1000,
                alpha=0.15, color='orange', label='IC 95% predicción')
ax.set_title('Costo de envío vs Distancia')
ax.set_xlabel('Distancia (km)')
ax.set_ylabel('Costo (miles CLP)')
ax.legend(fontsize=8)

# Panel derecho: diagnóstico de residuos
ax = axes[1]
ax.scatter(modelo_log.fittedvalues/1000, modelo_log.resid/1000,
           alpha=0.6, color='seagreen', edgecolor='white', s=50)
ax.axhline(0, color='red', lw=1.5, linestyle='--')
ax.set_title('Residuos vs Valores ajustados')
ax.set_xlabel('Valores ajustados (miles CLP)')
ax.set_ylabel('Residuos (miles CLP)')

plt.tight_layout()
plt.show()
# Visualización final
x_plot = pd.DataFrame({'distancia': np.linspace(10, 500, 300)})
pred_plot = modelo_log.get_prediction(x_plot).summary_frame(alpha=0.05)

fig, axes = plt.subplots(1, 2, figsize=(14, 5))

# Panel izquierdo: regresión con bandas
ax = axes[0]
ax.scatter(df_log['distancia'], df_log['costo']/1000, alpha=0.5, color='steelblue',
           edgecolor='white', s=50, label='Datos')
ax.plot(x_plot['distancia'], pred_plot['mean']/1000, 'black', lw=2, label='Ajuste')
ax.fill_between(x_plot['distancia'],
                pred_plot['mean_ci_lower']/1000, pred_plot['mean_ci_upper']/1000,
                alpha=0.4, color='steelblue', label='IC 95% media')
ax.fill_between(x_plot['distancia'],
                pred_plot['obs_ci_lower']/1000, pred_plot['obs_ci_upper']/1000,
                alpha=0.15, color='orange', label='IC 95% predicción')
ax.set_title('Costo de envío vs Distancia')
ax.set_xlabel('Distancia (km)')
ax.set_ylabel('Costo (miles CLP)')
ax.legend(fontsize=8)

# Panel derecho: diagnóstico de residuos
ax = axes[1]
ax.scatter(modelo_log.fittedvalues/1000, modelo_log.resid/1000,
           alpha=0.6, color='seagreen', edgecolor='white', s=50)
ax.axhline(0, color='red', lw=1.5, linestyle='--')
ax.set_title('Residuos vs Valores ajustados')
ax.set_xlabel('Valores ajustados (miles CLP)')
ax.set_ylabel('Residuos (miles CLP)')

plt.tight_layout()
plt.show()

Resumen¶

Concepto	Qué es	Cómo obtenerlo en statsmodels
Estimación OLS	$\hat{\beta}_0$, $\hat{\beta}_1$ que minimizan $\sum e_i^2$	`smf.ols().fit()`
Prueba t de $\beta_j$	$H_0: \beta_j = 0$ usando $t = \hat{\beta}_j / SE(\hat{\beta}_j)$	`modelo.tvalues`, `modelo.pvalues`
IC para coeficientes	Rango plausible del efecto real de $X$	`modelo.conf_int()`
$R^2$	Proporción de varianza explicada	`modelo.rsquared`
Prueba F global	¿El modelo en conjunto es significativo?	`modelo.fvalue`, `modelo.f_pvalue`
IC para la media	Rango para el valor promedio de $Y$ dado $X = x^*$	`get_prediction().summary_frame()` → `mean_ci_*`
IC para predicción	Rango para una observación individual dado $X = x^*$	`get_prediction().summary_frame()` → `obs_ci_*`
Diagnóstico residuos	Verificación de supuestos del modelo	`modelo.resid`, `modelo.fittedvalues`, `sm.qqplot()`

Ideas clave:

Un coeficiente significativo ($p < 0.05$) indica que el predictor tiene un efecto lineal detectable en los datos
El IC de predicción siempre es más ancho que el IC de la media — no confundirlos al comunicar resultados
Los gráficos de residuos son obligatorios antes de interpretar cualquier resultado inferencial
La significancia estadística no implica relevancia práctica: evaluar también la magnitud de $\hat{\beta}_1$

Referencias¶

James, G., Witten, D., Hastie, T., Tibshirani, R. (2021). An Introduction to Statistical Learning. Springer. Cap. 3.
Statsmodels OLS: https://www.statsmodels.org/stable/generated/statsmodels.regression.linear_model.OLS.html
Statsmodels get_prediction: https://www.statsmodels.org/stable/generated/statsmodels.regression.linear_model.RegressionResults.get_prediction.html

Regresión Lineal Simple¶

Librerías¶

1. El Modelo de Regresión Lineal Simple¶

2. Ajuste con statsmodels y Lectura del summary()¶

Cómo leer el summary()¶

3. Prueba t sobre los Coeficientes¶

4. Intervalos de Confianza para los Coeficientes¶

5. IC para Predicción vs IC para la Media¶

6. Visualización de Residuos¶

7. Ejemplo Integrador¶

Resumen¶

Referencias¶

2. Ajuste con statsmodels y Lectura del `summary()`¶

Cómo leer el `summary()`¶