CODEWITHSUNDEEP
rsurvival
gibsonmate
gibsonmate

Reputation: 25

How to calculate restricted mean survival time for covariates in R?

I am trying to calculate the restricted mean survival time grouped by the covariates but I get an error. Could you please have a look at the code example code below and let me know what the problem is?

library(survival)
library(flexsurv)

data <- survival::lung
formula <- Surv(data$time, data$status) ~data$sex
fsmodel_exponential<-flexsurvreg(formula,dist = "exp")

#Produces the expected results
rate_exponential<-fsmodel_exponential$res[1]
rmst_exp <- rmst_exp(t = 30, rate = rate_exponential, start = 0)
rmst_exp

#Gives error for the covariate
rate_exponential_sex<- fsmodel_exponential$res[2]
rmst_exp2 <- rmst_exp(t = 30, rate = rate_exponential_sex, start = 0)
rmst_exp2

Upvotes: 0

Views: 1242

Answers (1)

user12728748
user12728748

Reputation: 8506

Your fsmodel_exponential$res[2] is negative, which causes an error. In the exponential model covariates on these parameters represent an accelerated failure time (AFT) model. You do get those when you invert the factor levels of sex, which results in positive maximum likelihood estimates.

library(survival)
library(flexsurv)
data <- lung
data$sex <- relevel(factor(data$sex, labels = c("M", "F")), "F")
formula <- with(data, Surv(time, status) ~ sex)
fsmodel_exponential <- flexsurvreg(formula, dist = "exp")
rmst_exp(t = 30, rate = fsmodel_exponential$res[1], start = 0)
#> [1] 29.23162
rmst_exp(t = 30, rate = fsmodel_exponential$res[2], start = 0)
#> [1] 1.998406

plot(fsmodel_exponential, col = c("blue", "red"), 
    lwd.obs = 2, xlab = "months", ylab = "Recurrence-free survival")
legend("topright", levels(data$sex), col=c("red", "blue"), lty=1)

Edit: To get the RMST, simply run:

summary(fsmodel_exponential, type = "rmst", t= 30 )
#> sex=M 
#>   time     est      lcl      ucl
#> 1   30 28.7467 28.49528 28.94901
#> 
#> sex=F 
#>   time      est      lcl      ucl
#> 1   30 29.23162 28.98571 29.40346

Upvotes: 1

Related Questions