Esta nota pretende la introducción, de la forma más amena posible pero de forma extensa, de conceptos en relación con el manejo de los modelos de regresión logística incondicional, es decir, modelos basados en observaciones independientes. Existen modificaciones del modelo incondicional que nos permiten manejar datos dependientes, como ocurre para los estudios caso-control pareados.
A lo largo del texto general aparecerán las mínimas formulaciones necesarias, atendiendo principalmente a conceptos y estrategias de diseño. Aparte del texto general, en el apartado de Anexos, aparecerán siempre que se crea necesario, explicaciones breves, ejemplos o bien desarrollos matemáticos para que aquellos lectores que posean un conocimiento de cálculo suficiente puedan desarrollar sus aplicaciones informáticas propias.
Os recomiendo la siguiente bibliografía:
Los modelos de regresión son modelos estadísticos en los que se desea conocer la relación entre:
... siendo la ecuación inicial del modelo de tipo exponencial, si bien su transformación logarítmica (logit) permite su uso como una función lineal.
Como vemos, las covariables pueden ser cuantitativas o cualitativas. Las covariables cualitativas deben ser dicotómicas, tomando valores 0 para su ausencia y 1 para su presencia (esta codificación es importante, ya que cualquier otra codificación provocaría modificaciones en la interpretación del modelo). Pero si la covariable cualitativa tuviera más de dos categorías, para su inclusión en el modelo deberíamos realizar una transformación de la misma en varias covariables cualitativas dicotómicas ficticias o de diseño (las llamadas variables dummy), de forma que una de las categorías se tomaría como categoría de referencia. Con ello cada categoría entraría en el modelo de forma individual. En general, si la covariable cualitativa posee n categorías, habrá que realizar n-1 covariables ficticias. Un ejemplo con la covariable color de ojos:
| Categorías | F1 | F2 |
|---|---|---|
| Azules (categoría de referencia) | 0 | 0 |
| Verdes | 1 | 0 |
| Marrones | 0 | 1 |
En este ejemplo un sujeto de ojos azules (la categoría de referencia), entraría en el modelo con F1=0 y F2=0, mientras que un sujeto con ojos verdes entraría con F1=1 y F2=0.
Por sus características, los modelos de regresión logística permiten dos finalidades:
En esta nota me ocuparé únicamente de los modelos de regresión logística binaria. Resultan los de mayor interés ya que la mayor parte de las circunstancias analizadas en medicina responden a este modelo (presencia o no de enfermedad, éxito o fracaso, etc). Como hemos visto, la variable dependiente será una variable dicotómica que se codificará como 0 ó 1 (ausencia y presencia respectivamente). Este aspecto de la codificación de las variables no es vanal (influye en la forma en que se realizan los cálculos matemáticos), y habrá que tenerlo muy en cuenta si empleamos paquetes estadísticos que no recodifican automáticamente nuestras variables cuando éstas se encuentran codificadas de forma diferente (por ejemplo el frecuente uso de 1 para la presencia y -1 ó 2 para la ausencia).
La ecuación de partida en los modelos de regresión logística es como sigue:
... siendo P(y=1|X) la probabilidad de que y tome el valor 1 (presencia de la característica estudiada), en presencia de las covariables X (aquí X es un conjunto de n covariables x1, x2, ... , xn-1, xn). Los componentes de esta ecuación son:
Es lo que se denomina distribución logística. En la siguiente imagen vemos un ejemplo de esta distribución: la probabilidad de padecer enfermedad coronaria en función de la edad. Como puede verse, la relación entre la variable dependiente (cualitativa dicotómica), y la covariable (edad, cuantitativa continua en este caso), no es definida por una recta (lo que correspondería un modelo lineal), sino que decribe una forma sigmoidea (distribución logística).
Figura 1
Si dividimos la expresión anterior de la Ecuación 1 por su complementario, es decir, si construimos su odds (en el ejemplo de presencia o no de enfermedad, la probabilidad de estar enfermo entre la probabilidad de estar sano), obtenemos una expresión de de más fácil manejo matemático:
Pero esta expresión aún es difícil de interpretar. Su representación gráfica es como se ve en la Figura 2.
Figura 2
Si ahora realizamos su transformación logarítmica con el logaritmo natural, obtenemos una ecuación lineal que es lógicamente de manejo matemático aún más fácil y de mayor comprensión:
En la expresión de la Ecuación 3 vemos a la izquierda de la igualdad el llamado logit, es decir, el logaritmo natural de la odss de la variable dependiente (esto es, el logaritmo de la razón de proporciones de enfermar, de fallecer, de éxito, etc). El término a la derecha de la igualdad es la expresión de un recta, idéntica a la del modelo general de regresión lineal:
Ecuación 4
Siguiendo el ejemplo de las Figuras 1 y 2, podemos representar el logit frente a la edad de la forma que se observa en la Figura 3.
Pero la regresión lineal presenta una diferencia fundamental respecto al modelo de regresión logística. En el modelo de regresión lineal se asume que los errores estándar de cada coeficiente siguen una distribución normal de media 0 y varianza constante (homoscedasticidad). En el caso del modelo de regresión logística no pueden realizarse estas asunciones pues la variable dependiente no es continua (sólo puede tomar dos valores, 0 ó 1, pero ningún valor intermedio). Si llamamos ε al posible error de predicción para cada covariable xi,, tendremos que el error cometido dependerá del valor que llegue a tomar la variable dependiente y, tal como vemos en la Ecuación 5.
Ecuación 5
Esto implica que ε sigue una distribución binomial, con media y varianza proporcionales al tamaño muestral y a P(y=1|xi) (la probabilidad de que y=1 dada la presencia de xi).
Para la estimación de los coeficientes del modelo y de sus errores estándar se recurre al cálculo de estimaciones de máxima verosimilitud, es decir, estimaciones que maximicen la probabilidad de obtener los valores de la variable dependiente Y proporcionados por los datos de nuestra muestra. Estas estimaciones no son de cálculo directo, como ocurre en el caso de las estimaciones de los coeficientes de regresión de la regresión lineal múltiple por el método de los mínimos cuadrados. Para el cálculo de estimaciones máximo-verosímiles se recurre a métodos iterativos, como el método de Newton-Raphson. Dado que el cálculo es complejo, normalmente hay que recurrir al uso de rutinas de programación o a paquetes estadísticos. De estos métodos surgen no sólo las estimaciones de los coeficientes de regresión, sino también de sus errores estándar y de las covarianzas entre las covariables del modelo.
El siguiente paso será comprobar la significación estadística de cada uno de los coeficientes de regresión en el modelo. Para ello podemos emplear básicamente tres métodos:
El estadístico de Wald. Contrasta la hipótesis de que un coeficiente aislado es distinto de 0, y sigue una distribución normal de media 0 y varianza 1. Su valor para un coeficiente concreto viene dado por el cociente entre el valor del coeficiente y su correspondiente error estándar. La obtención de significación indica que dicho coeficiente es diferente de 0 y merece la pena su conservación en el modelo. En modelos con errores estándar grandes, el estadístico de Wald puede proporcional falsas ausencias de significación (es decir, se incrementa el error tipo II). Tampo es recomendable su uso si se están empleando variables de diseño.
El estadístico G de razón de verosimilitud. Se trata de ir contrastando cada modelo que surge de eliminar de forma aislada cada una de las covariables frente al modelo completo. En este caso cada estadístico G sigue una χ2 con 1 grado de libertad (no se asume normalidad). La ausencia de significación implica que el modelo sin la covariable no empeora respecto al modelo completo (es decir, da igual su presencia o su ausencia), por lo que según la estrategia de obtención del modelo más reducido (principio de parsimonia), dicha covariable debe ser eliminada del modelo ya que no aporta nada al mismo. Esta prueba no asume ninguna distribución concreta, por lo que es la más recomendada para estudiar la significación de los coeficientes.
La prueba Score. Su cálculo para el caso de una única variable viene dado por la Ecuación 6.
Ecuación 6
Al igual que en los casos anteriores, si alcanza significación nos indica que la covariable debería permanecer en el modelo. Su uso en algunos paquetes estadísticos ha quedado relegado a la selección de variables en métodos paso a paso (por la mayor rapidez de cálculo).
Cuando la covariable es cualitativa con n categorías (siendo n > 2), en el modelo se analizará la significación de cada una de sus n-1 variables ficticias, así como la significación global de la covariable comparando la presencia en bloque frente a la ausencia en bloque de sus n-1 covariables ficticias.
En el siguiente ejemplo, tomado de Hosmer y realizado con SPSS®, se analiza la variable edad (AGE) y la variable IVHX (usuario de drogas por vía parenteral); ésta segunda era una variable con tres categorías (nunca, previa y reciente), por lo que se crearon dos variables ficticias: IVHX(1) e IVHX(2); el resultado es una estimación de los β con sus errores estándar, la significación para IVHX(1) e IVHX(2), y la significación de IVHX considerada como la entrada frente a la salida en bloque del modelo de IVHX(1) e IVHX(2).
Una vez hemos estimado los coeficientes de regresión y sus correspondientes errores estándar debemos calcular los correspondientes intervalos de confianza para nuestras estimaciones. Cada intervalo de confianza se calculará bajo la hipótesis de que dichos coeficientes se distribuyen según respectivas distribuciones normales, por lo que para un determinado coeficiente su intervalo de confianza al 95% vendrá dado por la Ecuación 7.
Junto a la significación del estadístico que hayamos empleado para contrastar la significación de los coeficientes de regresión, la inclusión de la unidad en el intervalo de confianza es, lógicamente, indicativa de la ausencia de significación.
En ocasiones nos encontraremos con modelos que nos llaman la atención por la falta de sentido de sus estimaciones. Esta sorpresa suele venir dada por la presencia de estimaciones de grandes errores estándar, con frecuencia asociadas a estimaciones de coeficientes de regresión también anormalmente elevados. Las posibles causas de este hecho pueden ser:
Pero, del conjunto de variables que podemos tener en un estudio, ¿que variables deben introducirse en el modelo?
El modelo debe ser aquél más reducido que explique los datos (principio de parsimonia), y que además sea clínicamente congruente e interpretable. Hay que tener en cuenta que un mayor número de variables en el modelo implicará mayores errores estándar.
Deben incluirse todas aquellas variables que se consideren clínicamente importantes para el modelo, con independencia de si un análisis univariado previo se demostró o no su significación estadística. Por otro lado, no debería dejarse de incluir toda variable que en un análisis univariado previo demostrara una relación "suficiente" con la variable dependiente. Como puede verse no se habla de significación estadística (p<0.05), que sería un criterio excesivamente restrictivo, sino de un cierto grado de relación (por ejemplo p<0.25). La laxitud de esta recomendación se debe a que un criterio tan restrictivo como una p<0.05 puede llevarnos a dejar de incluir en el modelo covariables con una débil asociación a la variable dependiente en solitario pero que podrían demostrar ser fuertes predictores de la misma al tomarlas en conjunto con el resto de covariables.
Una cuestión importante a tener en cuenta es el correcto manejo de las variables cualitativas transformadas en varias variables ficticias.
Una vez tenemos nuestro modelo inicial debemos recurrir a su reducción hasta el modelo más reducido que siga explicando nuestros datos. Para ello podemos recurrir a métodos de selección paso a paso, bien mediante inclusión "hacia adelante" o por eliminación "hacia atrás", o a la selección de variables por mejores subconjuntos de covariables. Estos métodos se encuentran implementados en numerosos paquetes estadísticos, por lo que son muy populares. Dado que para la comprensión de los métodos de selección paso a paso se requiere un conocimiento previo acerca del ajuste del modelo, éste es un aspecto que debe ser tratado en otro momento; sugiero al lector que se introduzca en este aspecto una vez tenga conocimientos sobre el análisis del ajuste del modelo. No obstante hay que advertir que su uso nunca puede sustituir a la valoración juiciosa de los modelos que van surgiendo de forma seriada en cada paso y del modelo final. No hacerlo así puede llevarnos a dar por bueno un modelo surgido de forma automática (por criterios prestablecidos por el paquete estadístico muchas veces mal conocidos por el usuario del software), con escaso valor clínico.
Cada vez que nos encontremos ante un modelo de regresión logística (el inicial, cualquiera de los intermedios o el final), tendremos que contrastar su significación global. Esto se hace mediante las pruebas de ajuste global del modelo.
Una vez tenemos un modelo preliminar, podríamos incluir factores de interacción, es decir, estudiar cómo la asociación de dos o más covariables puede influir en la variable dependiente. Existen estrategias de desarrollo de modelos de regresión por las que se recomienda la inclusión en el modelo inicial de todas las covariables necesarias más las interacciones de las mismas, o por lo menos, las interacciones de primer orden (tomadas las covariables dos a dos), a los que se les llama modelos saturados. Interacciones de mayor orden suelen ser de difícil interpretación. En cualquier caso siempre hay que tener presente las limitaciones de tamaño muestral (que luego veremos), y de interpretación desde el punto de vista clínico (no se deberían incluir interacciones de significado incierto).
Otra estrategia en el desarrollo del modelo final es el diseño y ajuste de un modelo final preliminar sin interacciones, con el que luego ensayaríamos la inclusión, uno por uno, de términos de interacción que pudieran tener traducción clínica [Hosmer y Lemeshow], y valorar su significación respecto del modelo previo sin interacciones.
Una vez decidido que incluiremos un factor de interacción, tendremos en cuenta que siempre deberán estar incluidas también de forma aislada en el modelo las covariables que componen la interacción (principio jerárquico): si nuestra interacción es "HTA-diabetes", en el modelo se encontrarán como covariables HTA y diabetes:
logit = β0 + β1 · HTA + β2 · DM + β3 · HTA · DM + ...
Ecuación 8
Por otra parte, y en relación con la inclusión de interacciones, hay que tener en cuenta que la inclusión de las mismas puede generar multicolinealidad, tanto más probable cuanto mayor sea el número de interacciones.
Siempre debe considerarse la suficiencia del tamaño muestral para el número de covariables que se desea incluir en el modelo: modelos excesivamente grandes para muestras con tamaños muestrales relativamente pequeños implicarán errores estándar grandes o coeficientes estimados falsamente muy elevados (sobreajuste). En general se recomienda que por cada covariable se cuente con un mínimo de 10 individuos por cada evento de la variable dependiente con menor representación [Peduzzi]. Un ejemplo: si la variable dependiente Y es muerte y en nuestros datos hay 120 sujetos vivos y 36 sujetos muertos, el evento de Y menos representado es muerte, con 36 sujetos; de esta forma nuestro modelo no debería contener más de 36/10 ~ 3 covariables.
Lo anterior es válido siempre que se trate de covariables cuantitativas o cualitativas con distribuciones bien equilibradas. La situación se complica si una o más de las covariables cualitativas no tiene una distribución equilibrada (uno de sus dos valores tiene una mínima representación); en ese caso se recomienda que en su tabla de contingencia respecto a la variable dependiente, en cada celda haya un mínimo de 10 observaciones. Con un ejemplo lo entenderemos: En el siguiente ejemplo se debería disponer de suficiente tamaño muestral como para que en cada celda haya 10 ó más sujetos (es decir, que tanto a, b, c como d sean mayores de 10).
| Y|X | x=0 | x=1 |
|---|---|---|
| y=0 | a | b |
| y=1 | c | d |
Se dice que existe multicolinealidad cuando dos o más de las covariables del modelo mantienen una relación lineal.
Cuando la colinealidad es perfecta, es decir, cuando una covariable puede determinarse según una ecuación lineal de una o más de las restantes covariables, es posible estimar un único coeficiente de todas las covariables implicadas. En estos casos debe eliminarse la covariable que actua como dependiente.
Normalmente lo que encontraremos será una multicolinealidad moderada, es decir, encontraremos una mínima correlación entre covariables. Si esta correlación fuera de mayor importancia, su efecto sería, como ya vimos anteriormente, el incremento exagerado de los errores estándar, y en ocasiones, del valor estimado para los coeficientes de regresión, lo que hace las estimaciones poco creibles.
Un primer paso para analizar este aspecto puede ser examinar la matriz de coeficientes de correlación entre las covariables. Coeficientes de correlación muy elevados nos llevarán a investigar en mayor profundidad. Sin embargo este método, bueno para detectar colinealidad entre dos covariables, puede llevarnos a no poder detectar multicolinealidad entre más de dos de ellas.
Existen otros procedimientos analíticos para detectar multicolinealidad. Podemos desentendernos por el momento de la variable dependiente y realizar sendos modelos en los que una de las covariables actuará como variable dependiente y las restantes covariables como variables independientes de áquella. A cada uno de estos modelos podemos calcularle su R2 (o dispersión total, medida de ajuste que veremos más adelante). Se denomina tolerancia al complementario de R2 (1-R2), y factor de inflación de la varianza (FIV) al inverso de la tolerancia (1/(1-R2)). Cuando existe estrecha relación entre covariables la tolerancia tiende a ser 0, y por tanto FIV tiende al infinito. Como regla general nos deberían preocupar tolerancias menores de 0.1 y FIV mayores de 10. SPSS nos ofrece la matriz de correlaciones, pero no aporta índices de multicolinealidad para la regresión logística.
La solución a la multicolinealidad no es fácil:
Función de verosimilitud: Sabemos que cualquier variable dependiente de otra u otras variables, toma valores según los valores de las variables de las que depende. Por otra parte, esa variable dependiente irá tomando valores siguiendo o describiendo una determinada distribución de frecuencias; es decir, tomen los valores que tomen las variables independientes, si el experimento se repite múltiples veces, la variable dependiente tomará para esos valores de las independientes un determinado valor, y la probabilidad de ocurrencia de dicho valor vendrá dado por una distribución de frecuencias concreta: una distribución normal, una distribución binomial, una distribución hipergeoétrica, etc. En el caso de una variable dependiente dicotómica (como el caso que nos ocupa), la distribución de frecuencias que seguirá será la binomial, que depende de la tasa de éxitos (X sujetos de un total de N, que sería el elemento variable), para un determinado tamaño muestral (N) y probabilidad (p) de ocurrencia del evento valorado por la variable dependiente (parámetros constantes). La función de densidad de esta distribución de frecuencias vendrá dada por la siguiente expresión
Si en la expresión anterior introducimos los datos concretos de nuestra muestra de N sujetos (es decir, convertimos el elemento variable X en parámetro), y hacemos depender el resultado de la función de densidad del parámetro "probabilidad de ocurrencia" (p, que de esta forma lo convertimos en variable), estamos generando su función de verosimilitud, ƒ(p|x) (función dependiente de p dado el valor muestral de x) o L(p) (L del término likelihood), que nos ofrece como resultados las probabilidades de la función de densidad ajustada a nuestros datos.
Se deduce que, para una muestra concreta, esa probabilidad será diferente según qué valores tome el parámetro "probabilidad de ocurrencia".
Se demuestra que la mejor estimación de dicho parámetro (ê) es aquel valor que maximice esta función de verosimilitud, ya que son estimadores consistentes (conforme crece el tamaño muestral, la estimación se aproxima al parámetro desconocido), suficientes (aprovechan la información de toda la muestra), asintóticamente normales y asintóticamente eficientes (con mínima varianza), si bien no siempre son insesgados (no siempre la media de las estimaciones para diferentes muestras tenderá hacia el parámetro desconocido).
Método de Newton-Raphson: Se trata de un método iterativo, empleado en diversos problemas matemáticos, como en la determinación de las raíces de ecuaciones, y en nuestro caso, en la estimación de los coeficientes de regresión β por el procedimiento de máxima verosimilitud.
Por facilidad de cálculo toda la formulación la expresaremos en forma de matrices. Las particularidades del cálculo matricial escapan del ámbito de este documento. Tengamos en mente nuestra base de datos (una tabla con filas y columnas). Contaremos con:
El proceso se inicia construyendo la función de verosimilitud (likelihood function) de la ecuación de regresión logística,
... o mejor, su transformación logarítmica (LL, log likelihood):
Donde pi es la probabilidad de ocurrencia de y=1 con los valores muestrales de las covariables X ∈ {x1, ... , xm+1} , para el sujeto i ∈ {1, ... , N}. Llamamos devianza al valor -2·LL(β), y mide en qué grado el modelo se ajusta a los datos; cuanto menor sea su valor, mejor es el ajuste.
Se trata de conocer aquellos valores de β que hacen máxima la función de verosimilitud (o su logaritmo). Se sabe que si igualamos a 0 la derivada parcial de una función respecto a un parámetro, el resultado es unos valores de dicho parámetro que hacen llevar a la función a un valor máximo o un valor mínimo (un punto de inflexión de la curva). Para confirmar que se trata de un máximo y no de un mínimo, la segunda derivada de la función respecto a dicho parámetro debe ser menor de 0.
La primera derivada de LL(β) respecto de β (llamada función score) en su forma matricial es
siendo p una matriz de N filas y 1 columna que contiene las probabilidades de cada individuo de que tengan su correspondiente evento yi.
La segunda derivada, llamada matriz informativa o hessiana, es
siendo W una matriz diagonal (una matriz cuadrada en la que todos sus elementos son 0 excepto su diagonal) de N filas y N columnas, en la que los elementos de su diagonal vienen dados por los respectivos productos
de manera que W queda de la forma siguiente
y para cada fila su pi es
Una vez tenemos todos los elementos necesarios, pasamos a explicar como tal el método iterativo para la determinación de los coeficientes de regresión.
Entender esta formulación y el algoritmo de las iteraciones puede ser de gran utilidad, pues con conocimientos básicos de programación nos facilita el desarrollo de rutinas propias, por ejemplo en VisualBasic® dentro de una base de datos de Access®, que pueden librarnos de la dependencia de costosos paquetes estadísticos.
Odds ratio: Es un cociente de proporciones de enfermos por cada sano entre el grupo con un factor de riesgo y el grupo sin dicho factor de riesgo. Supongamos el siguiente ejemplo:
| Enfermedad | ||||
|---|---|---|---|---|
| sí | no | |||
| Factor de riesgo |
sí | 20 | 80 | 100 |
| no | 30 | 270 | 300 | |
| 50 | 350 | 400 | ||
En este caso, entre los que tienen el factor de riesgo hay 20 enfermos por cada 80 sanos (0.25), y entre los que no tienen el factor de riesgo hay 30 enfermos por cada 270 sanos (0.11), por lo que las personas con el factor de riesgo tienen un riesgo de enfermar 2.25 veces superior (0.25/0.11) que las personas sin el factor de riesgo.
Principio jerárquico: Siempre que se incluya en el modelo un término de interacción, el modelo debe incluir también todos los términos de orden inferior, y si el término de interacción resultase significativo y permaneciese en el modelo, también deberían permanecer lor términos de orden inferior, aunque no se lograra demostrar significación para ellos.
Modelo con interacción de primer orden:
y=a+b1x1+b2x2+b3x1x2
Modelo con interacción de segundo orden:
y=a+b1x1+b2x2+b3x3+b4x1x2+b5x1x3+b6x2x3+b7x1x2x3
Principio de parsimonia: En igualdad de condiciones la solución más sencilla que explique completamente un problema es probablemente la correcta [Guillermo de Ockham]. Según este principio, cuando más de un modelo se ajuste a nuestras observaciones, siempre deberíamos quedarnos con el modelo más simple que explique nuestras observaciones con un grado adecuado de precisión.
Tabla de contingencia: Una tabla de contingencia es una tabla de doble entrada para la comparación de variables cualitativas. En el texto, al tratarse de variables dicotómicas, se refiere a una tabla 2 x 2 como en la figura.
| Y|X | x=0 | x=1 |
|---|---|---|
| y=0 | a | b |
| y=1 | c | d |
Fco. Javier Caballero Granado
Servicio de Medicina Interna
Hospital Punta de Europa
Algeciras, Cádiz
