Methodological approaches to exploring the spatial variation in social impacts of protected areas: An intercomparison of Bayesian regression modeling approaches and potential implications

1.   Introduction

A large number of countries around the world have committed to protect 30% of land and 30% of water by 2030 ^[1]. This ambitious target has been set in order to both halt biodiversity loss and improve climate change adaptation via nature-based solution. Protected Areas (PAs) can be defined as specific geographical area, recognised, dedicated and managed, through legal or other effective means, to achieve the long-term conservation of nature with associated ecosystem services and cultural values ^[2]. This definition reveals two key characteristics of PAs. First, that PAs have specific boundaries within which certain regulations are enforced and second, that PAs have multidimensional targets focusing both on environmental and socio-economic aspects. Regarding the latter point scholars have highlighted the need to integrate local values into decision-making processes and put actions in place in order to maintain these values within PAs ^[3].

Due to the challenge of meeting these multidimensional targets for PAs, it is often the case that socio-economic aspects are not taken into consideration during PA designations or at least not to the same extent as ecological criteria. However, a growing body of literature has shown that PAs can have significant socio-economic benefits, such as the improvement of wellbeing, involvement in recreational activities and increase in personal income through tourism ^[4,5,6,7,8]. PAs can also have diverse negative impacts for local people ^[9]. This is because the designation of PAs often introduces key restrictions on human activity within certain geographical boundaries such as restrictions on fishing, logging and recreation ^[10,11]. Public perception of the social impacts of PAs, both positive and negative are important as they have been shown to influence public support for the PA ^[7,8,12]. It is therefore important to be able to identify and measure the type, extent and spatial distribution of such impacts in order to maximise PA management effectiveness both for nature and for people.

Despite the increase in research studies on the social impacts of PAs, there is very limited evidence regarding the distribution of social impacts across space within PAs, and relatedly, limited guidance on best practice for spatial assessment of impacts. Naidoo et al. ^[13] was one of the first studies that explored spatial aspects of social impacts and found that living near a PA has a positive effect on wellbeing levels in countries in the Global South. Jones et al. ^[9] also found that an individual's location in the PA is an explanatory parameter for people's wellbeing level, whilst McGinlay et al. ^[14] found that perceived impacts in the Eifel National Park, Germany influenced public support for the Park, and that impacts varied spatially in relation to settlement patterns and road infrastructure. Despite these contributions, our knowledge of how and why social impacts are distributed across space in PAs is very limited, with very few studies exploring whether social impacts differ between locations within a PA and whether the location of communities influence their perceptions.

The spatial aspect of social impacts is an important area of research, as defining a protected area is a spatial approach to biodiversity and landscape protection and enhancement. The fundamental aim is to set a specific boundary within which a new regulatory and governance regime will alter biophysical and social processes in a way that will benefit ecosystems, but which will also impact on people, potentially both negatively and positively. However, spatial aspects are often ignored when examining the factors that contribute to social impacts. Additionally, ignoring spatial autocorrelation, which refers to the correlation between nearby observations in space, can invalidate statistical modeling results if not accounted for ^[15].

Overall, given the scarcity of studies on the spatial distribution of social impacts in PAs, knowledge gaps are not only theoretical and conceptual, but also methodological, and so there is a need to consider what approaches to measuring, mapping and modeling impacts are likely to prove most helpful in supporting theory-building.

In this paper we compare different modeling approaches to explore links between perceived social impacts and the relative location of local communities in a PA. We additionally make recommendations on which approaches to take that can be used by researchers and practitioners when exploring the social impacts of new and existing PAs from a spatial perspective.

Our framework builds upon previous studies ^[13], while our approach is more comprehensive as we incorporate spatial information through kernel functions that consider proximity effects: i.e., the proximity of respondents to each other, as well as less detailed spatial information that only considers the neighboring proximity between municipalities of local residents.

2.   Materials and methods

2.1.   Data

Primary data were collected using structured questionnaires distributed in the Peak District National Park (United Kingdom) during the Summer of 2020 as part of the project FIDELIO (www.warwick.ac.uk/fidelio). The Peak District national park is located in central England and was established in 1951 (See Figure 1). The area is terrestrial (1, 438 km²) and socio-economic activities in the region are mainly tourism and recreation as well as agricultural activities. The specific research area was selected based on the following criteria that were suitable for our analysis: a) an area that has local communities living inside its boundaries; b) area that is large enough to allow the exploration of spatial distribution; c) area that is designated as National Park.

The National Park is designated as an IUCN category V acknowledging the interaction of local communities with the landscape and the contribution of this interaction to a number of aspects including the natural environment, cultural values and the local economy.

The research team sent 3100 postcards to a randomly selected sample of households in and immediately around the Peak District area inviting them to participate in the survey (for further details about this survey please see full report: ^[16]. This was estimated to be approximately 10% of the total population. The survey was also advertised online via social media and informal networks with the help of the Peak District National Park Authority. In total, 438 responses were received. Through the questionnaire, respondents were asked to provide information about their village of residence in order to capture their location which was essential for the spatial analysis in this paper. 227 respondents (51% of the sample) provided this information. The sampling frame included those who live inside the PA or within a 10km buffer around the PA boundary. The specific distance from the borders of the PAs was inspired by previous research ^[13,17]. A research protocol was designed and followed in order to ensure the validity and reliability of our analysis and derived results. Specifically, prior to the final collection of data a variety of proactive actions was performed, in terms of validity and respondent and non-respondent bias assessment. These actions followed include the pilot testing of the questionnaire to a small sample of respondents in all areas of the PA, to identify potential issues with question wording, terminology, response options etc. Potential selection bias was also considered prior to the data collection, by comparing the demographic characteristics of the target population of the PA in order to be aligned with the characteristics of survey respondents. Further, our sampling frame was distributed throughout all regions inside the PA and within the 10km outside buffer zone, in order to ensure the adequate regional representation of respondents.

Regarding perceptions of social impacts these were captured for 5 different aspects of personal well-being: personal income, quality of life, involvement in recreational activities, social relations and connectedness to nature (Table 1). Respondents were asked "How has the designation of the PA impacted you regarding the following issues in the past years?" and all questions were measured via a 5 point Likert Scale (1-Very negative impact, 2- Negative impact, 3-No Impact, 4-Positive Impact, 5-Very positive impact). These questions were part of a social capital and impact assessment questionnaire which was distributed in the area of the Peak District ^[18].

Furthermore, spatial information was utilized as explanatory variables for the variation in perceived social impacts of PAs, aiming to capture the role of location, and of mutual spatial proximity between local residents on their views of perceived impacts.

Data collected included the location in the form of geographical coordinates, in particular, the spatial data for the analysis involves coordinates that relate to the respondent's location (specifically the village of residence of respondents). To collect the coordinates, the names for each location in our survey were entered into an Excel worksheet. Once opened in the browser, an extension called 'Geocode' was used. To create the polygons, a layer from ArcGIS Online with UK wards was added to the map. Then, the wards were arbitrarily divided between the number of villages in each municipality. To ascertain the average distance of each respondent from the remaining respondents we calculated the Euclidean distance between each location and all other locations in the Peak District dataset using the haversine formula. This calculates the distance between two sets of coordinates. As a summary measure of distance, the average distance was utilized.

2.2.   Statistical analysis

2.2.1.   Bootstrap simulation

In order to understand the robustness of our analysis that uses a real dataset and validate the outcomes, we have simulated our data from Peak District PA using the Bootstrap simulation method ^[19,20]. Bootstrapping simulation approach involves resampling a single data set to create a multitude of simulated samples. Those samples are used to calculate standard errors, confidence intervals and for hypothesis testing. Bootstrap is also an appropriate way to control and check the stability of the results obtained by small datasets.

Simulating our research data prior to conducting the actual spatial analysis modeling can be an important step to help better understand the process, validate our methods, and gain more insights into potential outcomes. Specifically, we have simulated the main variables we are interested in for this study, such as the spatial proximity of local residents and the social impacts of respondents. We apply the non-parametric bootstrap technique for testing using the simulated samples if there are any spatial patterns present relating to the various social impacts. The procedure for bootstrapping essentially relies on resampling from an initial sample. Specifically, one creates B bootstrap samples by sampling with replacement from the original data.

The method of non-parametric resampling using bootstrap methodology for calculating estimates and confidence intervals for the parameter of interest is based on the following general scheme:

i. Sample n observations randomly with replacement from the initial sample of data, say vector y_obs.

ii. Calculate the bootstrap version of the statistic of interest, say $\hat \theta$.

iii. Repeat steps i and ii a large number of times, say B, to obtain an estimate of the bootstrap distribution.

The B samples are called bootstrap samples. In practice, the number of B samples that is chosen for our study is 10, 000.

2.2.2.   Statistical modeling

Due to the nature of the response variables for social impacts of PAs, that is, a 1-5 Likert scale, we fit suitable spatial regression models, such as ordinal logistic spatio-temporal regression models to the data. Within this category, the selection of the appropriate data distribution was investigated by comparing the former to the Poisson ^[21] and the Negative Binomial (NB) ^[22] distributional families. Hence, we test alternative spatial regression models, in terms of distributional specification of the responses, for example Poisson distribution and alternatives, such as the negative binomial and the generalized Poisson ^[23], both being suitable for data that are not equi-dispersed (i.e., the value of mean is different from the value of variance). Especially the negative binomial model is frequently considered as an alternative to the Poisson distribution in cases of over dispersed data.

Further, in the present paper a comparison is conducted, between candidate spatial approaches, in order to examine potential spatial effects of proximity between the local residents of the protected areas. Specifically, we compare stochastic spatial models (see, e.g., ^[24] that utilize precise spatial information in terms of exact co-ordinates of local residents to calculate the distances between respondents, to spatial regression models that use less precise information of proximity, based upon the conditional autoregressive spatial models (CAR) of Besag and Kooperberg ^[25].

This comparison is useful, since it can provide important insights into the spatial patterns of the social impacts of PAs and how spatial proximity between local residents can affect their views regarding these impacts.

We take a model comparison approach to assess the relative significance of the proposed spatial models, hence we chose to fit a series of four candidate distributional specifications for the response variables of the five impacts (Table 1), along with a series of candidate spatial models. For comparison reasons, we also fit a simple baseline regression model that includes only an intercept in order to assess the overall spatial effects by comparing this model to the spatial ones.

In the next sub-sections, we describe in detail the various specifications in terms of distributional assumptions and spatial approaches of the fitted regression models.

2.2.2.1.   Kernel-based spatial regression model(s)

As a measure of spatial proximity between the local residents, kernel-based spatial regression models utilized for our analysis use the average distance (in km) of each resident from all other local residents in the PA. The former measure of spatial proximity between local residents aims to capture patterns of spatial similarity (or dissimilarity) between local residents relative to one another regarding their views on perceived impacts as a result of their proximity to the PA.

The analytical expression of the Bayesian kernel-based statistical models used for the current study is mathematically formulated as follows:

Ordinal Logit model:

Poisson model:

Generalized Poisson (GP) model:

Negative binomial (NB) model:

where ${y_i}$ denotes the ith response value of 1-5 Likert scale $\left({i = 1, 2, ..., 227} \right)$, taking the values of 1, 2, 3, 4, 5, and in the ordinal logistic model, ${Y^*}$ denotes a continuous, unmeasured latent variable which is assumed to give rise to the observed categories 1, 2, 3, 4, 5 ^[26]. Also, ${\beta _0}$ is the intercept. In addition, ${\lambda _i}$ and $\lambda _i^*$ are the parameters of the Poisson and generalized Poisson distribution, respectively and $r, {\text{ }}{q_i}$ are the parameters of the negative binomial distribution. Finally, $\varepsilon$ is the error term.

The generalized Poisson regression model is a generalization of the standard Poisson regression model. When the dispersion parameter ω = 0, the probability function in Eq (3) reduces to the Poisson model (1). When ω > 0, the generalized Poisson model represents count data with overdispersion and when ω < 0, the generalized Poisson model represents count data with underdispersion.

The information collected on the perceived impacts of protected areas (PAs) on local residents is spatial in nature and varies across different locations. To address this, we have incorporated a set of distance measures in Eqs (1)-(4) as additional factors that account for the spatial differences in the prediction of perceived impacts of PAs on local residents. These distance measures will be used as covariates in the four distributional specifications, denoted as ${d_i}$. The average distances, ${d_i}$, were included in the regression models via a typical decaying function of distance.

A commonly used distance kernel function is the inverse distance power function, with:

for positive integer $\alpha$, with $\alpha$ often taken to be 1 ^[27,28]. Another typical function used is the negative exponential, given by:

with $\alpha > 0$, the parameter controlling for the rate of decay with distance ^[27,29].

In the current paper we utilize both kernel functions and compare their performance. In this way, we seek to capture the spatial pattern of density decay with average distance, ${d_i}$, of each local resident respondent from the remaining respondents of the Peak District protected area. The inverse distance power function has a fatter tail over long distances ^[30]. This is also true when compared to the negative exponential function.

2.2.2.2.   CAR spatial panel model(s)

Conditional autoregressive (CAR) models are regularly used to analyze data in a large range of disciplines, such as in demography, economy, epidemiology and geography ^[31]. A conditional autoregressive modeling approach following the Bayesian paradigm is employed in order to examine the potential spatial effects on the response variables of social impacts.

Let ${y_{is}}$ denote the ith response for the s^th municipality of Peak District region, UK. Hence, the CAR Bayesian regression model is subsequently written in its general form as:

for the CAR model with mean ${\mu _{is}}$, which is suitably adjusted for the various distributional specifications of ordinal logistic, Poisson and generalizations of the latter, as presented already in Eqs (1)-(4). In the model formulation of Eq (7), ${\beta _0}$ is the intercept.

CAR is a spatial regression method that incorporates spatial dependency into the data analysis through the inclusion of a spatial component, say, ${\eta _s}$ that is utilized for the modelling of the space random effects, and is represented as:

In the previous equation, $\backslash s$ denotes all municipalities excluding $s$, ${\bar \eta _s} = \frac{1}{{\sum\limits_j {{w_{sj}}} }} \cdot \sum\limits_{j = 1}^{{N_s}} {{w_{sj}} \cdot {n_j}}$ denotes the means of the neighbouring random effects for municipality $s$, and $\sigma _n^2$ is the corresponding variance.

A weight matrix (also called a spatial proximity matrix) is then utilized with its elements ${w_{ij}}$ being some type of measure of the spatial relationship between the ith and jth municipalities. A common definition for distance-based weights is the weights given by:

The assumption of "neighbours" typically refers to areas that share a common boundary ^[32,33]. We must also note that the property of symmetric matrices of weights ${\mathbf{W}}$, ${\mathbf{W'}}$ must hold for deriving valid results.

2.2.3.   Prior Specification and inference for the Bayesian regression models

Model implementation and inference in this paper is based upon a fully Bayesian approach. For inference by following the Bayesian paradigm, suitably vague priors for the parameters of the non-spatial, CAR and kernel-based spatial regression models are specified. In particular, the prior distributions for the intercept (${\beta _0}$) and the spatial parameters are specified through a Gaussian distribution, with zero-mean vector and a non-informative variance, i.e. $N\left({0, 1000} \right)$. Accordingly, the variance parameters $\sigma _{}^2$, were assigned weakly informative half-normal prior distributions. Finally, in function specifications 5 & 6, we do not assign a specific value on parameter α, instead allowing α to vary by assigning a weakly informative prior distribution.

In order to adopt the Bayesian paradigm for parameter estimation, the WinBUGS software has been utilized ^[34]. The model's parameters were estimated via Markov chain Monte Carlo (McMC) simulation. The convergence of the McMC chains was assessed through visual inspection of the posterior distributions. For selecting the optimal model among the candidate models fitted, suitable goodness-of-fit criteria are utilized, both from the Bayesian and frequentist perspective. The comparison of the models is done through the DIC criterion (Deviance Information Criterion, which is a generalization of the criterion AIC (Akaike Information Criterion) and which was proposed by Spiegenhalter et al. ^[35] and measures the variability of likelihood. Thus, the deviance information criterion (DIC) values reported, with smaller values of these measures indicating towards selection of the corresponding model. The DIC is based on the posterior mean deviance and the effective number of parameters, ${p_D}$, expressed as: $DIC = {p_D} + \bar D$ ($\bar D$ denotes the posterior mean deviance). The smaller the DIC value we have, the better fit the model makes to the data. The WinBUGS code of all spatial and non-spatial fitted models is available upon request by the corresponding author.

3.   Results

3.1.   Descriptive statistics

The following table (Table 2) shows descriptive statistics for the five variables of perceived social impacts reported in the survey at Peak District National Park. On average, local residents of the PA evaluate highly positively the effects of the national park on their lives, with the least positive effects observed on their personal income (average score 3.665 on the 5 point Likert scale). This is also the variable with the highest diversity in the responses of local residents (standard deviation: 1.052) and so the greatest variation from household to household across the Park.

The spatial distribution of the results depended on who responded to the survey out of the residents invited to participate, and who also gave their geographical location. Of greater interest is whether the magnitude of the reported impact scores of a household were related to, or independent of the geographical location in or around the Park, and in particular whether the score given by one household at one location was related to or independent of the scores given by other households nearby. In other words, firstly, did scores vary across the landscape and was location in the Park a significant factor affecting scores given, and secondly, might one household's perception influence that of others nearby. To address these questions, a range of models were used to predict how responses might vary across the Park, and to assess which modelling approach was most informative and useful.

Overall, the modeling results indicate that impact scores from a household were related to where they were in the Park and were not independent of the scores of other nearby households. For the 5 impacts studied, it appears that the impact of the PA on the social relations of local residents is the least affected by the proximity of other residents, and also varied least across the landscape (had the lowest standard deviation) (Table 2).

3.2.   Results of Bootstrap simulation

The Bootstrap subsampling simulation has been performed prior to the statistical spatial regression modeling by the utilization of a total sample of 10, 000 iterations to collect the Bootstrap sub-samples from the five samples of social impacts and the variable of average distances of each resident from all other local residents in the PA.

In Figure A1 in the Supplementary, the histograms and normal probability plots for the overall bootstrap samples, based upon the 10, 000 bootstrap sub-samples, are presented. As we observe, the bootstrap samples follow the normal distribution.

Upon the generation of the bootstrap samples, subsequently we examine potential spatial patterns in the data by calculating the Pearson's correlation coefficients between the five social impacts bootstrap samples and the average distance of respondents (Table A1).

As is seen by the correlations, the simulated values of the average distance between local residents are significantly correlated with all five simulated social impacts variables. This result is an indication that spatial dependence exists regarding the perceived social impacts, and subsequently that these variables are not randomly distributed but instead exhibit patterns or relationships based on their physical locations. Hence, these findings suggest that the spatial dependencies can be examined through techniques like spatial regression, which help uncover how the proximity or spatial arrangement influences social outcomes or impacts in the particular region.

3.3.   Results of spatial regression modelling and comparisons

In this sub-section, the results in terms of the deviance information criterion for all fitted regression-type Bayesian models are reported. Hence, Tables 3-6 show the DIC values for the simple baseline regression model, the CAR panel spatial model and the two alternative spatial regression kernel function models, for each one of the four alternative distributional specifications for the response variables of social impacts (ordinal logistic, Poisson, generalized Poisson and negative binomial).

The results clearly show that the fit of all models based on the ordinal logistic distribution (Table 3) outperformed all other distributional assumptions for the specific dataset, since the deviance information criterion values for the spatial ordinal logistic regression model are the lowest in comparison to all other modeling specifications.

The next best fit is observed for the generalized Poisson specification, as revealed by the DIC estimates in Table 5.

On the other hand, the worst model fit is seen for the negative binomial regression models (Table 6), with similar results observed for the Poisson models (Table 4).

After the selection of the best distributional specification (ordinal logistic, Table 3), we subsequently compare the various alternative spatial model formulations in comparison to the baseline regression model. When considering the alternative spatial specifications and their comparison to the non-spatial baseline regression model, it is observed that the best fit is achieved for the conditional autoregressive panel spatial model under an ordinal logistic distribution for the response variables (Table 3). CAR model outperforms the baseline regression model and the two alternative distance kernel function specifications in terms of model fit according to the DIC estimations (smallest DIC values in Table 3). On the other hand, the kernel-based spatial models have slightly better fit compared to CAR in the case of the Poisson and negative binomial specifications.

Overall, social impacts are best predicted by a spatial CAR regression model that takes into consideration the similarity of responses of residents in neighboring areas, and by considering an ordinal logistic distribution to link the responses to the spatial information.

As a final assessment of the best model performance, Figure 2 compares observed with predicted values of the best fitted spatial regression model (CAR spatial ordinal logistic model). Specifically, blue dots in the five impact graphs indicate a good predictive performance (observed social impact value falls within the 95% credible intervals of prediction), whereas red dots indicate poor performance (observed value falls outside the 95% credible intervals). Results indicate a good fit in most impacts, but with poorer fit observed for the impact on "Social relations" (Figure 2d).

4.   Discussion

Describing, mapping and explaining the spatial clustering of social impacts in and around PA, and of how they are perceived, is key to understanding their causes and to managing and mitigating them. In this paper we have examined the spatial variation in perceived social impacts and associated clustering effects in a large European Protected Area using survey data from local households. In particular, we have implemented generic Bayesian methodology, which offers flexibility in the ability to fit various distributional specifications, both typical and more complex as well as a variety of spatial regression models that are able to capture effectively the potential effects of spatial proximity between local residents of a protected area on their impacts assessment.

Our findings suggest that spatial proximity of local residents to each other in the Peak District National Park plays an important role in how people perceive social impacts, since all spatial Bayesian regression models showed a relatively improved fit compared to the baseline regression model that did not include any type of spatial information in the form of covariate.

Among the spatial regression models tested, the conditional autoregressive model that assigns a weight for adjacent municipalities and ignores spatial effects of non-adjacent areas gave the best fit. There are a few possible reasons why the conditional autoregressive (CAR) spatial model with binary weights outperformed the more detailed spatial regression model that uses exact coordinates to measure distance when analyzing data on local residents of a protected area.

Firstly, the CAR model is designed specifically to account for spatial dependence in the data ^[36]. This type of model assumes that the correlation between observations decreases as distance increases, but does not explicitly model the distance between observations. This can be appropriate when the spatial resolution of the data is low or when the research question is focused on identifying spatial patterns rather than estimating the magnitude of spatial effects. In the case of local residents of a protected area, the binary weights in the CAR model may have been able to capture the underlying spatial structure of the data, resulting in a better fit to the data compared to the more detailed spatial regression model.

Another reason why the CAR model may have outperformed the more detailed spatial regression model is that the latter may have suffered from overfitting ^[37]. When the spatial resolution of the data is high and the spatial regression model includes a large number of covariates or a complex functional form, it can be prone to overfitting. This means that the model fits the noise in the data rather than the true underlying patterns, which can result in poor out-of-sample performance. The binary weights in the CAR model may have provided a more parsimonious and robust approach for modelling the spatial dependence in the data, leading to better out-of-sample predictive performance.

It is worth noting that the choice of spatial model may also depend on the level of spatial resolution of the data ^[38]. If the data is highly aggregated (e.g., at the level of counties or postcodes/zip codes), then a CAR model with binary weights may be sufficient. However, if the data is more finely resolved (e.g., at the level of individual addresses), then a more detailed spatial regression model may be necessary to accurately capture the spatial variation in the data.

Overall, the choice of model will depend on the specific research question, the nature of the data being analyzed, and the level of spatial resolution of the data. Both types of models have their own strengths and weaknesses, and the choice between them should be based on careful consideration of these factors.

The key contribution of this paper lies in the implementation and intercomparison of different modeling strategies concerning the distributional specification of the dependent variable, as well as the focused implementation of spatial modeling in terms of accounting for spatial proximity effects in the estimation of social impacts in PAs.

In doing this, various alternative distributions have been utilized for the most accurate modeling of social impacts in terms of spatial effects assessment. The results of an ordinal logistic regression based on assigning a logistic distribution on the dependent variable have been compared with other suitable distributions for positive count data, such as the Poisson and the negative binomial, as well as to the generalized Poisson regression model, proposed as an alternative to the typical Poisson model when the response data suffers from underdispersion.

The Poisson regression model is typically used for count data, where the response variable represents the number of occurrences of an event in a fixed period of time or space. However, the Poisson distribution assumes that the mean and variance of the response variable are equal, which may not always be the case in real-world data. When the variance is greater than the mean, a negative binomial regression model may be more appropriate ^[39]. The generalized Poisson regression model, on the other hand, allows for the variance to be greater than or less than the mean, making it a more flexible option for count data. It achieves this by including an additional parameter, called the dispersion parameter, which allows for more flexibility in modeling the variance of the response variable ^[40].

The results indicated that the most suitable distribution for such type of response data, being a 5-point Likert scale, is the ordinal logistic regression modeling approach. A useful alternative has been shown to be the Generalized Poisson spatial regression model, which can be attributed to the fact that the specific dataset is an example of under-dispersed data. On the contrary, the negative binomial regression model did not performed well when compared to all other distributional specifications. The finding that the generalized Poisson regression model outperformed the typical Poisson and negative binomial regression models suggests that this particular model may be the most appropriate choice for analyzing the data at hand.

The modelling analysis of our survey findings suggests that the perceived social impacts reported by residents of a protected area are not independent of those reported by people living nearby. That is, their perceptions of impacts are influenced by the perceptions of others living in proximity, the effect diminishing with distance between neighbours. This may be due to a number of factors related to physical distance and social interactions.

Firstly, perceptions of impacts may be clustered because the causes of impacts themselves are clustered spatially, for example, in relation to settlements, infrastructure and processes (e.g., at traffic corridors, hotspots etc.). In the Eifel National Park, Germany, McGinlay et al. ^[14] found that many reported negative impacts were related to noise and disturbance from visitor traffic in the Park. Such disturbance (noise, air pollution, severance, traffic congestion) was spatially concentrated along road access corridors and around local honeypot sites with tourist attractions, local services and parking facilities. So positive perceptions are likely to cluster around sources of benefits, such as attractive areas for recreation and negative impacts around sources of disturbance.

Further to this, however, people's perceptions of the impacts may be amplified or dampened by social relational influences such as the norms, attitudes, opinions and experiences of their friends and neighbours ^[41,42]. People living closer together may have more frequent interactions with each other, leading to a greater sense of community and shared experiences related to living in the protected area. They may also have more opportunities to discuss and share their opinions about the social impacts of residing in the area, which could contribute to more nuanced and informed attitudes, modulating their individual impressions with a further relational aspect.

In contrast, residents living further apart from each other may have fewer opportunities to interact with their neighbours, and may be or feel more isolated and disconnected from the community in general. They may therefore be less likely to discuss and share their opinions about the social impacts of living in the protected area, which could contribute to more individualised perceptions and attitudes. Alternatively a given household may simply be more remote from the communities and locations most affected by a particular social impact, which will therefore have less relevance for them.

5.   Conclusions

We mapped the perceptions of a range of social impacts as reported by local people living in and around the Peak District National Park in an online survey, in order to investigate the spatial patterns of distribution of impacts. The resultant spatial patterns were modeled using a range of statistical modeling techniques in order to assess which approaches may prove most useful in assessing the spatial distribution of impacts.

Our statistical modeling and mapping work indicates that local people's perceptions of the social impacts of the Park tend to be clustered, that is to say, the reported impacts at one location are not independent of impacts at another location.

Among the tested spatial modeling techniques to capture clustering effects, the conditional autoregressive model (CAR) that assigns a weight for adjacent municipalities and ignores spatial effects of non-adjacent areas gave the best fit for the specific data. This could be attributed to technical reasons such as low spatial resolution or overfitting. This finding is also a strong indication that the effect of spatial dependence in this particular PA diminishes in a very rapid way, or in other words, the proximity effect is much stronger between residents living near to each other (e.g., between friends and close neighbors in same community) and has minimum effects between non-adjacent municipalities. The current analysis has generally shown that the choice of best modelling approach for assessing the spatial effects on social impacts may depend on the specific dataset, research question, and the level of spatial resolution of the data. The alternative types of models presented have their own strengths and weaknesses, and the choice between them should be based on careful consideration of the above factors. The findings of this study and the proposed modeling approach may assist in identifying localities inside or close to protected areas with similar social impacts on local residents.

Understanding the spatial distribution of perceived social impacts is key to understanding their causes, and to finding solutions to manage and mitigate them. Indeed, perceptions of social impacts may be clustered because impacts themselves are clustered in relation to settlements, infrastructure and processes (e.g., at traffic corridors, hotspots etc.), but perceptions may also be clustered because of social factors reinforcing perceptions.

Finally, in order to predict in advance the likely future impacts of a new PA, modeling and prediction tools will be needed to assist practitioners in designing PAs. This work should be of utility to policy makers and practitioners in understanding which modeling and prediction tools will be of most use in this regard.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgments

The project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research programme (Project FIDELIO, grant agreement no. 802605).

Conflict of interest

The authors have no conflict of interest to declare.

Appendix