When is spatial range > study extent a problem in INLA? How to diagnose?

[Courtney Taylor]

I'm fitting INLA spatial models across multiple species with multiple spatially-structured random effects. The models converge (mode.status = 0) and most show better WAIC than non-spatial models, but I'm getting some spatial range estimates that exceed my study extent (3-4 times the study extent) and I'm trying to understand if/when this is a problem.

[INLA help]

So the spde acts more like an intercept? 

—

Haavard Rue

Helpdesk

he...@r-inla.org

On 28 Jan 2026 at 09:28 +0300, Courtney Taylor <cou...@gmail.com>, wrote:

I'm fitting INLA spatial models across multiple species with multiple spatially-structured random effects. The models converge (mode.status = 0) and most show better WAIC than non-spatial models, but I'm getting some spatial range estimates that exceed my study extent (3-4 times the study extent) and I'm trying to understand if/when this is a problem.

--
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, visit https://groups.google.com/d/msgid/r-inla-discussion-group/5f0be926-7598-4d2a-9610-80495af7b7bdn%40googlegroups.com.

[Elias T. Krainski]

In that case, it is recommended to use 

  inla.spde2.pcmatern(..., constr = TRUE)

To view this discussion, visit https://groups.google.com/d/msgid/r-inla-discussion-group/9b64f708-472f-4732-a254-4f02525f6e41%40Spark.

[Finn Lindgren]

Whether it is a problem or not depends on what you're going to do with the estimates; if you only need spatial prediction/interpolation, and don't care about the range parameter itself, then things are usually ok, unless you have so many competing components that the resulting non-identifiability becomes a numerical problem. Otherwise, one can just look only at intercept+field, which is typically much more well-constrained by data than the field component itself.

Sometimes, adding an integrate-to-zero constraint can be helpful, but beware that it does change the model (and makes it more dependent on the specific mesh), and you still have to look at intecept+field to get an interpretable quantity.

Since you say you have multiple spatial structured random effects, you may need to be more careful with those posterior interpretations, and consider which combined quantities are well defined/identifiable.

If the range is very large, then you will run into boundary effect issues, even with the integrate-to-zero constraint, as the spde matern model imposes normal-derivative-zero at the boundary. To avoid this, make the outer extension of the mesh larger, so that it extends at least a full "range unit" out beyond the domain of interest.

Finn

On Wed, 28 Jan 2026 at 06:28, Courtney Taylor <cou...@gmail.com> wrote:

I'm fitting INLA spatial models across multiple species with multiple spatially-structured random effects. The models converge (mode.status = 0) and most show better WAIC than non-spatial models, but I'm getting some spatial range estimates that exceed my study extent (3-4 times the study extent) and I'm trying to understand if/when this is a problem.

--

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, visit https://groups.google.com/d/msgid/r-inla-discussion-group/5f0be926-7598-4d2a-9610-80495af7b7bdn%40googlegroups.com.

--

Finn Lindgren
email: finn.l...@gmail.com

[Courtney Taylor]

Thank you for your responses. 

I should clarify that the model structure does not contain a constant intercept, but a spatially variable one and this is the element for which the range is typically larger than the study extent (lambda). The components are defined as so: 

components <- ~ -1 +

lambda(... model = matern) +

conspecific(... model = matern)   +

heterospecific(... model = matern)  +

heteroneighbourhood_slo(idx, ...) # random slopes for heterospecific effect by neighbour type

With the estimates I want to be able to say firstly, either there is or isn't is spatial variation in that component at the observed scale. Secondly, I am less interested in the precise value of the estimated range but more in the relative differences between the components so I can say that the spatial variation in these processes operates at different scales. I realise there are a lot of competing components, what gives me some confidence is that all the spatial components have the same prior:

matern <- inla.spde2.pcmatern(
  mesh = mesh,
    #alpha = #what is alpha
  prior.range = c(range0 = 150, Prange = 0.5),
  prior.sigma = c(sigma0 = 1, Psigma = 0.05)
)

And although there is a lot of uncertainty around the median range estimates across all components, the median range is consistent amongst the parameters (see figure, each line/point is a separate model, lambdas range is always larger than conspecific and heterospecific). So the data rather than the prior seems to be driving these patterns. 

 

My next concern is that the range exceeding the study extent (plus the uncertainty) might be an indication of other problems with the models, e.g. that the range for lambda is not well estimated because there are too many competing spatial effects and these might be confounded as suggested by the correlations amongst the fitted values of the spatial components (see figure). I believe this becomes a problem if I want to do some kind of variance partitioning or further inference on the estimates from specific spatial parameters. As you say 

Since you say you have multiple spatial structured random effects, you may need to be more careful with those posterior interpretations, and consider which combined quantities are well defined/identifiable.

 

Based on posterior predictive checks I can see that the models do a good job at population level inference, so prediction is stable when the spatial parameters are combined. I also used adjusted PIT to look at model calibration and excluded models that had a D statistic greater than .20 (this threshold is probably a bit high but is perhaps ok for ecological data).

cpo <- res$cpo$cpo

pit <- res$cpo$pit

pit.adjusted <- pit - (0.5 * cpo) # for discrete responses

ks <- ks.test(pit.adjusted, "punif")

Again I think this is only useful for understanding how the combined parameters perform. My question then becomes, how can I diagnose the set of models where the individual components may be unstable and not suitable for further inference? 

[Courtney Taylor]

Further clarification: I am not interested in using the individual components to say what is the specific value of the heterospecific effect at some particular location for species_i, rather I want to say something about the range of values predicted from these individual components.