anova doesn't works with coxph function

[alexlaf314]

Good Morning,

I try your function anova.mitml.result with two Cox models (one empty and one with 1 variable), but the function doesn't work. Can you solve this bug

A simple exemple with the data lung from the package survival

library(Amelia)
library(mitml)
library(survival)
data(lung)

lung_mi<-amelia(lung,noms=c(1,3,5,6))
lung_mi2<-amelia2mitml.list(lung_mi)

fit0<-with(lung_mi2,coxph(Surv(time,status)~1))
fit1<-with(lung_mi2,coxph(Surv(time,status)~meal.cal))
anova.mitml.result(fit1,fit0)

Error in max(df) : invalid 'type' (list) of argument

Thank you for advance

Thank you for your work! it will be very usefull in my job !

[simongrund1]

Thanks for pointing this out. I am unable to reproduce this exact behavior. Instead, the problem seems to start at the call to with, which tells me:

Error in Surv(time, status) : Time variable is not numeric

Could you provide the exact code that produces this error message? Try initializing the random number generator with a fixed value to make the results reproducible, for example, with:

set.seed(1234)

[alexlaf314]

Sorry for this bad exemple. An other that I hope valid this time :

library(Amelia)
library(mitml)
library(survival)
data(lung)
set.seed(9)

lung_mi<-amelia(lung[,c(2,3,10)],noms=c(2))

-- Imputation 1 --

  1  2  3

-- Imputation 2 --

  1  2  3

-- Imputation 3 --

  1  2  3

-- Imputation 4 --

  1  2  3

-- Imputation 5 --

  1  2  3

lung_mi2<-amelia2mitml.list(lung_mi)

fit0<-with(lung_mi2,coxph(Surv(time,status)~1))
fit1<-with(lung_mi2,coxph(Surv(time,status)~wt.loss))
anova.mitml.result(fit1,fit0)
Error in max(df) : invalid 'type' (list) of argument

[simongrund1]

Thanks for the updated example! The error message is due to the fact that likelihood-based comparisons were not supported for coxph type models. However, I uploaded some changes that should add an initial support (d553679).

For example, the code above should now give:

Call:

anova.mitml.result(object = fit1, fit0)

Model comparison calculated from 5 imputed data sets.
Combination method: D2 (likelihood)

Model 1: Surv(time,status)~wt.loss
Model 2: Surv(time,status)~1

           F.value      df1      df2    P(>F)      RIV
  1 vs 2:    0.034        1 8578.430    0.855    0.022

Warnmeldung:
In anova.mitml.result(fit1, fit0) :
  The 'D3' method is currently not supported for models of class 'coxph'. Switching to 'D2'.

Note that method="D3" is still not supported, so anova switches to method="D2". This is intended because support for D3 is more difficult to implement and currently not a priority.

However, I also added a fix to testModels which should allow a better handling of "empty" models, including Wald-based comparisons of different models (2f04e8c). Because D1 tends to be more robust than D2, I recommend that you give it a try.

The respective command is simply (method="D1" is the default here):

testModels(fit1, fit0)

which gives:

Call:

testModels(model = fit1, null.model = fit0)

Model comparison calculated from 5 imputed data sets.
Combination method: D1

    F.value      df1      df2    P(>F)      RIV
      0.031        1 1569.748    0.861    0.053

Unadjusted hypothesis test as appropriate in larger samples.

I recommend that you stick with D1 whenever possible. If you decide to use D2, make sure that you specify a sufficient number of imputations (e.g., 100) because in my experience D2 tends to be more robust that way.

Let me know if the problem is solved for you. I rarely work with coxph, so suggestions are appreciated.

[alexlaf314]

Thanks you for your reactivity. In fact coxph object have not exactly the same format than a glm object for exemple. Unfortunately, I reuploaded mitml package and test the same code, and I again bugs.

library(Amelia)
library(mitml)
library(survival)
data(lung)
set.seed(9)
lung_mi<-amelia(lung[,c(2,3,10)],noms=c(2))

-- Imputation 1 --

  1  2  3

-- Imputation 2 --

  1  2  3

-- Imputation 3 --

  1  2  3

-- Imputation 4 --

  1  2  3

-- Imputation 5 --

  1  2  3

lung_mi2<-amelia2mitml.list(lung_mi)
fit0<-with(lung_mi2,coxph(Surv(time,status)~1))
fit1<-with(lung_mi2,coxph(Surv(time,status)~wt.loss))
anova(fit1,fit0)

Error in max(df) : invalid 'type' (list) of argument
In addition: Warning messages:
1: In anova.mitml.result(fit1, fit0) :
  The 'D3' method is currently not supported for models of class 'coxph'. Switching to 'D2'.
2: In df[mm] <- switch(lr.method, lmer = attr(logLik(modlist[[mm]][[1]]),  :
  number of items to replace is not a multiple of replacement length
3: In is.na(object$coefficients) :
  is.na() applied to non-(list or vector) of type 'NULL'
testModels(fit1,fit0)
Error in sapply(model, coef)[dpar, ] : incorrect number of dimensions

[simongrund1]

Your code runs fine on my machine with the current version from GitHub. Please make sure that you install the most recent development version from GitHub (not CRAN).

You can install the GitHub version as follows.

library(devtools)
install_github("simongrund1/mitml", force=TRUE)

Let me know if the errors persist.

[alexlaf314]

I'm sorry again, It's my first post on github. I try your modification and it's work. Unbelievable.

Thanks you again.

[simongrund1]

No problem, glad it works for you now! Let me know if anything isn't as it should be.