This paper develops a methodology for extending multilevel modelling to incorporate spatial interaction effects. The motivation is that classic multilevel models are not specifically spatial. Lower level units may be nested into higher level ones based on a geographical hierarchy (or a membership structure-for example, census zones into regions) but the actual locations of the units and the distances between them are not directly considered: what matters is the groupings but not how close together any two units are within those groupings. As a consequence, spatial interaction effects are neither modelled nor measured, confounding group effects (understood as some sort of contextual effect that acts 'top down' upon members of a group) with proximity effects (some sort of joint dependency that emerges between neighbours). To deal with this, we incorporate spatial simultaneous autoregressive processes into both the outcome variable and the higher level residuals. To assess the performance of the proposed method and the classic multilevel model, a series of Monte Carlo simulations are conducted. The results show that the proposed method performs well in retrieving the true model parameters whereas the classic multilevel model provides biased and inefficient parameter estimation in the presence of spatial interactions. An important implication of the study is to be cautious of an apparent neighbourhood effect in terms of both its magnitude and statistical significance if spatial interaction effects at a lower level are suspected. Applying the new approach to a two-level land price data set for Beijing, China, we find significant spatial interactions at both the land parcel and district levels.