16  Modelo de Fragilidad

16.1 Modelos de Fragilidad en Análisis de Sobrevivencia

16.1.1 Verosimilitud en datos independientes

Para \(n\) observaciones independientes \((t_i,\delta_i,z_i)\) bajo riesgos proporcionales:

\[ h(t;\beta)=h_0(t)\,e^{z_i\beta}, \]

la función de verosimilitud puede escribirse

\[ L(\beta;z_i) =\prod_{i=1}^n f(t_i)^{\delta_i}S(t_i)^{1-\delta_i} =\prod_{i=1}^n h(t_i)^{\delta_i}S(t_i) =\prod_{i=1}^n\bigl[h_0(t_i)e^{z_i\beta}\bigr]^{\delta_i} \exp\bigl(-H_0(t_i)\,e^{z_i\beta}\bigr), \]

donde \(H_0(t)=\int_0^t h_0(u),du\).

16.1.2 Introduciendo fragilidades para datos agrupados

Si las observaciones vienen en \(G\) clústeres (familias, hospitales, etc.) y cada clúster \(i\) comparte una “frailty” \(\omega_i\), asumimos

\[ h_{ij}(t_{ij})=h_0(t_{ij})\,\omega_i\,e^{z_{ij}\beta},\quad \omega_i\sim\Gamma\Bigl(\tfrac{1}{\theta},\theta\Bigr), \]

con densidad

\[ g(\omega;\theta) =\frac{\omega^{\frac1\theta-1}e^{-\omega/\theta}} {\Gamma(1/\theta)\,\theta^{1/\theta}}. \]

Si pudiéramos observar \(\omega_i\), el término de verosimilitud por sujeto sería

\[ L_{ij}(\beta,\theta;\omega_i) = g(\omega_i;\theta) \bigl[h_0(t_{ij})\,\omega_i\,e^{z_{ij}\beta}\bigr]^{\delta_{ij}} \exp\bigl(-H_0(t_{ij})\,\omega_i\,e^{z_{ij}\beta}\bigr), \]

y la verosimilitud total

\[ L(\beta,\theta) =\prod_{i=1}^G\prod_{j=1}^{n_i}L_{ij}(\beta,\theta;\omega_i). \]

Como las \(\omega_i\) son latentes, utilizamos un algoritmo EM que alterna:

  1. E-step: calcular \(E[\omega_i\mid\text{datos},\beta,\theta]\).
  2. M-step: maximizar la verosimilitud completa en \(\beta,\theta\) usando esas esperanzas.

En modelos paramétricos (p.ej. Weibull) es directo; en el modelo de Cox hay que también actualizar \(h_0(t)\) en cada iteración (véase Klein & Moeschberger, Therneau & Grambsch).

16.1.3 Alternativa normal de los efectos aleatorios

En lugar de usar una distribución gamma, podemos modelar

\[ \omega_i=e^{\sigma\,u_i},\quad u_i\sim N(0,1), \]

llevando el efecto aleatorio al mismo nivel que los fijos:

\[ h_{ij}(t_{ij})=h_0(t_{ij})\,\exp\bigl(z_{ij}\beta + u_i\sigma\bigr). \]

El EM es análogo, sustituyendo la distribución normal por la gamma.

16.2 Ejemplo: Fragilidad en el estudio “Washington Ashkenazi”

En el estudio “Washington Ashkenazi” (Struewing et al. [64]), se reclutaron 1 960 familias de ascendencia judía ashkenazi, cada una con al menos dos parientes femeninas de primer grado. A cada sujeto se le determinó su estado de mutación en el gen BRCA y se registró, para dos de sus parientes (edad de inicio de cáncer de mama o edad actual si no lo padecen), su edad de diagnóstico (“brcancer”: 1=cáncer, 0=sin cáncer). La variable clave es mutant (0/1 según la sujeto porte o no mutación). Como ejemplo, las familias 1, 9 y 94 muestran este formato:

famID brcancer age mutant
1 0 73 0
1 0 40 0
9 0 89 0
9 1 60 0
94 1 44 1
94 0 45 1

Dado que parientes de una misma familia comparten características no observables, no pueden tratarse como observaciones independientes.

16.2.1 Datos y objetivo

  • Familias: G = 1 960
  • Sujetos totales: N = 3 920 (dos parientes por familia)
  • Eventos de cáncer de mama: d = 473

Se busca estimar el efecto del estado de “portadora de mutación BRCA” de la probanda (mutant) sobre la edad de inicio de cáncer en sus parientes, teniendo en cuenta la posible dependencia entre miembros de una misma familia.

library(survival)
library(asaur)
data(ashkenazi)

16.2.2 1. Modelo de Cox “simple” (sin agrupar)

result.coxph <- coxph(Surv(age, brcancer) ~ mutant,
                     data=ashkenazi)
summary(result.coxph)
Call:
coxph(formula = Surv(age, brcancer) ~ mutant, data = ashkenazi)

  n= 3920, number of events= 473 

         coef exp(coef) se(coef)     z Pr(>|z|)    
mutant 1.1907    3.2895   0.1984 6.002 1.95e-09 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

       exp(coef) exp(-coef) lower .95 upper .95
mutant      3.29      0.304      2.23     4.853

Concordance= 0.523  (se = 0.006 )
Likelihood ratio test= 25.92  on 1 df,   p=4e-07
Wald test            = 36.03  on 1 df,   p=2e-09
Score (logrank) test = 40.49  on 1 df,   p=2e-10
  • Estimador: \(\hat\beta=1.1907\)
  • Hazard ratio: \(\exp(1.1907)=3.2895\)
  • Error estándar: se = 0.1984
  • Valor z: \(1.1907/0.1984\approx6.00\)
  • p-valor: \(1.95\times10^{-9}\)

Prueba de cociente de verosimilitudes La función coxph devuelve también los log-verosimilitudes:

result.coxph$loglik
[1] -3579.707 -3566.745
  • Sin covariables (null): \(\ell_0=-3579.707\)
  • Con mutant: \(\ell_1=-3566.745\)
  • Estadístico \(G^2 = 2(\ell_1 - \ell_0) = 2(-3566.745 + 3579.707) = 25.924\), que contra \(\chi^2_1\) da p ≪ 0.001.

16.2.3 2. Cox con frailty Gamma

result.coxph.frail <- coxph(Surv(age, brcancer) ~ mutant
                            + frailty(famID),
                            data=ashkenazi)
summary(result.coxph.frail)
Call:
coxph(formula = Surv(age, brcancer) ~ mutant + frailty(famID), 
    data = ashkenazi)

  n= 3920, number of events= 473 

               coef  se(coef) se2    Chisq  DF    p      
mutant         1.272 0.2317   0.2004  30.13   1.0 4.0e-08
frailty(famID)                       221.50 211.6 3.1e-01

       exp(coef) exp(-coef) lower .95 upper .95
mutant     3.567     0.2804     2.265     5.617

Iterations: 6 outer, 26 Newton-Raphson
     Variance of random effect= 0.4937594   I-likelihood = -3564.3 
Degrees of freedom for terms=   0.7 211.6 
Concordance= 0.938  (se = 0.004 )
Likelihood ratio test= 400.2  on 212.4 df,   p=1e-13
  • \(\hat\beta_{\rm mutant}=1.272\), se = 0.2317
  • Chisq(test único) = 30.13, p ≈ 4×10⁻⁸.
  • Frailty (varianza de ωₖ): estimada no significativa (p ≈ 0.31), lo que sugiere poca heterogeneidad adicional más allá de mutant.

Nota: la opción frailty() incorpora un efecto aleatorio Gamma \(\omega_i\) en el hazard: \[ h_{ij}(t)=h_0(t)\,\omega_i\,e^{\beta\,z_{ij}},\quad \omega_i\sim\Gamma\bigl(1/\theta,\theta\bigr). \]

16.2.4 3. Modelo mixto Cox con coxme

library(coxme)
Loading required package: bdsmatrix

Attaching package: 'bdsmatrix'
The following object is masked from 'package:base':

    backsolve
result.coxme <- coxme(Surv(age, brcancer) ~ mutant
                      + (1|famID),
                      data=ashkenazi)
summary(result.coxme)
Mixed effects coxme model
 Formula: Surv(age, brcancer) ~ mutant + (1 | famID) 
    Data: ashkenazi 

  events, n = 473, 3920

Random effects:
  group  variable        sd  variance
1 famID Intercept 0.5912135 0.3495334
                   Chisq    df        p   AIC     BIC
Integrated loglik  30.17   2.0 2.81e-07 26.17   17.85
 Penalized loglik 336.37 150.1 2.22e-16 36.16 -588.13

Fixed effects:
         coef exp(coef) se(coef)    z        p
mutant 1.2366    3.4439   0.2205 5.61 2.06e-08
  • Efecto fijo (mutant): \(\hat\beta=1.2366\), se = 0.2205 ⇒ HR=3.44, p ≈ 2×10⁻⁸.
  • Efecto aleatorio (variabilidad intrafamiliar): \(\widehat{\rm SD}(\text{famID})=0.591\) (Var=0.350).

16.2.4.1 Prueba de inclusión de frailty:

  • Tomamos la log-verosimilitud integrada \(\ell_{\rm int}=-3564.622\) y la comparamos con la verosimilitud del modelo sin frailty pero con mutant (\(\ell_1=-3566.745\) del paso 1):

    \[ G^2 = 2(\ell_1 - \ell_{\rm int}) = 2(-3566.745 + 3564.622) = 4.246, \]

    que contra \(\chi^2_1\) da p ≈ 0.039.

Interpretación: hay evidencia moderada de heterogeneidad adicional por familias.