I am trying to compare different Cox models to choose the “best” model. I am working with NHANES data and therefore use the survey package. With the svycoxph function the lrtest(), BIC(), AIC() and anova() don’t work. There is the function regTermTest() in the survey package. I am not really sure if I can use this function with the method=”LRT” to compare models or how it works.
library("survey")
NHANES_all <- svydesign(data=data_all, id=~Total.PSU, strata=~Total.stratum, weights=~Total.final.weight, nest=TRUE)
NHANES_subset <- subset(NHANES_all, cohortsubset==1)
M3<-svycoxph(Surv(PERMTH_INT, MORTSTAT==1)~relevel(group, ref = "2") + Age.at.interview..screener....qty + relevel(Sex, ref = "2") + Race.ethnicity + relevel(education, ref = "2") + relevel(pirgroup, ref = "2")+ relevel(Marital.status, ref = "3") + AlcFactCox + smoke + caffeinDec + LTPA + Body.mass.index + relevel(diabetes, ref = "2") + relevel(highBP, ref = "2") + eGFRgroup, subset(NHANES_subset, HS.primary.finding..recode==1))
regTermTest(M3, ~AlcFactCox + smoke, method="LRT", df=NULL)
> regTermTest(M3, ~AlcFactCox + smoke, method="LRT", df=NULL)
Working (Rao-Scott+F) LRT for AlcFactCox smoke
in svycoxph(formula = Surv(PERMTH_INT, MORTSTAT == 1) ~ relevel(group,
ref = "2") + Age.at.interview..screener....qty + relevel(Sex,
ref = "2") + Race.ethnicity + relevel(education, ref = "2") +
relevel(pirgroup, ref = "2") + relevel(Marital.status, ref = "3") +
AlcFactCox + smoke + caffeinDec + LTPA + Body.mass.index +
relevel(diabetes, ref = "2") + relevel(highBP, ref = "2") +
eGFRgroup,
design = subset(NHANES_subset, HS.primary.finding..recode ==
1))
Working 2logLR = 206.8198 p= 2.52e-07
(scale factors: 1.3 1 0.69 ); denominator df= 11
> regTermTest(M3, ~AlcFactCox + smoke, method="Wald", df=NULL)
Wald test for AlcFactCox smoke
in svycoxph(formula = Surv(PERMTH_INT, MORTSTAT == 1) ~ relevel(group,
ref = "2") + Age.at.interview..screener....qty + relevel(Sex,
ref = "2") + Race.ethnicity + relevel(education, ref = "2") +
relevel(pirgroup, ref = "2") + relevel(Marital.status, ref = "3") +
AlcFactCox + smoke + caffeinDec + LTPA + Body.mass.index +
relevel(diabetes, ref = "2") + relevel(highBP, ref = "2") +
eGFRgroup,
design = subset(NHANES_subset, HS.primary.finding..recode ==
1))
F = 42.20666 on 3 and 11 df: p= 2.5134e-06
In my example above with the method="LRT": do I compare the model M3 with its confounders with the model M3 without the confounders smoke and AlcfactCox? If I use the method=”Wald”, do I test if the estimated effect of both my tested confounders (smoke and AlcfactCox) is zero?
So basically my question is if the regTermTest() with method="LRT" is the right way to compare different svycoxph models or do you have any suggestions what is the best method or function to compare svycoxph models?
If I use the summary() function on my model, I get a concordance. But in my understanding, concordance is not sensitive enough for model comparison.
Thank you very much in advance!