On Jul 30, 2010, at 11:53 AM, Biau David wrote:
Thx for the answer.
I am using survival.
I didn't know that the Wald and score tests were the same for
individual variables in a coxph; I Thought the score test was the
"multivariate version" of the Log-rank.
However, say I have only one variable in the model, I should expect
the test for the full model and the one for a single variable to be
the same?
I don't understand what you are saying. What full model?
Then it seems to me that the default test is the Wald and that the
Wald and the Score are different.
Numerically it would appear.
> cox_lr_age <- coxph(Surv(tilr, ev_lr==1)~age, data=tam)
> summary(cox_lr_age)
Call:
coxph(formula = Surv(tilr, ev_lr == 1) ~ age, data = tam)
n=2156 (76 observations deleted due to missingness)
coef exp(coef) se(coef) z Pr(>|z|)
age 0.019504 1.019696 0.004651 4.193 2.75e-05 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
exp(coef) exp(-coef) lower .95 upper .95
age 1.020 0.9807 1.010 1.029
Rsquare= 0.008 (max possible= 0.669 )
Likelihood ratio test= 18.4 on 1 df, p=1.787e-05
Wald test = 17.58 on 1 df, p=2.751e-05
Score (logrank) test = 17.86 on 1 df, p=2.375e-05
You can find the methods used by executing:
#require(survival)
getAnywhere(summary.coxph)
You had not earlier indicated what functions and output you had
questions about. The Score test to which you refer above has the
additional parenthetical indication that it is a logrank statistic.
Looking in Therneau and Gramsch's, "Modeling Survival Data" in the
index, we find that an approximate intersection of cited pages for
score test and log rank-test includes page 53-54 where we see the
explanation that a) Wald test is deprecated as the least reliable and
b) that a Wald test based on one-iteration is closely related to the
"efficient score test" (but it is not clear that this is the test
referred to in this output. Reading between the lines it appears that
the score test above might not actually be from a Cox model but is
rather from a log-rank procedure. I hope Therneau can clarify.
--
David.
David Biau.
De : David Winsemius <dwinsem...@comcast.net>
À : Biau David <djmb...@yahoo.fr>
Cc : r help list <r-help@r-project.org>
Envoyé le : Ven 30 juillet 2010, 17h 34min 28s
Objet : Re: [R] COXPH: how to get the score test and likelihood
ratio test for a specific variable in a multivariate Coxph ?
On Jul 30, 2010, at 11:08 AM, Biau David wrote:
> Hello,
>
> I would like to get the likelihood ratio and score tests for
specific variables
> in a multivariate coxph model. The default is Wald, so the tests
for each
> separate variable is based on Wald's test. I have the other tests
for the full
> model but I don't know how to get them for each variable.
>
> Any idea?
>
The first idea would be to specify which function in which package
you are asking questions about. In the case of coxph in the survival
package, for instance, you do get a likelihood ratio test (==
differences in log-likelihoods) by default. A score test is, at
least as as I understand it for individual variables, equivalent to
a Wald test, so I don't really understand your question, since youa
re already getting all of that in the survival package.
(You can extract a "score" value and loglik values from a coxph
object by:
(with the first example in the coxph help page)
coxph(Surv(time, status) ~ x + strata(sex), test1)$score
xoxph(Surv(time, status) ~ x + strata(sex), test1)$loglik
But anova(coxph-object) would give you these values in a neater
bundle.
#Analysis of Deviance Table
# Cox model: response is Surv(time, status)
#Terms added sequentially (first to last)
# loglik Chisq Df Pr(>|Chi|)
# NULL -3.8712
# x -3.3277 1.0871 1 0.2971
The question about "getting them for each variable" does not make a
lot of sense to me, since likelihood tests are model comparisons.
You can only make such statements about the consequences of adding
or deleting a variable to/from an existing model.
--David Winsemius, MD
West Hartford, CT
David Winsemius, MD
West Hartford, CT
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.