Over the past four decades [1] there has been significant progress in the measurement and analysis of the social, economic and demographic determinants of subjective measures of happiness and well-being. It has now long and convincingly been argued that subjective measures of life satisfaction and happiness are reliable [2, 3, 4] . There has also been a considerable academic and policy-related interest in the analysis of what affects happiness and well-being, ranging from individual demographic (e.g. age, gender, household type) and socio-economic factors (e.g. income, education level) to social, environmental and spatial context (e.g. quality of the environment, natural and human-made amenities) [4, 5, 6, 7, 8] .
At the same time there has been a rapidly increasing and influential body of theoretical and empirical work highlighting the impact of inequality upon happiness and well-being [6, 9, 10, 11, 12, 13 ]. Yet, there is a relative paucity of research that explores the impact of interpersonal and inter-community social relations and of the impact of social and spatial inequalities. Notable exceptions include studies that considered social norms and attempted to model and quantify the effect of interactions between variables at different levels, such as the relationship between being unemployed and regional unemployment rate [14, 15, 16] or to use average regional income as a proxy to infer neighbourhood effects [17] . There have also been efforts to empirically consider social comparison effects by identifying and modelling reference groups [1, 18] .
There have been even more limited efforts to comprehensively consider social and spatial issues and possible spatial interdependencies and spillovers. This may be due to the significant data and methodological challenges regarding this type of socio-spatial analysis. In particular, there are very limited sources of data on subjective life satisfaction at the local level. In most cases the smallest area at which data is typically available and analysed is the region [19] or local authority district [16, 20] .
However, recent research published in this journal [21] addressed some of these challenges by combining suitable social survey and health survey microdata with small area data to build a new public use dataset for community-level life satisfaction in Canada, they identifying important inter-community differences. This new data set also opens up new possibilities for further comprehensive analysis that can provide insights regarding the importance of space, place and space for happiness and life satisfaction. As the researchers who created this data set point out (also inviting more researchers to use their novel new data set), possible 'next steps' and challenges include the examination of plausible sources of the substantial inter-community differences revealed in their original article as well as the analysis of social context variables. In this article we take up the challenge of examining these sources by focusing on possible spatial inter-dependencies and spillovers. In particular, we present a theoretical and methodological framework aimed at addressing the following questions:
Does the collective sense of well-being of an area affect that of a neighbouring area?
There is a large body of academic literature exploring the role of relative social status and inter-personal social comparisons pertaining to a wide range of disciplines. This includes theoretical work on conspicuous consumption highlighting the importance of relative social status [22] as well as the theory of "bandwagon effect" in which mob motivation and mass psychology are taken into account in demand theory [23] . Also of relevance is the "relative income hypothesis" proposing that individual consumption functions depend not only on an individual's own consumption but also on the consumption of other individuals [24] . Other influential work highlighted the importance of social comparisons and reference groups [25] and relative deprivation [26] .
More recently, there has also been very impactful work by social epidemiologists with compelling empirical evidence regarding 'status anxiety' [27] and of the detrimental impact of inequality upon social and individual well-being, resulting from psycho-social processes and social-evaluative threats [12, 13] .
The importance of social comparisons has also been extensively considered by the literature on the economics of happiness. Of particular relevance is the Easterlin paradox: many countries experienced high levels of growth in their GDPs but at the same time the levels of happiness for these countries remained constant [28, 29, 30] . The impact of inequality and social comparisons has been acknowledged as an important determinant of happiness and as a possible explanation for this paradox [1, 9, 10, 31, 32] .
Nevertheless, there has been a relatively limited number of empirical studies examining the impact of social comparisons on happiness. Most of these studies involved an analysis of comparisons of individual or household income or wealth or consumption to that of the rest of a country in which an individual resides in [1, 18, 33, 34] . There have also been efforts to perform similar analysis at the sub-national level with the inclusion of regional or subregional average income as an explanatory variable in happiness regression models as a proxy for social comparisons [17, 35, 36, 37] . Another approach involves the examination of interactions between individual level and regional or sub-regional level characteristics [14, 15, 16] finding that "unemployment hurts, but it hurts less when there are more unemployed people around" [14, p.346 ]. More recent work also considered the impact of wellbeing inequality at the regional level upon individual well-being finding that higher inequality is associated with lower levels of well-being [38, 39] . There have also been studies that used multilevel modelling approaches to implicitly capture the impact of social context upon (including social comparisons) upon happiness [16] .
Nevertheless, all of these studies considering social comparisons and happiness tend to only focus on individual geographical units (whether these are countries, regions or subregional units or community areas) without attempting to analyse possible spatial interactions, spillovers or interdependencies between geographical units. To the best of our knowledge, there is only one peer-reviewed study to date that explicitly involved analysis of spatial relationships and spillover between small areas and this was possible due the availability of high quality microdata sets made available via special license for the north of the Netherlands [40] . The paucity of studies considering geography, space, place, socioeconomic characteristics and happiness is to some extent due to the lack of publicly available data for geographical units smaller than the region. As also noted in the introduction, a significant positive development towards changing this has been a study published in this journal last year which we briefly discussed in the introduction. This study made available a new public use dataset for community-level life satisfaction in Canada, inviting other researchers to use it for further analysis. In this article we utilize the power of this new data set to make a contribution to the literature on social (and socio-spatial) comparisons and happiness by adopting a spatial econometrics approach.
To assess the role of socio-spatial comparisons in happiness levels, we enhance a standard utility model to include spatial spillover effects. The concept of utility has long been established in Economics in order to describe consumer behaviour and it has now long been suggested that happiness can be considered a suitable proxy to it [1, 18, 41] . Individual utility functions describe the elements that give utility to individuals. However, individuals consume many different goods that add to their utility. The concept of goods is expanded to include individual attributes and circumstances such as being in employment, in good health, owning a house etc. Given that there is an infinite number of goods that an individual can consume, we assume for simplicity that an individual's utility depends only on some personal characteristics. The variables often used to account for variation in utility include income, employment status, marital status, gender and age among others. Building on the idea of comparative life satisfaction which was discussed above, there is a need to include socio-spatial interdependencies in a standard utility model, building on previous pertinent theoretical work [33, 42, 43, 44, 45] .
When we empirically estimate the utility function, the latter is of the form:
where U is a happiness index usually proxied by subjective measures of life satisfaction, while the elements x, y and z represent the income, employment status, marital status etc.
For convenience, we denote all individual characteristics of individual i as Xi.
Previous theoretical and empirical work included the development of models that expand the textbook specification of utility and allow for externalities [14, 46] enabling the explicit modelling of the ways in which the actions of one agent affect the utility of another agent.
Under this scenario, the new utility function that accrues includes also characteristics of others that are likely to affect one's utility. Equation (1) is now of the form:
where Xj represents the characteristics of other individuals that affect individual i. It should also be noted that i≠j.
Equation (2) can be further expanded in order to allow for some new form of externalities.
Individuals interact with each other and based on this interaction they affect and they are affected. The influence that one individual has upon others is not limited to her own characteristics. It could be argued that an individual's utility might also help explain the utility of others. Hence, not only the income or the employment status of others may have a say upon her utility but also their utility per se, might play some role. Having said that, we need to modify equation (2) to account for this kind of externality. Hence, the new form is:
where Uj represents the utility of others and again i≠j.
An important research challenge is to identify a way to empirically measure the impact of Uj on Ui. In the following section, where we discuss the methodological framework, we provide the spatial econometric specification under which we can estimate this coefficient.
The main difference between the empirical part of the current article and the theory is that the latter is mostly focused on individual utility, whereas in our dataset we have aggregate (community level) measures of community's life satisfaction. Since the welfare of a community consists of the welfare of its individuals, we thus generalize the theory discussed above to describe the utility of communities. Therefore, we use the subscript i to denote a community instead of an individual. Accordingly the subscript j that used to denote the others, now denotes the other communities. At the community level however, the variables that might exert an influence towards the life satisfaction of another region should be modified accordingly. Hence, instead of having the personal income, employment status, age etc. as in the case of individual level, we now use the average income of the households, the unemployment rate of the community, and some geographical characteristics like whether the latter is rural or urban along with other community level variables. Having a dataset with the aforementioned variables and the geographical structure required for our purposes, we are able to empirically measure the interdependencies found in community level utilities. A crucial issue that needs to be considered is how 'others' are defined and specified from a modelling perspective. In theory, both at the individual and community area level, with others we include all those units that may influence one's utility. The next section details how we empirically establish the relevant others.
As also previously noted, we adopt a cross-sectional dataset for community life satisfaction in Canada which was made available as part of a study published in this journal [21] . The In the dataset for community life satisfaction in Canada, the measure of a communities' life satisfaction is based on aggregated individual measurements. As such, changes in life satisfaction can be attributed to changes in the community context as well as changes in the composition of the population following migration. In our empirical setup this is not expected to influence results as migration flows are relatively modest [47] . Also, we adopt a cross-sectional approach which is not sensitive to migration over time. It is also relevant to note here that 62.4% of the population in our sample lives for at least 5 years in the same house.
Regarding the Life Satisfaction variable per se, in both surveys, the question that
individuals had to answer in order to record their well-being was of the form " Taking However, in the absence of a better alternative and acknowledging the shortcomings of such measures, we argue that they can provide powerful insights into the debate about "What makes people better-off?" and particularly "Does the interaction between units (e.g.
individuals, firms, communities, countries etc) spur any spatial spillovers in subjective well-being?" The latter question is the underpinning of interdependencies both between individuals and communities. It should be noted that the likely mechanisms underpinning possible spatial interactions and interdependencies are very much dependent on the scale of spatial units. Our research focuses on the spatial level of community as defined above.
Apart from life satisfaction in each community, the dataset contains various community area variables. Table 1 Of particular interest is the standard deviation of Life Satisfaction, given the clear evidence (especially in highly developed economies) of a negative relationship between average well-being and well-being inequality [38, 39] . We use the standard deviation to capture the inequality of life satisfaction within communities. A high standard deviation (variance) of happiness in a community, means that there are both high and low levels in life satisfaction individuals in this community. We argue that the inequality of life satisfaction in a community exerts an (negative) effect upon the general level of happiness and can be considered a community area characteristic. For that reason, we consider it as a potential variable that explains the life satisfaction of neighbouring communities as well.
The approach we adopt originates from a regional economics perspective. The conceptualization of neighbours is based on a spatial weights matrix denoted by W. A weights matrix describes the spatial arrangement of units in a sample [53,54. It is always a squared N×N (WN×N) matrix that reflects the spatial connectivity among spatial units.
Specifically, the dimension of the matrix we are using is 1215x1215 as the number of the communities. Among the weights matrices that are most often used in empirical research are the p-order binary contiguity matrix and the inverse distance matrix. In the former case, assuming that p=1, if two communities share a common border they are assigned as neighbours while all other communities are not neighbours. In matrix terms, this will lead to the value of 1 in the cell that represents these two communities, while the value of 0 will be assigned to non-neighbours. On the other hand, the inverse distance matrix has as its elements the inverse of the distance between each pair of communities. It is based on what has been described as the first law of geography which states that the mutual influence of spatial units dissipates with distance [55] . Irrespective of the method used to construct the matrix, each element in the main diagonal of the matrix is always zero as it represents the connectivity of a community with itself. For all the empirical analysis, the weights matrices have been normalized and no cut-off point has been used for the inverse distance matrix.
Specifically, we have row normalized the binary contiguity weights matrices so as the sum of each row to equal one, while we have divided all elements of the inverse distance matrix with the largest eigenvalue of the matrix [56] . We are using both matrices and variations of them when we explore the data in order to examine the sensitivity of our results of the results to the formulation of the weight-matrix. We argue that a binary contiguity matrix can model better the underlying mechanisms that lead to the interdependencies in wellbeing compared to the inverse distance. It is reasonable to assume that a region can affect another region next to it rather than one at the other side of the country. 
Spatial econometrics is an expansion of the standard Ordinary Least Squares (OLS) regression in order to include spatial interactions using the weights matrix W. The general specification is as follows:
= + Equation (4) contains all possible spatial interactions: with the dependent variables (ρWYi)
, independent variables (WXθi) and with the error term (λWei) and is known as the General
Nesting Spatial model (GNS). Different combinations of spatial interactions give rise to different spatial econometric models. Table 2 , gives a brief description of the models. Following [56] , only models that include a spatial lag of X (WXθi) are able to produce flexible spatial spillover effects. Spatial spillovers are of interest in spatial econometrics as Likelihood. Equation (5) presents the reduced form for the SDM model:
Given that only models that include spatial interactions with the independent variables are able to produce flexible spatial spillovers effects we focus our analysis in such models.
Apart from that, since our intention is to estimate the interdependencies in utilities, we are particularly interested in the coefficient of Uj from Equation (3) which is the same as the coefficient rho (ρ) in Equation (4) . For that reason, we prefer models that include spatial interactions with the dependent variable as well. The significance threshold we follow is that of 5% (denoted in the tables by two asterisks), however, we also consider both the 10% and 1% significance levels (denoted with one and three asterisks respectively).
We also apply Moran's I test in order to estimate both the sign of the spatial autocorrelation and to detect whether there are clusters of high-high, low-low or mixed in life satisfaction communities. In addition, the same test is used along with different spatial weights matrices to examine whether adopting spatial econometric techniques is indeed the appropriate next stage of analysis. The software programs used for the current article are GeoDa 1.14.0, ArcMap 10.5.1 and Stata 16 [57, 58, 59] .
The main variable of interest is the Life Satisfaction of the regions. It has already been stated that we use this variable to proxy utility which in term represents the welfare of the communities. Hence, we start by presenting some descriptive statistics that provide insights for the regression analysis employed later in the article. We proceed by making a first exploration towards the geographical aspects of the communities under examination. It has often been suggested in the literature that rural areas are more satisfied with life than urban ones [61] and especially in more affluent countries [4, 39] . Our results are also consistent with this. Fig 6, This is consistent with relevant recent research presented in the World Happiness Report [39] . A possible explanation for this observation is the difference on how densely populated those two areas are. The mean value of the logarithm of population density in urban areas is 6.92 while for rural is 2.29. As a result, in urban areas inequality is more evident since more individuals reside in these areas and comparisons are more likely to take place. Life Satisfaction appears to be spatially dependent since the high satisfied communities appear to be clustered together. Particularly, it seems that there are many communities in the east of the country (Quebec) with high levels of life satisfaction. This observation suggests that cultural or institutional characteristics may be in place in these communities that explain this pattern. Spatial statistics can further explore the data and shed light to the underlying mechanisms governing these observations. The information provided by the correlation matrix is valuable since we get a first intuition about the relationship between our variables in all possible combinations. A general comment is that the majority of the correlations in Table 3 , are statistically significant. First of all, there is a positive and significant correlation at 1% level between Life Satisfaction and the average Household Income. In contrast, there is a negative correlation between
Life Satisfaction and Unemployment Rate as well as with standard deviation of Life Satisfaction. The former correlation is statistically significant at 10% level while the latter at 1% level. Apart from that, it is also interesting to note that there is a positive correlation between Urban areas, Commute Duration and Population Density as one might expect. The explanation for that is straightforward since the majority of economic activity is often concentrated in urban areas which leads to overpopulation resulting in traffic jams and increased commute duration. Finally, a thought-provoking correlation is observed between the community-level education attainment (Proportion of 4y degree) and community income (Household Income) on the one hand which is positive and significant and between the former and community level of Unemployment Rate (negative and significant) on the other.
24 Table 3 . Correlation Matrix.
(1) 
There are several diagnostic tests available that formally examine whether the residuals from the OLS regression are spatially correlated. If this is the case, we need to employ spatial regression models. In order to be able to examine spatial correlation, we first need to estimate a simple OLS model using our standard explanatory variables. Using then the weights matrices we discussed in the Weights matrix section, we apply those diagnostics. Table 4 , presents the OLS regression while Table 5 carries the test diagnostics on whether we should continue with error lag or spatial lag model The tests in Table 5 report that we cannot accept the hypothesis that the residuals based on the OLS specification are independent and identically distributed. Specifically, the test considered the alternative hypothesis that residuals are correlated with nearby residuals as defined by the different weights matrices. Namely all matrices point to the adoption of a model with spatial interaction in the error term while only for the first order binary contiguity matrix of a rook form and for the inverse distance matrix there is evidence for a model having spatial interaction with the dependent variable. However, based on the theory that we are interested in, we are dictated to use models that can estimate the coefficient rho (ρ) from Equation (4). It has already been stated that coefficient rho is the coefficient of Uj, the utility of others. Subsequently, in the following section, we present the results based on spatial econometric models.
Furthermore, using Moran's I test for spatial autocorrelation, we generate maps showing how the communities in the sample are clustered and whether there is any spatial dependence among them. Fig 12, shows the graphs for the spatial autocorrelation found in Life Satisfaction while Fig 13, carries the maps with clusters of communities for the two main weights matrices adopted. 
We start our analysis using the SLX model for some preliminary exploration of the flexible spatial spillover effects. Subsequently, we employ models with more spatial interactions that are able to estimate the coefficient of U * , for instance the Spatial Durbin Model (SDM) and the General Nested Model (GNS). We acknowledge the difficulty in obtaining the estimated coefficients from the latter model due to over-parameterization [62] . However, we use it for robustness checks. Tables 6 and 7 Table 7 , where we used the inverse distance matrix, the results are less conclusive. There is no significant indirect income effect while the rho coefficient is not significant either, meaning that there are no spillover effects of life satisfaction among the communities when the assumed spatial structure for possible interactions is increased to include a large part of Canada. Once again, only the unemployment rate has a significant indirect effect as in Table 6 . It is interesting that in both Tables 6 and 7, the indirect effect of life satisfaction inequality is not significant suggesting that apart from the direct effect, the inequality found in one community does not relate to the life satisfaction levels of neighbouring ones.
Given French speaking communities and hence there may be differences that could be attributed to cultural or linguistic effects [65] . By focusing in the regions of Ontario, we eliminate the cultural differences that might be responsible for results we found before. Furthermore,
given the different provinces in the entire sample, by focusing in just one province, we achieve a more homogenous sample in terms of institutions. The institutional differences that might have played a role before are now eliminated. Thus, two types of robustness take place in this exercise, the first concerns all the rural and urban communities of Ontario and the second the urban and rural communities around the city of Toronto.
In Tables 8 and 9 , we present the same regression as in Tables 6 and 7 but this time we are   using the new restricted sample with the 274 regions while in Tables 10 and 11 we are using the 355 communities that constitute Ontario, namely those that include the "islands" and the remote areas. Once again, we examine all the different spatial econometric models and we use the two main spatial weights matrices, hence the binary contiguity matrix of the rook form and the inverse distance matrix for sensitivity analysis. In Tables 8-11 
Although there have long been strong theoretical arguments for the analysis of social comparisons of happiness and well-being between people and areas, there are only few relevant empirical studies. However, a recent publication in this journal [21] paves the way for inter-community comparisons to become now possible. In this article, we modified the standard utility models that are typically developed and used by economists, allowing us to explore the spatial spillover effects of life satisfaction between community areas in Canada. We empirically estimated the model using spatial econometric techniques. By applying these techniques, we are able to explore the extent to which characteristics of one community (e.g. unemployment rate) can affect the life satisfaction of a neighbouring community. Furthermore, the impact of a community's life satisfaction on the life satisfaction of neighbouring communities can also be assessed. Thus, we have used life satisfaction data at the community area level in Canada in order to explore possible spatial interdependencies. Overall, we examined possible spatial spillover effects in communities'
happiness. Using data on life satisfaction, we proxy the utility of the communities and we empirically examine whether proximity plays any role.
Our results suggest that communities impact the life satisfaction levels of their neighbouring communities. The drivers of the interdependencies found are not limited to the characteristics of the neighbouring communities. Apart from the income level and unemployment rate, the level of satisfaction of a community affects the satisfaction level of its neighbours. These findings explain the clusters of high-high and low-low in wellbeing communities found in Canada. The results suggest that there is some evidence of spatial interdependencies between communities and that policy makers should take that into account when examining the well-being in rural and urban places. In particular, based on descriptive statistics, we observe that urban areas exhibit lower levels of life satisfaction, which is consistent with previous research [38, 39] . Moreover, the dispersion of life satisfaction in those areas is higher as higher life satisfaction inequality is found in urban areas. Life satisfaction inequality is measured by the standard deviation of Life Satisfaction variable.
Furthermore, we found positive spatially autocorrelated clusters of regions with similar levels of live satisfaction being neighbours (e.g. positive-positive or negative-negative).
The most persistent result from our spatial regressions is the significance of the indirect unemployment effect. In particular this finding suggests a possible 'could be worse' effect regarding unemployment: higher unemployment in neighbouring areas is associated with higher levels of community level subjective well-being. This highlights the possible important role of unemployment in shaping the happiness not only of individuals but also of communities and even of neighbouring communities. In addition, it illustrates the potential of our methodological approach to make a contribution to debates about sociospatial comparisons and well-being between people and places.
Nevertheless, our analysis would have been more powerful if we had access to small area microdata with individual addresses that could have enabled us to perform a similar analysis of spatial interdependencies between actual neighbours. In addition to the paucity of such suitable microdata an additional methodological challenge of such an exercise would have been the computational intensity and computer memory requirements and management, given that a spatial matrix of the total population of Canada would be needed.
Moreover, we acknowledge the limitations of using a cross-sectional dataset instead of panel data. Using panel data we could examine the dynamics in the relationships we are interested in. Having only contemporaneous relations may be misleading but they still provide the interdependencies we would have expected between communities' life satisfaction. We should be cautious when we interpret such results and further research is needed with more disaggregated regional and panel data.
Overall, the work reported in this article presents a new framework that can be used to explore and quantify the extent to which community area happiness measured by individuals' measures of subjective life satisfaction may be affected by neighbouring communities. Similar patterns could be explored regarding other variables and considering issues pertaining to other economic variables (e.g. economic growth in neighbouring regions) but also (and especially given the current Covid-19 crisis) health-related variables (e.g. how the high prevalence of Covid-19 in an area may affect subjective happiness or sense of anxiety in neighbouring areas). The methodological framework specified in this article can be used and adapted to explore a wide range of regional and sub-regional sociospatial interdependencies