Hi @SDBsjr,
Welcome to the Stan forum. I can’t say how it might impact Bayes factors, but here’s some general information on the topic.
Say x is a set of samples taken across some defined space. If the observations are correlated in space, the variance of the observations will tend to be inflated, and you’ll effectively have less information about the population than you would if you had the same number of uncorrelated observations (as if you had …somehow partially…re-sampled the same point/person). If you also have some spatially correlated samples of another variable y, those properties will similarly impact the covariance of the two variables: in repeated sampling, you see more extreme values.
Say for each variable you’re able to estimate a spatial trend, the non-stationary mean around which values fluctuate. If x and y share a common spatial pattern—say, they both exhibit an east-west gradient—then the covariance will tend to be inflated; whereas if they have orthogonal patterns, the covariance will tend to be deflated. A few people out there will tell you that you should not ‘adjust’ for spatial trends, on the grounds that the common pattern is the covariation we’re interested in; others see it more like time-series data and say that the common trend is probably caused by some other factors which logically increases the chances of ‘nonsense correlations’ appearing in the data. (I certainly favor the latter view for what its worth.)
This is why you shouldn’t rely on non-spatial hypothesis tests with spatial data, and that is equally true for Bayesian and non-Bayesian models. Similarly, regression coefficient estimates will tend to have larger errors when you don’t estimate the spatial trend.
Usually people only focus on SA in the outcome variable. However, when you have independent variables with strong SA, this will independently impact your model results. Because the SA increases the variance of x but includes essentially duplicate information, the standard deviation of the posterior will tend to be overconfident. You can see this if you just simulate a bunch of data with SA; you’ll see that as the degree of SA in the covariates increases, the posterior standard deviation for its coefficient \beta (on average) decreases and the coverage rate (proportion of say 90% quantile intervals for the posterior distribution of \beta containing the correct value) also goes way down. Again, Bayesian and frequentist models respond identically.
You can see these results in this paper on Bayesian spatial filtering https://osf.io/fah3z/ (also in Spatial Statistics)particularly figures 2, 3, and 4 show what happens to regression results when you vary the degree of SA in the independent and dependent variables. Those spatial filtering models, BYM, and intrinsic autoregressive models are in the geostan R package; though its still being developed.