Regresión lineal en Python para epidemiólogos en 6 pasos
En este tutorial cubriremos los siguientes pasos
1. Abre el conjunto de datos
2. Explora los datos
3. Haz una pregunta de investigación (que se pueda responder usando un modelo de regresión lineal)
4. Ejecutar y leer una regresión lineal simple
5. Ejecutar y leer una regresión lineal múltiple
6. Respondiendo a nuestra pregunta de investigación
¡Empecemos!
1. Abre el conjunto de datos
Vamos a extraer una versión disponible públicamente del Framingham heart study; https://www.kaggle.com/amanajmera1/framingham-heart-study-dataset/version/1
En esta web también se puede consultar la codificación de las diferentes variables. Por ejemplo, con respecto al género (variable masculino; masculino=1 significa masculino y masculino=2 significa femenino)
Primero descarga el archivo a tu computadora.
Aquí la Versión de Jupyter Notebook (Conjunto de datos descargado de Kaggle).
Aquí está el Kaggle Notebook para que pueda ejecutar todo en tu navegador
Primero importa todos los paquetes relevantes
import numpy as np
import pandas as pd
import statsmodels.formula.api as smf
Si es la primera vez que usa StatsModels, es posible que deba usar esto antes (instalación desde Anaconda, distribución de Python).
conda install -c conda-forge statsmodels
Luego abra la base de datos de Framingham. Lo llamaremos df para facilitar las cosas.
df = pd.read_csv('framingham.csv')
A partir de ahora, cada vez que queramos "llamar" a nuestro marco de datos, usaremos df en nuestro código.
Echemos un vistazo a las primeras 5 filas y todas nuestras columnas.
df.head()
Veamos cuántas columnas y filas tenemos…
df.shape()
(4240, 16)
El resultado significa que la base de datos tiene 4240 observaciones (filas) y 16 variables (columnas). Recuerda .head() solo nos dio un vistazo de las 4240 observaciones/filas.
Ahora describamos nuestras variables
df.describe()
El paréntesis vacío en df. describe()le dice a Python que muestre todas las variables disponibles.
Aquí podemos ver el número de observaciones por columna (recuento), seguido de la media y la desviación estándar (sd). Ten en cuenta que para las variables binarias discretas, la interpretación es algo diferente y la media será igual a la proporción relativa (0 a 1, en lugar de 0–100%).
3. Haz una pregunta de investigación
¿Existe una asociación entre la presión arterial sistólica y la edad, después de ajustar los factores de confusión relevantes?
¿Cómo podemos responder a esto?
Aquí tenemos un resultado continuo (presión arterial sistólica en mmHg) y una covariable continua independiente (edad en años).
Dado que la presión arterial sistólica (PAS) (variable denominada sysbp) tiene una distribución normal y queremos probar una asociación, debemos optar por regresión lineal .
Consideraremos las siguientes variables disponibles como posibles factores de confusión relevantes para la asociación entre PA y edad; género (variable binaria: masculino), educación (variable categórica: educación), tabaquismo (variable continua: cigsPerDay para cigarrillos por día), ingesta de medicamentos para la hipertensión (variable binaria: BPMeds) y colesterol total ( variable continua: totChol).
Antes de ejecutar los números tenemos que definir una hipótesis nula y una hipótesis alternativa:
- Hipótesis nula (H0): NO existe asociación entre la PAS y la edad (considerando los factores de confusión elegidos);
- Alternative hypothesis (H1): There is an association between sBP and age (considering the chosen confounders).
4. Ejecutar y leer una regresión lineal simple
Vamos para una regresión lineal simple.
Describamos el modelo. Una alternativa clásica en ciencia de datos sería usar el paquete Scikilt-learn . En su lugar, estamos utilizando el paquete StatsModels , que es más fácil de usar para aquellos con experiencia en epidemiología que ya están familiarizados con otros lenguajes y programas, como Stata y R /RStudio.
model1 = smf.ols(formula='sysBP ~ age', data=df).fit()
- smf llama al paquete StatsModels
- ols le dice a Python que estamos usando una regresión de mínimos cuadrados ordinarios (OLS) (un tipo de regresión lineal)
- formula= usada para escribir la(s) variable(s) dependiente(s) y todas las independientes
- primera variable dentro del paréntesis antes de “~”/variable dependiente/resultado: la primera variable es nuestra única variable dependiente. Este es nuestro resultado, la variable que determina qué tipo de regresión usar y la que se asociará con todas las demás covariables;
- ~ dentro de paréntesis: marca el límite entre el resultado (variable dependiente) a la izquierda y las covariables (variables independientes) a la derecha;
- covariables independientes/variables independientes: todas las demás variables después de "~", entre paréntesis;
- signo + entre paréntesis: el signo + se usa para separar diferentes variables independientes dentro del mismo modelo (útil para modelos multivariables, también conocidos como muchas variables independientes)
- ,data=: This marks the name of the data frame
- .fit() le dice a Python que queremos llenar nuestra función ("ejecutar la función")
model1 = smf.ols(formula='sysBP ~ age', data=df).fit()
print(model1.summary())
En la mitad superior de la salida (dentro del cuadrado rojo) podemos ver algunos datos generales importantes sobre el modelo.
En la columna de la izquierda vemos
- La variable dependiente es sysBP
- Estamos ejecutando un modelo de regresión mínimo cuadrado ordinario (OLS) (una forma de regresión lineal)
- Se incluyeron 4240 observaciones (pacientes).
En la columna de la derecha, podemos ver información importante sobre el diagnóstico del modelo:
- R-cuadrado = 0,155, lo que significa el modelo explica el 16 % de la variabilidad de los resultados. El R-cuadrado varía de 0 a 1 y mide qué tan cerca están los datos de la línea de regresión ajustada. En otras palabras, indica con %, la variabilidad del resultado (datos de respuesta) es explicada por el modelo;
- El valor F-statistic se utiliza para calcular el valor p del modelo, Prob (F-statistic), que aquí es <0,05. Esto también nos dice que Python está usando una prueba ANOVA, lo que implica una distribución F.
Dado que estamos realizando un estudio epidemiológico para observar/probar una asociación, en lugar de construir un modelo predictivo, no dividiremos nuestros datos en entrenamiento y prueba. fuerte>. Por lo tanto, no se sorprenda si este enfoque parece diferente de lo que leerá en la mayoría de los tutoriales de ciencia de datos.
Con respecto ao output en la mitad inferior
Vamos a leer el output, columna por columna. Aquí estamos evaluando la asociación entre nuestro resultado (sysbp) y nuestra variable independiente, la edad (primera columna), que tiene una asociación positiva de 1,0128 (o coeficiente, o segunda columna), con un error estándar de 0,036 (tercera columna), con una relación t (para un distribución t) de 27,911 (cuarta columna) que da un valor p <0,001 (quinta columna) , correspondiente a un intervalo de confianza del 95 % para el coeficiente beta entre 0,94 2 (límite inferior, la sexta columna) y 1,084 (límite superior, séptima columna).
Para hacer las cosas más simples…
Podemos rechazar la hipótesis nula de no asociación entre la edad y la PAS. Existe una asociación positiva (coeficiente beta > 1,0) entre la edad y la presión arterial sistólica, que es estadísticamente significativa (p < 0,01 y el IC del 95 % no incluye un coeficiente = 0,0).
Recuerda que el coeficiente beta puede tomar un valor entre -∞ y +∞.
- Un coeficiente positivo representa una asociación lineal positiva;
- Un coeficiente negativo representa una asociación lineal negativa.
Con respecto al intervalo de confianza del 95 % (IC del 95 %) para el coeficiente beta
- Si el IC del 95 % incluye 0, la asociación no es estadísticamente significativa (p≥0,05), no se puede rechazar la hipótesis nula de no asociación;
- Si el IC del 95 % no incluye el 0, la asociación es estadísticamente significativa (p<0,05), la hipótesis nula de asociación puede rechazarse.
¿Qué pasa con la “(intercepción)”?
En este caso (Intercepto)= 82,14 representa el intercepto, o el β0. Esto significa que, cuando la edad = 0, sBP = 82 mmHg. En este contexto, el β0 es irrelevante para nuestra pregunta de investigación, porque NO estamos construyendo un modelo predictivo. Estamos evaluando una asociación.
5. Ejecutar y leer una regresión lineal múltiple
El proceso aquí es bastante simple, simplemente agregue todas las covariables (variables independientes) después de su variable independiente (en este caso, después de sysbp).
model2 = smf.ols(formula='sysBP ~ age + male + education + cigsPerDay + BPMeds + totChol', data=df).fit()print(model2.summary())
Antes de leer el resultado, vemos que este es un modelo válido pero... no es exactamente lo que queríamos... Las variables discretas con múltiples categorías (por ejemplo, educación) no se estratifican por grupos. Podemos especificar el tipo de covariables.
model3 = smf.ols(formula='sysBP ~ age + male + C(education) + cigsPerDay + BPMeds + totChol', data=df).fit()print(model3.summary())
Genial, ahora podemos pasar a la interpretación de nuestro output
Echemos un vistazo a las covariables continuas.
En primer lugar, considera que hay un valor ligeramente diferente para el coeficiente y el IC del 95 % para edad, en comparación con la regresión lineal simple anterior. Este es un fenómeno normal y esperado. Como se observó antes, hay una regresión lineal positiva estadísticamente significativa.
Solo tenemos 4007 observaciones. Esto se debe a que 233 pacientes (4240–4007) tenían al menos una observación faltante para cualquiera de las variables del modelo.
Recuerda el OLS clásico solo puede ejecutarse con datos completos, lo que significa que cada vez que una observación (fila) tenga datos faltantes para al menos una covariable, esta fila no ingresará al modelo.
How do we read the output; taking as example age:
- The null hypothesis (H0) is: In the population, there is no linear association between age and sBP, when controlling for gender, education, smoking, anti-hypertensors intake, and total cholesterol (in other words, when controlling for all other covariates present in the model);
- Since p-value<0.05, we reject H0, meaning there is a linear association between age and sBP, when controlling for gender, education, smoking, anti-hypertensors intake, and total cholesterol;
- Interpretation of the slope: For each 1 unit increase in age, mean sBP increases +0.84mmHg.
Notice an important pattern regarding p-values and 95%CI
- When p<0.05, the confidence interval does NOT include 0.0 (eg: age, beta coefficient (95%CI): 0.84 [0.76;0.91]);
- When p≥0.05 the confidence interval includes 0.0 (eg; male, beta coefficient (95% CI): -0.087 [-1.40;+1.22]).
Echemos un vistazo a las covariables discretas binarias
Veamos el género (male= 0 significa paciente femenino; male = 1 significa paciente masculino):
- La hipótesis nula (H0) es: En la población, no hay diferencia entre hombres y mujeres en su valor promedio de PAS, cuando se controla por edad, educación, tabaquismo, consumo de antihipertensores y colesterol total;
- Aquí coeficiente beta (95% IC) = -0.09 (-1.4;+1.23). Dado que el IC del 95 % incluye 0,0 (y p≥0,05) no podemos rechazar H0.
Echemos un vistazo a las covariables discretas no binarias
Ahora llegamos a la parte más desconcertante del modelo... La covariable no binaria discreta; Education (Educación).
Pero, ¿qué falta? Category one (education=1)
¿Por qué? Porque todas las demás categorías se comparan en relación con la categoría de referencia (educación=1)
df.education.value_counts()
Como vemos, hay 4 categorías (1, 2, 3, 4). En la base de datos original
- 1 = algo de escuela secundaria
- 2= escuela secundaria o equivalente
- 3= Algún tiempo de universidad o escuela profesional
- 4= Educación Superior
Por convención, el modelo utiliza la categoría con el índice más bajo (en este caso 1, algún instituto) como referencia. Todas las categorías restantes en el modelo se comparan con la categoría de referencia.
Pero, ¿cómo leo los coeficientes de educación?
En el caso de educación, mirando primero la primera fila, cuando educación=2
- La hipótesis nula (H0) es: En la población, no hay diferencia entre aquellos con “algo de secundaria” (education=1) y aquellos con “bachillerato” (education=2 ) en su valor medio de PAS, controlando por edad, sexo, tabaquismo, ingesta de antihipertensores y colesterol total;
- Dado que p-value≥0.05, rechazamos H0, no podemos rechazar la hipótesis nula para esta categoría.
¿Qué pasa con las diferentes categorías de educación?
Vemos que para la 3ra categoría de Educación (“algún Colegio”), el coeficiente beta (95% IC): -3.57 (-5.35;-1.79), lo que significa que; “En comparación con los pacientes con algo de educación secundaria (categoría de educación = 1), los pacientes con algo de educación universitaria (educación = 3) tienen una PAS media más baja, y esta diferencia es estadísticamente significativa”.
Considera que en este modelo no podemos comparar la asociación entre los pacientes de la categoría Education=2 y los pacientes de la categoría Education=3. Para ello, tendríamos que cambiar la referencia.
6. Respondiendo a nuestra pregunta de investigación
¿Existe una asociación entre la presión arterial sistólica y la edad?
Sí, existe una asociación positiva entre la edad y la presión arterial sistólica (coeficiente de Betta [95%IC]: 0,84[0,76;0,91]), que es estadísticamente significativa (p<0,05 e intervalo de confianza no incluyen 0), incluso cuando se considera el efecto de posibles factores de confusión relevantes (covariables en el modelo).
ES 
EN 
PT