Nature’s Contributions to People
We included global maps of ten NCP including vulnerable ecosystem carbon storage, coastal risk reduction, flood regulation, sediment retention, nitrogen retention, crop pollination, fodder production for livestock (including grazing), fuel wood production, timber production, and access to nature1. Vulnerable ecosystem carbon storage is mapped as the above-ground and below-ground ecosystem carbon lost in a ‘typical’ disturbance event2. This includes terrestrial and coastal (mangrove, salt marsh, seagrass) ecosystem carbon pools (aboveground, belowground, and soils), based on what carbon is likely to be released if the ecosystem were converted. Coastal protection, sediment retention, nitrogen retention, and crop pollination were modeled using InVEST models, adapted to be run at global scales1,3,4. Fodder production for livestock, timber production, fuelwood production were modeled using Version 3 of Co$ting Nature5,6. Flood regulation was modeled using Version 2 of WaterWorld7–9. Access to nature is expressed as the number of urban and rural people10 within one hour travel of natural and semi-natural habitat, taking the least-cost path (by foot, road, rail or boat) across a friction layer11. Data sources, units of measurement, and the original spatial resolution of each modeled NCP are summarized in Table S1. Full model details are available in the supplementary materials from Chaplin-Kramer et al.1
All NCP are “realized,” either as an end use or benefit (e.g., timber harvest or livestock fodder production per unit area of land), or, where possible given current data, weighted by number of beneficiaries. Beneficiaries include people downstream of habitats providing flood regulation or water quality benefits, coastal populations protected from coastal storm surge, or people within a certain travel time of natural habitats. We attributed all NCP to the natural and semi-natural land cover classes providing the benefit, excluding developed lands (croplands and urban areas) and unvegetated areas (Table S3). We excluded Antarctica due to lack of data on NCP from that continent.
Biodiversity (Area of Habitat, AOH)
As a measure of biodiversity, we used species “Area of Habitat” (AOH) for 26,709 terrestrial vertebrate species12. AOH are based on species range maps from IUCN but refined using habitat preferences and elevational limits from IUCN Red List data12. AOH are more specific than extent-of-occurrence (EOO) which can overestimate species range sizes13. AOH areas “exclude areas of unsuitable habitat from each species’ range, which reduces commission errors and more closely approximates the actual occurrence of the species”14.
Species AOH ranges were produced for 10,774 species of birds, 5,219 mammals, 4,462 reptiles and 6,254 amphibians with available IUCN range polygon data12. Species range polygons obtained from the IUCN Red List spatial data portal15 and the Birdlife International spatial data zone16 were first filtered for ‘extant’ range then rasterized to a global one km grid in the Eckert IV equal area projection. Individual species range rasters were then modified to only include land cover classes that match the habitat associations for each species. Habitat associations were obtained from the IUCN Red List species habitat classification scheme and were matched to ESA land cover classes for the year 201817. ESA land cover classification data was aggregated from its native 300 m resolution to match the global ten km grid using a majority rule. Species ranges were additionally filtered so that only areas within a species accepted elevational range were included. Global elevation data derived from SRTM was obtained from WorldClim v. 218. For bird species, seasonal range codes 1-3 (1=year-round; 2=breeding range; 3=non-breeding range) were processed individually and stored as separate range files where applicable.
Spatial optimization
We used integer linear programming techniques19 to estimate how much land area is needed to provide different levels of NCP and/or achieve species representation targets. We identified areas that provide the highest value across all ten NCP using spatial prioritization procedures. Specifically, we generated prioritizations using the minimum set formulation of the reserve selection problem, and completing optimization procedures using integer linear programming techniques20. These procedures were completed using the prioritizr R package20,21 and Gurobi22.
The objective function (equation 1) is to minimize the cost of selected planning units. Constraints (equation 2) are used to ensure that the representation targets are met.
To explore the land area required to maintain different levels of NCP provision, we ran the optimization using 19 different targets ranging from 5% to 95% of total NCP value, across all ten NCP, at 5% increments. To explore how much additional area is required to conserve biodiversity, we generated an additional 19 scenarios which added species representation targets, using data on extent of suitable habitat, or Area of Habitat (AOH) data for mammals, reptiles, amphibians, and birds12. We followed previous studies which established targets based on species’ habitat size, with the goal of ensuring that both restricted-range and wide-ranging species are represented12,20,23,24. We assigned a 100% threshold to species with less than 1,000 km2 of suitable habitat, a 10% threshold to species with more than 250,000 km2 of suitable habitat, and log-linearly interpolated thresholds for species with intermediate amounts of suitable habitat. We also assigned a cap of 1,000,000 km2 for species with a large amount of suitable habitat (>10,000,000 km2). These targets should be considered minimum representation targets as they do not account for habitat connectivity, ecological intactness25, species traits26, evolutionary processes, ecosystem representation 27, genetic diversity, or other important dimensions of biodiversity.
Separately, we also ran scenarios in which protected areas from the World Database on Protected Areas28 were “locked in” to the spatial prioritization (that is, protected areas were required as part of the solution in each scenario). This allowed us to estimate how much additional land area would be required to achieve NCP and species representation targets, beyond the current system of protected areas (Figs. S5-S6, Table S4).
For the scenarios that included NCP (but not species), we ran the optimizations globally and at a spatial resolution of 2 km. Due to computational constraints and the large number of species (26,709), 10 km was the finest resolution at which we were able to run optimizations for scenarios that contained both NCP and species targets. NCP models assume low or no provision of NCP in sparsely vegetated areas (e.g. deserts) and human-modified habitats (e.g. croplands.) However, many species rely on deserts and modified landscapes. After running the optimizations, we masked the 10 km optimization results using natural and semi-natural landcover data29 at 2 km. This allowed us to more precisely map and quantify the natural and semi-natural habitats providing NCP, and to align our results with NCP-only scenarios run at a higher spatial resolution1. As a second step, we re-included the prioritized areas for species (at 10 km) to ensure that species that rely on deserts and human-modified habitats were included.
In the present study, solutions which achieved all targets in the least amount of land area (minimum land area) are used in place of a minimum cost objective. Current globally available data on costs (e.g. Naidoo and Iwamura30) have limitations and were considered unsuitable for this analysis. Proxy indicators such as gridded GDP were also considered a poor measure of cost, since costs of safeguarding NCP and biodiversity may be poorly correlated with national or local estimates of economic productivity. Instead of a measure of cost, therefore, we separately include data on development potential across 14 major sectors, as areas with high value for NCP and biodiversity that also have high suitability for development may have high opportunity costs for alternative land uses.
For subsequent analyses, we focused on “prioritized areas”, defined as areas providing 90% of current levels of all ten NCP which also meet all species targets. We overlaid prioritized areas with data on development potential to examine overall conversion risk as well as risks from major sectors.
Development Potential Index
We also integrated spatially-explicit estimates of the potential for habitat conversion by development across several economic sectors31. To create an development pressure map (Fig. S7), we used published Development Potential Indices (DPIs)31 for renewable energy (concentrated and photovoltaic solar power, wind power, and hydropower), oil and gas (conventional and unconventional oil and gas), mining (coal, metallic and non-metallic mining), and commercial agriculture (crop and biofuels expansion). For urban expansion pressure, we created an Urban Pressure Index (UPI) based on global urban growth projections from 2020 to 205032 (see Supplementary Materials for details on the UPI).
DPIs are global, spatially-explicit 1-km resolution maps that depict the suitability of land for potential expansion development expansion by agriculture, renewable energy, oil and gas, mining, and urbanization. Each DPI has standardized 0–1 values that account for sector-specific land constraints that restrict development (e.g., suitable land cover, slope); land suitability for sector expansion based on resource availability (sector-specific yields); and siting feasibility of new development (e.g., ability to transport resources or materials, access to demand centers, existing development, and other economic costs associated with resource siting). For each of the 14 DPIs, we binned the range of values represented into six categories based on standardized z-score ranges to characterize development pressure as “very low” (≤10th percentile), “low” (>10th – 25th percentile), “medium-low” (>25th – 50th percentile), “medium-high” (>50th – 75th percentile), “high” (>75th – 90th percentile), and “very high” (>90th percentile). We calculated z-scores by mean-standardizing values per country to capture national-level domestic demand coupled with global-level demand likely to drive national-level resource extraction to occur within each countries’ highest development suitability for that resource. To identify regions of high development pressure (HDP), we retained the highest value within the 14 classified DPIs (Fig. S7A) and then selected the “high” and “very high” classified cells (i.e., values 5 and 6) (Fig. S7B).
Protected areas
We evaluated the extent to which currently protected areas and other effective area-based conservation measures (OECM) might achieve targets and maintain NCP, and how much additional land area might require conservation or stewardship. For this step, we compiled spatial data to delineate the boundaries of protected areas and other conservation areas worldwide. To achieve this, we obtained the World Database on Protected Areas (WDPA) and the World Database on Other Effective Area-Based Conservation Measures (WDOECM)28. We prepared these data for analysis following standard practices (using the wdpar R package)33. Briefly, we (i) excluded sites within an unknown or proposed designation, (ii) excluded UNESCO Biosphere Reserves34, (iii) transformed the site boundaries to avoid numerical issues associated with geometries that cross the dateline, (iv) reprojected the data to an equal-area coordinate reference system (World Behrmann; ESRI:54017), (v) replaced sites represented as point localities with circular reserves matching their reported area35, (vi) removed slivers and (vii) excluded marine protected areas. Additionally, throughout this process, we implemented routines to detect and repair invalid geometries.
We included WDPA and OECM areas in our analysis in two different ways. First, to calculate the percentage of prioritized areas that are currently protected or effectively conserved, we overlaid prioritized areas with the combined WDPA and OECM areas (Fig. S3). Second, to calculate how much additional land area would be required to achieve NCP and species representation targets, beyond the current system of protected areas, we ran separate spatial optimization scenarios in which WDPA and OECM areas were “locked in” (Fig. S4). Because protected areas and OECM do not necessarily overlap with prioritized areas, locking them in to the prioritization solutions results in larger total land areas to achieve targets (Fig S5).
Methods References
1. Chaplin-Kramer, R. et al. Mapping the planet’s critical natural assets. Nat Ecol Evol 7, 51–61 (2022).
2. Noon, M. L. et al. Mapping the irrecoverable carbon in Earth’s ecosystems. Nat Sustain 5, 37–46 (2022).
3. Chaplin-Kramer, R. et al. Global modeling of nature’s contributions to people. Science 366, (2019).
4. Sharp, R. et al. InVEST 3.8.0 User’s Guide. http://releases.naturalcapitalproject.org/invest-userguide/latest/ (2020).
5. Mulligan, M. CostingNature v3.x Model Documentation. https://docs.google.com/document/d/136OvAO6PSyVBp0gNl9f0_pIAIg-h-4JZGArnf2V6-Us/edit (2022).
6. Mulligan, M. et al. Mapping nature’s contribution to SDG 6 and implications for other SDGs at policy relevant scales. Remote Sensing of Environment 239, 111671 (2020).
7. Gunnell, K., Mulligan, M., Francis, R. A. & Hole, D. G. Evaluating natural infrastructure for flood management within the watersheds of selected global cities. Science of The Total Environment 670, 411–424 (2019).
8. Mulligan, M. WaterWorld: a self-parameterising, physically based model for application in data-poor but problem-rich environments globally. Hydrology Research 44, 748 (2013).
9. Mulligan, M. WaterWorld (AguAAndes) v2.x Model Documentation. https://docs.google.com/document/d/1GKheQFp5_rsZyazwCJxCzeEStQF2jh04x-Dl_oG5yoY/edit (2022).
10. Cattaneo, A., Nelson, A. & McMenomy, T. Global mapping of urban–rural catchment areas reveals unequal access to services. Proceedings of the National Academy of Sciences 118, e2011990118 (2021).
11. Weiss, D. J. et al. A global map of travel time to cities to assess inequalities in accessibility in 2015. Nature 553, 333–336 (2018).
12. Brooks, T. M. et al. Measuring terrestrial Area of Habitat (AOH) and its utility for the IUCN Red List. Trends in Ecology & Evolution 34, 977–986 (2019).
13. Hurlbert, A. H. & Jetz, W. Species richness, hotspots, and the scale dependence of range maps in ecology and conservation. PNAS 104, 13384–13389 (2007).
14. Soto-Navarro, C. et al. Mapping co-benefits for carbon storage and biodiversity to inform conservation policy and action. Philosophical Transactions of the Royal Society B: Biological Sciences 375, 20190128 (2020).
15. IUCN Red List & UNEP WCMC. Species Range Polygons. https://www.iucnredlist.org/resources/other-spatial-downloads (2020).
16. BirdLife International and Handbook of the Birds of the World. Bird species distribution maps of the world. (2019).
17. Santini, L. et al. Applying habitat and population-density models to land-cover time series to inform IUCN Red List assessments. Conservation Biology 33, 1084–1093 (2019).
18. Fick, S. E. & Hijmans, R. J. WorldClim 2: new 1-km spatial resolution climate surfaces for global land areas. International Journal of Climatology 37, 4302–4315 (2017).
19. Rodrigues, A. S., Orestes Cerdeira, J. & Gaston, K. J. Flexibility, efficiency, and accountability: adapting reserve selection algorithms to more complex conservation problems. Ecography 23, 565–574 (2000).
20. Schuster, R. et al. Optimizing the conservation of migratory species over their full annual cycle. Nature Communications 10, 1754 (2019).
21. Hanson, J. O. et al. prioritizr: Systematic Conservation Prioritization in R. (2020).
22. Gurobi Optimization, LLC. Gurobi Optimizer Reference Manual. (2022).
23. Hanson, J. O. et al. Global conservation of species’ niches. Nature 580, 232–234 (2020).
24. Rodrigues, A. S. L. et al. Global gap analysis: Priority regions for expanding the global protected-area network. BioScience 54, 1092–1100 (2004).
25. Riggio, J. et al. Global human influence maps reveal clear opportunities in conserving Earth’s remaining intact terrestrial ecosystems. Global Change Biology 26, 4344–4356 (2020).
26. Harris, T., Mulligan, M. & Brummitt, N. Opportunities and challenges for herbaria in studying the spatial variation in plant functional diversity. Systematics and Biodiversity 19, 322–332 (2021).
27. Sayre, R. et al. An assessment of the representation of ecosystems in global protected areas using new maps of World Climate Regions and World Ecosystems. Global Ecology and Conservation 21, e00860 (2020).
28. UNEP-WCMC and IUCN. Protected Planet: The World Database on Protected Areas (WDPA). https://www.protectedplanet.net (2021).
29. ESA Climate Change Initiative - Land Cover project. Land Cover CCI. https://maps.elie.ucl.ac.be/CCI/viewer/download.php (2017).
30. Naidoo, R. & Iwamura, T. Global-scale mapping of economic benefits from agricultural lands: Implications for conservation priorities. Biological Conservation 140, 40–49 (2007).
31. Oakleaf, J. R. et al. Mapping global development potential for renewable energy, fossil fuels, mining and agriculture sectors. Scientific Data 6, 101 (2019).
32. Zhou, Y., Varquez, A. C. G. & Kanda, M. High-resolution global urban growth projection based on multiple applications of the SLEUTH urban growth model. Scientific Data 6, 34 (2019).
33. Hanson, J. O. wdpar: Interface to the World Database on Protected Areas. (2022).
34. Coetzer, K. L., Witkowski, E. T. F. & Erasmus, B. F. N. Reviewing Biosphere Reserves globally: effective conservation action or bureaucratic label?: Reviewing Biosphere Reserves globally. Biol Rev 89, 82–104 (2014).
35. Visconti, P. et al. Effects of errors and gaps in spatial data sets on assessment of conservation progress: Errors and gaps in spatial data sets. Conservation Biology 27, 1000–1010 (2013).