Moran I test under randomisation
data: eng324$realNetPre
weights: lw
Moran I statistic standard deviate = 6.9357, p-value = 4.042e-12
alternative hypothesis: two.sided
sample estimates:
Moran I statistic Expectation Variance
0.244871308 -0.003095975 0.001278221
I = \frac{n \sum_{(2)} w_{ij} z_i z_j}{S_0 \sum_{i=1}^{n} z_i^2}
where x_i, i=1, \ldots, n are n observations on the numeric variable of interest, z_i = x_i - \bar{x}, \bar{x} = \sum_{i=1}^{n} x_i / n, \sum_{(2)} = \stackrel{\sum_{i=1}^{n} \sum_{j=1}^{n}}{i \neq j}, w_{ij} are the spatial weights, and S_0 = \sum_{(2)} w_{ij}. A comparison of implementations of measures of spatial autocorrelation shows that a wide range of measures is available in R in a number of packages, chiefly in the spdep package, and that differences from other implementations can be attributed to design decisions (Bivand and Wong 2018). The spdep package also includes the only implementations of exact and Saddlepoint approximations to global Moran’s I (Tiefelsdorf 2002; Bivand, Müller, and Reder 2009).
Global measures of spatial autocorrelation using spatial weights objects based on graphs of neighbours are, as we have seen, rather blunt tools, which for interpretation depend critically on a reasoned mean model of the variable in question. If the mean model is just the intercept, the global measures will respond to all kinds of mis-specification, not only spatial autocorrelation (Schabenberger and Gotway 2005; McMillen 2003).
7.2 Tests of residual spatial autocorrelation
When the model is the null model, containing only the intercept, the univariate test moran.test on the residuals with randomisation=FALSE and the test on residuals lm.morantest are the same:
Moran I test under normality
data: residuals(lm_null)
weights: lw
Moran I statistic standard deviate = 6.7433, p-value = 1.548e-11
alternative hypothesis: two.sided
sample estimates:
Moran I statistic Expectation Variance
0.244871308 -0.003095975 0.001352204
Moran I test under normality
data: residuals(lm_obj_pre)
weights: lw
Moran I statistic standard deviate = 4.521, p-value = 6.156e-06
alternative hypothesis: two.sided
sample estimates:
Moran I statistic Expectation Variance
0.163150227 -0.003095975 0.001352204
Global Moran I for regression residuals
data:
model: lm(formula = form_pre, data = eng324)
weights: lw
Moran I statistic standard deviate = 4.7671, p-value = 1.869e-06
alternative hypothesis: two.sided
sample estimates:
Observed Moran I Expectation Variance
0.163150227 -0.010433124 0.001325885
This is because lm.morantest takes into account the structure of the regression model, as do lm.morantest.sad and lm.morantest.exact:
Global Moran I statistic with exact p-value
data:
model:lm(formula = form_pre, data = eng324)
weights: lw
Exact standard deviate = 4.4649, p-value = 8.01e-06
alternative hypothesis: two.sided
sample estimates:
[1] 0.1631502
In all these cases, the measured value of Moran’s I is the same, but the inferences differ somewhat.
7.3 Model specification
A recent review of spatial regression in a spatial econometrics setting is given by Kelejian and Piras (2017); note that their usage is to call the spatial coefficient of the lagged response \lambda and that of the lagged residuals \rho, the reverse of other usage (Anselin 1988; J. P. LeSage and Pace 2009); here we use \rho_{\mathrm{Lag}} for the spatial coefficient in the spatial lag model, and \rho_{\mathrm{Err}} for the spatial error model. One interesting finding is that relatively dense spatial weights matrices may downweight model estimates, suggesting that sparser weights are preferable (Smith 2009). Another useful finding is that the presence of residual spatial autocorrelation need not bias the estimates of variance of regression coefficients, provided that the covariates themselves do not exhibit spatial autocorrelation (Smith and Lee 2012).
In trying to model these spatial processes, we may choose to model the spatial autocorrelation in the residual with a spatial error model (SEM).
where {\mathbf y} is an (N \times 1) vector of observations on a response variable taken at each of N locations, {\mathbf X} is an (N \times k) matrix of covariates, {\mathbf \beta} is a (k \times 1) vector of parameters, {\mathbf u} is an (N \times 1) spatially autocorrelated disturbance vector, {\mathbf \varepsilon} is an (N \times 1) vector of independent and identically distributed disturbances and \rho_{\mathrm{Err}} is a scalar spatial parameter.
If the processes in the covariates and the response match, we should find little difference between the coefficients of a least squares and a SEM, but very often they diverge, suggesting that a Hausman test for this condition should be employed (Pace and LeSage 2008). This may be related to earlier discussions of a spatial equivalent to the unit root and cointegration where spatial processes match (Fingleton 1999).
A model with a spatial process in the response only is termed a spatial lag model (SLM, often SAR - spatial autoregressive) (J. P. LeSage and Pace 2009).
{\mathbf y} = \rho_{\mathrm{Lag}} {\mathbf W}{\mathbf y} + {\mathbf X}{\mathbf \beta} + {\mathbf \varepsilon},
where \rho_{\mathrm{Lag}} is a scalar spatial parameter.
Work reviewed by Mur and Angulo (2006) on the Durbin model; the Durbin model adds the spatially lagged covariates to the covariates included in the spatial lag model, giving a spatial Durbin model (SDM) with different processes in the response and covariates:
where {\mathbf \gamma} is a (k' \times 1) vector of parameters. k' defines the subset of the intercept and covariates, often k' = k-1 when using row standardised spatial weights and omitting the spatially lagged intercept.
This permits the spatial processes to be viewed and tested for as a Common Factor (Burridge 1981; Bivand 1984). The inclusion of spatially lagged covariates lets us check whether the same spatial process is manifest in the response and the covariates (SEM), whether they are different processes, or whether no process is detected. The Common Factor is present when {\mathbf \gamma} = - \rho_{\mathrm{Lag}} {\mathbf \beta}:
{\mathbf y} = \rho_{\mathrm{Lag}} {\mathbf W}{\mathbf y} + {\mathbf X}{\mathbf \beta} - \rho_{\mathrm{Lag}} {\mathbf W}{\mathbf X} {\mathbf \beta} + {\mathbf \varepsilon},
\qquad ({\mathbf I} - \rho_{\mathrm{Lag}} {\mathbf W}){\mathbf y} = ({\mathbf I} - \rho_{\mathrm{Lag}} {\mathbf W}) {\mathbf X}{\mathbf \beta} + {\mathbf \varepsilon},
where {\mathbf I} is the N \times N identity matrix, and {\mathbf I} - \rho_{\mathrm{Lag}} {\mathbf W} is the Common Factor:
If we extend this family with processes in the covariates and the residual, we get a spatial error Durbin model (SDEM). If it is chosen to admit a spatial process in the residuals in addition to a spatial process in the response, again two models are formed, a general nested model (GNM) nesting all the others, and a model without spatially lagged covariates (SAC, also known as SARAR - Spatial AutoRegressive-AutoRegressive model). If neither the residuals nor the residual are modelled with spatial processes, spatially lagged covariates may be added to a linear model, as a spatially lagged X model (SLX) (Elhorst 2010; Bivand 2012; J. LeSage 2014; Halleck Vega and Elhorst 2015).
This may be constrained to the double spatial coefficient model SAC/SARAR by setting {\mathbf \gamma} = 0, to the spatial Durbin (SDM) by setting \rho_{\mathrm{Err}} = 0, and to the error Durbin model (SDEM) by setting \rho_{\mathrm{Lag}} = 0. Imposing more conditions gives the spatial lag model (SLM) with {\mathbf \gamma} = 0 and \rho_{\mathrm{Err}} = 0, the spatial error model (SEM) with {\mathbf \gamma} = 0 and \rho_{\mathrm{Lag}} = 0, and the spatially lagged X model (SLX) with \rho_{\mathrm{Lag}} = 0 and \rho_{\mathrm{Err}} = 0.
The specification problem is then to choose between (from most complex to simplest): GNM (three sets of spatial parameters: \rho_{\mathrm{Lag}}, \rho_{\mathrm{Err}}, {\mathbf \gamma}), SAC, SDM or SDEM (two sets of spatial parameters among the three), SEM, SLM, SLX (one set of spatial parameters among the three) or ordinary least squares (no spatial parameters) - assuming that the dependent variable is Gaussian. Effectively the same specification choice is present with limited dependent variables.
7.4 Rao’s score (Lagrange multiplier) tests
When working over 25 years ago on Bivand and Szymanski (2000), the only Rao’s score (Lagrange multiplier) tests available were those newly proposed by Anselin et al. (1996). The underlying concern is to test whether the spatial parameter of interest could be taken as zero, without estimating the parameter. One works from the information matrix of the spatial model; in the case of the first tests, a mixed regressive-spatial autoregressive model with a spatial autoregressive disturbance (SAC) was considered, as this had also been used for residual autocorrelation in an SLM model. Five tests are proposed, a test for missing \rho_{\mathrm{Lag}} {\mathbf W}{\mathbf y}, "RSlag" and pointing to SLM, another for missing \rho_{\mathrm{Err}} {\mathbf W} {\mathbf u}, "RSerr" pointing to SEM, two tests matching the first two tests but adjusted for the presence of the other ("adjRSlag", "adjRSerr"), and a "SARMA" test combining both alternatives, equal to "RSlag" plus "adjRSerr" and "RSerr" plus "adjRSlag".
Here, for the independent variables initially included in the model to account for the logarithm of real net expenditure on garbage collection in each district, the logarithm of the number of pick-up points (units), the proportion of dwellings among the pick-up points (house), the logarithm of pick-up point density (dens), a categorical variable Metrop with levels: "district", "london", "metrop", and the logarithm of the real regional wage rate for regions. In this case, the tests point to a SEM model rather than an SLM model or falling back to OLS.
Koley and Bera (2022) propose a similar set of tests to distinguish between either an SLM or an SLX model based on an SDM model, falling back to ordinary least squares if neither alternative is chosen. Because the interpretation of spatially lagged categorical variables is unclear, we drop Metrop from the tests
In this case, an SLX model appears to be indicated rather than SDM, SLM or OLS. However, we also need to look at SDEM as the tests against SAC indicated SEM rather than SLM. Koley (2024) proposes approriate tests:
which indicate SEM rather than SDEM. Clearly, the underlying difficulty is that the GNM is poorly identified, and as Koley (2024) finds:
… identification conditions have not been established for the GNS model in the spatial literature that also suffers from an identification problem. In fact, to the best of our knowledge, the non-identification of the GNS model has not yet been analytically studied in the literature (Koley 2024, 3).
However, it has been found that the lm.RStests tests can be used on a fitted SLX model using the lmSLX function in spatialreg(Koley 2024, 18), again omitting Metrop from the {\mathbf W}{\mathbf X} term in addition to omitting the intercept:
In Bivand and Szymanski (2000), it was found that the inclusion of the categorical variable Majority, with levels: "coalition", "conservative" and "labour", as anecdotal reports indicated that neighbouring districts with the same political party in control shared more information than otherwise, so the yardstick used would be weighted towards the same political party:
The tests here indicate some support for the inclusion of \rho_{\mathrm{Lag}} {\mathbf W}{\mathbf y}. On turning to the new tests for SDM compared to SLM or SLX:
Anselin, Luc, Anil K. Bera, R. Florax, and M. J. Yoon. 1996. “Simple Diagnostic Tests for Spatial Dependence.”Regional Science and Urban Economics 26: 77–104. https://doi.org/10.1016/0166-0462(95)02111-6.
Bivand, Roger. 1984. “Regression Modeling with Spatial Dependence: An Application of Some Class Selection and Estimation Methods.”Geographical Analysis 16: 25–37.
———. 2012. “After ‘Raising the Bar’: Applied Maximum Likelihood Estimation of Families of Models in Spatial Econometrics.”Estadística Española 54: 71–88.
Bivand, Roger, W. Müller, and M. Reder. 2009. “Power Calculations for Global and Local Moran’s I.”Computational Statistics and Data Analysis 53: 2859–72.
Bivand, Roger, and Stefan Szymanski. 2000. “Modelling the Spatial Impact of the Introduction of Compulsory Competitive Tendering.”Regional Science and Urban Economics 30: 203–19.
Bivand, Roger, and David W. S. Wong. 2018. “Comparing Implementations of Global and Local Indicators of Spatial Association.”TEST 27 (3): 716–48. https://doi.org/10.1007/s11749-018-0599-x.
Burridge, P. 1981. “Testing for a Common Factor in a Spatial Autoregression Model.”Environment and Planning A 13: 795–800.
Cliff, Andrew D., and J. Keith Ord. 1981. Spatial Processes. London: Pion.
Elhorst, J. Paul. 2010. “Applied Spatial Econometrics: Raising the Bar.”Spatial Economic Analysis 5: 9–28.
Fingleton, B. 1999. “Spurious spatial regression: Some Monte Carlo results with a spatial unit root and spatial cointegration.”Journal of Regional Science39: 1–19.
Halleck Vega, Solmaria, and J. Paul Elhorst. 2015. “The SLX Model.”Journal of Regional Science 55 (3): 339–63. https://doi.org/10.1111/jors.12188.
Koley, Malabika, and Anil K. Bera. 2022. “Testing for Spatial Dependence in a Spatial Autoregressive (SAR) Model in the Presence of Endogenous Regressors.”Journal of Spatial Econometrics 3 (11): 1–46. https://doi.org/10.1007/s43071-022-00026-7.
LeSage, James P., and R. Kelley Pace. 2009. Introduction to Spatial Econometrics. Boca Raton FL: Chapman; Hall/CRC.
McMillen, Daniel P. 2003. “Spatial Autocorrelation or Model Misspecification?”International Regional Science Review 26: 208–17.
Mur, Jesús, and Ana Angulo. 2006. “The Spatial Durbin Model and the Common Factor Tests.”Spatial Economic Analysis 1 (2): 207–26. https://doi.org/10.1080/17421770601009841.
Smith, Tony E., and K. L. Lee. 2012. “The effects of spatial autoregressive dependencies on inference in ordinary least squares: a geometric approach.”Journal of Geographical Systems 14 (January): 91–124. https://doi.org/10.1007/s10109-011-0152-x.
Tiefelsdorf, Michael. 2002. “The Saddlepoint Approximation of Moran’s I and Local Moran’s {I}_i Reference Distributions and Their Numerical Evaluation.”Geographical Analysis 34: 187–206.