Elsevier

Spatial Statistics

Volume 21, Part B, August 2017, Pages 440-459
Spatial Statistics

A comparison of estimation methods for multilevel models of spatially structured data

https://doi.org/10.1016/j.spasta.2017.01.002Get rights and content

Abstract

Two recent contributions (Dong et al., 2015; Osland et al., 2016) point to the relevance of multilevel models for spatially structured data. In Osland et al. (2016) these models are used to examine the importance of district-level covariates for house prices in Stavanger, Norway, in Dong et al. (2015) similarly for land parcel prices in Beijing; we use these data sets in our comparison. In Osland et al. (2016), a district-level spatial random effect was fitted using an intrinsic CAR model estimated using WinBUGS. Dong et al. (2015) used R code provided in supplementary materials to their article, and subsequently improved in an R package (Dong et al., 2016a); computation there used custom MCMC C++ code to fit a SAR district-level spatial random effect. This article compares approaches to estimating models of this kind, using the R packages R2WinBUGS, HSAR, INLA, R2BayesX, hglm and the new package mclcar for Monte Carlo maximum likelihood estimation (Sha, 2016b). We show that multilevel models of spatially structured data may be estimated readily using a variety of approaches, not only the intrinsic CAR model more typically found in the existing literature. We also point to a range of issues for further research in situations in which data acquired at different levels of spatial resolution are combined.

Introduction

Spatial data are often organized in multilevel structures, and can be characterized by a mixture of two important features: spatial autocorrelation and spatial heterogeneity (Anselin, 1988). For example, when analysing housing markets, individual houses are located in neighbourhoods which are located within census tracts, which are situated in even higher or more aggregated levels such as municipalities, counties and regions. Typical feature of multilevel data are correlations between observations located within one of the spatial zones or defined neighbourhoods. Houses and, hence, housing prices within a neighbourhood may share similar features, inter alia because they are built at the same time, or because there may be a relatively large proportion of wealthy people living there, or because they are located close to environmental amenities or disamenities. The result could be within-zonal correlations at various levels. So when using geographical data, it could be unrealistic to assume that observations within delimited spatial zones are independent. Following, Corrado and Fingleton (2012), analysts should consider using multilevel models in these situations.

Similarly to spatial regression models, but in contrast with the ordinary least squares estimator, multilevel estimators are not based on the assumption that observations are independent. The multilevel models account for unexplained local correlations, by introducing for instance a random intercept, one for each local zone. Frequently, these random intercepts are not the main focus of the research, given that they capture the impact of excluded factors found in the zones. However, by studying the significance, sign and magnitude of the random effects, they may provide important or useful information in the modelling process. This is explained in for instance Rabe-Hesketh and Skrondal (2008) or Osland et al. (2016). We choose not to examine an alternative random slope approach to these models here (Rabe-Hesketh and Skrondal, 2008), but note that a very recent paper by Dong et al. (2016b) takes up the equivalent spatial random slope approach.

In addition to inducing correlations within zones, variants of multilevel models as described in e.g. Rabe-Hesketh and Skrondal (2008), capture unobserved spatial heterogeneity which is a common feature in spatial data. By way of example, coefficients related to housing market structures may vary between the most urbanized area and semi-urban areas. In this way, classical multilevel models are useful in that they via estimation of the variance of random effects also account for the interzonal variations in the data.

One important assumption in the conventional multilevel models is that there is independence between the zonal random effects (Rabe-Hesketh and Skrondal, 2008, p. 61). For many types of data and analyses, this may be a reasonable assumption. For instance, when studying what determines birthweight of children, one may find that birthweights of children of the same mother are correlated. Birthweight of children of different mothers which represents the group level can probably be assumed to be uncorrelated. Another frequently used example is related to academic results in schools: it could be reasonable to assume that there are correlations or more similar results for pupils within the same school. It may also be a reasonable assumption that correlations of educational results between different schools could be ignored. So in these situations the classical multilevel models could be a useful approach, and may provide improved efficiency of estimates at the lower observational level (Dong et al., 2015).

The conventional random effect estimator may not be efficient when analysing spatial data, however. For this type of data it could be highly relevant to account for correlations between random effects located in close geographic proximity (Osland et al., 2016). When studying regional housing prices, for instance, urban residential neighbourhoods may have many similar features. As we move further away from these neighbourhoods, the zones may more or less gradually change character, from being highly urban to being more semi-urban and even rural. Zones located close to each other, will thus be more similar compared to zones located further away. In this way, there may exist correlations between the higher level zones located nearby, similar to what we may find at the individual observational data level in the traditional spatial regression models.

The new class of models termed Hierarchical Spatial Autoregressive (HSAR) is introduced and developed in contributions by Dong and Harris (2015) and Dong et al. (2015). Their proposed model adds an observation level spatial autoregressive component to the model, the spatially lagged response Wy, where y is a vector of n observations on the response variable, and W is a spatial weights matrix. The HSAR model may also be estimated without this component. In the following, Wy will be omitted to maintain comparability with Osland et al. (2016), although HSAR fits a simultaneous autoregressive (SAR) spatial random effect rather than a conditional autoregressive (CAR) spatial random effect. Both Dong et al. (2015) and Osland et al. (2016) analyse housing data in a very similar format, with some covariates observed for each individual-level observation, and some others observed only at the aggregated, district level. Both employ district level spatial random effects.

It is worth noting that similar approaches to spatial multilevel models have been used by Savitz and Raudenbush (2009) in studying what they term as collective efficacy–the fusion of social cohesion and informal social control–in Chicago neighbourhoods, and by Pierewan and Tampubolon (2014) in examining well-being across regions in Europe. Spatial multilevel approaches in regional and housing economics are inviting, but do not seem to have been applied commonly before Dong et al. (2015) and Osland et al. (2016), apart from Gelfand et al. (2007) and Liu and Roberts (2012); other exceptions may exist in a very dispersed literature. In other research domains, work is documented in, for example, the Munich (Fahrmeir and Lang, 2001, Rue and Held, 2005) and Zambia (Kandala et al., 2001) data sets included in several R packages (Umlauf et al., 2015b, Rue et al., 2015). Both of these examples are for non-Gaussian responses, and many of the estimation methods discussed below extend to discrete responses (see Table 1). In an associated but separate discussion, Bivand (2016) uses weighted aggregation to the district level in the case of a data set in which the estimate of the coefficient of a district level covariate is of central concern.

In this contribution, we examine the software implementations available for estimating district-level multilevel spatially structured random effects models of the kind described by Dong et al. (2015) and Osland et al. (2016); we use the data sets and variables from these articles below. The intrinsic CAR model, used by Osland et al. (2016) and common in disease mapping, may be considered as a standard CAR model with a constrained process coefficient set to its maximum. We will not only look at implementations of the intrinsic CAR model, but also at alternatives such as HSAR and multilevel CAR models in which the autoregressive coefficient is estimated. Comparisons will be carried out using scatterplot matrices of the estimated spatial random effects, tables of the estimated fixed effects and “caterpillar” plots (Spiegelhalter et al., 2003, p. 28). These plots of random effects sorted by mean value and shown with 95% error bars are useful for presenting condensed displays of estimates of the posterior marginal variances; Blangiardo and Cameletti (2015) recommend that these model outputs be shown in detail.

Our approach is first to apply the chosen estimation methods and implementations to the subset of data from the Beijing land parcel data set used in Dong et al. (2015) and published in the R HSAR package (Dong et al., 2016a), as this data set is available for readers to replicate using a script provided as supplementary material (see Appendix A). Next, we apply the same methods to the Stavanger data set used in Osland et al. (2016), but where the house price data have confidential status. Most of the methods used will be described briefly here, as they are well-documented in the literature, but special attention will be given to the Monte Carlo likelihood approach to estimating multilevel CAR models (Sha, 2016a, Sha, 2016b).

Section snippets

Estimation methods

Our approach in this comparison is to advance the methodological proposals made by Osland et al. (2016) and Dong et al. (2015). If we have a multilevel data structure, and the upper level entities (districts, zones) contain at least one lower level observation, we can attempt to apply a multilevel model in some specification. We will first describe the general class of models that we will be considering, before going on to examine implementation issues which will mostly be covered by references

Model estimation results

We move now to consider the use of software implementations listed above and in Table 1 to obtain estimators for multilevel spatial models. The implementations we have used include iid random effects at the district level, and SAR or CAR spatially structured random effects at the district level. Finally two implementations of the Besag–York–Mollie model with both an iid random effect and an intrinsic CAR spatially structured random effect will also be used. Results from a selection of these

Concluding remarks

Recent research has shown that multilevel models have the potential to be useful when analysing spatial data. In this paper we have shown that estimation of multilevel models, both CAR and SAR, may be carried out by using several implementations included in R packages. We have also tried out a new package mclcar for Monte Carlo maximum likelihood estimation. This method was long neglected, and has been retrospectively evaluated only recently, including its implementation in prototype form in R.

Acknowledgements

The authors would like to thank Guanpeng Dong and Lars Rönnegård for helpful clarifications, and the editor and two referees for constructive comments. Any remaining confusion is the responsibility of the authors alone. This paper is based on funding from the Research Council of Norway through an Industrial Ph.D. for Ingrid Sandvig Thorsen.

References (55)

  • J.E. Besag et al.

    On conditional and intrinsic autoregression

    Biometrika

    (1995)
  • J. Besag et al.

    Bayesian image restoration, with two applications in spatial statistics

    Ann. Inst. Statist. Math.

    (1991)
  • R. Bivand

    Revisiting the boston data set (harrison and rubinfeld, 1978c): changing the units of observation affects estimated willingness to pay for clean air

    Region

    (2016)
  • R. Bivand et al.

    Spatial data analysis with R-INLA with some extensions

    Journal of Statistical Software

    (2015)
  • M. Blangiardo et al.

    Spatial and Spatio-temporal Bayesian Models with R-INLA

    (2015)
  • O. Cappé et al.

    On the convergence of the Monte Carlo maximum likelihood method for latent variable models

    Scand. J. Stat.

    (2002)
  • Y. Chung et al.

    A nondegenerate penalized likelihood estimator for variance parameters in multilevel models

    Psychometrika

    (2013)
  • L. Corrado et al.

    Where is the economics in spatial econometrics?

    J. Reg. Sci.

    (2012)
  • G. Dong et al.

    Spatial autorgressive models for geographically hierarchical data structures

    Geogr. Anal.

    (2015)
  • G. Dong et al.

    Multilevel modeling with spatial interaction effects with application to an emerging land market in Beijing, China

    PLOS One

    (2015)
  • Dong, G., Harris, R., Mimis, A., 2016a. HSAR: Hierarchical Spatial Autoregressive Model (HSAR). URL:...
  • G. Dong et al.

    Spatial random slope multilevel modeling using multivariate conditional autoregressive models: A case study of subjective travel satisfaction in beijing

    Ann. Am. Assoc. Geogr.

    (2016)
  • Dorie, V., 2015. blme: Bayesian linear mixed-effects models. URL: https://CRAN.R-project.org/package=blme  R package...
  • M. Evans et al.

    Approximating Integrals via Monte Carlo and Deterministic Methods

    (2000)
  • L. Fahrmeir et al.

    Bayesian inference for generalized additive mixed models based on Markov random field priors

    J. R. Stat. Soc. Ser. C

    (2001)
  • A.E. Gelfand

    Misaligned spatial data: The change of support problem

  • A. Gelman

    Inference and monitoring convergence

  • Cited by (27)

    • The (in)stability of Bayesian model selection criteria in disease mapping

      2021, Spatial Statistics
      Citation Excerpt :

      For example, CARBayes (Lee, 2013), OpenBUGS (Lunn et al., 2009), NIMBLE (de Valpine et al., 2017), and Stan (Carpenter et al., 2017), (Stan Development Team, 2018a) make use of MCMC methods, each with their own specifications, whereby Stan also offers the possibility to use an approximation method, namely variational inference. Most comparisons between Bayesian software applications focus on parameter estimation and prediction (Bivand et al., 2017; Vranckx et al., 2019). Model selection, which aims at appointing a representative model from a set of candidate models, given the data, is usually not considered in the comparison of the software tools.

    • A spatial multilevel analysis of the impacts of housing conditions on county-level life expectancy at birth in China

      2020, Applied Geography
      Citation Excerpt :

      The significant provincial-level effects on county-level LEB, the significant spatial autocorrelation of county-level LEB, and the comparisons of the spatial autocorrelation of residuals suggested the better fit of spatial multilevel models, compared to traditional regression models. Similarly, to avoid the potential biased estimation of standard errors caused by ignoring spatial homogeneity or spatial heterogeneity (Pierewan & Tampubolon, 2014; Steele, 2008), the application of spatial dependence multilevel models in studies linking aggregated health data and environmental data was discussed and recommended in some studies (Bivand et al., 2017; Ren et al., 2013). Therefore, this empirical study further underscores the importance of considering both spatial homogeneity and spatial heterogeneity in studies of spatially structured health data.

    • Comparison of different software implementations for spatial disease mapping

      2019, Spatial and Spatio-temporal Epidemiology
      Citation Excerpt :

      On the other hand, Neyens et al. (2018) concluded that estimates of INLA in its default state and WinBUGS are very similar, based on a simulation study that focused on a negative binomial CAR model for count data. Bivand et al. (2017) recently made a comparison of several software implementations of the spatial intrinsic CAR and SAR (simultaneous autoregressive) model on two Gaussian examples in the context of spatial econometrics. They show that spatially structured data may indeed be estimated readily using a variety of approaches and conclude that robust results across methods of estimation are obtained from INLA and WinBUGS implementation of the intrinsic CAR model for Gaussian responses.

    • Spatiotemporal model for estimating electric vehicles adopters

      2019, Energy
      Citation Excerpt :

      Unlike other studies that use a large set of socioeconomic variables to find the profile of EV buyers [53], we sought to minimize the set of these variables to reduce the influence of socioeconomic variables and to better characterize the stochastic nature. Other influential measures could be also evaluated in future studies such as EV utilization costs or benefits through the locational vicinity of charging infrastructure [48]. Such a criterion may be adopted based on the experience of planners and recommendations from policy makers.

    • Inferring neighbourhood quality with property transaction records by using a locally adaptive spatial multi-level model

      2019, Computers, Environment and Urban Systems
      Citation Excerpt :

      ϵ and ζ are model residual terms at the pupil and school levels, assumed to follow independent Normal distributions with variances of σϵ2 and σζ2, respectively. Directly applying standard multi-level models to spatial data using spatial groups as the “second (or higher) level” has been shown to be problematic (Arcaya et al., 2012; Bivand, Sha, Osland, & Thorsen, 2017; Dong & Harris, 2015; Dong, Ma, Harris, & Pryce, 2016). A key concern is that geography exists both “vertically,” in the sense that observations have regional membership, and “horizontally,” in that spillovers between adjacent observations often exist regardless of grouping (Haining, 2003; Owen, Harris, & Jones, 2016).

    View all citing articles on Scopus
    View full text