example for calculating the cure fraction

[markdanese]

Apologies if I am missing this, but I am trying to estimate the cure proportion. I see that this can only be estimated using the predict method. But I am having trouble and example syntax would be helpful to make sure I am doing this properly.

fst <- 
    stpm2(
       Surv(days_to_death, event_flag) ~ tx_2L_01,
       tvc = list(tx_2L_01 = 3),
       df = 4,
       data = dtw, 
       cure = TRUE,
       weights = iptw,
       link.type = "PH")

predict(fst, type = "probcure")
Error in .subset2(x, i, exact = exact) :  attempt to select less than one element in get1index

[mclements]

The predict function for stpm2 needs newdata to be specified. Can you try, for example,

predict(fst, type = "probcure",

newdata=data.frame(tx_2L_01=1))

I would be interested to hear how you go.

Sincerely, Mark.

On 8 Jun 2018 8:07 pm, Mark Danese notifications@github.com wrote:

Apologies if I am missing this, but I am trying to estimate the cure proportion. I see that this can only be estimated using the predict method. But I am having trouble and example syntax would be helpful to make sure I am doing this properly.

fst <- stpm2( Surv(days_to_death, event_flag) ~ tx_2L_01, tvc = list(tx_2L_01 = 3), df = 4, data = dtw, cure = TRUE, weights = iptw, link.type = "PH")

predict(fst, type = "probcure") Error in .subset2(x, i, exact = exact) : attempt to select less than one element in get1index

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHubhttps://github.com/mclements/rstpm2/issues/9, or mute the threadhttps://github.com/notifications/unsubscribe-auth/AAiBCcCL14NByrYAJP23ymfFIsLWmLAhks5t6r14gaJpZM4Ugu31.

[markdanese]

I just ran this: predict(fst, type = "probcure", newdata = data.frame(tx_2L_01 = 1)) and this was the error. Error in nsx(log(days_to_death), knots = c(4.1190041246092, 4.83628190695148, : object 'days_to_death' not found

I tried adding "days_to_death = 365" to the data.frame but that returned a different error.

[mclements]

Apologies: it's late Friday evening here.

You did either specify the time variable or use grid=TRUE.

Can I ask: what was the error when you included days_to_death?

--- Mark

On 8 Jun 2018 9:15 pm, Mark Danese notifications@github.com wrote:

I just ran this: predict(fst, type = "probcure", newdata = data.frame(tx_2L_01 = 1)) and this was the error. Error in nsx(log(days_to_death), knots = c(4.1190041246092, 4.83628190695148, : object 'days_to_death' not found

I tried adding "days_to_death = 365" to the data.frame but that returned a different error.

— You are receiving this because you commented. Reply to this email directly, view it on GitHubhttps://github.com/mclements/rstpm2/issues/9#issuecomment-395861608, or mute the threadhttps://github.com/notifications/unsubscribe-auth/AAiBCR1w9ZHljonl9ciIjk_Bulqw9BsBks5t6s1ZgaJpZM4Ugu31.

[LasseHjort]

Excuse me for interfering.

type = "probcure" will not provide the cure proportion, but the conditional probability of cure given survival until some time point. To obtain the cure proportion, you simply use predict(fst, newdata = data.frame(tx_2L_01 = 1, days_to_death = last_knot_of_splines).

Does this help?

Best, Lasse

[markdanese]

Thanks to both of you for your help. It is not a rush. @mclements the error was Error in .subset2(x, i, exact = exact) : attempt to select less than one element in get1index

@LasseHjort thank you. I would be interested in estimating the function over a variety of times.

Full run below:

> fst <- stpm2(
+     Surv(days_to_death, event_flag) ~ 
+         tx_2L_01,
+     tvc = list(tx_2L_01 = 3),
+     df = 4,
+     data = dtw, 
+     cure = TRUE,
+     weights = iptw,
+     link.type = "PH")
> summary(fst)
Maximum likelihood estimation

Call:
mle2(minuslogl = negll, start = coef, eval.only = TRUE, vecpar = TRUE, 
    gr = function (beta) 
    {
        localargs <- args
        localargs$init <- beta
        localargs$return_type <- "gradient"
        return(.Call("model_output", localargs, PACKAGE = "rstpm2"))
    }, control = list(parscale = c(1, 1, 1, 1, 1, 1, 1, 1, 1), 
        maxit = 300), lower = -Inf, upper = Inf)

Coefficients:
                                                       Estimate Std. Error z value               Pr(z)
(Intercept)                                            -9.79046    1.23379 -7.9353 0.00000000000000210
tx_2L_01                                                4.04954    1.26358  3.2048           0.0013516
nsx(log(days_to_death), df = 4, cure = TRUE)1           9.30146    1.22476  7.5945 0.00000000000003089
nsx(log(days_to_death), df = 4, cure = TRUE)2           3.00399    0.25280 11.8830           < 2.2e-16
nsx(log(days_to_death), df = 4, cure = TRUE)3          19.20863    2.72692  7.0441 0.00000000000186688
nsx(log(days_to_death), df = 4, cure = TRUE)4           8.02233    0.78228 10.2550           < 2.2e-16
tx_2L_01:nsx(log(days_to_death), cure = TRUE, df = 3)1 -0.97615    0.29785 -3.2773           0.0010481
tx_2L_01:nsx(log(days_to_death), cure = TRUE, df = 3)2 -9.10907    2.83778 -3.2099           0.0013277
tx_2L_01:nsx(log(days_to_death), cure = TRUE, df = 3)3 -2.81481    0.74911 -3.7575           0.0001716

-2 log L: 6176.666 
> cure <- predict(fst, type = "probcure", newdata = data.frame(tx_2L_01 = 1, days_to_death = 500))
Error in .subset2(x, i, exact = exact) : 
  attempt to select less than one element in get1index

[mclements]

Mark: can you call traceback() after the error, please?

On 8 Jun 2018 10:34 pm, Mark Danese notifications@github.com wrote:

Thanks to both of you for your help. It is not a rush. @mclementshttps://github.com/mclements the error was Error in .subset2(x, i, exact = exact) : attempt to select less than one element in get1index

@LasseHjorthttps://github.com/LasseHjort thank you. I would be interested in estimating the function over a variety of times.

Full run below:

fst <- stpm2(

• Surv(days_to_death, event_flag) ~
• tx_2L_01,
• tvc = list(tx_2L_01 = 3),
• df = 4,
• data = dtw,
• cure = TRUE,
• weights = iptw,
• link.type = "PH") summary(fst) Maximum likelihood estimation

Call: mle2(minuslogl = negll, start = coef, eval.only = TRUE, vecpar = TRUE, gr = function (beta) { localargs <- args localargs$init <- beta localargs$return_type <- "gradient" return(.Call("model_output", localargs, PACKAGE = "rstpm2")) }, control = list(parscale = c(1, 1, 1, 1, 1, 1, 1, 1, 1), maxit = 300), lower = -Inf, upper = Inf)

Coefficients: Estimate Std. Error z value Pr(z) (Intercept) -9.79046 1.23379 -7.9353 0.00000000000000210 tx_2L_01 4.04954 1.26358 3.2048 0.0013516 nsx(log(days_to_death), df = 4, cure = TRUE)1 9.30146 1.22476 7.5945 0.00000000000003089 nsx(log(days_to_death), df = 4, cure = TRUE)2 3.00399 0.25280 11.8830 < 2.2e-16 nsx(log(days_to_death), df = 4, cure = TRUE)3 19.20863 2.72692 7.0441 0.00000000000186688 nsx(log(days_to_death), df = 4, cure = TRUE)4 8.02233 0.78228 10.2550 < 2.2e-16 tx_2L_01:nsx(log(days_to_death), cure = TRUE, df = 3)1 -0.97615 0.29785 -3.2773 0.0010481 tx_2L_01:nsx(log(days_to_death), cure = TRUE, df = 3)2 -9.10907 2.83778 -3.2099 0.0013277 tx_2L_01:nsx(log(days_to_death), cure = TRUE, df = 3)3 -2.81481 0.74911 -3.7575 0.0001716

-2 log L: 6176.666

cure <- predict(fst, type = "probcure", newdata = data.frame(tx_2L_01 = 1, days_to_death = 500)) Error in .subset2(x, i, exact = exact) : attempt to select less than one element in get1index

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHubhttps://github.com/mclements/rstpm2/issues/9#issuecomment-395882428, or mute the threadhttps://github.com/notifications/unsubscribe-auth/AAiBCcLkK211zEvsT44b4Bs-sOUqJWTJks5t6t-6gaJpZM4Ugu31.

[markdanese]

with days_to_death:

8: (function(x, i, exact) if (is.matrix(i)) as.matrix(x)[[i]] else .subset2(x, 
       i, exact = exact))(x, ..., exact = exact)
7: `[[.data.frame`(data, var[i])
6: data[[var[i]]]
5: exposed(newdata)
4: predict.stpm2.base(object = object, newdata = newdata, type = type, 
       grid = grid, seqLength = seqLength, se.fit = se.fit, link = link, 
       exposed = exposed, var = var, keep.attributes = keep.attributes, 
       use.gr = use.gr, level = level, type.relsurv = type.relsurv, 
       scale = scale, rmap = substitute(rmap), ratetable = ratetable, 
       ...)
3: .local(object, ...)
2: predict(fst, type = "probcure", newdata = data.frame(tx_2L_01 = 1, 
       days_to_death = 300))
1: predict(fst, type = "probcure", newdata = data.frame(tx_2L_01 = 1, 
       days_to_death = 300))

with grid=TRUE:

8: (function(x, i, exact) if (is.matrix(i)) as.matrix(x)[[i]] else .subset2(x, 
       i, exact = exact))(x, ..., exact = exact)
7: `[[.data.frame`(data, var[i])
6: data[[var[i]]]
5: exposed(newdata)
4: predict.stpm2.base(object = object, newdata = newdata, type = type, 
       grid = grid, seqLength = seqLength, se.fit = se.fit, link = link, 
       exposed = exposed, var = var, keep.attributes = keep.attributes, 
       use.gr = use.gr, level = level, type.relsurv = type.relsurv, 
       scale = scale, rmap = substitute(rmap), ratetable = ratetable, 
       ...)
3: .local(object, ...)
2: predict(fst, type = "probcure", newdata = data.frame(tx_2L_01 = 1), 
       grid = TRUE)
1: predict(fst, type = "probcure", newdata = data.frame(tx_2L_01 = 1), 
       grid = TRUE)```

[LasseHjort]

You can obtain the conditional cure probability given survival until specific time points, say 10 and 50 days, using the following code:

exposed <- function(newdata){
  newdata$days_to_death <- last_knot_of_splines
  newdata
}

cure <- predict(fst, type = "probcure", 
                newdata = data.frame(tx_2L_01 = 1, days_to_death = c(10,  50)), 
                exposed = exposed)

I just updated updated the code on the develop branche, so you will have to reinstall the package to make it work.

Lasse

[markdanese]

I installed your commit and it looks like it works. For some reason, I can't seem to find the documentation for the predict function to understand what "exposed" is doing. Maybe it isn't in this branch.

I assume you mean the boundary knot when you say "last knot of spline" or do you mean the last non-boundary knot?

If I interpret this correctly, the population "cure fraction" would be the value at about day 1. This would be the expected proportion cured given treatment initiation. This probability increases over time, because it is the conditional probability given survival until a given time.

> exposed <- function(newdata){
+   newdata$days_to_death <- exp(6.81563999007433)
+   newdata
+ }
> predict(fst, type = "probcure", 
+                 newdata = data.frame(tx_2L_01 = 1, days_to_death = c(10,  50)), 
+                 exposed = exposed)
[1] 0.3560970 0.4015381
> predict(fst, type = "probcure", 
+                 newdata = data.frame(tx_2L_01 = 1, days_to_death = c(10,  50, 500)), 
+                 exposed = exposed)
[1] 0.3560970 0.4015381 0.9344874