The use of mediation analysis is common in epidemiology to decompose an exposure-outcome total effect into direct and indirect effects via a given set of mediators. The counterfactual framework in causal inference allowed for great improvements in mediation analysis to define identification assumptions and clarify sources of bias in different settings. Multiple mediation methods have been proposed to consider confounders of various types and for different settings such as longitudinal designs. However, in some settings we may reasonably expect spatially heterogeneous estimates of both natural direct and indirect effects. In various research areas, especially environmental and social epidemiology, we would expect there to be spatial dependence that could be leveraged. In this paper, we describe and compare a set of approaches to perform spatial mediation analyses and conduct an extensive set of simulations to quantify their performance. We also present a case study assessing the spatially heterogeneous mediating role of PM2.5 in the association between income and cardiovascular diseases. We provide all code via GitHub for replication and reproducibility purposes. We find that using a framework combining geographically weighted regression (GWR) and other spatially varying coefficient (SVC) models and the product method of mediation analysis is accurate and precise at estimating spatially heterogeneous effect estimates and handling spatial confounding in a grid setting. Such spatial approaches can be extended to any mediation approach and various confounding structures.