Protected areas have remarkable spillover effects on forest conservation on the Qinghai‐Tibet Plateau

Although protected areas (PAs) are assumed to reduce natural threats within boundaries, their spillover effects remain equivocal. It is necessary to determine whether PAs truly achieve conservation targets and whether they promote or inhibit natural habitat degradation in adjacent areas by blockage or leakage spillover. This study aims to choose 54 nature reserves (NRs) focusing on forest protection as a case study to assess PA conservation effectiveness and spillover prevalence.


| INTRODUC TI ON
Covering approximately 30% of the land surface (Keenan et al., 2015), forest ecosystems are primary habitats for the majority of terrestrial species and provide essential ecosystem services such as carbon sequestration, hydrologic cycling, energy balance and wood supply (Bonan, 2008;Joppa et al., 2008;Lui & Coomes, 2016;Naughton-Treves et al., 2005). However, with an estimated annual deforestation rate of 0.5% since the 1990s (Taubert et al., 2018), forests are under tremendous pressure from human disturbances and climate change, severely threatening global biodiversity and ecological security (Curtis et al., 2018;Keenan et al., 2015;Lui & Coomes, 2016;Taubert et al., 2018). Designed to restrain forest degradation and conserve biodiversity, protected areas (PAs) are perceived as a vital tool for natural ecosystem conservation Geldmann et al., 2013;Joppa et al., 2008;Schleicher et al., 2019). In recent decades, global protected area networks have rapidly expanded, and the coverage of terrestrial PAs has increased to approximately 15.4% of Earth's land surface (UNEP-WCMC, 2021), whereas their conservation impact and effectiveness in safeguarding forests remain disputed (Pfaff et al., 2014;Spracklen et al., 2015). Accordingly, comprehensive assessments of PA performance on forest conservation are required to evaluate whether these areas successfully achieve their conservation targets.
The effectiveness of PAs can be measured through multiple approaches that focus on extent, strategic location, management adequacy and ecology outcomes such as habitat loss abatement or biodiversity conservation (Joppa & Pfaff, 2011;Schleicher et al., 2019). Outcomes such as forest conservation are not directly measurable because we must compare them with counterfactual scenarios, that is, scenarios of what would have happened if PAs had not been established (Adams et al., 2019;Andam et al., 2008;Joppa & Pfaff, 2010). PAs are sometimes assigned nonrandomly with less human disturbances, deeper slopes, higher altitudes and lower agricultural suitability across landscapes (Joppa & Pfaff, 2010). Therefore, simply comparing protected and unprotected areas is not a sufficient counterfactual analysis and can easily lead to biased results (Andam et al., 2008;Bowker et al., 2017). Recently, statistical methods matching sites according to environmental variables have been used to account for protected area location bias (Ament & Cumming, 2016;Andam et al., 2008;Carranza et al., 2014;Ren et al., 2015). These methods can balance the pre-existing covariate heterogeneity by comparing habitat degradation such as deforestation inside PAs with that of matched control areas (CAs) (Ren et al., 2015).
Several studies have shown the positive conservation effect of PAs (Ament & Cumming, 2016;Andam et al., 2008;Geldmann et al., 2013;Joppa & Pfaff, 2011;Zhao et al., 2019). However, these studies rarely quantified influences or spillovers of PAs on adjacent unprotected areas that sometimes compromise conservation values (Ewers & Rodrigues, 2008;Ford et al., 2020;Fuller et al., 2019;Spracklen et al., 2015). Negative spillover, which is also called 'leakage', usually indicates that land use inside protected areas is displaced to nontarget unprotected areas, which would potentially diminish or entirely offset the net impacts of conservation while accelerating PA isolation and landscape fragmentation (Bode et al., 2015;DeFries et al., 2005;Fuller et al., 2019Fuller et al., , 2020. Conversely, positive spatial spillovers termed 'blockages' may occur if the natural cover loss is restrained beyond PA boundaries. In forest ecosystems, blockages manifest lower deforestation or increased forest cover in the periphery, perhaps due to ecotourism incentives for adjacent land users to maintain forest cover or low investment in transportation infrastructure and industry (Ament & Cumming, 2016;Fuller et al., 2019). Aichi Target 11 emphasizes the importance of buffer zones of PAs and states that area-based conservation should integrate into the wilder landscape (CBD, 2010). Nevertheless, the actual extent of spillovers varies by region. For example, PAs in the Peruvian Amazon have been assessed for leakage (Oliveira et al., 2007), while PAs in the Costa Rica or Brazilian Amazon have been found either inappreciable leakage or blockage (Andam et al., 2008;Fuller et al., 2020;Pfaff & Robalino, 2017). Given the complexity of regional differences and the effect of further global change on biodiversity (Naughton-Treves et al., 2005), lucubrating the direction, magnitude, prevalence and drivers of spillovers are vital for comprehensive management of PAs.
China harbours more than 30,000 types of vascular plants and approximately 2340 species of terrestrial vertebrates Liu et al., 2003). Conserving this flora and fauna is highly dependent on the PA network. Since the establishment of the first PA, the Dinghushan Nature Reserve, China has established more than 11,800 PAs with multifarious types, including national parks, natural reserves (NRs), forest parks, wetland parks, scenic spots, and other special reserves, covering 18% of its land area and 4.6% of its sea area . Among these PAs, a large proportion takes forest ecosystems or arboreal wildlife as their primary conservation targets. Previous studies have evaluated the effectiveness of PA networks in controlling deforestation at the national (Ren et al., 2015;Yang et al., 2019), regional (Wang et al., 2013;Zhao et al., 2019) or single-PA (Gao et al., 2020) scale in China, yet few studies have taken buffer zones into account for systematic assessments of PAs. Furthermore, NRs always serve as the mainstay of PAs (Wu et al., 2018) due to their stricter protection level and perfect management system; thus, habitat loss or deforestation across both a reserve and its surroundings should be identified and quantified. As the 'Third Pole' (in addition to the North Pole and the South Pole) of the world (Qiu, 2008), the Qinghai-Tibet Plateau (QTP) plays an important role in the ecological security of China and Asia (Sun et al., 2012). It provides water, pasture, timber as well as recreational K E Y W O R D S conservation effectiveness, forest ecosystem, matching methods, protected areas (PAs), Qinghai-Tibet Plateau (QTP), spillover opportunities to billions of people inhabiting the plateau and the surrounding regions . In addition, the QTP lies across two global biodiversity hotspots (the Mountains of Southwest China and the Himalaya) (Mittermeier et al., 2011;Myers et al., 2000), and forest ecosystems in its southeast region nourish rich species, especially endangered and endemic species such as the giant panda (Ailuropoda melanoleuca) (Hou et al., 2021) and black snub-nosed monkey (Rhinopithecus bieti) (Su et al., 2019). To date, Chinese governments have established over 150 NRs on the QTP to protect its sensitive ecosystems and biodiversity, while the conservation status and effectiveness of these NRs remain poorly understood. Thus, integrated NR effectiveness with spillovers on inhibiting deforestation is a matter of great urgency.
In this study, rather than simply focusing on impact assessments of PAs or NRs, we applied a matching method to jointly examine the NR effectiveness and the spillover influence on forest protection on the QTP. We identified forest cover changes within NRs and investigated whether the conservation efficiency varied with NR management levels, ages and areas. Especially, we investigated the evidence of NR nontarget impacts (spillover) with distance and sought to determine possible drivers of the direction and magnitude of spillovers.
Based on the results, we further provided suggestions for PA management in the future. cn/), we chose NRs on the QTP designated primarily for protecting forest ecosystems or wild animals whose habitat was the forest.

| Site selection
Because our forest cover data baseline was 2000, we excluded NRs that were established after 2001 (2 natural reserves established in 2001 were included in this study). We also removed proposed reserves that did not have a specific spatial boundary. Through these processes, we finally selected 54 NRs, including 30 national NRs and 24 provincial NRs, as the study area (Figure 1, Figure S1 and Table   S2) to evaluate conservation effectiveness and spillover effect on forest protection.

| Forest cover change data
Forest cover and loss data were sourced from the Global Forest Change dataset (GFCD) 2000-2019 (version 1.7) (Table S1) (Hansen et al., 2013). The GFCD showed raster files of forest cover We did not use the forest gain data in our study because of the unmatched time scale and relatively short duration, which could not represent the forest regeneration in NRs. With more than 2 million users of this dataset by 2018, it is now widely applied both globally and locally (Curtis et al., 2018;Spracklen et al., 2015;Zhao et al., 2019). For local calibration, former research has approved a relatively robust performance of this dataset on forest F I G U R E 1 Spatial distributions of the 54 NRs on the QTP cover and loss in China with local reference dataset and Google Earth visual verifications Zhang et al., 2020).

| Environmental variables
Due to the location bias of NRs, landscape characteristics needed to be considered when we assessed the impact and spillover of selected NRs. We focused on six environmental variables (Table S5) that were potentially associated with forest loss: annual precipitation, annual mean temperature, slope, elevation, distance to major roads and distance to settlements (Ament & Cumming, 2016;Andam et al., 2008;Ford et al., 2020;Joppa & Pfaff, 2010;Ren et al., 2015).
The annual precipitation and annual mean temperature data were obtained from World Clim version 2 with a spatial resolution of 1 km (Fick & Hijmans, 2017). These data were the average annual data for the years 1970-2000. Elevation and slope were derived from the DEM data with a spatial resolution of 90 m downloaded from the United States Geological Survey (USGS) (http://earth explo rer.usgs. gov/). Digital vector data (scale 1:1 million) of major roads (lines) and major settlements (points) were downloaded from the National Catalogue Service for Geographic Information (http://www.webmap.cn/). Major roads included national roads and provincial roads, and major settlements consisted of cities, towns and villages. All spatial data were projected to a WGS 1984 Albers projection in ArcGIS software (version 10.5.1, ESRI, Redlands, CA, USA). The detailed information of the data sources is shown in Table S1.

| Protection effectiveness assessment in NRs
For each NR, we considered three zones: the NR area, a buffer zone (in a 20-km buffer around the NR, where we expect spillover effects to occur), and a control area (CA, in a 20-120-km buffer around the NR). We generated 1000 random sample points (pixels) within each NR (n = 54) and 80,000 random sample points in the CA of each NR. The distance between points was set to at least 100 m to limit spatial autocorrelation while still ensuring the quantity of samples with smaller NRs. Only points that held forest cover in 2000 were selected for matching analysis. Values of forest cover, forest loss and environmental variables (annual precipitation, annual mean temperature, slope and elevation) were extracted at the location of each point. The distances to roads and settlements were calculated for each point.
The genetic matching method was adopted to balance landscape heterogeneity and made an 'apple to apple' comparison between each NR sample and sample in its CA. Different from propensity score matching or Mahalanobis distance matching, genetic matching used a multivariate matching method to find the optimal balance where an evolutionary search algorithm determined the weight each covariate was given and maximizes the balance of covariates across matched treated and control samples (Diamond & Sekhon, 2013;Sekhon, 2011). We divided initial points into four groups with initial forest cover 1-25%, 25-50%, 50-75% and 75-100% for each NR and CA, matching points in each NR with pixels in CA with similar forest cover. We used a calliper of 0.25 standard deviations for each covariate during the matching process.
The NR effectiveness or avoided deforestation was assessed with the following formula: where D ca and D nr are the deforestation rates of CA and NR, that is, the percentage of CA (NR) sample points where forest loss occurs over the total matched CA (NR) points, and E stands for effectiveness of NRs.
A positive E value indicates that the NR achieves a beneficial target impact on curbing forest loss.
We ran a Wilcoxon rank-sum test to compare the median E value of NR with zero. Additionally, we investigated whether NR properties affected forest loss by dividing 54 NRs into different groups based on their size (area < medium area of NRs 479.35 km 2 , area ≥ 479.35 km 2 ), age (establish year 1963-1993, 1994-2001) and management level (national or provincial). Wilcoxon rank-sum tests were then conducted to compare the median E value of NRs in each group with zero and with each other respectively. All statistical analyses were performed in R software version 4.0.4 (The R Project for Statistical Computing, http://www.r-proje ct.org/). Genetic matching was performed individually for each NR and CA using MatchIt Package (Ho et al., 2011) in R.

| Spillover assessment of NRs
We explored NR spillover effects by assessing whether NRs drove forest loss to their circumjacent areas. Previous studies constantly analysed spillover or PA effectiveness using a single buffer radius such as 10 km buffer zones (Carranza et al., 2014;Fuller et al., 2019;Poor et al., 2019). However, spatial heterogeneity was de facto scale-dependent (Lui & Coomes, 2016), which meant that the extent of spillover possibly varied with NR sizes. To better analyse the spatial variation of spillovers, we delineated 4-km-wide concentric buffer zones around each NR at five levels: 4 km, 8 km, 12 km, 16 km and 20 km. We then erased areas of buffer zones that overlapped with other NRs to ensure that all the buffer zones were set fully within unprotected areas.
Analogous sample and matching procedures were implied as those detailed for NR effectiveness analysis to each buffer zone.
Similarly, the spillover value was calculated as the difference between the deforestation in matched CA and in buffer zones. Consequently, we then treated the 20-km buffer zone as a whole area and calculated the final spillover value for each NR by taking into account samples of all five buffer zones and their corresponding CAs.

| Drivers of spillover
We set a generalized linear model (GLM) to estimate potential drivers of the direction and magnitude of spillovers. Seven possible drivers were selected to predict spillover values, three of which were characteristics of NRs: area, age and management level. The remaining four potential predictors were natural or socio-economic factors of spillover zones related to forest loss: mean elevation, mean slope, mean gross domestic product (GDP) and mean population density, which were calculated within 20-km buffer zones using Raster Calculator and Zonal Statistic tools in ArcGIS. Where possible, we intended to calculate the average value for the period 2000-2019 to reflect the whole deforestation period in this study. However, due to the data limitation, GDP spatial distribution data from 2016 to 2019 were not included. Driver details (Table   S2) and data sources (Table S1) could be checked in the Supporting Information.
The leaps library (Lumley, 2020) in R was used to test all possible combinations of potential predictors, and k-fold (k = 10) crossvalidation (James et al., 2013) was applied to select the model with the best predictive performance. The major steps of data analyses in this study were summarized in Figure 2.

| Matching analysis
For most NRs and their buffer zones, we could successfully sample and match 1000 random points. However, inside the buffer zones of Erhai-Zechai NR and Qomolangma NR, the number of successfully matched pairs was drastically lower with 37, 64, 148, 68 and 88 pairs around Erhai-Zechai and 5, 3, 3, 3 and 3 pairs around Qomolangma NR respectively. The first two buffer zones of Wanglang NR had 939 and 443 pairs, yet the remaining three buffers of this NR held 0 (zero) sample points. Therefore, we could not calculate the deforestation rate and estimate the spillovers of Eerhai-Zecha, Qomolangma and Wanglang NRs, leading to 51 NRs with 255 buffer zones for the final spillover calculation. Sample sizes and the numbers of matched pairs in each NR and its corresponding buffer zones were provided in the Supporting Information (Table S3).

| Spillover effect of NRs on forest protection
Spillover effects of NRs on deforestation on the QTP were heterogeneous over different distances (Figures 5a-e and 6). Buffers in the first two distance gradients (0-12 km) held a more significant number of leakages than blockages, while buffers 12-20 km from the NR boundary had more blockages (Figure 6). In total, 122 of 255 buffer zones (51 NRs excluded Eerhai-Zecha NR, Qomolangma NR and Wanglang NR where we were unable to collect sufficient samples, 5 zones per NR) had less forest loss than their matched control samples with positive spillover values ranging from 0.1% to 5.3%. 119 buffer zones had more forest loss than matched CAs, with negative values from −8.84% to −0.1%. Hence, blockages slightly outnumbered leakages at different distances outside NR boundaries ( Figure 6). Furthermore, when examining NR individually, we found both blockage and leakage in five buffer zones for 35 NRs  (Table S4) for more results.
We treated buffer zones outside NRs as a whole spillover area and calculated one spillover outcome (Figure 5f). Spillover value ranged from −3.59% to 2.52% (n = 51). Leakages occurred more frequently than blockages.

| Drivers of NR spillover
Leakages were slightly more prevalent in the national NR group while the number of blockages was nearly the same (Figure 7a). The percentage of leakages increased with NR age ranging from 42.85% to 63.63% (Figure 7c). The spillover direction (blockage, leakage and no obvious effect) was distributed without apparent patterns across NR area quartiles (Figure 7d). The lack of forest-conservation NRs in Qinghai province prohibited the calculation of different types of spillover outcomes and made it hard to classify NRs by province ( Figure 7b). For 51 NR spillover values, the best-selected model had 30.8% DE with R 2 = .42 and included three indicators (Table S6 and Figure S2): NR age (coefficient = −0.0223, p = .0537), mean elevation of NR buffer zones (coefficient = 0.0002, p = .1652) and population density of buffer zones (coefficient = −0.004, p < .001), among which only population density was significant. The mean elevation of buffer zones was positively correlated with spillover values, that is, NR surroundings located at higher altitudes were associated with blockages. NR age and population density of buffer zones were negatively correlated with spillover value, that is, buffer zones of older NRs and with higher population density were associated with leakages. Management level, NR area, mean slope of buffer zones and GDP of buffer zones were not main predictors of spillover value.

| Protection effectiveness of NRs
Our results indicated that the majority (75.9%) of NRs on the QTP were effective for forest conservation within their boundaries compared to matched unprotected CAs. This finding is congruent with previous studies about PA conservation in China or other regions, which reveals that establishing PAs is a practical method to mitigate deforestation and sustaining habitats for flora and fauna (Ament & Cumming, 2016;Andam et al., 2008;Bowker et al., 2017;Carranza et al., 2014;Feng et al., 2021;Wang et al., 2013;Zhao et al., 2019).
Due to the absence of spatial data, we only analysed national and provincial NRs on the QTP, regardless of NRs at the municipal or county level; thus, our results may be conservative to some degree.
The forest loss data we applied in our study helped us reveal the overall trend of NR performance on forest conservation on the QTP. Forest loss results from both natural and anthropogenic disturbances such as logging, cultivation, wildfire, climate change and catastrophic events (Curtis et al., 2018;Dietze & Moorcroft, 2011;Garwood et al., 1979). However, the detailed drivers of forest loss were not a focus of this research.
In our study, the ineffective NRs which failed to avoid forest loss (Figure 3) are primarily geographically aggregated around the Longmen Shan mountain range at the eastern margin of the QTP (Qiu et al., 2015;Zeng et al., 2016), where the Wenchuan earthquake in 2008 occurred, causing substantial damage to the local ecosystem and biodiversity. To some extent, this devastating seismic earthquake provides a possible explanation for the large percentage of forest loss inside these NRs compared to matched CAs in addition to intensive human activities (Huang & Fan, 2013). Therefore, we suggest that further studies conduct finer-scale analyses and consider forest details such as forest type or biomass. PA properties such as management level, area and age interact in complex ways to drive PA effectiveness (Bowker et al., 2017). Our study shows no significant difference between national and provincial NRs on forest conservation, indicating that stricter protected areas are not always more effective at forest conservation (Butsic et al., 2017;Nelson & Chomitz, 2011;Pfaff et al., 2014). One potential explanation is that most national NRs assessed in our study are upgraded from provincial NRs, some of which have just been upgraded for less than 10 years. Even with stricter policies and management, some national NRs cannot make remarkable progress in forest conservation within a relatively short period.

| Spillover and potential drivers
An increasing number of scholars have realized spillover impacts on biodiversity conservation and calls for precise assessment of PA spillovers (Bowker et al., 2017;Ewers & Rodrigues, 2008;Pfaff & Robalino, 2017;Ren et al., 2015;Spracklen et al., 2015), but only a few studies have attempted to quantify them (Ford et al., 2020;Fuller et al., 2019;Lui & Coomes, 2016). Our assessment provides a preliminary exploration of the PA spillover on the QTP of China despite including only forest-related NRs. Results from 255 buffer zones at five distance gradients showed that blockages were slightly F I G U R E 3 Map of 54 NR impacts on forest conservation from 2001 to 2019 after matching process F I G U R E 4 NR effectiveness in avoiding deforestation from 2000 to 2019, including all NRs and three groups with different NR properties: management levels (national NRs, provincial NRs), NR ages (established 1963-1993, established 1994-2001) and NR sizes (area < 479 km 2 , area ≥ 479 km 2 ). All NRs: N = 54, National NRs: N = 30, Provincial NRs: N = 24, NRs established 1963-1993: N = 33, NRs established 1994-2001: N = 21, NRs with area < 479 km 2 : N = 27, NRs with area ≥ 479 km 2 : N = 27. *p < .05, **p < .01, ***p < .001 more prevalent at different distances outside NRs and the positive effect fell off as the distance increased. Although at the regional level, the average positive spillover value decreased with distance, the individual NR did not always exhibit this pattern (Figure 5a-e).
We further treated five separate buffer zones of each NR as an entire area and tried to identify the overall spillover situation. In contrast, leakages became the primary spillover outcome occurring around 26 NRs. One factor that may explain the inconsistency is that spillovers are offset among buffer zones with different distances.
Blockage of a certain distance could easily be neutralized if larger leakages happen nearby. Furthermore, we noticed that the spatial distribution of severe leakages was somehow consistent with NRs suffering high deforestation (Figures 3 and 5f). The earthquake disturbances mentioned before in Sichuan province have also brought extensive forest loss, which influenced deforestation irrespective of protection status. Even so, the results of our study suggest a higher percentage of leakage outcomes from PAs than previous studies (Ament & Cumming, 2016;Andam et al., 2008;Fuller et al., 2019).
These high deforestation rates surrounding NRs in our study underscore the possible tendency of PA isolation. Since PAs are parts of larger ecosystems, the isolation phenomenon is likely to rescale ecosystems, restrict material and gene flows and destroy habitat integrity (Arturo Sánchez-Azofeifa et al., 2003;DeFries et al., 2005;Hansen & DeFries, 2007;Seiferling et al., 2012). Given the negative influence of adjacent leakages, we recommend that periodic spillover monitoring should start when PAs are established, and PA effectiveness should be assessed at a wider range beyond boundaries in future PA management.
Clarifying the drivers of PA spillover is complex because it involves both natural and anthropogenic factors. Our results indicated that spillover-prone NRs did not show precise patterns when classified by management level, NR age, area or location. Model selection F I G U R E 5 Maps of NR spillover impacts on forest conservation from 2001 to 2019 after the matching process. The spillover outcomes were calculated within (a) 0-4 km buffer, (b) 4-8 km buffer, (c) 8-12 km buffer, (d) 12-16 km buffer, (e) 16-20 km buffer and (f) 0-20 km buffer as a whole area F I G U R E 6 Number of blockages/leakages/no effect for buffer zones at different distances from NR boundaries analysis indicates that NR age and population density of NR surroundings are the most relevant predictors of spillover value. Older NRs with large population densities in surrounding areas are more readily associated with leakages. Unsurprisingly, older NRs with longer periods of land-use limitation may transfer resource demand pressure to the surrounding landscape, especially when high population communities live nearby.

| Implications of further PAs management
We chose national and provincial NRs with clear boundaries as our samples in this study due to the spatial PA border limitation. It is important to acknowledge that except for these NRs, the performance of PAs with lower management levels or other types also requires systematic assessments. Although we have assessed forest conservation within and beyond PA boundaries, two challenges remain for further research. First, buffer zones could have equal or higher conservation values or irreplaceability than PAs (Ford et al., 2020;Pressey et al., 1993). Irreplaceability and species richness, especially the richness of threatened species, should be evaluated beyond borders for PAs with high magnitude leakages to prioritize areas for management. Second, in this study, we used deforestation rates to quantify NR impacts and spillover. Still, non-forest land covers such as grasslands, wetlands or lakes are significant habitats for much of the biodiversity and threatened species. With intensified ecosystem degradation and land use expansion (Bond & Parr, 2010;Fuller et al., 2019), PA evaluation and spillover analyses of these habitats should become a focus, especially in biodiversity hotspots.
Our results showed that most NRs in this study exhibited less deforestation, but leakages detected around NRs undermined the conservation achievements. It is crucial to ensure that PAs are developed ecologically, economically and politically sustainable (Cumming & Allen, 2017) and that their conservation achievements are not offset by leakages. In our study, a large number of NRs with leakages are geographically close together in the east of Sichuan province. Livelihoods of households surrounding these NRs highlight the importance of sustainable development for PA neighbourhoods. Incorporating the social-ecological framework at this point may be an essential step for recognizing the plurality of PAs before contrasting conservation successes and failures (Cumming, 2016;Cumming & Allen, 2017;Ostrom, 2009;Ostrom, 2009). Therefore, we suggest that further PA management consider both biodiversity conservation and local socio-economic development, applying improved policies such as payment for ecosystem services or ecological migration. With increasing in PA coverage planned in the future years (CBD, 2021), it is critical to shift our focus more centred on the connectivity and quality of PAs (Barnes et al., 2018) to ensure that they are able to achieve targets of biodiversity protection .

| CON CLUS ION
PAs are considered as a pivotal tool for the conservation of biodiversity and ecosystem integrity. Given intensified land use expansion, it is necessary to evaluate the performance of PAs as well as their spillover effects to achieve desired conservation targets truly. In this study, we used the genetic matching method to give the first evaluation of NR effectiveness and its spillover on forest conservation on the QTP in the recent 20 years. Considering counterfactual impacts, our study revealed that 75.9% of the NRs F I G U R E 7 Proportional distribution and number of blockages/leakages/ no effect, generated across (a) NR management levels, (b) NR provinces, (c) NR age quartiles and (d) NR area quartiles effectively prevented deforestation within their boundaries on the QTP. Spillover effects are spatially heterogeneous. Both positive and negative spillover were detected in NR adjacent areas. NRs age and population density of buffer zones were the most relevant factors to predict spillover value. PAs could easily displace deforestation pressure to close regions reducing connectivity or diminishing the net benefits from conservation. Therefore, we recommend a new framework that involves both PAs and their adjacent areas when assessing conservation successes and failures of PAs to promote PA management in the future.

CO N FLI C T O F I NTE R E S T
The authors declare there is no conflict of interest.

PEER R E V I E W
The peer review history for this article is available at https://publo ns.com/publo n/10.1111/ddi.13466.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings of this study are openly available in Dryad at https://doi.org/10.5061/dryad.k0p2n gf9c.