Beyond Conditional Independence: Modeling Observation-Level Spatial Dependence in Hierarchical Models

For normally distributed data, modeling dependence is straightforward—we specify a covariance matrix. For non-Gaussian data, hierarchical models can introduce dependence on a latent scale through spatially correlated Gaussian fields, but observations remain conditionally independent given the latent parameters. This ignores dependence in the data beyond what the latent structure captures. In extreme rainfall modeling, for instance, nearby locations may produce observations with similar quantiles even after conditioning on their location, scale and shape parameters—when rainfall at one site exceeds its median, neighboring sites likely do too. I introduce a computationally tractable approximate Bayesian method that models this observation-level dependence using copulas. Simulations show that ignoring such dependence biases parameter estimates, while modeling it improves point estimation. However, an interesting trade-off emerges: the copula model can produce overly narrow parameter intervals when observation correlation is high, even as observation-level intervals maintain proper coverage