Desempaquetando la Cox Un Oscuro Secreto Oculto de la Regresión de Cox
Unveiling the Hidden Dark Secret of Cox Regression
¿Por qué los predictores perfectos dan como resultado un valor p de 0.93?
Profundizando en los predictores perfectos
Si ha seguido mis publicaciones de blog anteriores, recordará que la regresión logística encuentra un problema al tratar de ajustar datos perfectamente separados, lo que lleva a una razón de probabilidades infinita. En la regresión de Cox, donde el riesgo reemplaza las probabilidades, puede preguntarse si surge un problema similar con los predictores perfectos. Ocurre, pero a diferencia de la regresión logística, es mucho menos evidente cómo ocurre aquí e incluso lo que constituye “predictores perfectos”. Como se hará más claro más adelante, los predictores perfectos se definen como predictores x cuyas posiciones coinciden exactamente con las posiciones de los tiempos de eventos (su correlación de Spearman es uno).
Anteriormente, en “Desempaquetar a Cox”:
Desempaquetar a Cox: Guía intuitiva para las regresiones de Cox
¿Cómo predicen los riesgos y las estimaciones de máxima verosimilitud la clasificación de eventos?
towardsdatascience.com
Explicamos la estimación de máxima verosimilitud e introdujimos un conjunto de datos ficticio con cinco individuos, donde un solo predictor, x, representaba la dosis de un medicamento para prolongar la vida. Para convertir x en un predictor perfecto de los tiempos de eventos, aquí intercambiamos los tiempos de eventos de los sujetos C y D:
import numpy as npimport pandas as pdimport plotnine as p9from cox.plots import ( plot_subject_event_times, animate_subject_event_times_and_mark_at_risk, plot_cost_vs_beta,)perfect_df = pd.DataFrame({ 'subject': ['A', 'B', 'C', 'D', 'E'], 'time': [1, 3, 4, 5, 6], 'event': [1, 1, 1, 1, 0], 'x': [-1.7, -0.4, 0.0, 0.9, 1.2],})plot_subject_event_times(perfect_df, color_map='x')
Para entender por qué estos “predictores perfectos” pueden ser problemáticos, retomemos donde lo dejamos y verifiquemos el costo del logaritmo negativo de la verosimilitud trazado frente a β:
- ¿Pandas 2.0 ¿Un cambio de juego para los científicos de datos?
- Evolución del panorama de los datos.
- Inteligencia Artificial Verde Métodos y Soluciones para Mejorar la ...
negloglik_sweep_betas_perfect_df = neg_log_likelihood_all_subjects_sweep_betas( perfect_df, betas=np.arange(-5, 5, 0.1))plot_cost_vs_beta(negloglik_sweep_betas_perfect_df, width=0.1)
De inmediato se puede ver que ya no hay un valor mínimo de β: si usamos valores negativos de β muy grandes, terminamos con ajustes de verosimilitud que son casi perfectos para todos los eventos.
Ahora, profundicemos en la matemática detrás de esto y echemos un vistazo a la verosimilitud del evento A. Analizaremos cómo cambian el numerador y el denominador a medida que ajustamos β:
Cuando β es alto o un número positivo grande, el último término en el denominador (con el valor x más grande de 1.2), que representa el riesgo del sujeto E, domina todo el denominador y se vuelve extremadamente grande. Por lo tanto, la verosimilitud se vuelve pequeña y se acerca a cero:
Esto resulta en un gran logaritmo negativo de verosimilitud. Lo mismo ocurre para cada verosimilitud individual porque el último riesgo de sujeto E, siempre superará cualquier riesgo en el numerador. Como resultado, el logaritmo negativo de verosimilitud aumenta para los sujetos A a D. En este caso, cuando tenemos un alto β, hace bajar todas las verosimilitudes, dando como resultado pobres ajustes para todos los eventos.
Ahora, cuando β es bajo o un número negativo grande, el primer término en el denominador, que representa el riesgo de sujeto A, domina ya que tiene el valor x más bajo. Como el mismo riesgo de sujeto A también aparece en el numerador, la verosimilitud L(A) puede ser arbitrariamente cercana a 1 haciendo β cada vez más negativo, creando así un ajuste casi perfecto:
Lo mismo ocurre para todas las otras verosimilitudes individuales: β negativos ahora aumentan las verosimilitudes de todos los eventos al mismo tiempo. Básicamente, tener un β negativo no tiene desventajas. Al mismo tiempo, ciertos riesgos individuales aumentan (sujetos A y B con x negativo), algunos permanecen iguales (sujeto C con x = 0), y otros disminuyen (sujeto D con x positivo). Pero recuerde, lo que realmente importa aquí son las relaciones entre riesgos. Podemos verificar esta matemática aproximada de manera visual representando los riesgos individuales:
def plot_likelihoods(df, ylim=[-20, 20]): betas = np.arange(ylim[0], ylim[1], 0.5) subjects = df.query("event == 1")['subject'].tolist() likelihoods_per_subject = [] for subject in subjects: likelihoods = [ np.exp(log_likelihood(df, subject, beta)) for beta in betas ] likelihoods_per_subject.append( pd.DataFrame({ 'beta': betas, 'likelihood': likelihoods, 'subject': [subject] * len(betas), }) ) lik_df = pd.concat(likelihoods_per_subject) return ( p9.ggplot(lik_df, p9.aes('beta', 'likelihood', color='subject')) + p9.geom_line(size=2) + p9.theme_classic() ) plot_likelihoods(perfect_df)
La forma en que se combinan las verosimilitudes, como una relación entre un riesgo y la suma de riesgos de todos los sujetos que aún están en riesgo, significa que los valores β negativos hacen un ajuste perfecto para la verosimilitud de cada sujeto cuyo rango de tiempo de evento es mayor o igual que el rango del predictor! Como nota adicional, si x tuviera una correlación negativa de Spearman perfecta con los tiempos de evento, las cosas se invertirían: los valores β arbitrariamente positivos nos darían ajustes arbitrariamente buenos.
Rangos de predictor y tiempo desalineados
Podemos ver esto y mostrar lo que sucede cuando los rangos de tiempo de evento y los rangos de predictor no coinciden usando otro ejemplo ficticio:
sample_df = pd.DataFrame({ 'subject': ['A', 'B', 'C', 'D', 'E', 'F', 'G', 'H'], 'x': [-1.7, -0.4, 0.0, 0.5, 0.9, 1.2, 1.3, 1.45], 'time': [1, 2, 4, 3, 5, 7, 6, 8], 'rank_x': [1, 2, 3, 4, 5, 6, 7, 8], 'event': [1, 1, 1, 1, 1, 1, 1, 0],})sample_df
En este ejemplo particular, la columna time varía de 1 a 8, donde cada valor representa su propio rango. También tenemos una columna x_rank, que clasifica los predictores x. Ahora, aquí está la observación clave: para los sujetos D y G, su x_rank es en realidad más alto que su rango de time. Como resultado, las verosimilitudes de D y G no experimentarán el efecto de cancelación entre el numerador y el denominador cuando tenemos valores β negativos grandes:
Ahora sus probabilidades son máximas en algunos valores finitos intermedios de β . Echemos un vistazo a una gráfica de probabilidades individuales para ver esto:
plot_likelihoods(sample_df)
Estos rangos “desalineados” entre un tiempo y un predictor juegan un papel crucial: detienen todas las probablidades de colapsar esencialmente en una cuando tenemos valores de β significativamente negativos.
En conclusión, en la regresión de Cox, para obtener un coeficiente β finito para un predictor x , necesitamos al menos una instancia donde el rango del predictor x sea menor que el rango del tiempo del evento.
Perfecto es enemigo de lo bueno (valor-p)
Entonces, ¿cómo se comportan estos predictores perfectos en escenarios de la vida real? Bueno, para descubrirlo, volvamos una vez más a la biblioteca de lifelines para investigar:
from lifelines import CoxPHFitterperfect_cox_model = CoxPHFitter()perfect_cox_model.fit( perfect_df, duration_col='time', event_col='event', formula='x')perfect_cox_model.print_summary()
#> /.../coxph_fitter.py:1586: ConvergenceWarning:#> The log-likelihood is getting suspiciously close to 0 and the delta is still large.#> There may be complete separation in the dataset.#> This may result in incorrect inference of coefficients.#> See https://stats.stackexchange.com/q/11109/11867 for more.#> /.../__init__.py:1165: ConvergenceWarning:#> Column x has high sample correlation with the duration column.#> This may harm convergence.#> This could be a form of 'complete separation'.#> See https://stats.stackexchange.com/questions/11109/how-to-deal-with-perfect-separation-in-logistic-regression#> /.../coxph_fitter.py:1611:#> ConvergenceWarning: Newton-Rhaphson failed to converge sufficiently.#> Please see the following tips in the lifelines documentation:#> https://lifelines.readthedocs.io/en/latest/Examples.html#problems-with-convergence-in-the-cox-proportional-hazard-model
Al igual que en la regresión logística, encontramos advertencias de convergencia y obtenemos intervalos de confianza extremadamente amplios para nuestro coeficiente predictor β . Como resultado, obtenemos un valor-p de 0.93 !
Si simplemente filtramos modelos basados en valores-p sin tener en cuenta este problema o realizar una investigación adicional, pasaríamos por alto estos predictores perfectos.
Para abordar este problema de convergencia, la documentación de la biblioteca de lifelines y algunas útiles preguntas de StackOverflow sugieren una solución potencial: incorporar un término de regularización en la función de costo. Este término aumenta efectivamente el costo para valores grandes del coeficiente, y puede activar la regularización L2 estableciendo el argumento penalizer en un valor mayor que cero:
perfect_pen_cox_model = CoxPHFitter(penalizer=0.01, l1_ratio=0)perfect_pen_cox_model.fit(perfect_df, duration_col='time', event_col='event', formula='x')perfect_pen_cox_model.print_summary()
Este enfoque soluciona la advertencia de convergencia, pero realmente no reduce significativamente ese molesto valor-p. Incluso con este truco de regularización, el valor-p para un predictor perfecto todavía ronda un valor algo grande de 0.11.
El tiempo es relativo: solo importan los rangos
Por último, verificaremos que los valores absolutos de los tiempos de evento no tienen impacto en un ajuste de regresión de Cox, usando nuestro ejemplo anterior. Para hacer esto, introduciremos una nueva columna llamada time2, que contendrá números aleatorios en el mismo orden que la columna time:
sample_df = pd.DataFrame({ 'subject': ['A', 'B', 'C', 'D', 'E', 'F', 'G', 'H'], 'x': [-1.7, -0.4, 0.0, 0.5, 0.9, 1.2, 1.3, 1.45], 'time': [1, 2, 4, 3, 5, 7, 6, 8], 'rank_x': [1, 2, 3, 4, 5, 6, 7, 8], 'event': [1, 1, 1, 1, 1, 1, 1, 0],}).sort_values('time')np.random.seed(42)sample_df['time2'] = sorted(np.random.randint(low=-42, high=888, size=8))sample_df
Sus ajustes son de hecho idénticos:
sample_cox_model = CoxPHFitter()sample_cox_model.fit( sample_df, duration_col='time', event_col='event', formula='x')sample_cox_model.print_summary()
sample_cox_model = CoxPHFitter()sample_cox_model.fit( sample_df, duration_col='time2', event_col='event', formula='x')sample_cox_model.print_summary()
Conclusión
¿Qué aprendimos de todo esto?
- Los predictores perfectos en los modelos de supervivencia son aquellos cuyas clasificaciones coinciden perfectamente con las clasificaciones de los tiempos de evento.
- La regresión de Cox no puede ajustar estos predictores perfectos con un coeficiente finito β, lo que conduce a intervalos de confianza amplios y valores de p grandes.
- Los valores reales de los tiempos de evento no importan realmente, todo se trata de sus clasificaciones.
- Cuando las clasificaciones de los tiempos de evento y los predictores no se alinean, no obtenemos ese efecto de cancelación útil para valores grandes de β en las verosimilitudes. Por lo tanto, necesitamos al menos un caso en el que las clasificaciones no coincidan para tener un ajuste con un coeficiente finito.
- Incluso si intentamos algunas técnicas de regularización elegantes, los predictores perfectos todavía pueden darnos esos intervalos de confianza molesto amplios y valores de p altos en situaciones de la vida real.
- Al igual que en la regresión logística, si realmente no nos importan esos valores de p, el uso de un método de regularización todavía puede proporcionarnos un ajuste de modelo útil que acierte en la predicción.
Si desea ejecutar el código usted mismo, no dude en utilizar los cuadernos de IPython en mi Github: https://github.com/igor-sb/blog/blob/main/posts/cox_perfect.ipynb
¡Hasta la próxima publicación! 👋
