inlabru-org / inlabru

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

Accessing log-intensities at explicit points #69

Closed shinneuro closed 2 years ago

shinneuro commented 4 years ago

Hello inlabru,

I am working with point data on sphere with radius 1, which is not on the globe. I could fit the data successfully, and I want to access posterior estimates at some points, but I have no idea how to do it. For example, I want to know posterior mean at point (0, 0, 1).

Below is my code, which is just a direct application of example.

cmp.L <- coordinates ~ mySmooth(map = coordinates, model = spde.L) + Intercept fit.L <- lgcp(cmp.L, coords.L) lambda.L <- predict(fit.L, pixels(mesh.L), ~ exp(mySmooth + Intercept))

Thank you in advance for your answers. And I really appreciate your commitment of building this useful tool.

Yours Sincerely, Minho

finnlindgren commented 4 years ago

Hi Minho, The pixels(mesh.L) call you use for the predict call sets up a spatial object with locations in a grid covering the mesh. To get the value at some specific point, one can either extract the value by accessing the value corresponding to the closest point in that grid, or instead construct an alternative input data object that just contains that specific location.

However, the current inlabru version doesn’t fully support working directly on the sphere; the internal code has several places that transform the domain to a flat projection, which breaks models/code that should operate on the entire sphere, so I wouldn’t expect a spherical mesh to work at the moment (but INLA does support it, and modifying inlabru to support it is high on my TODO list). How are you currently building you mesh? Do you get any warnings/errors when running lgcp()"?

shinneuro commented 4 years ago

Dear finn,

Thank you for the quick reply. I may briefly explain what I am doing now.

I just modified a bit from https://www.tandfonline.com/doi/full/10.1080/01621459.2019.1611582 and I just realized that you are one of the authors! So I assume you are quite familiar with the topic.

The sphere that I am working is actually brain cortical surface mapped to a sphere, having 162k vertices. So I just used all the coordinates of vertices (162k x 3) to generate the mesh.

Also, my data is a set of coordinates, each matching one of the vertices above. So the point data that I have exactly match to part of the vertices that I used to generate the mesh.

So far, I didn't encounter any errors. Do you think everything worked fine then?

vert.L = fs.lh$vertices
vert.R = fs.rh$vertices
# normalize to a sphere with radius 1
loc.norm.L = vert.L/sqrt(rowSums(vert.L^2)) #left
loc.norm.R = vert.R/sqrt(rowSums(vert.R^2)) #right
mesh.L <- inla.mesh.2d(loc.norm.L, cutoff = .022, max.edge = .2) #left
mesh.R <- inla.mesh.2d(loc.norm.R, cutoff = .022, max.edge = .2) #right
shinneuro commented 4 years ago

Hi, for your information, everything just went fine. I also could generate posteriors using predict providing coordinates that I want to access.

finnlindgren commented 2 years ago

Spherical mesh support was added at some point since this report, so I'm closing this Issue.