Introducción

En el análisis de supervivencia se han desarrollado dos enfoques principales: los métodos no paramétricos y los semiparamétricos. En el primer caso, el estimador de Kaplan–Meier permite estimar la función de supervivencia sin asumir una forma funcional para el tiempo hasta el evento. En el segundo, el modelo de riesgos proporcionales de Cox incorpora covariables sin especificar la función de riesgo basal.

Estos enfoques se utilizan ampliamente en investigación biomédica debido a su baja especificación estructural. Sin embargo, al no modelar explícitamente una distribución paramétrica para el tiempo de supervivencia, su inferencia suele limitarse al rango observado de los datos, restringiendo la extrapolación y la predicción a largo plazo.

En este contexto, los modelos paramétricos asumen que el tiempo de supervivencia sigue una distribución perteneciente a una familia conocida, caracterizada por una forma funcional específica y parámetros desconocidos a estimar. Esta especificación permite describir completamente el comportamiento temporal del evento y obtener inferencias basadas en toda la distribución de supervivencia.

En contraste con el modelo de Cox, los modelos paramétricos determinan completamente la distribución del tiempo de supervivencia mediante familias como la exponencial, Weibull, log-normal, log-logística o gamma generalizada.

En este documento se revisan los principales modelos paramétricos de supervivencia y sus supuestos. En particular, se analiza el enfoque de tiempo de fallo acelerado (AFT) y su relación con el modelo de riesgos proporcionales (PH), dos formulaciones que pueden ser equivalentes en ciertos casos pero que difieren en su interpretación.

Finalmente, se examinan los modelos paramétricos más utilizados en la práctica, junto con sus supuestos, propiedades y criterios de selección basados en ajuste, plausibilidad teórica y capacidad predictiva.

Riesgos proporcionales (PH) vs tiempo de falla acelerado (AFT)

El modelo de riesgos proporcionales (PH) se expresa como:

\[h(t\mid X)=h_0(t)\exp(\beta^{\top}X)\]

donde las covariables actúan multiplicativamente sobre la función de riesgo. En este enfoque, la medida de asociación corresponde al hazard ratio (HR).

Por su parte, el modelo de tiempo de falla acelerado (Accelerated Failure Time, AFT) asume que las covariables modifican directamente la escala temporal del proceso de supervivencia.

Siguiendo la notación de Kleinbaum:

\[\log(T)=\mu+\beta^{\top}X+\sigma W\]

donde:

  • \(T\) representa el tiempo hasta el evento,
  • \(\mu\) es el intercepto,
  • \(\beta\) corresponde a los coeficientes de regresión,
  • \(\sigma\) es un parámetro de escala,
  • y \(W\) es un término aleatorio cuya distribución determina el modelo paramétrico.

Según la distribución asumida para \(W\):

  • distribución de valor extremo → modelo Weibull,
  • distribución normal → modelo log-normal,
  • distribución logística → modelo log-logístico.

De forma equivalente:

\[T=\exp(\mu+\beta^{\top}X)\exp(\sigma W)\]

lo que muestra que las covariables aceleran o desaceleran el tiempo hasta el evento.

En este contexto, los coeficientes se interpretan mediante razones de tiempo (time ratios, TR) o factores de aceleración (acceleration factors, AF). Valores positivos de \(\beta\) indican una prolongación del tiempo de supervivencia, mientras que valores negativos indican una reducción.

Conceptualmente, el modelo PH describe efectos relativos sobre el riesgo, mientras que el modelo AFT describe efectos relativos sobre el tiempo de supervivencia.

Relación entre la Función de Densidad de Probabilidad con la Función de Riesgo y con la de Supervivencia

Para los modelos paramétricos de supervivencia, se asume que el tiempo sigue una distribución cuya función de densidad de probabilidad \(f(t)\) puede expresarse en términos de parámetros desconocidos. Una vez que se especifica una función de densidad de probabilidad para el tiempo de supervivencia, se pueden determinar las funciones de supervivencia y de riesgo correspondientes.

Como vimos previamente, la función de supervivencia \(S(t)\) se puede obtener a través de la función de densidad y la función de riesgo se puede encontrar a través de la función de supervivencia. Finalmente, la función de densidad de probabilidad puede expresarse como el producto de la función de riesgo y la función de supervivencia: \(f(t) = h(t) \cdot S(t)\).

Recordemos el gráfico:

  • \(S(t)\) es la función de supervivencia
  • \(f(t)\) es la función de densidad del tiempo hasta el evento
  • \(h(t)\) es la función de riesgo instantáneo (hazard)
  • \(H(t)\) es la función de riesgo acumulado

Distribuciones paramétricas principales

Una vez especificada la estructura del modelo, es necesario definir la distribución del tiempo de supervivencia, ya que esta determina la forma de la función de riesgo en el tiempo.

En los modelos paramétricos de supervivencia, las distribuciones más utilizadas son la exponencial, Weibull, log-logística, log-normal y gamma generalizada. Cada una de ellas impone diferentes comportamientos sobre la función de riesgo, permitiendo representar patrones temporales diversos.

  • Distribución exponencial: asume un riesgo constante en el tiempo.

  • Distribución Weibull: permite riesgos monótonos crecientes o decrecientes según el parámetro de forma \(p\). Cuando \(p=1\), se reduce al modelo exponencial.

  • Distribución log-logística: admite funciones de riesgo no monótonas, con incrementos iniciales seguidos de descensos.

  • Distribución log-normal: genera riesgos unimodales, adecuados cuando el riesgo alcanza un máximo en tiempos intermedios.

  • Distribución gamma generalizada: constituye una familia flexible que incluye varios modelos clásicos como casos particulares.

Estas distribuciones proporcionan un marco flexible para describir distintos patrones de supervivencia y riesgo a lo largo del tiempo. La elección del modelo paramétrico más apropiado debe considerar tanto el contexto clínico y biológico del fenómeno estudiado como el grado de ajuste que la distribución presenta frente a los datos empíricos.

Familias paramétricas comunes en supervivencia

A continuación se presenta un resumen de las principales distribuciones paramétricas utilizadas en análisis de supervivencia, junto con sus funciones de supervivencia, riesgo y características principales. Cada una de ellas se presentará en detalle, más adelante.

Distribución Función de supervivencia \[S(t)\] Función de riesgo \[h(t)\] Riesgo acumulado \[H(t)\] Característica principal
Exponencial \[S(t)=\exp(-\lambda t)\] \[h(t)=\lambda\] \[H(t)=\lambda t\] Riesgo constante en el tiempo
Weibull \[S(t)=\exp(-\lambda t^p)\] \[h(t)=\lambda p t^{p-1}\] \[H(t)=\lambda t^p\] Riesgo monótono (creciente o decreciente)
Log-logística \[S(t)=\frac{1}{1+\lambda t^p}\] \[h(t)=\frac{\lambda p t^{p-1}}{1+\lambda t^p}\] \[H(t)=\log(1+\lambda t^p)\] Riesgo no monótono
Log-normal \[S(t)=1-\Phi\left(\frac{\log(t)-\mu}{\sigma}\right)\] \(f(t)/S(t)\)* \(-\log S(t)\)* Riesgo unimodal
Gamma generalizada \(1-F(t;\lambda,p,k)\)** \(f(t)/S(t)\)* \(-\log S(t)\)* Alta flexibilidad funcional

Nota. Para las distribuciones Log-normal y Gamma generalizada, los valores de \(h(t)\) y \(H(t)\) no poseen formas cerradas simples. La función de riesgo se obtiene como \(h(t) = f(t)/S(t)\), donde \(f(t)\) es la función de densidad correspondiente y \(S(t)\) es la función de supervivencia. El riesgo acumulado se calcula como \(H(t) = -\log S(t)\) mediante integración numérica. En la distribución log-normal, \(\Phi\) es la función de distribución acumulada de la normal estándar. En la distribución gamma generalizada, \(F(t;\lambda,p,k)\) es la función de distribución acumulada parametrizada por \(\lambda\), \(p\) y \(k\).

Principales modelos paramétricos de supervivencia

Modelo Exponencial

El modelo exponencial representa la estructura más simple y fundamental en el análisis de supervivencia ya que asume un riesgo constante en el tiempo. Esta definido por la siguiente función de riesgo instantáneo:

\[h(t) = \lambda\]

con \(\lambda > 0\). Bajo este supuesto, la probabilidad de que el evento ocurra en el siguiente instante es siempre la misma, independientemente de cuánto tiempo haya sobrevivido el sujeto. Esta simplicidad se refleja en sus funciones derivadas.

El riesgo acumulado crece de forma lineal:

\[H(t) = \lambda t\],

mientras que la función de supervivencia decae exponencialmente:

\[S(t) = \exp(-\lambda t)\].

En términos de sus parámetros, la esperanza o media de vida se define como \(E[T] = 1/\lambda\), con una varianza de \(1/\lambda^2\). Un detalle matemático relevante es que su desviación estándar es igual a su media, lo que implica un coeficiente de variación de 1 y sugiere que los tiempos de supervivencia presentan una variabilidad intrínseca muy alta.

Una consecuencia única de este modelo es la lo que se conoce como falta de memoria, la cual puede sr expresada formalmente como:

\[P(T > t + x \mid T > t) = P(T > x)\]

Esto implica que la probabilidad de ocurrencia de un evento futuro no depende del tiempo previamente transcurrido, es decir, el proceso no presenta memoria del pasado. Aunque este supuesto suele ser poco realista en contextos biomédicos, donde por ejemplo el envejecimiento, la progresión de la enfermedad o el daño acumulado si modifican el riesgo con el tiempo, el modelo exponencial resulta útil para describir eventos agudos o impredecibles, así como para servir como aproximación simplificada en estudios de corta duración donde las variaciones temporales del riesgo son mínimas o despreciables.

Ejemplo con modelo exponencial

Para ilustrar su implementación, consideremos el estudio de los 42 pacientes con leucemia dividido en dos grupos:

  • tratamiento experimental (\(\text{TRT}=1\)) y
  • placebo (\(\text{TRT}=0\)).

Al incorporar esta covariable, reparametrizamos la tasa de riesgo mediante una función exponencial para asegurar valores positivos:

\[\lambda(X) = \exp(\beta_0 + \beta_1 \text{TRT})\]

Bajo esta formulación, el riesgo basal (grupo placebo) queda definido por \(h(t) = \exp(\beta_0)\), mientras que el riesgo del grupo tratado se expresa como \(h(t) = \exp(\beta_0 + \beta_1)\). Al calcular la relación entre ambos, obtenemos el HR:

\[\text{HR} = \frac{h(t \mid \text{TRT}=1)}{h(t \mid \text{TRT}=0)} = \exp(\beta_1)\]

Podemos ver como en el modelo exponencial, tanto los riesgos individuales de cada grupo como su razón permanecen constantes en el tiempo. Esta es una forma restrictiva pero poderosa de los riesgos proporcionales.

Resultados obtenidos previamente con el modelo de Cox:

## Call:
## coxph(formula = Surv(t, d) ~ trat_f, data = base)
## 
##   n= 42, number of events= 30 
## 
##                      coef exp(coef) se(coef)      z Pr(>|z|)    
## trat_fTratamiento -1.5721    0.2076   0.4124 -3.812 0.000138 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                   exp(coef) exp(-coef) lower .95 upper .95
## trat_fTratamiento    0.2076      4.817   0.09251    0.4659
## 
## Concordance= 0.69  (se = 0.041 )
## Likelihood ratio test= 16.35  on 1 df,   p=5e-05
## Wald test            = 14.53  on 1 df,   p=1e-04
## Score (logrank) test = 17.25  on 1 df,   p=3e-05

El hazard ratio estimado por Cox es \(\text{HR} = \exp(-1.5721) \approx 0.208\). Esto indica que el tratamiento reduce el riesgo instantáneo en un 79% respecto al placebo.

Resultados aplicando el modelo exponencial (forma AFT):

## 
## Call:
## survreg(formula = Surv(t, d) ~ trat_f, data = base, dist = "exponential")
##                   Value Std. Error    z       p
## (Intercept)       2.159      0.218 9.90 < 2e-16
## trat_fTratamiento 1.527      0.398 3.83 0.00013
## 
## Scale fixed at 1 
## 
## Exponential distribution
## Loglik(model)= -108.5   Loglik(intercept only)= -116.8
##  Chisq= 16.49 on 1 degrees of freedom, p= 4.9e-05 
## Number of Newton-Raphson Iterations: 4 
## n= 42

Observamos que \(\hat{\beta}_1 = 1.527\). En el contexto de Tiempo de Falla Acelerado (AFT), el efecto se interpreta mediante el factor de aceleración (o Time Ratio (TR)):

\[\text{TR} = \exp(1.527) \approx 4.60\]

Esto sugiere que el tratamiento prolonga en 4.6 veces el tiempo libre de progresión en comparación con el grupo control, lo cual es consistente con un efecto protector.

Nóte que el signo del coeficiente es opuesto al de Cox: mientras que Cox modela el riesgo (signo negativo = menos riesgo), ATP modela el tiempo (signo positivo = más tiempo), por lo que en este modelo se cumple que \(\beta_{\text{PH}} \approx -\beta_{\text{AFT}}\).

Relación entre la forma PH y AFT en el modelo exponencial

La distribución exponencial posee la propiedad particular de permitir una representación equivalente tanto bajo el enfoque PH como bajo el enfoque de AFT.

En el enfoque de los PH, el parámetro actúa sobre el riesgo así:

\[h(t \mid x) = h_0(t)\exp(\beta_{\text{PH}} x)\]

Y en el enfoque AFT, el parámetro actúa sobre la escala del tiempo así:

\[T = T_0 \exp(\beta_{\text{AFT}} x)\]

Aunque este modelo es un excelente punto de partida conceptual, el supuesto de riesgo constante rara vez se ajusta a la realidad clínica a largo plazo. Para capturar comportamientos donde el riesgo aumenta o disminuye con el tiempo, es necesario recurrir a modelos más flexibles como la distribución Weibull.

Modelo Weibull

El modelo Weibull es una de las distribuciones paramétricas más utilizadas en análisis de supervivencia debido a su gran flexibilidad. A diferencia del modelo exponencial, que asume un riesgo constante en el tiempo, el modelo Weibull incorpora un parámetro de forma \(p\) que permite que la función de riesgo cambie de manera monótona a lo largo del tiempo.

Su estructura matemática se define mediante las siguientes funciones fundamentales:

  • Función de riesgo: \(h(t) = \lambda p t^{p-1}, \quad \lambda > 0, \, p > 0\)

  • Supervivencia: \(S(t) = \exp(-\lambda t^p)\)

  • Media (\(E[T]\)): \(\frac{1}{\lambda^{1/p}} \Gamma\left(1 + \frac{1}{p}\right)\)

El comportamiento del riesgo depende enteramente del valor de \(p\):

  1. \(p = 1\) (Riesgo constante): La distribución se reduce al modelo exponencial (ej. Riesgo de infección hospitalaria en pacientes con exposición constante).

  2. \(p > 1\) (Riesgo creciente): Modela procesos de envejecimiento o acumulación de daño (ej. mortalidad en pacientes con insuficiencia cardíaca crónica avanzada, enfermedades degenerativas o cáncer).

  3. \(p < 1\) (Riesgo decreciente): Describe situaciones de “selección de sobrevivientes” o eventos tempranos (ej. complicaciones agudas post-cirugía, riesgo de infección durante los primeros días de hospitalización).

Por esta razón, el modelo Weibull resulta especialmente útil en contextos clínicos y epidemiológicos donde el riesgo puede incrementarse por progresión de la enfermedad o disminuir después de una fase crítica del tratamiento.

Al igual que el modelo exponencial, el modelo Weibull puede formularse tanto bajo el enfoque de riesgos proporcionales como bajo el enfoque de tiempo de falla acelerado. Sin embargo, Weibull ofrece mayor flexibilidad al permitir riesgos crecientes o decrecientes a través de su parámetro de forma (p).

En este contexto, los coeficientes de ambas parametrizaciones se relacionan mediante:

\[\beta_{PH}=-p\alpha_{AFT}\]

donde (p) corresponde al parámetro de forma de la distribución Weibull.

Esta propiedad permite interpretar el efecto de las covariables tanto en términos de riesgo relativo como de aceleración o desaceleración del tiempo hasta el evento.

Ejemplo con distribución Weibull

Supongamos una Weibull con \(p=2\) (riesgo creciente) y un tratamiento que reduce el riesgo a la mitad (\(\text{HR} = 0.5\)). Bajo el enfoque de riesgos proporcionales, el tratamiento se modela como una reducción constante del riesgo en el tiempo. En cambio, bajo el enfoque de tiempo de falla acelerado, el mismo efecto se interpreta como una prolongación del tiempo hasta el evento.

Bajo el enfoque de riesgos proporcionales (PH), esto significa que el tratamiento reduce en un 50% el riesgo instantáneo del evento respecto al grupo control.

Sin embargo, debido a la dualidad PH/AFT del modelo Weibull, el mismo efecto puede interpretarse desde el enfoque de tiempo de falla acelerado (AFT). La relación entre ambas parametrizaciones viene dada por:

{PH}=-p{AFT}

Como:

\[ \beta_{PH}=\log(0.5)\approx -0.693 \]

y (p=2), entonces:

\[ \alpha_{AFT}=-\frac{\beta_{PH}}{p} =\frac{0.693}{2} \approx 0.347 \]

Por tanto, el factor de aceleración o razón de tiempo es:

\[ TR=\exp(0.347)\approx 1.41 \]

Esto indica que, bajo la formulación AFT, el tratamiento prolonga aproximadamente 1.4 veces el tiempo hasta el evento en comparación con el grupo control.

En consecuencia, ambos enfoques describen el mismo fenómeno clínico: el modelo PH lo expresa como una reducción del riesgo instantáneo, mientras que el modelo AFT lo interpreta como una prolongación del tiempo de supervivencia.

En el gráfico anterior se observa que las curvas obtenidas bajo las formulaciones PH y AFT se superponen completamente. Esto ilustra la equivalencia matemática de ambas parametrizaciones en el modelo Weibull.

La diferencia entre ambos enfoques no se refleja en las predicciones del modelo, sino en la interpretación de los coeficientes: PH expresa el efecto de las covariables sobre el riesgo instantáneo, mientras que AFT lo expresa sobre el tiempo hasta el evento.

Veamos ahora el ejemplo real con los datos de leucemia, aplicando el modelo Weibull:

## 
## Call:
## survreg(formula = Surv(t, d) ~ trat_f, data = base, dist = "weibull")
##                    Value Std. Error     z       p
## (Intercept)        2.248      0.166 13.55 < 2e-16
## trat_fTratamiento  1.267      0.311  4.08 4.5e-05
## Log(scale)        -0.312      0.147 -2.12   0.034
## 
## Scale= 0.732 
## 
## Weibull distribution
## Loglik(model)= -106.6   Loglik(intercept only)= -116.4
##  Chisq= 19.65 on 1 degrees of freedom, p= 9.3e-06 
## Number of Newton-Raphson Iterations: 5 
## n= 42

En la salida del modelo Weibull, el coeficiente estimado para el tratamiento corresponde a la parametrización AFT es:

\[ \hat{\alpha}_{AFT}=1.267 \]

Por lo tanto, la interpretación la debemos hacer calculando la razón de tiempo (TR), así:

\[ TR=\exp(1.267)\approx 3.55 \]

Podemos decir, que los pacientes del grupo tratamiento presentan un tiempo libre de progresión aproximadamente 3.6 veces mayor que los pacientes del grupo control, resultado similar, pero no igual al obtenido con el modelo exponencial (TR=4.6).

La diferencia se debe a que el modelo Weibull ajusta un riesgo creciente, lo que implica que el efecto del tratamiento se amplifica a medida que transcurre el tiempo, mientras que el modelo exponencial asume un efecto constante.

Nota: En R, para obtener el parámetro de forma, debe utilizarse:

\[p = \frac{1}{\text{scale}}\]

En R, el parámetro de forma para el modelo Weibull se obtiene a partir del parámetro scale:

\[ p=\frac{1}{0.732}\approx 1.37 \]

Dado que p>1, el riesgo aumenta progresivamente con el tiempo. En este estudio de leucemia, esto sugiere que los pacientes presentan una probabilidad instantánea creciente de recaída o muerte conforme transcurre el seguimiento.

La relación entre ambas parametrizaciones viene dada por:

{PH}=-p{AFT}

Sustituyendo los valores estimados:

\[ \beta_{PH}=-(1.37)(1.267)\approx -1.74 \]

Por consiguiente:

\[ HR=\exp(-1.74)\approx 0.18 \]

Esto indica que el tratamiento reduce el riesgo instantáneo aproximadamente en un 82% en comparación con el grupo control.

Como vemos, las interpretaciones bajo los enfoques PH y AFT son consistentes entre sí: mientras el modelo PH describe el efecto del tratamiento como una reducción del riesgo instantáneo, el modelo AFT lo expresa como una prolongación del tiempo hasta el evento.

Modelo Log-Normal

A diferencia de los modelos exponencial y Weibull, el modelo log-normal permite describir patrones de riesgo no monótonos. En este caso, el riesgo instantáneo puede aumentar inicialmente, alcanzar un máximo y posteriormente disminuir, generando una función de riesgo unimodal o con forma de “joroba”.

El modelo asume que el logaritmo del tiempo de supervivencia sigue una distribución normal:

\[\log(T)\sim N(\mu,\sigma^2)\]

donde () representa el parámetro de localización y () el parámetro de escala.

Equivalentemente:

\[T\sim \text{Lognormal}(\mu,\sigma^2)\]

Esta distribución produce tiempos de supervivencia positivos y asimétricos hacia la derecha, patrón frecuente en datos clínicos y epidemiológicos. Su flexibilidad resulta útil para modelar procesos con un período intermedio de mayor vulnerabilidad, como el riesgo de rechazo tras un trasplante o la recaída de ciertos cánceres después del tratamiento.

Utilizando la función de distribución acumulada normal estándar (\(\Phi\)) y su densidad (\(\phi\)), las principales funciones del modelo log-normal se expresan como:

  • Función de supervivencia:

\[S(t)=1-\Phi\left(\frac{\log(t)-\mu}{\sigma}\right)\]

  • Función de riesgo:

\[h(t)=\frac{\phi\left(\frac{\log(t)-\mu}{\sigma}\right)}{\sigma t\left[1-\Phi\left(\frac{\log(t)-\mu}{\sigma}\right)\right]}\]

  • Mediana del tiempo de supervivencia:

\[\text{Mediana}(T)=\exp(\mu)\]

En este modelo, el parámetro \(\mu\) controla la localización de la distribución, desplazando los tiempos de supervivencia hacia valores mayores o menores, mientras que \(\sigma\) determina la dispersión y la forma de la función de riesgo.

Diferencias con el modelo Weibull

Aunque los modelos Weibull y log-normal pertenecen a la familia de modelos paramétricos de supervivencia, presentan diferencias importantes en la forma del riesgo y en sus supuestos estadísticos.

La principal diferencia radica en la función de riesgo: mientras el modelo Weibull permite riesgos monótonos (constantes, crecientes o decrecientes), el log-normal describe riesgos unimodales, donde el riesgo aumenta inicialmente y luego disminuye.

Además, la distribución log-normal presenta colas más pesadas, lo que implica una mayor probabilidad de eventos tardíos en comparación con el Weibull.

Finalmente, a diferencia del modelo Weibull, la distribución log-normal no satisface el supuesto de riesgos proporcionales. Por ello, este modelo se interpreta exclusivamente bajo el enfoque de tiempo de falla acelerado, y no mediante hazard ratios constantes en el tiempo.

Ejemplo con el modelo log-normal

Utilizando nuevamente los datos de leucemia (aunque en este caso un modelo Weibull monótono podría resultar más apropiado), ajustamos un modelo log-normal para ilustrar la interpretación de sus parámetros:

## 
## Call:
## survreg(formula = Surv(t, d) ~ trat_f, data = base, dist = "lognormal")
##                     Value Std. Error     z       p
## (Intercept)        1.8251     0.2016  9.05 < 2e-16
## trat_fTratamiento  1.3468     0.3165  4.26 2.1e-05
## Log(scale)        -0.0792     0.1323 -0.60    0.55
## 
## Scale= 0.924 
## 
## Log Normal distribution
## Loglik(model)= -106.7   Loglik(intercept only)= -115.4
##  Chisq= 17.38 on 1 degrees of freedom, p= 3.1e-05 
## Number of Newton-Raphson Iterations: 4 
## n= 42
## Parámetro de ubicación (μ - Intercepto): 1.825
## Parámetro de escala (σ): 0.924
## Factor de aceleración (TR): 3.845

El coeficiente asociado al tratamiento fue:

\[ \hat{\beta}=1.3468 \]

En el modelo AFT, la interpretación se realiza mediante la razón de tiempo (TR):

\[ TR=\exp(1.3468)\approx 3.845 \]

Esto indica que los pacientes del grupo tratamiento presentan un tiempo libre de progresión aproximadamente 3.85 veces mayor que los pacientes del grupo control.

Obsérvese que el efecto estimado del tratamiento es consistente con el obtenido en los modelos anteriores, aunque la interpretación cambia debido a la estructura no monótona del riesgo en el modelo log-normal.

El parámetro de localización \((\mu)\) estimado fue:

\[ \hat{\mu}=1.825 \]

mientras que el parámetro de escala \((\sigma)\) fue:

\[ \hat{\sigma}=0.924 \]

En el modelo log-normal, \(\mu\) controla la posición de la distribución de supervivencia, de modo que valores mayores se asocian con tiempos de supervivencia más prolongados. Por su parte, \(\sigma\) describe la dispersión de los tiempos hasta el evento; valores más altos indican mayor variabilidad entre individuos.

A partir del intercepto, la mediana del tiempo de supervivencia para el grupo control puede estimarse como:

\[\text{Mediana}(T)=\exp(\hat{\mu})\]

Por tanto:

\[ \text{Mediana}(T)=\exp(1.825)\approx 6.20 \]

Esto sugiere que aproximadamente el 50% de los pacientes del grupo control experimentan el evento antes de los 6.2 meses de seguimiento, mientras que el tratamiento prolonga sustancialmente este tiempo.

Modelo Log-Logístico

El modelo log-logístico constituye otra alternativa flexible dentro de los modelos paramétricos de supervivencia. Al igual que el modelo log-normal, permite describir funciones de riesgo no monótonas, donde el riesgo puede aumentar inicialmente, alcanzar un máximo y posteriormente disminuir con el tiempo.

Este modelo asume que el logaritmo del tiempo de supervivencia sigue una distribución logística:

\[\log(T)\sim \text{Logística}(\mu,\sigma)\]

Equivalentemente:

\[T\sim \text{Log-Logística}(\mu,\sigma)\]

Una característica distintiva del modelo log-logístico es que presenta colas más pesadas que la distribución log-normal. En consecuencia, asigna mayor probabilidad a tiempos de supervivencia muy prolongados, lo que resulta útil en estudios donde algunos individuos presentan supervivencias excepcionalmente largas.

Además, este modelo posee una interpretación basada en odds proporcionales. En particular, la razón de odds acumuladas del evento entre grupos permanece constante en el tiempo.

Una ventaja importante del modelo log-logístico es que su función de supervivencia posee una expresión algebraica cerrada y relativamente simple:

\[S(t)=\frac{1}{1+(\lambda t)^p}\]

donde:

  • \(\lambda\) representa un parámetro de escala relacionado con la mediana de supervivencia,
  • y (p) corresponde al parámetro de forma.

La mediana del tiempo de supervivencia puede expresarse como:

\[\text{Mediana}(T)=\exp(\mu)\]

Al igual que en el modelo Weibull, el parámetro de forma se relaciona con \(\sigma\) mediante:

\[p=\frac{1}{\sigma}\]

El comportamiento de la función de riesgo depende directamente del valor de (p):

  • Si \(p>1\), la función de riesgo es unimodal.
  • Si \(p\leq1\), el riesgo disminuye continuamente con el tiempo.

Diferencias con los modelos Weibull y log-normal

El modelo log-logístico comparte con el modelo log-normal la capacidad de representar riesgos unimodales. Sin embargo, sus colas más pesadas permiten modelar con mayor flexibilidad eventos muy tardíos o supervivientes de larga duración.

A diferencia del modelo Weibull, el modelo log-logístico no satisface el supuesto de riesgos proporcionales. Por ello, su interpretación se realiza principalmente bajo el enfoque de tiempo de falla acelerado, utilizando razones de tiempo en lugar de HR constantes.

Además, la forma cerrada de su función de supervivencia facilita algunos cálculos analíticos e interpretaciones probabilísticas en comparación con el modelo log-normal.

Ejemplo práctico con el modelo log-logístico

Ajustamos un modelo log-logístico para evaluar la interpretación de sus parámetros bajo el enfoque de tiempo de falla acelerado (AFT) al ejemplo de leucemia.

## 
## Call:
## survreg(formula = Surv(t, d) ~ trat_f, data = base, dist = "loglogistic")
##                    Value Std. Error     z       p
## (Intercept)        1.893      0.208  9.12 < 2e-16
## trat_fTratamiento  1.265      0.326  3.89   1e-04
## Log(scale)        -0.604      0.150 -4.02 5.7e-05
## 
## Scale= 0.547 
## 
## Log logistic distribution
## Loglik(model)= -107.7   Loglik(intercept only)= -115.4
##  Chisq= 15.38 on 1 degrees of freedom, p= 8.8e-05 
## Number of Newton-Raphson Iterations: 4 
## n= 42
## Parámetro de ubicación (μ): 1.893
## Parámetro de escala (σ): 0.547
## Factor de aceleración (TR): 3.545

El coeficiente asociado al tratamiento fue:

\[ \hat{\beta}=1.265 \]

Calculando la razón de tiempo:

\[ TR=\exp(1.265)\approx 3.545 \]

Nuevamente, el tratamiento se asocia con un incremento importante en el tiempo hasta la progresión de la enfermedad, aproximadamente 3.5 veces superior al observado en el grupo control.

El parámetro de localización (\(\mu\)) estimado fue:

\[ \hat{\mu}=1.893 \]

mientras que el parámetro de escala (\(\sigma\)) fue:

\[ \hat{\sigma}=0.547 \]

A partir del parámetro de escala, el parámetro de forma se obtiene como:

\[p=\frac{1}{\sigma}\]

Por tanto:

\[ p=\frac{1}{0.547}\approx 1.83 \]

Dado que \(p>1\), el modelo es compatible con una función de riesgo que aumenta durante una etapa inicial del seguimiento y posteriormente disminuye.

En el contexto del estudio de leucemia, este patrón podría reflejar un período de mayor susceptibilidad a recaída tras el tratamiento, seguido de una reducción gradual del riesgo en pacientes que permanecen libres de progresión durante tiempos más prolongados.

Además, el valor relativamente bajo de \(\sigma\) sugiere una menor dispersión de los tiempos de supervivencia y una concentración más marcada del riesgo alrededor de su punto máximo.

Modelo Gamma Generalizado

El modelo Gamma Generalizado constituye una de las familias paramétricas más flexibles en análisis de supervivencia. A diferencia de distribuciones más específicas como Weibull o log-normal, este modelo puede adaptarse a una amplia variedad de formas de riesgo y patrones de supervivencia.

Su función de densidad puede expresarse como:

\[f(t)=\frac{|\beta|\theta^{\beta k}t^{\beta k-1}\exp(-(\theta t)^\beta)}{\Gamma(k)}\]

donde:

  • \(k\) y \(\beta\) son parámetros de forma,
  • \(\theta\) es un parámetro de escala,
  • y \(\Gamma(k)\) corresponde a la función Gamma.

Una característica importante del modelo Gamma Generalizado es que incluye como casos particulares varios de los modelos paramétricos más utilizados:

  • Si \(k=1\), el modelo se reduce a una distribución Weibull.
  • Si \(k \to \infty\), converge hacia una distribución log-normal.
  • Si \(k=1\) y \(\beta=1\), se obtiene el modelo exponencial.

En consecuencia, esta familia permite representar riesgos constantes, monótonos o unimodales dentro de un mismo marco matemático.

Ventajas y limitaciones

La principal ventaja de este modelo es el de tener una elevada flexibilidad, ya que se puede aproximar a diferentes distribuciones paramétricas mediante variaciones de sus parámetros.

No obstante, esta flexibilidad también implica mayor complejidad estadística. A diferencia de modelos como Weibull o log-normal, el Gamma Generalizado requiere estimar múltiples parámetros de forma, lo que puede aumentar la incertidumbre de las estimaciones y dificultar la interpretación clínica de los resultados.

Además, este modelo no posee una interpretación PH/AFT tan directa como Weibull, por lo que suele emplearse principalmente cuando los modelos paramétricos más simples no logran describir adecuadamente los datos observados.

En la práctica, este modelo resulta especialmente útil como familia de referencia para comparar distintas distribuciones paramétricas y evaluar cuál de ellas ofrece el mejor ajuste a los datos observados, sin imponer restricciones a priori sobre la forma del riesgo.

Ejemplo con el modelo Gamma Generalizado

Usando los datos de leucemia, podemos ajustar un modelo Gamma Generalizado:

## Call:
## flexsurvreg(formula = Surv(t, d) ~ trat_f, data = base, dist = "gengamma")
## 
## Estimates: 
##                    data mean  est     L95%    U95%    se      exp(est)  L95%  
## mu                     NA      2.063   1.386   2.741   0.346      NA        NA
## sigma                  NA      0.830   0.540   1.278   0.183      NA        NA
## Q                      NA      0.551  -0.807   1.910   0.693      NA        NA
## trat_fTratamiento   0.500      1.306   0.680   1.932   0.319   3.691     1.973
##                    U95%  
## mu                     NA
## sigma                  NA
## Q                      NA
## trat_fTratamiento   6.903
## 
## N = 42,  Events: 30,  Censored: 12
## Total time at risk: 541
## Log-likelihood = -106.3895, df = 4
## AIC = 220.779

El efecto del tratamiento estimado fue:

\[β_1=1.306\]

En el contexto del modelo AFT, la interpretación también se realiza mediante la razón de tiempo:

\[TR=exp(1.306)≈3.69\]

En los pacientes del grupo tratamiento presentan un tiempo libre de progresión aproximadamente 3.7 veces mayor que los pacientes del grupo control, consistente con los resultados obtenidos en los modelos paramétricos anteriores.

El modelo Gamma generalizado introduce parámetros adicionales de forma que permiten mayor flexibilidad en la estructura del riesgo:

El parámetro μ = 2.063 representa el nivel de localización del tiempo de supervivencia en escala logarítmica. Este parámetro refleja el tiempo central esperado hasta el evento (progresión o muerte, según la definición del estudio) en la población analizada. Al transformarlo a la escala original, se obtiene un tiempo típico aproximado de supervivencia de \(e^{2.063} \approx 7.87\) unidades de tiempo, lo cual proporciona una referencia del horizonte temporal central de la enfermedad bajo el esquema del modelo ajustado.

El parámetro σ = 0.830 corresponde a la dispersión de los tiempos de supervivencia en la escala logarítmica y describe la heterogeneidad entre los pacientes con leucemia. En términos clínicos, este valor indica que existe una variabilidad moderada a alta en la evolución de la enfermedad, lo cual es esperable en este tipo de patologías, donde factores biológicos, respuesta al tratamiento y características individuales generan trayectorias de supervivencia bastante diferenciadas.

El parámetro Q controla la forma de la función de riesgo (hazard) y permite caracterizar distintos patrones de evolución del riesgo en el tiempo. En términos generales:

Cuando Q = 1, el modelo se reduce al caso de Weibull, lo que implica un hazard monótono, ya sea creciente o decreciente.

Cuando Q tiende a 0, el modelo se aproxima a una distribución log-normal, asociada típicamente a un hazard no monótono con forma de “campana”, es decir, con aumento inicial del riesgo seguido de una disminución.

Cuando Q < 1, el modelo permite que el riesgo de ocurrencia del evento no siga una tendencia única y constante en el tiempo. En lugar de aumentar o disminuir de manera uniforme, el riesgo puede cambiar de comportamiento a lo largo del seguimiento, presentando distintas fases, como un incremento inicial seguido de estabilización o incluso una disminución posterior. Esto refleja situaciones clínicas en las que la evolución del riesgo no es homogénea, sino que depende del momento del proceso de la enfermedad.

Cuando Q > 1, el hazard suele adoptar formas más estructuradas y monótonas, asociadas a procesos de deterioro progresivo más sostenido.

En el contexto del ejemplo, el valor de Q = 0.551 indica que el riesgo de progresión o muerte sigue una dinámica más flexible en el tiempo, coherente con la heterogeneidad clínica de los pacientes y la variabilidad en la evolución de la enfermedad. En particular, este valor sugiere que la función de riesgo no se ajusta a un patrón monótono simple, como el asumido por modelos exponenciales o Weibull, sino que puede presentar cambios en su intensidad a lo largo del tiempo. Esto implica que el riesgo puede ser mayor en ciertos periodos del seguimiento y disminuir en otros, reflejando posibles fases clínicas diferenciadas en la evolución de la leucemia.

Finalmente, los intervalos de confianza y errores estándar asociados a los parámetros reflejan la precisión de las estimaciones obtenidas mediante máxima verosimilitud. En conjunto, los resultados indican que el modelo Gamma Generalizado proporciona una representación flexible y adecuada del comportamiento de los tiempos de supervivencia en pacientes con leucemia, capturando tanto la heterogeneidad individual como la complejidad en la dinámica del riesgo a lo largo del tiempo.

Elección del modelo paramétrico

En la práctica, la elección del modelo paramétrico debe ser un análisis robusto que requiere balancear algunos aspectos fundamentales:

Plausibilidad Biológica y Teórica

Es el criterio más importante. Antes de mirar los datos, debemos preguntarnos evaluar si esa forma de riesgo tiene sentido para la enfermedad en estudio.

  • Riesgo Constante (Exponencial): ¿Es razonable asumir que la probabilidad de morir no cambia con el tiempo? (Ej: muerte accidental).

  • Riesgo Monótono (Weibull): ¿El riesgo debería aumentar por acumulación de daño (envejecimiento) o disminuir por selección de sobrevivientes (recuperación post-quirúrgica)?

  • Riesgo Unimodal (Log-Normal/Log-Logística): ¿Existe un periodo crítico de vulnerabilidad seguido de una etapa de estabilización? (Ej: rechazo de trasplante o recaída tras quimioterapia).

Si un modelo ajusta bien estadísticamente pero propone una curva de riesgo que contradice la naturaleza de la patología, sus predicciones a largo plazo no serán confiables.

Principio de Parsimonia

El principio de parsimonia establece que, al construir modelos estadísticos, se debe preferir aquel que sea lo más simple posible, siempre que mantenga un adecuado poder explicativo en comparación con modelos más complejos.

En este sentido, los modelos con un número reducido de parámetros, como el modelo exponencial, suelen ser más robustos y presentan menor riesgo de sobreajuste, ya que imponen una estructura más restrictiva sobre los datos.

Por el contrario, modelos más flexibles y con mayor número de parámetros, como la distribución gamma generalizada, pueden lograr un ajuste muy preciso en una muestra específica, pero esa misma flexibilidad puede convertirse en una desventaja al momento de generalizar, ya que tienden a capturar ruido propio de los datos observados en lugar de patrones estructurales reales.

En consecuencia, aunque estos modelos complejos puedan ofrecer un mejor ajuste aparente, su capacidad predictiva fuera de la muestra o en otras poblaciones puede ser considerablemente más limitada.

3. Criterios Estadísticos de Información (AIC y BIC)

Como apoyo para la selección de modelos cuando la teoría no es suficiente, se recurre a criterios de comparación como el AIC y el BIC. En el caso del AIC, este tiende a ser más permisivo con la complejidad del modelo, por lo que resulta especialmente útil cuando el objetivo principal es la predicción, mientras que el BIC impone una penalización más fuerte sobre el número de parámetros, especialmente a medida que aumenta el tamaño de la muestra, lo que lo hace más adecuado cuando se busca un modelo más parsimonioso y con mayor capacidad interpretativa.

4. Validación Visual y Diagnóstico

Finalmente, la inspección gráfica constituye una herramienta indispensable en la validación del modelo. La comparación de las curvas de supervivencia estimadas mediante el modelo paramétrico con aquellas obtenidas a través del estimador de Kaplan–Meier permite evaluar de forma visual si el modelo tiende a subestimar o sobrestimar la supervivencia en distintos puntos del seguimiento. Cuando la curva paramétrica se desvía de manera sistemática respecto a la de Kaplan–Meier, debe cuestionarse la adecuación del supuesto de distribución y, en consecuencia, la validez del modelo propuesto

La Función de Verosimilitud

A diferencia del modelo de Cox, que utiliza una verosimilitud parcial (basada solo en el orden de los eventos y no en el tiempo exacto), los modelos paramétricos se basan en la verosimilitud completa. Esto significa que aprovechan toda la información disponible, incorporando explícitamente la distribución de probabilidad asumida para el tiempo de supervivencia.

La verosimilitud es, esencialmente, la probabilidad de observar los datos que tenemos dado un conjunto de parámetros. Dependiendo de lo que le ocurra a cada paciente durante el estudio, su contribución a la función de verosimilitud (\(L_i\)) varía:

  • Caso 1: Evento observado (Exacto). Conocemos el momento preciso del evento. Su contribución es la densidad de probabilidad en ese instante: \(L_i = f(t_i \mid \mathbf{x}_i)\).
  • Caso 2: Censura a la derecha. Solo sabemos que el paciente estaba vivo hasta el tiempo \(t\). Su contribución es la probabilidad de sobrevivir más allá de ese punto: \(L_i = S(t_i \mid \mathbf{x}_i)\).
  • Caso 3: Censura por intervalos. Sabemos que el evento ocurrió entre dos visitas (tiempo \(L\) y tiempo \(U\)), pero no el momento exacto. Su contribución es la diferencia de supervivencia entre ambos puntos: \(L_i = S(L_i) - S(U_i)\).

La Log-Verosimilitud Conjunta

Para un estudio con \(n\) individuos independientes, la verosimilitud total es el producto de todas las contribuciones. Sin embargo, matemáticamente es mucho más sencillo trabajar con la suma de los logaritmos, conocida como log-verosimilitud (\(\ell\)):

\[\ell(\theta) = \sum_{i=1}^{n} \left\{\delta_i \log f(t_i \mid \mathbf{x}_i) + (1-\delta_i) \log S(t_i \mid \mathbf{x}_i)\right\}\]

Donde \(\delta_i\) es un indicador que vale 1 si el evento ocurrió y 0 si fue censurado. Esta estructura es elegante y eficiente: los eventos observados aportan información sobre la forma de la densidad, mientras que los censurados aportan información sobre la probabilidad de supervivencia.

Los parámetros del modelo (\(\theta\)), que incluyen tanto los coeficientes de las covariables como los parámetros de forma y escala de la distribución, se obtienen mediante el método de Máxima Verosimilitud (MLE). Este proceso busca los valores que maximizan la probabilidad de haber observado nuestros datos.

Bajo una especificación correcta del modelo, estos estimadores poseen propiedades estadísticas óptimas como la de convergencia, eficiencia y normalidad adintótica.

Finalmente, los modelos paramétricos ofrecen ventajas insuperables cuando la distribución elegida es adecuada:

  • Mayor Eficiencia: Al usar la verosimilitud completa, proporcionan estimaciones más precisas y mayor potencia estadística.

  • Capacidad Predictiva Extrema: Permiten extrapolar resultados y predecir la supervivencia en tiempos que van más allá del último evento observado en los datos.

  • Estimación de Cuantiles: Facilitan la obtención directa de métricas clínicas intuitivas, como la mediana o los cuartiles de tiempo de vida, algo que Cox solo puede aproximar de forma indirecta.

Ejemplo práctico

Contexto clínico del estudio

La insuficiencia cardíaca crónica es una enfermedad progresiva en la que el corazón pierde la capacidad de mantener un gasto adecuado, lo que conduce a episodios recurrentes de descompensación y hospitalización. Las rehospitalizaciones cardiovasculares son un desenlace clínico clave, ya que reflejan progresión de la enfermedad y mayor riesgo de mortalidad.

En este ejemplo se estudia el tiempo hasta la primera rehospitalización en pacientes seguidos durante varios meses, con el objetivo de evaluar si un programa intensivo de seguimiento prolonga el tiempo libre de evento frente al manejo convencional. El evento de interés es la rehospitalización cardiovascular, aunque no todos los pacientes lo presentan durante el seguimiento, ya sea por finalización del estudio o pérdida en el seguimiento, lo que genera observaciones censuradas.

Paso 1. Abriendo la base de datos

library(readxl)
datos <- read_excel("parametricos.xlsx")
summary(datos)
##       time           status       tratamiento             edad      
##  Min.   : 0.60   Min.   :0.0000   Length:60          Min.   :50.00  
##  1st Qu.: 6.15   1st Qu.:0.0000   Class :character   1st Qu.:59.75  
##  Median :11.70   Median :1.0000   Mode  :character   Median :66.00  
##  Mean   :12.25   Mean   :0.7333                      Mean   :66.75  
##  3rd Qu.:17.40   3rd Qu.:1.0000                      3rd Qu.:73.00  
##  Max.   :31.70   Max.   :1.0000                      Max.   :88.00

Como podemos ver, esta base de datos contiene cuatro variables: time, status, tratamiento y edad.

La variable time representa el tiempo observado hasta que ocurre la rehospitalización o hasta que el paciente es censurado. La variable status indica si realmente ocurrió el evento:

  • 1 = rehospitalización
  • 0 = censura

La variable tratamiento identifica el grupo terapéutico al cual pertenece el paciente (Convencional o Intensivo).

Finalmente, la edad se incluye como una covariable clínica importante, debido a que pacientes mayores suelen presentar peor pronóstico cardiovascular.

Paso 2. Exploración inicial de los datos

Antes de construir modelos, es importante hacer un análisis exploratorio de la información con el fin de detectar errores, verificar estructuras de variables y obtener una comprensión preliminar de la población estudiada.

datos$tratamiento <- factor(datos$tratamiento,
                             levels = c("Convencional", "Intensivo"))

# Crear variable factor para status solo para reportes, mantener status numérico para Surv()
datos$status_label <- factor(datos$status,
                           levels = c(0, 1),
                           labels = c("Censurado", "Rehospitalizado"))

summary(datos)
##       time           status             tratamiento      edad      
##  Min.   : 0.60   Min.   :0.0000   Convencional:27   Min.   :50.00  
##  1st Qu.: 6.15   1st Qu.:0.0000   Intensivo   :33   1st Qu.:59.75  
##  Median :11.70   Median :1.0000                     Median :66.00  
##  Mean   :12.25   Mean   :0.7333                     Mean   :66.75  
##  3rd Qu.:17.40   3rd Qu.:1.0000                     3rd Qu.:73.00  
##  Max.   :31.70   Max.   :1.0000                     Max.   :88.00  
##           status_label
##  Censurado      :16   
##  Rehospitalizado:44   
##                       
##                       
##                       
## 
table(datos$tratamiento)
## 
## Convencional    Intensivo 
##           27           33
library(tableone)
tabla <- CreateTableOne(vars = c("time", "edad"), 
                        strata = "tratamiento",
                         data = datos)
print(tabla, showAllLevels = TRUE)
##                   Stratified by tratamiento
##                    level Convencional Intensivo    p      test
##   n                         27           33                   
##   time (mean (SD))       11.91 (7.64) 12.52 (7.78)  0.762     
##   edad (mean (SD))       66.81 (9.43) 66.70 (9.26)  0.961
# Tabla de status por grupo
table_status <- table(datos$tratamiento, datos$status_label)
print(table_status)
##               
##                Censurado Rehospitalizado
##   Convencional         7              20
##   Intensivo            9              24

Observamos una mediana de seguimiento cercana a 12 meses, lo que indicaría que aproximadamente la mitad de los pacientes presentó el evento o fue censurado antes de ese tiempo.

Al observar las frecuencias de cada grupo, notamos que el número de pacientes en el grupo intensivo es ligeramente mayor que en el convencional. Esto es importante para tener en cuenta al interpretar los resultados, ya que un desequilibrio en el tamaño de los grupos puede influir en la potencia estadística y la precisión de las estimaciones.

Adicionalmente, al comparar el status entre grupos, se observa que el grupo intensivo tiene una proporción ligeramente menor de rehospitalizaciones en comparación con el grupo convencional. Esto podría ser un indicio preliminar de que el programa intensivo tiene un efecto protector, aunque es necesario realizar análisis formales para confirmar esta hipótesis.

En términos de edad, ambos grupos presentan una distribución similar, con una media alrededor de 67 años. Aunque pareciera que la edad no es un factor de confusión importante, se debe evaluar si es necesario ajustar por ella en modelos posteriores para mayor precisión.

Paso 3. Curvas de Kaplan-Meier

Estas curvas permiten estimar la probabilidad acumulada de permanecer libre del evento a través del tiempo. En este caso, la curva representa la probabilidad de continuar sin rehospitalización cardiovascular conforme avanzan los meses de seguimiento.

fit_km <- survfit(
  Surv(time, status) ~ tratamiento,
  data = datos
)

plot(fit_km, col  = c("red", "blue"), lwd  = 2, xlab = "Meses",
  ylab = "Probabilidad libre de rehospitalización",
  main = "Curvas Kaplan-Meier"
)

legend("topright", legend = levels(datos$tratamiento),
  col = c("red", "blue"), lwd = 2
)

El gráfico muestra las curvas de Kaplan–Meier para el tiempo libre de rehospitalización en dos grupos: manejo convencional (rojo) e intensivo (azul). En los primeros meses del seguimiento, ambos grupos presentan una supervivencia muy similar, e incluso el grupo intensivo muestra una ligera desventaja inicial, lo que sugiere rehospitalizaciones relativamente más tempranas en ese periodo.

Sin embargo, a partir de aproximadamente los 6–8 meses las curvas comienzan a separarse e incluso se observa un cruce, lo que indica que el efecto del tratamiento no es constante en el tiempo. Después de este punto, el grupo intensivo tiende a mostrar una mayor probabilidad de permanecer libre de rehospitalización durante gran parte del seguimiento, especialmente entre los meses 15 y 25.

Este patrón de cruce sugiere una posible violación del supuesto de riesgos proporcionales, ya que la ventaja relativa entre grupos cambia de dirección en el tiempo.

Paso 4. Prueba log-rank

Aunque la inspección visual de las curvas es muy útil, necesitamos la prueba de log-rank para determinar si las diferencias observadas son compatibles con el azar o representan un efecto real del tratamiento.

library(survival)
survdiff(Surv(time, status) ~ tratamiento, data = datos)
## Call:
## survdiff(formula = Surv(time, status) ~ tratamiento, data = datos)
## 
##                           N Observed Expected (O-E)^2/E (O-E)^2/V
## tratamiento=Convencional 27       20     19.4    0.0157    0.0295
## tratamiento=Intensivo    33       24     24.6    0.0124    0.0295
## 
##  Chisq= 0  on 1 degrees of freedom, p= 0.9

Dado que \(p = 0.9\), se concluye que no existe evidencia estadísticamente significativa de diferencia entre los tratamientos. Esto significa que el programa intensivo no modifica significativamente el tiempo hasta rehospitalización.

Paso 5. Evaluación gráfica de la distribución

Una de las decisiones más importantes en modelos paramétricos de supervivencia consiste en identificar la forma funcional más adecuada para el tiempo hasta el evento. A diferencia de los métodos no paramétricos, estos modelos requieren asumir una distribución específica (exponencial, Weibull, log-normal, gamma, entre otras), por lo que una inspección exploratoria es fundamental antes del ajuste formal.

Una forma útil de explorar la estructura del tiempo de falla es mediante representaciones de la función de riesgo acumulado o transformaciones logarítmicas de la función de supervivencia, ya que estas permiten evaluar si el patrón de los datos es consistente con riesgos constantes, crecientes o decrecientes.

library(survival)

fit_haz <- survfit(Surv(time, status) ~ tratamiento, data = datos)

plot(fit_haz, fun = "cumhaz",
     xlab = "Tiempo (meses)",
     ylab = "Riesgo acumulado",
     main = "Función de riesgo acumulado (Nelson-Aalen)",
     lwd = 2)

Se observa que el riesgo acumulado aumenta de forma monótona a lo largo del tiempo, sin embargo, el patrón de crecimiento no es completamente uniforme:

  • En los primeros meses (0–10 meses), el incremento del riesgo es relativamente lento y estable.

  • Entre los meses 10 y 20 se observa una aceleración progresiva del riesgo, con incrementos más frecuentes.

  • A partir de aproximadamente los 20 meses, el riesgo acumulado presenta una pendiente más pronunciada, lo que sugiere mayor incidencia de eventos en etapas tardías del seguimiento.

Este comportamiento indica que el riesgo de rehospitalización no es constante en el tiempo, sino que tiende a incrementarse conforme avanza la enfermedad.

Paso 6. Elección de la distribución paramétrica

Una vez explorado el comportamiento de la función de riesgo acumulado, el siguiente paso en el modelado paramétrico consiste en identificar qué distribución podría representar adecuadamente el tiempo hasta la rehospitalización. Esta decisión es crucial, ya que cada distribución impone una forma específica sobre el riesgo en el tiempo.

En términos generales, recordemos que las distribuciones más utilizadas en análisis de supervivencia se diferencian por la forma de su función de riesgo:

  • La distribución exponencial asume un riesgo constante en el tiempo.

  • La distribución Weibull permite riesgos crecientes o decrecientes dependiendo del parámetro de forma.

  • La distribución log-normal o log-logística permite riesgos no monótonos, con un pico en el tiempo.

  • La distribución gamma generalizada ofrece mayor flexibilidad, incluyendo casos que se aproximan a las anteriores.

Dado el comportamiento observado en la función de riesgo acumulado, donde el riesgo aumenta progresivamente con el tiempo, es poco probable que el modelo exponencial sea adecuado, ya que este implicaría una pendiente aproximadamente constante en el riesgo acumulado.

La forma observada es consistente con un escenario de riesgo creciente, lo que sugiere que la distribución Weibull podría ser una alternativa apropiada. Recordemos que en el modelo Weibull la forma del riesgo depende del parámetro de forma \(p\):

\[h(t)=\lambda p t^{p-1}\]

donde:

  • Si \(p = 1\) → riesgo constante (equivalente a exponencial)

  • Si \(p > 1\) → riesgo creciente en el tiempo

  • Si \(p < 1\) → riesgo decreciente en el tiempo

En este contexto, la evidencia gráfica sugiere un comportamiento compatible con \(p > 1\). Esto es coherente con la progresión clínica de la insuficiencia cardíaca, donde la probabilidad de rehospitalización aumenta en etapas más avanzadas del seguimiento.

Sin embargo, la elección de la distribución no debe basarse únicamente en inspección visual. Es necesario confirmar esta hipótesis mediante ajuste de modelos paramétricos y comparación de criterios de bondad de ajuste como AIC, BIC o log-verosimilitud.

Paso 7. Verificación formal mediante ajuste paramétrico

Para validar la selección de la distribución, se ajustan los distintos modelos paramétricos y se comparan sus criterios de ajuste. Esto permite cuantificar cuál distribución describe mejor los datos observados.

library(survival)

# =========================
# Modelo exponencial (AFT)
# =========================
fit_exp <- survreg(
  Surv(time, status) ~ tratamiento,
  data = datos,
  dist = "exponential"
)

# =========================
# Modelo Weibull (AFT)
# =========================
fit_weib <- survreg(
  Surv(time, status) ~ tratamiento,
  data = datos,
  dist = "weibull"
)

# =========================
# Modelo log-normal
# =========================
fit_lnorm <- survreg(
  Surv(time, status) ~ tratamiento,
  data = datos,
  dist = "lognormal"
)

# =========================
# Modelo log-logístico (alternativa común)
# =========================
fit_llogis <- survreg(
  Surv(time, status) ~ tratamiento,
  data = datos,
  dist = "loglogistic"
)

AIC(fit_exp,fit_weib,fit_lnorm,fit_llogis)
##            df      AIC
## fit_exp     2 339.7292
## fit_weib    3 324.4731
## fit_lnorm   3 328.6836
## fit_llogis  3 328.9708

Al comparar los valores de AIC, vemos que el modelo weibull es el que presenta menor valor, por lo que se considera el mejor, confirmando que el supuesto de riesgo constante es inadecuado y que la dinámica del evento es dependiente del tiempo.

Paso 7. Interpretación detallada del modelo seleccionado

El siguiente paso consiste en revisar cuidadosamente los coeficientes del modelo Weibull.

summary(fit_weib)
## 
## Call:
## survreg(formula = Surv(time, status) ~ tratamiento, data = datos, 
##     dist = "weibull")
##                        Value Std. Error     z       p
## (Intercept)           2.7816     0.1279 21.74 < 2e-16
## tratamientoIntensivo  0.0519     0.1727  0.30    0.76
## Log(scale)           -0.5614     0.1192 -4.71 2.5e-06
## 
## Scale= 0.57 
## 
## Weibull distribution
## Loglik(model)= -159.2   Loglik(intercept only)= -159.3
##  Chisq= 0.09 on 1 degrees of freedom, p= 0.76 
## Number of Newton-Raphson Iterations: 7 
## n= 60
Interpretación del intercepto y escala

El intercepto representa el logaritmo del tiempo base esperado hasta el evento en el grupo de referencia (tratamiento convencional). Indica que, en promedio, los pacientes presentan tiempos relativamente largos hasta la rehospitalización.

Recordemos la relación entre la escala del modelo AFT y el parámetro de forma \(p\) en Weibull:

\(p = \frac{1}{\text{scale}} \approx 1 / 0.57 \approx 1.75\)

Esto implica que \(p > 1\), lo cual confirma formalmente lo observado en la exploración gráfica: el riesgo de rehospitalización es creciente en el tiempo.

Interpretación del efecto del tratamiento

Como vimos previamente, en modelos AFT, la interpretación se hace mediante razones de tiempo.

TR <- exp(
  coef(fit_weib)["tratamientoIntensivo"]
)
TR
## tratamientoIntensivo 
##             1.053272

Este resultado indica que, en promedio, los pacientes en el grupo intensivo presentan un tiempo hasta la rehospitalización cerca de un 5% mayor en comparación con el grupo convencional. Sin embargo, este efecto no es estadísticamente significativo, lo que sugiere que la evidencia disponible no es suficiente para concluir un beneficio real del programa intensivo sobre la prolongación del tiempo libre de evento.

Desde una perspectiva clínica, estos hallazgos sugieren que, aunque descriptivamente podría observarse cierta separación tardía entre los grupos, el análisis paramétrico ajustado no respalda una diferencia consistente en la evolución temporal de la rehospitalización. En consecuencia, el programa intensivo evaluado no modifica de manera relevante la historia natural del evento en este conjunto de datos.

Finalmente, podemos concluir que el tiempo hasta la rehospitalización sigue una dinámica de riesgo creciente, consistente con la fisiopatología progresiva de la insuficiencia cardíaca, mientras que el efecto del tratamiento no muestra impacto significativo sobre dicha trayectoria temporal.

Bibliografía

  • Survival analysis: A self-learning text. David G. Kleinbaum and Mitchel Klein. Springer, New York, 2012.

  • Función de riesgo instantáneo: Cox, D. R. (1972). Regression models and life-tables. Journal of the Royal Statistical Society: Series B (Methodological), 34(2), 187-202.

  • Función de riesgo acumulado: Nelson, W. (1972). Theory and applications of hazard plotting for censored failure data. Technometrics, 14(4), 945-966.

  • R: Therneau, T. M. (2026). A package for survival analysis in R. R package version 3.8-6, https://CRAN.R-project.org/package=survival.

  • klein, J. P., & Moeschberger, M. L. (2003). Survival analysis: techniques for censored and truncated data. Springer Science & Business Media.

  • Cox, D. R., & Oakes, D. (1984). Analysis of Survival Data. London: Chapman & Hall.

  • Hosmer, D. W., Lemeshow, S., & May, S. (2008). Applied Survival Analysis: Regression Modeling of Time-to-Event Data (2nd ed.). Hoboken, NJ: Wiley.

  • R Core Team. (2023). R: A language and environment for statistical computing. Vienna, Austria. Retrieved from https://www.R-project.org/