Managing Forests for Biodiversity Conservation and Climate Change Mitigation

We include biodiversity impacts in forest management decision making by incorporating the countryside species area relationship model into the partial equilibrium model GLOBIOM-Forest. We tested three forest management intensities (low, medium, and high) and limited biodiversity loss via an additional constraint on regional species loss. We analyzed two scenarios for climate change mitigation. RCP1.9, the higher mitigation scenario, has more biodiversity loss than the reference RCP7.0, suggesting a trade-off between climate change mitigation, with increased bioenergy use, and biodiversity conservation in forests. This trade-off can be alleviated with biodiversity-conscious forest management by (1) shifting biomass production destined to bioenergy from forests to energy crops, (2) increasing areas under unmanaged secondary forest, (3) reducing forest management intensity, and (4) reallocating biomass production between and within regions. With these mechanisms, it is possible to reduce potential global biodiversity loss by 10% with minor changes in economic outcomes. The global aggregated reduction in biodiversity impacts does not imply that biodiversity impacts are reduced in each ecoregion. We exemplify how to connect an ecologic and an economic model to identify trade-offs, challenges, and possibilities for improved decisions. We acknowledge the limitations of this approach, especially of measuring and projecting biodiversity loss.

(using the simulation units (SimU) that correspond to the intersection of countries, grid cells, altitude, slope and soil class), the demand is represented on a regional basis and bilateral trade between regions is included.Since it maximizes total surplus, 3 land allocation decisions are based on the profitability of the activities by land use type.It incorporates information from the EPIC agricultural model and is coupled with the G4M forest management model.The model is run recursively with 10 year time steps, from 2000 to 2100.GLOBIOM is used as the land module of two integrated assessment models, WITCH (from the RFF-CMCC European Institute on Economics and the Environment) and MESSAGEix (from the International Institute of Applied System Analysis (IIASA)).
GLOBIOM-Forest is a version of GLOBIOM that focuses on a more detailed representation of the forest sector, while simplifying the representation of the agricultural and bioenergy sectors.
It is a bottom-up partial equilibrium model in which total economic surplus (see equation S14) is maximized.In it, spatially explicit supply-related decisions are made on spatial units that correspond to the intersection of a grid, that can be 200km ⇥ 200km or 50km ⇥ 50km resolution, with country boundaries.On the other hand, demand is represented on a regional basis (up to 58 regions -Table S4 4 ).The model includes: (a) transportation costs of woody biomass from forest to mill gate within each region, (b) harvest costs, (c) process costs, (d) investment costs, (e) trade costs, and (f) land use change costs.a,b,f are spatially explicit and c,d,e are on a regional basis.
The model includes a representation of both forestry and the forest industry.Biomass production is via the primary harvested products (pulplogs, sawlogs, industrial plantations biomass, other industrial roundwood, fuelwood, logging residues) and one non-harvested product (deadwood).
Forest industry is represented via the by-products (sawdust, woodchips, bark, black liquor, recycled wood, recycled paper, recycled pulp), the intermediate products (chemical pulp, mechanical pulp), and final products (sawnwood, plywood, fiberboard, other industrial roundwood, newsprint, paper for printing or writing, packaging, other paper, fuelwood, energy wood).
Regarding the production of biomass, GLOBIOM-Forest decides on (1) the area of forest to be harvested during the rotation period for each spatial unit and forest management type and (2) 3 Which is a monetary measure of welfare. 4For the results presented here, the model was run with 58 regions that include 180 countries.
the harvested quantities of a particular primary product in each spatial unit and under each forest management type.In comparison to GLOBIOM, GLOBIOM-Forest includes more than one forest management intensity (Low, Medium and High), has detail on tree species (distinguishing between coniferous (softwood) and non-coniferous (hardwood)) and includes details on age-class dynamics.
Management intensities are defined as a combination of assumptions on (1) the percentage of the increment that is harvested, (2) the limit on logging residues that can be obtained and (3) a minimum amount of the increment that has to be left as deadwood.
Regarding the production of intermediate and final demand products, the model decides on (1)   the quantity of final products to produce by processing primary products and (2) the processing capacity of the main final products.It also decides on the level of investment in each region and each product, which will increase the production capacity.
GLOBIOM-Forest also includes bilateral trade of forest products between regions, deciding the quantities of each product that are exported or imported from region to region.
Energy crops are incorporated through short rotation plantations (SRP).These are represented separately from the previously mentioned forest management types because by a sustainability assumption energy crops are not located in forestland 5 .Instead, the model decides the amount of area devoted to these industrial plantations by transforming from natural land 6 , grasslands or cropland.The representation of these land use changes is simplified in this version compared to GLOBIOM.The biodiversity impact of both forest management and land use transformations to SRPs are represented in this paper.
Biomass for energy can be produced both in forestland and in SRPs.However, the rest of the biomass for forest products can only be produced in forestland.
GLOBIOM-Forest is a recursive optimization model that is calibrated from years 2000 to 2020, and runs in ten-year intervals up to 2100.Data sources for calibration include the Global Forest Resource Assessment (FRA) which provides regional data on forest types and harvest potential between coniferous and non-coniferous tree species; the World Database on Protected Areas (WDPA) for grid level data on forest management and Nature Map Explorer.These three databases were used to improve the allocation of forest and forest management areas during the calibration period.Additionally, the FAOSTAT database was used for reference volumes for demand functions, forest industry production capacities, the separation between coniferous and non-coniferous final products, harvest volumes and net trade quantities.Finally, BACI trade data was used for the bilateral trade quantities.Other sources of data, besides those for calibration purposes, include G4M which provides increments, harvest costs, and total forest area as a result of deforestation and a↵orestation decisions.
The model is run under the Shared Socioeconomic Pathway 2 (SSP2) -the middle of the road scenario for global socioeconomic development to 2100 [3,4] .Within this pathway, two RCPs are explored: RCP1.9, representing a scenario in which global warming is kept below 1.5°C by keeping the increase in radiative forcing to less than 1.9 W/m2 in 2100 [5] , and RCP7.0, the reference scenario for MESSAGE in which increased radiative forcing reaches 7 W/m2 by 2100, here called RCPref.In the model, the SSP a↵ects the GDP and population data which then a↵ect the demand functions.The RCPs a↵ect the bioenergy and wood pellet demand as in [6] .

S2 Data Mappings
Three data components have to be considered, the spatial units, the forest management types and the time periods.

Spatial units
As described in biodiversity model description section 2.2, the biodiversity impact indicator is calculated on an ecoregion basis.On the other hand, the data for a nities in the biodiversity model (from [7] ), is on a continental level (excluding Antarctica).This requires a mapping between continents and ecoregions.The assignation of an ecoregion to a continent is through countries, therefore an initial mapping between ecoregions and countries is required.This was done in ArcGIS using the UIA World Country Boundaries Layer and the WWF ecoregions layer.It is assumed that the a nity value h for the ecoregion will correspond to the a nity value assigned to the continent to which the ecoregion belongs.If an ecoregion has area in more than one country, a weighted average of the a nity factors of the continents was calculated, based on the proportion of the ecoregion area in each continent.
The GLOBIOM-Forest model was run on a 200km x 200km (2º) grid resolution for all countries.
The spatial units in GLOBIOM-Forest correspond to the intersection between country boundaries and this grid.To connect the ecoregion level and these spatial units, an intersection of the two layers was done in ArcGIS.From this, a mapping was created.It includes (1) if the combination between the ecoregion and the spatial unit exists, and (2) the weight mW s,j used in equations 3, S1, S2 and S3.
Management Types [7] contains information on the response ratios, and therefore a nities, for ten management types, whereas GLOBIOM-Forest includes three management intensity types.The mapping used is presented in Table S1.See Table S3 for the definition of each forest management type in GLOBIOM-Forest.
Table S1: Mapping between [7] and GLOBIOM-Forest management types [7]  It is assumed that a nity factors will not change between tree species (coniferous and nonconiferous) and the a nity factor for a GLOBIOM-Forest management intensity will correspond to the average of the [7] management types according to mapping in Table S1.
Land use Types [8] contains information to estimate the a nities for seven land use types.The following table show the mapping used to connect to GLOBIOM-Forest land use types for SRP.

Time periods
The biodiversity cSAR model estimates potential regional species loss by comparing a reference scenario to a future scenario.To combine this with the 10-year periods (from 2000-2100) of GLOBIOM-Forest, the "pristine" reference scenario remains the same for all GLOBIOM-Forest time periods, while the future scenario corresponds to each of the GLOBIOM-Forest runs.
This means that the potential species loss estimated is always with respect to the reference scenario.In tableS3, the maximum increment available for harvest comes from a G4M estimate of the mean annual increment (MAI) based on Net Primary Production (NPP) maps and the normal forest assumption [9] .

S7
S4 Bioenergy demands for the two climate change mitigation scenarios

S5 Mathematical formulation changes to GLOBIOM-Forest
The required changes to incorporate biodiversity to the partial equilibrium model mathematical formulation via the constraint approach [10]  This variable (SRP s,i 0 ) represents the area changed from suitable land use type i 0 to SRP in spatial unit s, measured in 1000 ha.

Auxiliary variables:
A F M i,l : Area under each forest management type i in each ecoregion l in the scenario being analyzed.
A LU i 0 ,l : Suitable area for SRP under each land use type i 0 in each ecoregion l in the scenario being analyzed.
A SRP l : Area under SRP in each ecoregion l in the scenario being analyzed.
Slost Regional g,l : corresponds to the calculation from the countryside SAR model for taxon g and ecoregion l.It represents the potential regional species loss due to habitat loss. Constraints: 1. Defines A F M i,l .Adds over the spatial units that belong to each ecoregion.
Defines A LU i 0 ,l .Adds over the spatial units that belong to each ecoregion.
3. Defines A SRP l .Adds over the spatial units that belong to each ecoregion.
5. Constraint methodology Defines an upper limit for biodiversity loss of current forest management allocation for each taxon About Bmax g For the results presented here Bmax g is a function of the total biodiversity loss per taxa, calculated on an ex-post basis, for the baseline scenario (without the incorporation of biodiversity).First, the baseline model is run.Second, using the results for the harvest areas under each type of management (HARV EST V AR) and the amount of area transformed to SRP (SRP V AR), the biodiversity impact for baseline model (Bmax 0g ) is estimated using the same cSAR model representation that is then used in the model that includes biodiversity.Third, Bmax g is defined according to equation S6.
where % corresponds to the percentage of reduction desired.Here 10%, 20%, 30% and 40% were tested.For the regional-per taxa implementation only 10% was mathematically feasible.The baseline scenario, without the biodiversity constraint, is represented in this document by a 0% desired reduction.1. Defines cSAR1.cSAR1 is defined over interval [0, 1].

Piecewise linear approximation for cSAR
Note that with constraint S7, the countryside SAR model (constraint S4) can be rewritten as, 2. Implement linearization.Replacing cSAR1 z l g,l in equation S8 with a piecewise linear approximation L(cSAR1 z l g,l ).
where the linearization L(cSAR1 z l g,l ) is defined as: 3. Define cSAR1 g,l as the convex combination of break points.
4. Define the lower bound on t k,g,l As mentioned in [10] , since y = cSAR1 z l is a concave function when 0  h  1, then it is guaranteed that no more than two consecutive weights can be nonzero.The decision on m implies a trade-o↵ between model accuracy and the computational time to run.As expected, with increased breakpoints, the linear approximation becomes a better representation of the non-linear cSAR function, but because of the increased number of t k,g,l decision variables, computational time increases significantly when m increases.We defined m = 6 and the domain of cSAR1 as [a 1 = 0, a 6 = 1].

S6 GLOBIOM-Forest regions
In GLOBIOM-Forest, the model includes 180 countries that could be grouped into 58 regions to facilitate the analysis of the results.Table S4.

S7 Why don't economic outcomes change much?
There are small changes in economic outcomes, measured by biomass production in forests, for two reasons.First, to adapt to the imposed biodiversity constraint, the model uses three main strategies: (1) shifting biomass production to SRP, (2) reducing management intensity in managed forests, and (3) reallocating biomass production between and within regions.These strategies allow the model to supply almost the same levels of biomass while reducing the biodiversity impacts.
Second, demanded quantities do not change much.This is because endogenous demand for forest products for material use is highly inelastic [ 0.3, 0.1], and exogenous demand for bioenergy biomass must be satisfied.For endogeneous demands, on a global basis, the largest reduction in demanded quantities is of -4.53% for fuelwood biomass (FW biomass C) for the RCPref, 10% reduction scenario in 2030.
Similar results are obtained for the total global economic surplus, for more details see supplementary S8.1.

S8 Total global surplus
The objective function of GLOBIOM-Forest is to maximize the total global economic surplus CSP S according to where indices a, b correspond to economic regions (See S6), k are products, s are spatial units, i, j are forest management types, f are forest industry production activities, and i 0 are the land use types than can be transformed to SRP.The following are the decision variables: x a,k are demanded quantities, y a,k are supplied quantities per region, y p a,f are the levels of production activities execution, y f s,i,k are the supplied quantities per spatial unit under each management type, e a,b,k are traded quantities, L s,i are harvested areas (HARV EST V AR) in forests, SRP s,i 0 are the harvested areas in SRP (SRP VAR), z f a,i,j is the forest management area changed from i to j, z a,i 0 is the land use change area from i 0 to SRP, and I a,k are the investment quantities.
S8.1 Why does the total global economic surplus have small changes?
On a global basis, total benefits and total costs do not change significantly.Total benefits depend on endogenous demands and endogenous demands are highly inelastic, as mentioned in a previous section S7, resulting in unchanged demanded quantities.
Total costs, as shown in equation S14 depend on several sub costs.Figure S2 shows how global cost components, change when the biodiversity constraint is introduced in 2030 for the scenario with 10% reduction for both climate change mitigation scenarios.
The graph shows how some costs increase while others decrease, resulting in a net small change in total costs.To exemplify this, note that SRP costs, forest management costs and trading costs increase, whereas harvest costs and transportation costs decrease.This is a result of producing less biomass in forests, more biomass in SRP, reallocating production, and changing forest management toward less intensification.

S10 Carbon storage in forests outcomes
Age class dynamics and carbon storage GLOBIOM-forest includes spatial explicit ageclass dynamics which determines available biomass for harvests and biomass stock of forests for each period.Each grid cell of the model includes di↵erent management systems with management specific growth curves and initial age-class distribution.The S-shaped growth curves are based on a Chapman-Richards biomass growth model similar to [11] and [12] .Initial age-class areas are based on [13] .At the grid level, growth curves and initial age-class distributions are calibrated to match biomass data from the G4M model [9,14] .In addition to this, the model is calibrated to match FRA's [15] country level data on forest areas and biomass stocks.After 2020, the age-class dynamics develops endogenously based on periodic harvest volumes, growth curves, and mortality.
Figure S11: Comparison between three ways to restrict biodiversity loss, for the RCP1.9 and 10% reduction scenario for 2030.Regional per taxa is used in this study.The biodiversity impact is presented as total potential regional species loss (total sum of extirpations) as a percentage of the reference number of species included in the model for each taxon.
In this figure, when comparing between Global per taxa and Regional per taxa, we can see the e↵ect of including the level of endemism of the ecoregion in the constraint.Because of how the level of endemism is included in equation S16, i.e. by multiplying V S l , it was expected that the model would not have the same incentives to protect biodiversity in regions with lower levels of endemism (those with V S l very close to zero) in comparison to using the Regional per taxa implementation.
For those ecoregions with very low levels of endemism, for the model, the regional biodiversity impact becomes irrelevant to satisfy the constraint under the Global per taxa implementation.
However, figure S11 shows that on an aggregated level, the total sum of extirpations remains relatively the same.
Between including the constraint per taxa or aggregated (Regional all versus Regional per taxa), there are some di↵erences, with the aggregated constraint resulting in a larger number of global extirpations.Protection of each taxon may force the model to incur in additional mechanisms to protect biodiversity for a particular taxon, for example, amphibians, which as we saw earlier are the taxon with higher shadow prices of satisfying the biodiversity constraint.
For the same scenario, RCP1.9 with 10% reduction in 2030, the Global per taxa implementation is the one with the largest decreases in RW biomass production (9.3%), closely followed by Regional per taxa (8.4%) and Regional all (2.6%).Again, there are no major di↵erences between Global per taxa and Regional per taxa.When comparing between aggregating taxa or constraining each taxon separately, constraining taxa separately results in less biomass production in 2030 and more biodiversity protection with respect to the baseline.
Shadow prices vary in magnitude for the di↵erent ways of incorporating biodiversity as can be seen in Figure S12.This reflects how the di↵erent ways of incorporating biodiversity have an e↵ect on how the model prioritizes biodiversity impacts among taxa and ecoregions and therefore how hard it is to satisfy the constraints.For example, when using Global per taxa, amphibians are no longer the taxon with the highest shadow price and when using Regional all shadow prices decrease significantly (to 576 (2020 USD)) because now the model can interchange taxon protection to satisfy the constraint.This may be the reason why Regional all is the only implementation feasible for a 20% reduction, among the three tested.
Figure S12: Comparison between shadow prices of the three ways to restrict biodiversity loss, for the RCP1.9 and 10% reduction scenario for 2030.Regional per taxa is used in this study.

Figure S1 :
Figure S1: Global bioenergy demand for each climate mitigation scenario.
For computational e ciency, GLOBIOM-Forest is modeled as an LP.With constraint S4 the model becomes non-linear.The following are the changes used to linearize the problem by using a piecewise linear approximation for the countryside SAR model.Add the following set, parameter, auxiliary variables, and constraints: Set: K: {1, 2, ..., m} Number of breakpoints for linearization.Parameter: a k : Breakpoints for linearization of y = cSAR1 z l g,l where a 1 < a 2 < ... < a m .a 1 = 0 and a m = 1.Auxiliary variables: cSAR1 g,l : Ratio between areas available for species relationship on SAR model.t k,g,l : Weight given to break point k for each combination of ecoregion l and taxon g.Constraints:

Figure S2 :
Figure S2: How cost components change for the mitigation scenarios and the introduction of the biodiversity loss constraint (for a 10% reduction) for 2030.SRP costs (SRPC) include both land use change cost and harvest cost in SRP.

Figure S3 :
Figure S3: Potential regional species loss for mammals in 2030 under the biodiversity loss reduction scenario assessed on a regional basis.

Figure S4 :
Figure S4: Land use change resulting from increased SRP, with its origin between agricultural land or grasslands, for the scenarios assessed, on a global basis, for 2030.

Figure S5 :
Figure S5: Land use change resulting from increased SRP, with its origin between agricultural land or grasslands, for the scenarios assessed, on a global basis, for 2100.

Figure S6 :
Figure S6: Forest areas under each management type, on a global basis, for the scenarios assessed in 2030.Unmanaged forests are Primary forests and Secondary forests.Managed forests are shown under low, medium and high intensity management.

Figure S7 :
Figure S7: Forest area under each management type, for 2030 and 2100, under the scenarios analyzed.Primary forest and secondary forest correspond to unmanaged forests and Low, Medium and High correspond to the three levels of intensity management in managed forests.

Figure S8 :
Figure S8: Forest areas under each management type, on a world basis, for the scenarios assessed in 2100.Unmanaged forests are Primary forests and Secondary forests.Managed forests are the ones under low, medium and high intensity

Figure S13 :
Figure S13: Reference species richness Sorg g,l for all ecoregions for each taxon.

Figure S14 :
FigureS14: Box plot of a nities h F M g,i,l (rows 5-9), h LU g,i 0 ,l (rows 1,2,4) and h SRP g,l (row 3) for all taxon and ecoregions land use type GLOBIOMf land use type

Table S3 :
Description of forest management types in GLOBIOM forest when deciding on forest land types harvest for the rest of the world (ROW).Its a mix of several low intensity management types: Retention forestry, "Nature" management, Selective logging with uneven age management.
Sorg g,l : Number of species of taxa g present in ecoregion l in the reference scenario.Aorg l : Natural habitat area in the reference scenario in ecoregion l.See assumptions in the biodiversity model description section 2.2.This variable (L s,i ) represents the area of forest that will be harvested in each spatial unit under each forest management type during the rotation time, measured in 1000 ha.SRP s,i 0 =SRP VAR(COUNTRY,ALLCOLROW,AltiClass,SLPCLASS, SoilClass, AEZCLASS, LC TYPE SRP) are explained next.It is important to note that the constraint on biodiversity was activated after the GLOBIOM-Forest calibration period, i.e., for 2030-2100.Sets:Defining the following sets according to GLOBIOM-Forest, S: set of spatial units in GLOBIOM-Forest indexed in s *.F : set of forest management types in GLOBIOM-Forest indexed in i. Includes {Primary forest, Secondary forest, Low intensity management, G: set of taxonomic groups indexed in g.Includes {Mammals, Birds, Amphibians, Plants}.L: set of ecoregions of the world indexed in l.*In GLOBIOM-Forest, the spatial units correspond to the intersection of the indices that correspond to the sets COUNTRY, ALLCOLROW, AltiClass, SLPCLASS, SOILCLASS and AEZ-CLASS.ALLCOLROW represents the gridcell which can be defined and used under two di↵erent resolutions, 200km x 200km (2°) or 50km x 50km (0.5°).The other sets represent Altitude, Slope, Soil and Agro-ecological zones (AEZ), respectively.hLUg,i 0,l : SAR model parameter that reflects the a nity of taxonomic group g to land use type i 0 in ecoregion l. h SRP g,l : SAR model parameter that reflects the a nity of taxonomic group g to SRP in ecoregion l.Bmax g : Maximum number of species allowed to disappear regionally across the ecoregions, due to habitat loss caused by forest management decisions, for each taxon g. mW s,l : Weight of spatial unit s in ecoregion l.Based on area.Decision variables from GLOBIOM-Forest:The decision variables from GLOBIOM-Forest that will be connected with the biodiversity model Countryside SAR are: L s,i =HARVEST VAR(COUNTRY, ALLCOLROW, AltiClass, SLPCLASS, SOILCLASS, AEZ-CLASS, ForMngType)

Table S4 :
Regions and associated countries in GLOBIOM forest