Examples using both 'block' and 'group' in inla.spde.make.A?

252 views
Skip to first unread message

Sam Royston

unread,
Jul 8, 2022, 11:31:13 AM7/8/22
to R-inla discussion group
Hi

I am attempting to run some R-INLA code but please note I a not a statistician so apologies if this is a trivial enquiry. Does anyone know of example vignettes / code where the inla.spde.make.A function uses both "block" and "group" terms?

I'll outline my issue below, but in summary we require a spatio-temporal model with some of the observations being regional (area-averages over a polygon area) and so we would like to integrate the process value over the corresponding polygon. We were hoping to use "block" to define the blocks of spatial data and "group" for time. 

Thanks in advance,
Sam 

-------------------------------------------------------------------------------------------------------------------------------------------

We are running a spatio temporal model using the "group" for the temporal element. We have data that is areal-average over a polygon and were hoping to use "block" in the A matrix to sample and area-weight integrate our data. We have done each separately (using "block" on 2D time-invariant problem and using "group" with point data in 3D). 

When we only set "groups" in the A matrix and assume the lon-lat centres of polygons as point locations, the A matrix is (n loc * n times) in length. When I attempt to use both together, I find I have to manipulate / replicate the inputs to make the sizes work with inla.spde.make.A and then when passed to inla.stack it gives an error that length(group.index) != length(index), with n.A being the number time steps (group) smaller than when only "group" is used and "block" is omitted. 

Attached is a toy example that fails at the call to inla.stack()
............................................................................................................................................


block_and_group_example.R

Finn Lindgren

unread,
Jul 12, 2022, 4:20:11 PM7/12/22
to R-inla discussion group
Hi Sam,

it's not entirely clear from the code you supplied what the data
format is, so I'll describe what you need to give as input to
inla.spde.make.A to make it work, with your current call as starting
point.

This is your code:
A_eg = inla.spde.make.A(mesh, loc = df_xyz_by_yr, index = df_block_by_yr,
block = df_block_by_yr, block.rescale = "count",
group=df_obsdat$time,n.group=n_years)

The block feature works so that when index is supplied, you must have
index, group, and block to be vectors of the length of the number of
total desired integration/summation points.
The contents of block should go from 1 through the total number of
observations, and after constructing the initial evaluation matrix,
A_raw, for the integration points, then it aggregates rows of A_raw so
that the resulting matrix A has
A[i,] <- colSums(A_raw[block == i, ])
(apart from some weighting defined by the block.rescale argument.

The fact that you have both space and time isn't really relevant; it
does not care whether you have observations for all time points for
each block or not. You must supply the complete information, so that
index[j], group[j], block[j] tells it that integration point #j should
influence observation block[j], at spatial location determined by
loc[index[j],] and "time"(=group) point group[j]

Does the above make sense?
You didn't say what actual values you expected where the error message
said it expected something other than what it was given.

We're also working on helpers for aggregated data models for the
inlabru interface, where setting up these things will be a lot easier,
just as it is for point processes now.
That also provides
ipoints(samplers, mesh)
to construct proper integration weights that doesn't rely on gridding
and inside-counting, but rather constructs a properly weighted
integration scheme for the mesh.
When the samplers is SpatialPolygonsDataFrame with one row for each
observation, and an observation ID column, this can be used to help
construct the A-matrix.
It should go something like this:
ips$ID <- seq_len(nrow(samplers))
ips <- ipoints(samplers, mesh, group = "ID")
ips$time <- samplers$time[ips$ID] # not sure if this step is needed;
check the ips contents first
A <- inla.spde.make.A(mesh, ips, group = ips$time, block = ips$ID,
weights = ips$weight, block.rescale = "none", n.group=n_years)

Finn
> --
> 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/bf142093-7edb-4e73-a8a4-0712a422c85dn%40googlegroups.com.



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

Finn Lindgren

unread,
Jul 12, 2022, 4:22:25 PM7/12/22
to R-inla discussion group
For area averaging, change the block.rescale argument in my final example:

ips$ID <- seq_len(nrow(samplers))
ips <- ipoints(samplers, mesh, group = "ID")
ips$time <- samplers$time[ips$ID] # not sure if this step is needed;
A <- inla.spde.make.A(mesh, ips, group = ips$time, block = ips$ID,
weights = ips$weight, block.rescale = "sum", n.group=n_years)

Can you try this with your data? It should generate the correct
integration matrix for your problem.

Finn

Finn Lindgren

unread,
Jul 12, 2022, 4:27:37 PM7/12/22
to R-inla discussion group
(install and load the inlabru package first, for the ipoints function...)

Sam Royston

unread,
Jul 18, 2022, 5:32:05 AM7/18/22
to R-inla discussion group
Sincere thanks for the assistance. I have resolved the issue in my problem. The inlabru ipoints() function looks really useful so we will try to implement that. For completeness I have attached the 'toy' problem with code that should work using inlabru. 

Thanks again
Sam 

block_and_group_example_inlabru_sampling.R
Reply all
Reply to author
Forward
0 new messages