Abstract Protected areas (PAs) are fundamental for biodiversity conservation, yet their impacts on nearby residents are contested. We synthesized environmental and socioeconomic conditions of >87,000 children in >60,000 households situated either near or far from >600 PAs within 34 developing countries. We used quasi-experimental hierarchical regression to isolate the impact of living near a PA on several aspects of human well-being. Households near PAs with tourism also had higher wealth levels (by 17%) and a lower likelihood of poverty (by 16%) than similar households living far from PAs. Children under 5 years old living near multiple-use PAs with tourism also had higher height-for-age scores (by 10%) and were less likely to be stunted (by 13%) than similar children living far from PAs. For the largest and most comprehensive socioeconomic-environmental dataset yet assembled, we found no evidence of negative PA impacts and consistent statistical evidence to suggest PAs can positively affect human well-being.

INTRODUCTION The world has committed, through the Aichi Biodiversity Targets and the Sustainable Development Goals (SDG), to halt biodiversity loss and increase protected area (PA) coverage (Aichi Target 11 and SDG 15) and to reduce multidimensional poverty by half by 2030 (SDG 1.2) (1, 2). It is crucial to determine whether these goals are synergistic or antagonistic. Recent calls to evaluate interactions between SDGs have highlighted that achieving one goal in isolation may actually have negative consequences for sustainable development foci of other goals (3). Therefore, is the expansion of the world’s PA network—a cornerstone of biodiversity conservation strategies (4–6)—likely to enhance the prospects of achieving global goals around poverty alleviation and human health or to hamper them? Whether conservation activities benefit or harm people living near PAs has been debated extensively (7, 8). The empirical foundation for the debate has been shaped by research using different methodologies across varying temporal and spatial scales (9–13), making it difficult to derive general insights. A recent meta-analysis of 1043 studies concluded that empirical evidence for impacts of PAs on human well-being remains thin: Only 8% of studies examining impacts on material living standards and 1% of studies analyzing impacts on health used rigorous, quantitative methods and data (14, 15). In addition, a separate systematic review found that the few studies that used comparable, quantitative approaches produced a mix of positive and negative outcomes that were highly dependent on context and methodology, making it virtually impossible to detect any global patterns in PA impacts on human well-being (16). To detect these patterns, we need data on PAs, environmental conditions, and indicators of well-being that are sufficiently fine-grained to reflect complex dynamics at local scales but that are consistent and comprehensive enough to enable analyses at global scales. We also require an analytical approach that can disentangle the many, complex factors that shape multidimensional human well-being, allowing the independent impacts of PAs to be revealed. To address these challenges, we developed a georeferenced database comprising information on ~300,000 children and ~190,000 households across 34 countries in the developing world (Fig. 1) (17). We merged household Demographic and Health Surveys (DHS; table S1) on maternal and reproductive health, childhood growth, and household assets with spatial data layers characterizing the biophysical environment and the world’s PAs (18). While human well-being includes multiple dimensions that can be measured in many ways (19), our database allowed us to select proxies for two important aspects of well-being: health and material living standards (15). For both, we examined average PA impacts and whether there was evidence of “pro-poor” impacts [i.e., differential impacts of PAs on the least well-off people (20)]. Our outcome variables for health were early childhood (age 6 to 60 months) height-for-age growth scores relative to internationally consistent World Health Organization (WHO) standards and whether a child is stunted (stunting affects more than 160 million children, often limiting physical and cognitive growth for life, and is defined as whether a height-for-age score is more than two SDs lower than WHO benchmarks) (21). For material living standards, the outcome variables were an internationally standardized household wealth score (derived from the presence or absence in households of a variety of durable goods and assets related to living standards) and whether a household is poor (defined as a household wealth score of less than 1000 international dollars) (22). Rather than construct a multidimensional index of well-being or poverty (23), we analyzed each of these outcome variables separately to allow for possible differential PA impacts on each metric and to avoid any perception that our well-being indicators are comprehensive enough to warrant their own multidimensional index. Fig. 1 Geographic distribution of developing country household surveys. (A) Global distribution of surveys. (B) Inset of Nepal. Dots represent sampling clusters (blue, further than 10 km from a PA; red, within 10 km) in relation to International Union for Conservation of Nature (IUCN) categories I to VI PAs (green polygons) in countries with surveys. To analyze the impact of proximity to a PA on these dimensions of well-being, we identified all households within the database that were located within 10 km of a PA of International Union for Conservation of Nature (IUCN) classes I to VI (24). This distance conforms with previous thresholds at which PAs are thought to exert ecological and socioeconomic impacts (10), although we tested the sensitivity of our results to this threshold (see Materials and Methods). Since PAs are not situated randomly in landscapes but rather tend to occur in more isolated, less productive areas (25), people living near PAs may also systematically differ in socioeconomic attributes that may confound any assessment of well-being (26). We therefore used a conceptual model (Fig. 2) along with quasi-experimental matching techniques (27, 28) to create a “control” group located further than 10 km from PAs that were, on average and in relevant ways, similar to people living close to PAs. We aggregated all children and households living near and far from PAs across all countries and then used Bayesian regression modeling techniques (29) to estimate the impact of PA proximity on our four outcome variables while accounting for the hierarchical, nonindependent nature of our data (see Materials and Methods). We also assessed whether PAs with different characteristics—age, size, IUCN categorization, and the documented presence of tourism—exerted differential impacts on the health and wealth of nearby households (table S2) and examined how sensitive our results were to possible hidden bias due to unobserved confounding variables, using Rosenbaum bounds (table S3). Fig. 2 Conceptualizing PA impacts. Possible mechanisms of PA impacts on the health and wealth of nearby people. Individual pathways can be combined to conceptualize an impact mechanism; e.g., pathway ADG suggests how PAs can lead to better health outcomes via income gains from PA-related tourism employment that are then spent on improving children’s health.

RESULTS After matching (fig. S1 and tables S2 and S4), the best impact estimation regression models showed strong effects, in the expected direction, of factors typically associated with human well-being gains. A mother’s education level was the strongest predictor of height-for-age scores and likelihood of stunting, while living in an urban (versus rural) area was the strongest predictor of increased wealth and decreased likelihood of being poor (Fig. 3 and fig. S2). We also observed strong effects for breastfeeding (children not breastfed had lower height-for-age scores and a higher likelihood of being stunted) and for distance to the nearest road (households closer to roads had higher household wealth and a lower likelihood of being poor). For all well-being outcomes, there were also strong effects of survey year; height-for-age and household wealth scores increased, and the likelihood of being stunted or poor decreased, over the 15 years of DHS surveys. This result reflects the general advances in development seen around the world during this period (30). That our statistical models demonstrate the same well-being associations that have been extensively documented elsewhere provides confidence in our results for PA proximity on health and wealth outcomes. However, they also demonstrate the value of integrating environmental and socioeconomic data: Elevation and tree cover (negative) and human population density and rainfall (positive) had similar effects on health and wealth outcomes as several of the socioeconomic variables (Fig. 3 and fig. S2). Fig. 3 Postmatching regression model results. Regression coefficients and 95% credible intervals from Bayesian hierarchical models for the impacts of proximity to PA, as well as additional matching covariates and interactions (e.g., “Within 10 km × PA tourism”), on height-for-age growth scores (A) and likelihood of poverty (B). For (A), positive regression coefficients indicate variables that are associated with increased height-for-age scores in children under 5 years old. For (B), negative regression coefficients indicate variables that are associated with a reduction in the likelihood of household poverty. See fig. S2 for regression results for likelihood of stunting and household wealth scores. Colored symbols represent different categories of predictor variables: green, PAs; blue, environmental conditions; brown, socioeconomic information. Table S2 provides a detailed description of the matching covariates. We used our statistical models to simulate predictions for how proximity to PAs of different types affects the health and wealth of people (Fig. 4). We find that all else equal, a hypothetical move of rural households to within 10 km of PAs with documented tourist visits would result in significantly higher wealth scores (by 16.7% on average) and a lower likelihood of poverty (by 16.1%) compared to similar rural households living further than 10 km from a PA. These impacts rise to 20.1 and 25.7% for wealth and poverty likelihood, respectively, for a scenario where households shift to living close to multiple-use PAs (IUCN categories V and VI), rather than those under stricter protection (IUCN categories I to IV), where tourism has been documented. Similarly, a hypothetical shift to living near multiple-use PAs where tourism has been documented would, all else equal, increase children’s height-for-age growth scores by 9.8% and reduce the likelihood of stunting by 13.4%, compared to similar children living further than 10 km from a PA. The likelihood of poverty would also be 8.8% lower for households that shift to live near multiple-use PAs, even with no documented tourism at these PAs. In contrast, no early childhood growth gains were observed for scenarios where children hypothetically move close to PAs where no tourism has been documented, nor would wealth scores be higher in households moving adjacent to PAs without such tourism. There was also no evidence for any negative impacts of PAs on human well-being in any of our scenarios. Fig. 4 Simulated well-being impacts of PA proximity. Predicted impacts (%) of proximity to PAs of various types, as well as impacts of changes in socioeconomic condition, relative to a baseline scenario, for height-for-age scores and likelihood of stunting of young children and household wealth scores and likelihood of poverty. Baseline = rural household located greater than 10 km from a strict (IUCN categories I to IV) PA having no tourism, with children that are breastfed. x axis is broken because of high percentage impacts of urban residency on household wealth and likelihood of poverty. Context for these PA impacts can be generated by using our models to simulate well-being impacts for variables whose human development effects are more commonly studied (Fig. 4). For example, a hypothetical switch from a rural to an urban household, holding everything else constant, results in a 14.7% increase in height-for-age growth scores and a 20.1% reduction in stunting likelihood, while ceasing breastfeeding for children would result, all else equal, in a 15.3% greater chance of being stunted and a 15.7% reduction in height-for-age growth scores. For wealth, a simulated increase in the number of years of education (from a median of 4 to 7) results in household wealth scores that are 34% higher and a likelihood of being poor that is 34% lower. These examples underscore the fact that the PA impacts we describe are not only statistically significant but also of comparable magnitude to changes in socioeconomic conditions that are typically associated with improved well-being or reduced poverty in the developing world (30). The exception to this comparability was the impact on wealth for a rural-to-urban switch of households [a dominant driver of improvements in multidimensional poverty (30)], which results in a greater than doubling of household wealth scores and an 84% reduction in the likelihood of being poor (Fig. 4).

DISCUSSION Our results demonstrate that for a truly widespread dataset, going far beyond the spatial scope of previous studies, there is empirical evidence that PAs can positively affect human well-being in developing countries. We suggest that there are at least four possible pathways or mechanisms (31) through which this could be occurring. PAs with documented tourist visits (~15% of all PAs in our dataset), regardless of management class, had strong positive impacts on household wealth outcomes. This suggests firstly that such PAs may improve household wealth by generating income or other material benefits via tourism-related employment or affiliated markets that can then be spent on household assets (Fig. 2, pathway ADH). Second, multiple-use PAs where tourism has been documented also resulted in increased height-for-age scores and reduced likelihood of stunting among children. The tourism component of this impact may reflect increased household income that is being spent, in part, on additional food, medicine, or medical clinic visits that improve children’s health (Fig. 2, pathway ADG). These tourism-related pathways for PA impacts provide further evidence that the impacts of nature-based tourism can be positive for people and for wildlife (32–35). The third pathway through which PAs affect human well-being was unrelated to tourism. The likelihood of being poor was reduced in households living near multiple-use PAs (IUCN categories V and VI, ~1/ 3 of all PAs in our database), as compared to similar households living further than 10 km from a PA. This suggests that multiple-use PAs lead to improved environmental conditions experienced by nearby households and that their accessibility—unlike categories I to IV PAs—then allows people to benefit from a greater abundance of useful plants and animals via harvest and sales at markets, resulting in income that can be spent on household assets (Fig. 2, pathway BFH). Last, tourism alone did not improve children’s health outcomes; improvements were seen only in combination with proximity to multiple-use PAs. This suggests a role for improved environmental conditions to positively affect health via pathways BE and BFG (Fig. 2), as has been documented elsewhere (36), although the fact that benefits are seen only at multiple-use PAs suggests that an increased availability of natural resources, rather than enhanced air or water quality, drives the positive impacts. More generally, the synthetic relationships described here across multiple countries can motivate further field studies that test mechanisms for PA impacts in specific countries and PAs; data from such additional empirical research will also strengthen the global evidence base used to assess PA–human well-being relationships. We recognize that a more complete evaluation of well-being impacts of PAs would include additional aspects we could not capture here because of limitations in the availability of high-resolution, global datasets. These include additional components of multidimensional well-being (19), social equity (37), historical displacements and exclusions (38), the opportunity costs of PAs (39), environmental governance (40), and less-tangible PA benefits that are difficult to quantify (19, 41). In addition, our results are based on the location of current PAs, meaning that there is no guarantee that they will hold if PA expansion occurs in areas that are systematically different from existing locations. Last, given DHS survey limitations, our analysis is largely cross-sectional. While preprocessing data by matching and then estimating impacts via regression can perform as well as or better than difference-in-difference estimators (42), it would nevertheless be desirable, where data exist for particular PAs and/or countries, to assess how well-being impacts evolve over time after PA establishment. Despite these caveats, our results suggest that rather than displaying any negative effects, several types of PAs across the developing world have positive impacts on important aspects of human well-being. That multiple-use PAs and PAs with documented tourism improve health and wealth outcomes for the least well-off people in the developing world suggests that the expansion of appropriately managed PAs can make an important contribution to SDGs that target poverty reduction, food security, health, and livelihoods (SDGs 1, 2, 3, and 8). Advancing this area of research will be critical to further inform how targeted investment in PAs can support global goals around both biodiversity conservation and human development.

MATERIALS AND METHODS Database—Children, households, environmental conditions, and PAs The database from which health, wealth, and some of the matching covariates were extracted is more fully described in (17). Briefly, data from all DHS surveys (18) that have been conducted in 39 developing countries since 2000 were included in an initial database. The DHS program conducts nationally and subnationally representative surveys, implemented using a stratified two-stage cluster sampling design, across the developing world. These publicly accessible surveys contain detailed demographic and socioeconomic data at both the individual and household level, obtained by interviewing women and men aged 15 to 49 on a variety of issues related to livelihoods, household assets, reproductive health, family planning, and child health. After eliminating instances where relevant DHS data for our analyses were missing, this resulted in 312,727 observations across 33 countries for early childhood growth and 190,794 observations across 34 countries for household wealth. Key blocks of missing data occurred in Indonesia, Peru, and the Philippines, where questions on stunting were not asked in some or all DHS survey years. We used global spatial data layers on elevation, annual precipitation, tree cover, roads, anthropogenic land transformation, and human population density to characterize the biophysical environment of DHS sampling clusters within which the households in our database were contained. For elevation and annual precipitation, we extracted the corresponding value at each cluster, while for roads, we calculated the distance to the nearest road. For tree cover, anthropogenic land transformation, and human population density, we calculated the average value in a 10-km buffer around each DHS cluster. We used the 2013 version of the World Database on Protected Areas (22) to assess whether DHS sampling clusters were located within or outside a 10-km radius to a PA. We restricted our analysis to PAs in IUCN categories I to VI, removing those PAs that were unclassified since we could not be sure of their management objectives. We also merged previously developed databases on the prevalence of tourism at PAs, classifying PAs with any demonstrable level of visitation as having associated tourism (33, 34). Matching methods Within countries and DHS survey years, we used a genetic matching algorithm that weighted the matching covariates to achieve optimal balance (43). The set of matching covariates reflected our understanding of conditions that (i) might differ among near/far from PA respondents and (ii) were likely to influence the health and wealth outcome variables. Controlling for variables that impact the likelihood of being in treated versus control groups while simultaneously affecting treatment outcomes is a key principle for the assessment of causal impacts from treatment (28). Description of the matching covariates and the rationale for their inclusion are found in table S2. We assessed the resulting covariate balance across one to four matches, ultimately using two matches, sampled without replacement, for each near-PA respondent as a compromise between final sample size (more matches = greater power to detect impacts), bias (larger samples = lower balance), and country-level sampling constraints (some countries did not have enough observations that were further than 10 km from a PA to reach n = 4 or n = 3 matching without replacement). Our country-by-country matching yielded a final dataset of 87,033 children and 60,041 households across 34 countries, including 28,913 children and 20,022 households situated within 10 km of one of 603 PAs (Fig. 2). For the analysis of early childhood growth, matching to the most similar controls resulted in 91% of cases having a resulting mean standardized difference below the 0.25 threshold that is recommended for subsequent regression estimation (28). For the final dataset on which we ultimately assessed PA-proximity impacts on early childhood growth, mean standardized differences for covariates were all below 0.25 (table S4 and fig. S1). Matching to the most similar control households for the wealth analyses produced similar values, with mean standardized differences below 0.25 for 89% of the 333 country-variable combinations and for all variables in the final wealth dataset. Impact estimation regressions After matching, we used Bayesian hierarchical regression models to quantify the impact of PA proximity on the health/wealth outcomes, while accounting for the nonindependence of respondents at three levels: (i) individuals/households within DHS clusters, (ii) DHS clusters associated with individual PAs, and (iii) PAs within countries. We constructed a set of candidate models that varied in (i) hierarchical structure, including random intercepts for sampling cluster and country, and random slopes and intercepts for PAs and countries; and (ii) whether the original matching variables were included in the final regression, since including matching variables in the estimation regression may further reduce bias in estimating the treatment effect, as well as allow any remaining impacts of these variables to be explicitly quantified within the treated-matched set (table S2) (44). Continuous variables were standardized by centering and dividing by two SDs so that their coefficients were directly comparable and on the same scale as binary variables, including the proximity-to-PA treatment variable (45). Support for these impact estimation models was calculated using leave-one-out cross-validation (loo); resulting values and interpretation are analogous to Akaike information criterion values in frequentist model-selection procedures, and 95% SEs for the difference in loo scores were calculated and used to assess whether best models were statistically better than next-best models (46). The best impact estimation model included all of the matching variables in the final regression, as well as the full hierarchical structure (i.e., random intercepts for sampling cluster and random slope plus intercept for PA and country). Support for this model was dominant with respect to the other candidate models for all four outcomes (table S5). All analyses were conducted in the statistical computing software R, particularly using the packages cobalt (47), Matching (48), arm (49), brms (50), and loo (51). Simulating impacts of PA proximity To assess the well-being impacts of proximity to PAs of various types, we simulated a baseline scenario in which predicted values of our outcome measures of well-being were derived from posterior simulations of the Bayesian hierarchical models. For this baseline scenario, we generated predictions holding all continuous predictor variables at their mean values, while binary variables were assigned either a 1 or 0, depending on which value was most frequently encountered in the data. The baseline scenario assumed a rural household living greater than 10 km from a strictly protected (IUCN classes I to IV) PA that had no associated tourism and whose mothers breastfed their children. We subsequently simulated a variety of scenarios that reflected changes in these conditions (Fig. 4), collecting median, lower 95%, and upper 95% values over 1000 simulations for the following impact measure I = W S 1 − W S 0 / W S 0 × 100 where I represents the percentage change in the well-being outcome under the alternative scenario relative to the baseline scenario, W S1 is the predicted value for the well-being indicator of interest (i.e., height-for-age score, probability of stunting, household wealth score, and probability of being poor) under the alternate scenario, and W S0 is the equivalent measure under the baseline scenario. Assessing sensitivity of our model results We assumed that people living within 10 km of a PA were near enough to be affected by its presence, while those further than 10 km were unlikely to be. While there is strong empirical and theoretical evidence for the validity of this threshold (10), we also tested for impacts at two alternative distance thresholds. We reran all of our matching models and subsequent impact estimation regressions using 5 and 20 km, respectively, as alternate distance-to-PA thresholds. The resulting impact estimation models showed few differences in regression coefficient values (fig. S3). The 95% Bayesian credible intervals for coefficients of all predictor variables overlapped with one another regardless of the proximity-to-PA distance threshold used. When examining impacts of PAs on well-being outcomes using these different thresholds (fig. S4), we saw that despite a greater sample size and therefore greater power to detect impacts, all well-being impacts at a distance of 20 km are no longer statistically greater than zero except for those associated with household wealth scores, and even these have declined in absolute value relative to estimates at 10 km. This suggests that a distance of 20 km from a PA may be too great for PA impacts to extend out to. On the other hand, at a distance of 5 km, household wealth and poverty impacts increase, but height-for-age and stunting impacts decrease and are no longer statistically greater than zero. This may be a function of reduced sample size and power; 95% Bayesian credible intervals are larger at the 5-km versus 10-km threshold, reflecting increased variability of estimates. However, it also suggests that potential impacts at a 5-km threshold are being dampened because of the presence in the control group of households that are between 5 and 10 km from a PA and yet are seeing well-being improvements, as per the analysis presented in Results. Note that the well-being impacts associated with changes in socioeconomic variables (Fig. 4) remain relatively robust to PA-proximity threshold changes. We also assessed how sensitive our model results were to the presence of hidden bias via an unobserved covariate that might strongly affect selection into the treatment. We used the Rosenbaum bounds approach as implemented in the “sensitivitymult” package in R (52), which calculates whether differences in outcomes between treated and untreated observations remain statistically significant as the value Γ, which represents the odds of an observed covariate affecting differences between treated and untreated, increases. Lower values of Γ (i.e., close to 1) indicate models that are highly sensitive to the presence of hidden bias, with greater values of Γ indicating models that are more robust to such bias. In our case, the values of Γ at which treatment differences are no longer significant due to hidden bias (table S3) are similar to those from other studies that have looked at PA impacts on poverty (13, 53, 54) and can be characterized as emanating from models that are moderately sensitive to possible hidden bias. Note that this test does not imply that a powerful and unobserved confounding variable does exist; it merely assesses the sensitivity of matching models to hidden bias if such a variable were in existence.

SUPPLEMENTARY MATERIALS Supplementary material for this article is available at http://advances.sciencemag.org/cgi/content/full/5/4/eaav3006/DC1 Fig. S1. Assessment of matching effectiveness. Fig. S2. Postmatching regression model results. Fig. S3. Sensitivity of impact regression model results to PA-proximity threshold. Fig. S4. Sensitivity of scenario simulations to PA-proximity threshold. Table S1. Countries with DHS and associated number of observations used to assess impacts of PAs on human well-being. Table S2. Summary of matching covariates and PA characteristics used in quasi-experimental evaluation of the impacts of PA proximity on human well-being. Table S3. Critical P values from sensitivity tests (Rosenbaum bounds) to hidden bias, showing Γ values at which significant differences between observations within versus beyond 10 km from a PA disappear. Table S4. Absolute values of the mean standardized differences for unmatched versus matched comparison groups of children and households within and beyond 10 km of a PA. Table S5. Evaluation of candidate models for estimating impact of proximity to PA on growth scores, stunting, household wealth, and poverty.

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.