Modelling minimally-sampled spatial data

Modeling
1

I have data collected spatially at the 4 corners of a square, with lots of observations of a continuous measurement at each corner. My underlying generative model is that there is a single spatially-located latent signal that drives the measurements, so I’ve been modelling a latent signal with an x & y location plus a scale parameter encoding the diminishing influence of the signal at increasing distances from its location (just using a isotropic exponentiated quadradic). But with a naive parameterization of x & y being unbounded (with Gaussian priors to assign most of the mass inside the square; the latent signal occuring outside the square is possible but unlikely in my case) and scale being positive-bound, I feel like I need a prior on the scale that’s a bit more nuanced than usual as I would like to steer the sampling away from solutions where the latent signal scale is too small to manifest substantial influence at any of the measurement locations, and what counts as “too small” will obviously be influenced by the x & y locations. That is, a very small signal scale would be fine if it is very close to one of the measurement locations, but not fine if directly in the center of the square or way outside the square. Thoughts on solutions? Here’s what I have so far:

data{
	int<lower=2> n_obs_per_loc ;
	int<lower=4,upper=4> n_loc ;
	matrix[n_obs_per_loc,n_loc] obs ;
transformed data{
	row_vector[n_loc] loc_x = [-1,1,-1,1] ;
	row_vector[n_loc] loc_y = [-1,-1,1,1] ;
}
parameters{
	real x ;
	real y ;
	real<lower=0> scale ;
	real<lower=0> amplitude ;
	vector<lower=0> loc_noise_scale ;

}
transformed parameters{
	vector[n_loc] loc_mean = (
		amplitude
		.* exp(
			-(
				pow( (loc_x-x), 2 )
				+ pow( (loc_y-y), 2 )
			)
			./ ( 2 .* pow(scale,2) )
		)
	) ;
}
model{
	x ~ std_normal() ;
	y ~ std_normal() ;
	scale ~ weibull(2,1) ;
	amplitude ~ weibull(2,1) ;
	loc_noise_scale ~ weibull(2,1) ;
	for(i_loc in 1:n_loc){
		obs[,i_loc] ~ normal(loc_mean[i_loc],loc_noise_scale[i_loc]) ;
	}
}
2

Why? Is this for identifiability reasons, or is there some underlying domain knowledge at play here? Given that there’s some signal at one or more sampling location, can you just let the likelihood enforce this constraint?

Tangentially, note that with enough data this likelihood will likely contain some multimodalities, because a source near-ish to one of the edges of the square will likely induce a secondary mode reflected across that edge.

3

Yes, purely to avoid the identifiability issue of location becoming uninformed when the scale parameter is small. The likelihood should constrain under conditions of low measurement noise, and I guess I just need to evaluate whether I’m in that domain for my particular case (maybe run some prior predictive stuff to discern at what levels of noise the pathology I’m worried about starts to become an issue).

Ah, good point. This would occur when the source has a scale small enough to influence the nodes on the edge but have negligible influence the parallel edge’s nodes.

4

I don’t think you’ll have a serious computational issue here as long as you have some weak priors that prevent the location from flying off to infinity when the likelihood is uninformative.

In the case where noise is such that one cannot distinguish presence of signal from absence of signal, I think you’d want to be very uncertain about the location and scale of the signal, rather than reasonably certain that the signal is substantially present on at least one of the four corners.

1 Like
5

Oh, I belatedly discovered another identifiability issue: as the latent signal approaches the center, the scale and amplitude parameters become non-identified with respect to each other. I think the solution is to have a mechanism for smoothly transitioning from the latent-signal-with-scale model and common-amplitude model as the former approaches the center:

data{
	int<lower=2> n_obs_per_loc ;
	int<lower=4,upper=4> n_loc ;
	matrix[n_obs_per_loc,n_loc] obs ;
transformed data{
	row_vector[n_loc] loc_x = [-1,1,-1,1] ;
	row_vector[n_loc] loc_y = [-1,-1,1,1] ;
}
parameters{
	real x ;
	real y ;
	real<lower=0> scale ;
	real<lower=0> amplitude ;
	vector<lower=0> loc_noise_scale ;

}
transformed parameters{
	vector[n_loc] dist_xy_loc = (
		exp(
			-(
				pow( (loc_x-x), 2 )
				+ pow( (loc_y-y), 2 )
			)
			./ ( 2 .* pow(scale,2) )
		)
	) ;

	vector[n_loc] dist_xy_center = (
		exp(
			-(
				pow( x, 2 )
				+ pow( y, 2 )
			)
			./ ( 2 .* pow(scale,2) )
		)
	) ;

	vector[n_loc] loc_mean = (
		amplitude
		.* (
			(dist_xy_loc .* (1-dist_xy_center) // 0 when at center
			+ dist_xy_center // 1 when at center
		)
	) ;

}
model{
	x ~ std_normal() ;
	y ~ std_normal() ;
	scale ~ weibull(2,1) ;
	amplitude ~ weibull(2,1) ;
	loc_noise_scale ~ weibull(2,1) ;
	for(i_loc in 1:n_loc){
		obs[,i_loc] ~ normal(loc_mean[i_loc],loc_noise_scale[i_loc]) ;
	}
}

Edit: oops, removed unnecessary exponentiations.

1 Like