# [R] Efficient way to update a survival model

Frank S. |_j_rod @end|ng |rom hotm@||@com
Fri Aug 30 11:54:16 CEST 2019

```Hi everyone,

Vito, perhaps my previous mail was not clear.  It is true that I used a loop, but the key point is that such a loop
cannot compute the desired result. For example, for k = 3 the following loop

Cox <- list()
Cox[[1]] <- coxph(Surv(time,status == 2) ~ v + cos(v), data =  pbc)
for (k in 2:10) {
Cox[[k]] <- update(Cox[[k-1]], . ~ . + cos(k * v), data =  pbc)
}

leads to a model Cox[[3]] which accounts for terms {v, cos(v), cos(3*v)}, but does not include the term cos(2*v).
I think that this could be one way to solve my question:

library("survival")
set.seed(1)
v <- runif(nrow(pbc), min = 0, max = 2)
Cox0 <- coxph(Surv(pbc\$time,pbc\$status == 2) ~ v, data =  pbc)
k.max <- 9
Z <- outer(v, 1:k.max, function (x, y) {sin(x * y)})  # Matrix with the outer product of the two arrays

Cox <- list()
for (k in 1:k.max){
Cox[[k]] <-
update(Cox0, substitute(. ~ . + Z[, 1:k]), data =  pbc)
attr(Cox[[k]]\$coefficients, "names")[2:(k+1)] <- paste0("sin(", 1:k, "* v)")
}
Cox

Best,

Frank

________________________________
De: Frank S. <f_j_rod using hotmail.com>
Enviado: jueves, 29 de agosto de 2019 12:38
Para: Vito Michele Rosario Muggeo <vito.muggeo using unipa.it>
Cc: r-help using r-project.org <r-help using r-project.org>
Asunto: RE: [R] Efficient way to update a survival model

Hi Vito,

Cox[[3]] <- update(Cox[[2]], . ~ . + cos(3 * v), init=c(coef(Cox[[1]]), 0, 0), data =  pbc)
Cox[[3]] <- update(Cox[[2]], . ~ . + cos(3 * v), data =  pbc)

and both expressions lead to the same result. Is that OK?

Additionally, in my original question I wondered about the possibility of reducing the
10 lines of code to one general expression or some  loop. Is it possible?

Best,

Frank
________________________________
De: Vito Michele Rosario Muggeo <vito.muggeo using unipa.it>
Enviado: jueves, 29 de agosto de 2019 8:54
Para: Frank S. <f_j_rod using hotmail.com>
Cc: r-help using r-project.org <r-help using r-project.org>
Asunto: Re: [R] Efficient way to update a survival model

dear Frank,

update() does not update actually.. It just builds a new call which is
evaluated. To speed up the procedure you could try to supply starting
values via argument 'init'. The first values come from the previous
fit, and the last one referring to new coefficients is set to zero (or
any other appropriate value).

Something like (untested), for instance

update(Cox[[2]], . ~ . + cos(3 * v), init=c(coef(Cox[[1]]),0), data =  pbc)

Hope this helps,
best,
vito

"Frank S." <f_j_rod using hotmail.com> ha scritto:

> Hello everybody, I come with a question which I do not know how to
> conduct in an efficient way. In order to
> provide a toy example, consider the dataset "pbc" from the package
> "survival". First, I fit the Cox model "Cox0":
>
> library("survival")
> set.seed(1)
> v <- runif(nrow(pbc), min = 0, max = 2)
> Cox0 <- coxph(Surv(pbc\$time,pbc\$status == 2) ~ v, data =  pbc)
>
> Then, from the above model, I can fit recursively 10 additional models as:
>
> Cox <- list()
>
> Cox[[1]] <- update(Cox0, . ~ . + cos(1 * v), data =  pbc)
> Cox[[2]] <- update(Cox[[1]], . ~ . + cos(2 * v), data =  pbc)
> Cox[[3]] <- update(Cox[[2]], . ~ . + cos(3 * v), data =  pbc)
> Cox[[4]] <- update(Cox[[3]], . ~ . + cos(4 * v), data =  pbc)
> ...
> Cox[[10]] <- update(Cox[[9]], . ~ . + cos(10* v), data =  pbc)
>
> Since in practice I have to repeat above step until Cox[[100]], say,
> do you know an efficient way to
> wrap this code chunk in a loop or similar?
>
>
> set.seed(1)
> v <- runif(nrow(pbc), min = 0, max = 2)
> Cox0 <- coxph(Surv(pbc\$time,pbc\$status == 2) ~ v, data =  pbc)
>
> Cox <- list()
> Cox[[1]] <- update(Cox0, . ~ . + cos(1 * v), data =  pbc)
> for (k in 1:10) {
>   Cox[[k + 1]] <- update(Cox[[k]], . ~ . + cos((k + 1) * v), data =  pbc)
> }
>
> However, from Cox[[3]] onwards, the intermediate values of integer k
> are not included here (for
>  instance, the model Cox[[10]] would only include the cosinus terms
> for cos(1*v) and cos(10*v)).
>
[[elided Hotmail spam]]
>
> Frank
>
>        [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help using r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help