The supporting data contains two types of datasets:

1. Social cost of carbon (SCC) values, extrapolated from
   UK government carbon values (2010-2100) to cover the range of 2000-2500. The dataset is stored in a csv file ("scc_extended.csv") and contains four columns:

1. year: numeric
2. low: numeric, central series minus 50% sensitivity range
3. central: numeric, central estimated time series of modelled monetary value that society places on one tonne of carbon dioxide equivalent (£/tCO2e)
4. high: numeric, central series plus 50% sensitivity range

The detailed method that the UK government used to model carbon value estimates can be found at: https://www.gov.uk/government/publications/valuing-greenhouse-gas-emissions-in-policy-appraisal/valuation-of-greenhouse-gas-emissions-for-policy-appraisal-and-evaluation#annex-1-carbon-values-in-2020-prices-per-tonne-of-co2

2. sampled yearly carbon loss (tCO2e) in project area and in counterfactual scenario, and yearly additionality (Mg CO2e) in four ongoing REDD+ projects, estimated using a combination of JRC-TMF Landsat-based annual time series of land use cover and GEDI L4A footprint-level aboveground biomass density estimates to track forest cover and carbon stock through time. The datasets are stored in csv files (Gola: "Gola_country.csv"; Alto Mayo: "CIF_Alto_Mayo.csv"; RPA: "VCS_1396.csv"; Mai Ndombe: "VCS_934.csv") in long format, containing five columns:

1. year: numeric
2. var: character, either "project", "counterfactual" or "additionality"
3. val: numeric, total carbon flux (or additionality) from the previous year to this year (Mg CO2e)
4. n_sim: numeric, number of repetition
5. boolean, whether the year is larger than project start (t0)

To quantify how forest cover changes over time, we used the annual change collection in the JRC-TMF dataset (Vancutsem et al. 2021), which provides the spatial extent and the annual change of the tropical moist forest (TMF) biome at the 0.09-hectare (30 m × 30 m pixels) resolution from 1990 to 2022, derived from the L1T archive imagery (orthorectified top of atmosphere reflectance). The six following land cover classes were mapped: 1) undisturbed forest, 2) degraded forest, 3) deforested land, 4) forest regrowth, 5) permanent and seasonal water, and 6) other land cover. To use information on forest cover to quantify how carbon stock changes over time in NBS projects, we assume that for each project, we can calculate a reference carbon density value for each land cover class that is stable over time. For this, we used the GEDI Level 4A dataset, which contains footprint-level aboveground biomass density (AGBD) estimates (Mg ha-1) for each 25-m GEDI shot (Dubayah et al. 2020, Duncanson et al. 2022). The AGBD estimates are generated from models linking GEDI waveform-derived canopy height metrics with field AGBD estimates for multiple regions and plant functional types.

We selected GEDI shots occurring from 1st January 2020 to 1st January 2021, and which falling within the project area plus a 30-km buffer around it. The inclusion of a 30-km buffer around the project area is to ensure that enough GEDI shots can be found for each land cover class. For each land cover class, we selected the subset of GEDI shots associated with it as shots that overlap with a JRC-TMF pixel 1) that belongs to the land cover class in question and 2) whose eight neighboring pixels also belong to the land cover class in question. The second condition was included to account for the potential geolocation error up to 10 m of GEDI shots. For each land cover class, we calculated the median AGBD value of all the GEDI shots associated with it. We then estimated belowground biomass and deadwood biomass to be 20% and 11% of AGB, respectively, calculated the total biomass as the sum of aboveground, belowground and deadwood biomass, and converted total biomass to total carbon density by multiplying it by the average carbon density of biomass, taken to be 0.47 for this study (Cairns et al. 1997, Penman et al. 2003, Martin & Thomas 2011).

We adopted a pixel-based matching approach to find the counterfactual scenario for each project, following the PACT Tropical Moist Forest Accreditation Methodology (doi:10.33774/coe-2023-g584d-v5). We sampled pixels in the project area at a density of 0.25 points/ha for smaller projects (≤ 250k ha) and 0.05 points/ha for large projects (> 250k ha). We then sampled candidate matching pixels from the match destination to the amount of ten times the number of sampled project pixels. The match destination is defined as the area of a 2000-km buffer around the project that falls within the project’s country boundary (from the LSIB dataset) and the RESOLVE ecoregion boundaries for all the ecoregions that lie within the project (Dinerstein et al. 2017), excluding all other REDD+ project areas and a 5-km leakage buffer around each of the REDD+ projects (including the project being matched).

For each project pixel in a 10% sample of the sampled project pixel set, we matched it to one candidate matching pixel which has the exact same value for the following categorical variables: (1) Land cover class at t−10, t−5, and t0 (where t is the project start year), (2) Country, and (3) Ecoregion, and which has the minimum Mahalanobis distance (Mahalanobis 2018) across the following continuous variables: (1) Elevation from the SRTM data (Jarvis et al. 2008), (2) SRTM-derived slope, (3) Accessibility (Weiss et al. 2018), and (4) Coarsened proportional cover of undisturbed forest and deforested land, at 1200 m × 1200 m resolution, within a 1-km radius buffer around the pixels at t−10, t−5, and t0.

We deemed the matching results valid if all the standardized mean differences (SMD) of each continuous matching variable between the sampled project pixels and the matched pixels is smaller than 0.2 (unless if a continuous matching variable with an SMD > 0.2 is distributed in the range [0, 1], and the value in one of the pixel sets is close to 0 or 1: this is because (as near those values SMD becomes misleading). We performed 20 repetitions of the matching process, each time using an independent sample of the project pixels as input and producing a set of matched pixels as output. The matched pixel sets are the ”counterfactual scenarios” of the project area, and can be considered to be representative of the trajectory of forest cover change and carbon flux in the project area if the project had not existed.

We then evaluated the carbon losses of both the pixels in the project area and the pixels in the counterfactual scenarios, at a yearly interval. For each year within the JRC-TMF time series (1990-2021), for both the project area and the counterfactual scenarios, we calculated the proportion of pixels in each JRC-TMF land cover class, and used the GEDI L4A-derived estimates of total carbon density for each land cover class, described in the previous section, to calculate the total carbon stock (Mg CO2), and calculated the mean total carbon stock value of all 100 counterfactual repetitions. We then calculated carbon losses (lt) of each year t in both the project area (p) and the counterfactual (c) as the difference between the carbon stock of that year (bt) and that of the previous year (bt−1): lt = bt−1 − bt. Finally, we calculated annual carbon drawdown (at) as the difference between the project carbon loss (ltp) and counterfactual carbon loss (ltc): at = ltc − ltp.