Error specifying Cox PH survival model
Michael Smith
I've been trying to get a spatial term into a survival analysis using the following code, which probably has a beginner’s error somewhere:
# load spatial points dataframe with columns grid8@data
grid8 <-
readOGR(paste(getwd(), "/grid_10km.shp", sep = ""))
surv <- inla.surv(grid8@data$def_mean, grid8@data$status) # year and
censoring status
# this seems to work OK: 3 predictor variables
coxinla
<- inla(surv ~ 1 + c_mean + e_mean + p_mean, data = grid9,
family = "coxph",
control.hazard = list(hyper = list(prec = list(param = c(0.001, 0.001)))))
#set up mesh and spde using a spatial polygons dataframe named ‘boundary’
mesh1 <- inla.mesh.2d(boundary=boundary, max.edge = c(25, 250), cutoff = 3)
matern1 <- inla.spde2.pcmatern(mesh1,
prior.sigma = c(1, 0.01),
prior.range = c(25, 0.01))
A.index <- inla.spde.make.A(mesh = mesh1, loc = coordinates(grid8))
s.index <- inla.spde.make.index(name = "spatial.field",
n.spde = matern1$n.spde)
### inla.stack throws an error here: ###
c.stack <- inla.stack(data = list(surv = surv),
A = list(A.index, 1),
effects = list(c(s.index, list(Intercept = 1)),
list(c_mean = grid9$c_mean,
e_mean = grid9$e_mean,
p_mean = grid9$p_mean
)),
)
cox2 <- inla(surv ~ 1 + c_mean + e_mean + p_mean + f(spatial.field, model = matern1),
data = grid9, family = "coxph",
control.hazard = list(hyper = list(prec = list(param = c(0.001, 0.001)))))
...the error is:
Error in lapply(X = X, FUN = FUN, ...) :
argument is missing, with no default
> traceback()
3: lapply(X = X, FUN = FUN, ...)
2: sapply(list(...), function(x) inherits(x, "inla.data.stack"))
1: inla.stack(data = list(surv), A = list(A.index, 1, 1, 1), effects = list(c(effect = list(spatial = 1:matern1$n.spde),
list(Intercept = 1)), list(cost_mean = grid9$cost_mean, elev_mean = grid9$elev_mean,
popn_mean = grid9$popn_mean)), )
....can anyone help? I’m not great with stacks....
Helpdesk
can you attach the complete code and data used (at least something that
have this error). thx.
On Thu, 2021-06-03 at 08:26 -0700, 'Michael Smith' via R-inla discussion
> You received this message because you are subscribed to the Google
> Groups "R-inla discussion group" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to r-inla-discussion...@googlegroups.com.
> To view this discussion on the web, visit
> https://groups.google.com/d/msgid/r-inla-discussion-group/60b336b8-f21d-47cf-a663-72583e8139bbn%40googlegroups.com
> .
--
Håvard Rue
he...@r-inla.org
Michael Smith
Helpdesk
On Thu, 2021-06-03 at 10:18 -0700, 'Michael Smith' via R-inla discussion
> I'm having trouble posting the files so will try emailing the help @r-
> inla? Thank you....
>
> You received this message because you are subscribed to the Google
> Groups "R-inla discussion group" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to r-inla-discussion...@googlegroups.com.
> To view this discussion on the web, visit
Finn Lindgren
On 3 Jun 2021, at 16:26, 'Michael Smith' via R-inla discussion group <r-inla-disc...@googlegroups.com> wrote:
Hi everyone,
--
You received this message because you are subscribed to the Google Groups "R-inla discussion group" group.
To unsubscribe from this group and stop receiving emails from it, send an email to r-inla-discussion...@googlegroups.com.
To view this discussion on the web, visit https://groups.google.com/d/msgid/r-inla-discussion-group/60b336b8-f21d-47cf-a663-72583e8139bbn%40googlegroups.com.
Helpdesk
alternative, is to use
inla.coxph
which will transform the the coxph model into poisson regression, and
then you're back to a normal case.
see ?inla.coxph for details. this is what happens internally.
H
On Thu, 2021-06-03 at 21:54 +0100, Finn Lindgren wrote:
> Hi,
> Inla.stack doesn’t now how to handle the special response objects,
> like inla.surv() output. I think what you’d need to do (unless
> skipping inla.stack, and constructing the inputs manually instead) is
> one of these two:
> 1) call inla.stack using the inputs to inla.surv instead, and then
> call inla.surv using the inla.stack.LHS() contents, or
> 2) call inla.surv() first, and then extract its contents and use as
> input to inla.stack, and then extract the resulting vectors from the
> stack and replace in the surv object.
>
> I’m not sure which of these two methods is the correct one.
>
> When I’m back from annual leave (July) my plan is to work on support
> for these response objects in inlabru, and if the objects are
> sufficiently predictable I may be able to add direct support for them
> in inla.stack as well, or write the inlabru support in such a way that
> the support is encapsulated in a new function that can go in the inla
> package.
>
> Finn
>
> > On 3 Jun 2021, at 16:26, 'Michael Smith' via R-inla discussion group
> > <r-inla-disc...@googlegroups.com> wrote:
> >
he...@r-inla.org
Michael Smith
grid8@data$id2 <- 1:nrow(grid8@data) # set up ID same as rows
lattice_temp <- poly2nb(grid8, row.names = grid8@data$id2)
nb2INLA(paste(getwd(), "/Lattice.graph", sep = ""), lattice_temp)
lattice.adj <- paste(getwd(), "/Lattice.graph", sep = "")
intercept1 = rep(1, nrow(grid8@data))
grid9 <- grid8@data
coxinla = inla.coxph(surv ~ -1 + intercept1 + c_mean + e_mean + p_mean +
f(id2, model = "bym", graph = lattice.adj),
list(surv = surv, c_mean=grid9$c_mean, e_mean=grid9$el_mean,
p_mean=grid9$p_mean, id2=grid9$id2, intercept1 = intercept1))
r = inla(coxinla$formula,
family = coxinla$family,
data=c(as.list(coxinla$data), coxinla$data.list),
E = coxinla$E)