inlabru-org / inlabru

inlabru
https://inlabru-org.github.io/inlabru/
76 stars 21 forks source link

Add a proper inlabru rerun method (was: inla.rerun() impact on predict()) #78

Closed ghost closed 3 years ago

ghost commented 4 years ago

I have been experimenting with specifying temporal effects using the 1D SPDE based approach within a lgcp() model. The model has 3 fixed effects and 3 random effects (Spatial SPDE, Temporal 1D SPDE, Seasonal IID)

Whilst fitting this model, the iteration status messages return some Hessian warnings and one of the fixed coefficients had a clearly wild/erroneous estimate.

I believe using inla.rerun() may help in the case of the Hessian warnings?

After running inla.rerun(model.fit), coefficient estimates seem more reasonable.

However, now predict() does not function as before:

The predict output is split into a multi list of class “prediction” and “list”, named “ID” and “value”

predict(fit.model, ppxl.data.frame, ~ exp( Intercept + spat.eff + time.eff + season.eff + cov.A + cov.B + cov.C), n.samples = 200) returns the following error message: Error in Ops.data.frame(Intercept + spat.eff, time.eff) : ‘+’ only defined for equally-sized data frames

I get a single prediction if I only include fixed effects , e.g. mean sd q0.025 median q0.975 smin smax cv var 1 1.269469 0.01517862 1.248578 1.276495 1.287422 1.248445 1.288849 0.01195667 0.0002303906

I get a number of predictions equal to the number of levels within that random effect if I include either the temporal or seasonal effects. Including both returns the equally-sized data frame error.

If I include the spatial effect, I get a number of predictions, but I can’t see how that number relates to anything.

I was wondering if you could clarify what exactly is going on here? And whether there is an inlabru specific solution to running inla.rerun()?

Package versions: inlabru_2.1.13.999 INLA_20.07.12

finnlindgren commented 4 years ago

Unfortunately inla.rerun() forgets the inlabru-specific information, and the inlabru::predict.inla() method that tries to convert a pure inla result to a bru-object that inlabru::predict.bru() works for was a mistake that has never properly worked. I plan on adding an inla.rerun wrapper to inlabru to solve this problem.

However, the wrapper might be as simple as copying the old inlabru-specific information over to the new object. Can you try the following approach and report here whether it works?

result <- lgcp(...)
result2 <- inla.rerun(result)
result2$sppa <- result$sppa # Copying the internal bru-information
class(result2) <- class(result)
predict(result2, ...)

(First double-check that names(result) contains the "sppa"; at some point the name of that entry will change, but I think your version still uses that one.)

Also note that the default value for control.inla=list(h = ...) is being lowered from 0.01 to 0.005, to help avoid the root cause of the problem that is improved by inla.rerun().

finnlindgren commented 4 years ago

A proper inlabru version of inla.rerun would also restart the non-linear iterations; inla.rerun will only rerun the final approximation step. Lowering the h-setting will make the whole nonlinear inlabru iteration more robust. Fortunately, if your LGCP model only uses a pure log-linear predictor, the model doesn't need non-linear iterations. (It's when parameterised distance sampling models are included the non-linear iterations are needed.)

ghost commented 4 years ago

My version contains the "sppa" information.

Copying this over to the output of inla.rerun() and redefining the class appears to have worked as intended!

I havent had a distance sampling component in this model but i will keep it in mind. I will setup another run lowering the h-value nonetheless.

Cheers, H

finnlindgren commented 3 years ago

Note: In a recent devel branch update, the $sppa object changed name from $bru_info, and also changed some internal storage structure; the methods that use the estimation object will attempt to detect and upgrade stored objects from the old structure to the new, but the manual method for inla.rerun() in the thread above will need to use the bru_info name for newly generated results.

finnlindgren commented 3 years ago

inlabru now has a bru_rerun() method (added in 0227028ffa9e41cffcc1f4bee44230d3b268b56f) that takes a previous bru() or lgcp() result as starting point, as well as options, that can be used to modify algorithm parameters for the rerun (e.g. changing inla int.strategy, or setting bru_max_iter=1 to just do one more inla iteration.