[R] Score test for a subset of parameters estimated with coxph

[Denise b]

An update regarding my question to compute the score test for a subset of
parameters with coxph,
when the option frailty is used:

After further explorations, based on my previous examples, I realized that
in presence of frailty term,
coxph returns only the results from the likelihood ratio test, whereas when
there is no frailty term in the model specified, it returns the likelihood
ratio, wald and score test as illustrated in the following example.

It could be due to the penalized-likelihood approach used by coxph for
frailty model in coxph.
Is there any other potential explanation?

# R code
library(survival)

fit1<-coxph(Surv(time,status)~sex+ ph.karno + pat.karno , data=lung)
summary(fit1)
Call:
coxph(formula = Surv(time, status) ~ sex + ph.karno + pat.karno,
    data = lung)

  n= 224, number of events= 161
   (4 observations deleted due to missingness)

               coef exp(coef)  se(coef)      z Pr(>|z|)
sex       -0.511878  0.599369  0.169275 -3.024  0.00249 **
ph.karno  -0.006155  0.993864  0.006822 -0.902  0.36695
pat.karno -0.017020  0.983124  0.006535 -2.604  0.00920 **
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

          exp(coef) exp(-coef) lower .95 upper .95
sex          0.5994      1.668    0.4301    0.8352
ph.karno     0.9939      1.006    0.9807    1.0072
pat.karno    0.9831      1.017    0.9706    0.9958

Concordance= 0.646  (se = 0.026 )
Rsquare= 0.097   (max possible= 0.999 )
Likelihood ratio test= 22.91  on 3 df,   p=4.208e-05
Wald test            = 22.64  on 3 df,   p=4.789e-05
Score (logrank) test = 23.2  on 3 df,   p=3.673e-05

frailtyfit1<-coxph(Surv(time,status)~sex+ ph.karno + pat.karno
+frailty(inst), data=lung)
summary(frailtyfit1)
Call:
coxph(formula = Surv(time, status) ~ sex + ph.karno + pat.karno +
    frailty(inst), data = lung)

  n= 223, number of events= 160
   (5 observations deleted due to missingness)

              coef     se(coef) se2      Chisq DF   p
sex           -0.50984 0.169534 0.169498 9.04  1.00 0.0026
ph.karno      -0.00620 0.006850 0.006842 0.82  1.00 0.3700
pat.karno     -0.01703 0.006544 0.006541 6.77  1.00 0.0093
frailty(inst)                            0.25  0.21 0.3600

          exp(coef) exp(-coef) lower .95 upper .95
sex          0.6006      1.665    0.4308    0.8373
ph.karno     0.9938      1.006    0.9806    1.0073
pat.karno    0.9831      1.017    0.9706    0.9958

Iterations: 5 outer, 21 Newton-Raphson
     Variance of random effect= 0.001494078   I-likelihood = -711.9
Degrees of freedom for terms= 1.0 1.0 1.0 0.2
Concordance= 0.649  (se = 0.026 )
Likelihood ratio test= 23.22  on 3.21 df,   p=4.703e-05

Le jeu. 15 août 2019 à 18:11, Denise b <denise.jeane.b using gmail.com> a écrit :

> Dear R help list,
>
> I would like to perform a score test for a subset of the parameters
> estimated with coxph using the frailty() option.
> As illustrated in the following reproducible example, I am able to perform
> the score test with the standard coxph() but not in presence of frailty()
> argument. See Examples 1 and 2
>
> library(survival)
>
> # Example 1: score test using standard coxph()
> # score test for ph.karno & pat.karno coefficients
> fit1<-coxph(Surv(time,status)~sex+ ph.karno + pat.karno , data=lung)
> fit0<-coxph(Surv(time,status)~pat.karno , data=lung) # fit under null
> scorestat <- coxph(Surv(time,status)~sex+ ph.karno + pat.karno,
> init=c(0,0,coefficients(fit0)) ,iter=0, data=lung)
> scorestat$score
> [1] 10.50409
>
> # Example 2: coxph with frailty term
> # score test for ph.karno & pat.karno coefficients
> frailtyfit1<-coxph(Surv(time,status)~sex+ ph.karno + pat.karno
> +frailty(inst), data=lung)
> frailtyfit0<-coxph(Surv(time,status)~sex+ frailty(inst), data=lung)
> scorefrailtyfit1<-coxph(Surv(time,status)~sex+ ph.karno + pat.karno
> +frailty(inst), data=lung,
> init=c(0,0,coefficients(frailtyfit0)) ,iter=0)
>
> I get the following error message:
> Error in coxph(Surv(time, status) ~ sex + ph.karno + pat.karno +
> frailty(inst),  :
>   wrong length for init argument
>
>
> When I checked in the code I found
> if (length(init) != ncol(X))
>             stop("wrong length for init argument")
> Which seems that the length of the init argument doesnt match the number of
> predictors in the model.
>
> It seems that the code is expecting that I specify another parameter for
> frailty term.
> So, I tried to specify 2 as an inital value for frailty but got another
> error message later:
>
> scorefrailtyfit1<-coxph(Surv(time,status)~sex+ ph.karno + pat.karno
> +frailty(inst), data=lung,
> + init=c(0,0,coefficients(frailtyfit0),4.5),iter=0)
> Error in coxpenal.fit(X, Y, strats, offset, init = init, control, weights =
> weights,  :
>   Wrong length for inital values
>
> I have also checked coxme() but it seems that there is no option
> implemented
> for the score test.
> Do you have an idea of what could be the problem?
>
> Thanks,
> Denise
>
>
>

	[[alternative HTML version deleted]]