Forest management practice for enhancing carbon sequestration in national forests of Korea

Abstract This paper examines forest management planning and its possible outcomes using linear programing (LP). More specifically, the most appropriate forest harvesting schedule was selected that can maximize the carbon sequestration in the current forest areas considering forest manager’s income. The LP model allows the managers to segment forests into cutting units under rotation basis logging activities. Through harvest prescription from LP, we derived the balanced age-class distribution that constitutes improved conditions for sustainable use of forest resource. However, the solutions from LP did not achieve normal forests with perfectly even aged distribution. Instead, it produced a left-skewed age-class distribution due to the cost restriction of management ruling out the achievement of a normal forest as an optimal solution. The results from our LP model also confirm that the forest management activities will enhance yearly carbon sequestration in forests for all scenarios compared to baseline, and the shorter rotation ages tend to call for more carbon sequestration and economic profit. However, it is difficult to ensure that 50 years rotation is the optimal rotation age for the target forests, since we do not consider the benefit of biodiversity conservation.


Introduction
The forest carbon pool is the largest terrestrial reservoir, holding more than 3/4 of the entire above-ground terrestrial carbon (Watson et al. 2000). In the past few decades, the world's forests have sequestrated 30% of annual global anthropogenic CO 2 emission (Bellassen and Luyssaert 2014. Forests in Korea store approximately 55 million tons of CO 2, which is about 0.08% of all the carbon stored in Asia's forests (Reich 2013).
The forests have multiple roles to play in climate change, and forest management can help to optimize their roles (Bowyer et al. 2011). Many recent studies show that a policy of active and responsible forest management is effective for enhancing carbon sequestration in forests (Perez-Garcia et al. 2007;Fahey et al. 2010;Lippke et al. 2011) which is accomplished by allocating optimum age distribution in forest land (Fahey et al. 2010). Net forest carbon sink, inferred from changes in storage, sharply decreases with standing age (Bradford and Kastendick 2010). Thus, conserving forests without harvesting on time only works in short-term but call for missing opportunities for greater carbon mitigation in the long term by losing net carbon sequestration capacity (American Hardwood Expert Council 2012).
The objective of this study is to bring the concept of normal forest into forest carbon management in Korea, in order to alleviate the problems that the forests in the nation are now facing. More specifically, we are developing harvest regulation to make Korean forests move toward the desired forests conditions that can guarantee the sustainable revenues and carbon sequestration. To achieve these goals, we would like to introduce the concept of "normal forest". One of the distinguished characteristics of normal forests is that there is an even distribution of age classes in the forest from 0 to rotation age. However, due to the intensive forestation in the 1970s, the current age structure of forests in Korea is highly unbalanced. Without alleviating this problem, the carbon sequestrations in forest will sharply decrease after the major age group class reach the mature stage. Under the basic concept of normal forest theory, this study would like to answer the following research questions: sequestration. The goal of this study suggests 'sustainable forest management' over one rotation. Sustainable forest management will provide integrated benefits such as providing economic advantage to local livelihood and mitigating some of the effects of climate change (LEDS 2016).

Application of the concept of normal Forest to Korea
To seek for forest management strategies, the concept of a normal forest provides a good starting point. The term "normal forest" is a concept of an ideal forest that provides an even flow of timber volume harvested in each period (Amacher et al. 2009), which meet the needs of forest industries supplying wood through ensuring the wood-based revenues to local communities (Amacher et al. 2009). In classical approaches, attaining the normal forest was the direct or indirect goal of forest management. However, the limitation was that it only focused on timber-oriented silviculture. Forest management problems are complex due to complicated interaction between different components of the forest and diversity of values associated with natural resources such as carbon sequestration and biodiversity (Buongiorno and Gilless 2003). We should consider the environmental benefits from forests when we build forest management strategies since the environmental value from forests is much greater than their economic values. One of the greatest environmental benefits of forests is carbon sequestration, and net carbon sequestrations in trees largely depend on their standing age. Although there is variability among species, net carbon storage decline in trees generally occur between ages 20 and 30 for major Korean tree species (see Figure 1). After this age, the sequestration rate decreases gradually until maturity at forest age 80 to 100 (Johnson and Coburn 2010). Thus, conserving forests without harvesting on time only works in shortterm but gives for missing opportunities for greater carbon mitigation in the long term (American Hardwood Expert Council 2012).
To suggest the harvest scheduling that maximizes the discounted profits from the forest considering carbon sequestration in trees, we formulate a linear program (LP). The LP is suitable for establishing forest management planning because it can incorporate multiple decision criteria in the model. Therefore, many agencies such as the Forest Service, U. S. Department of Agricultureoften use the LP package for their timber management planning to calculate the potential yield under forest regulation (Field 1978). Carbon sequestration has become a crucial ecological service that forests provide because of increasing attention to global climate change. The forest owners can gain profits from harvest; however, loss of carbon storage caused by forest cutting activities can generate additional harvesting costs. We hereby introduce a profit maximization formula, which incorporates the benefit and cost of harvest considering cutting-related loss of living trees. The research target of this study is the total area of national forests in Korea. We are only regarding national forests to avoid problems related to complicated ownership, but the method can be expanded to private forests in future researches. The target specie is red pine (Pinus densiflora), which is one of the major forest types in Korea and covers about 40% of the total forest area (Table 1).
3. The profit maximizing harvest scheduling model

Harvest planning period
A profit maximizing harvest scheduling model using linear programing developed by McDill (1999) is useful to derive the optimum harvest schedule; as such, this study applied McDill's approach for maximizing the harvest profit. However, we modified several equations for our objective function. We assume some finite period of forest management schedules including 50 years, 60 years and 70 years planning horizon, and 10 years are one planning period. The general rule to set the length of the planning horizon is that the planning horizon should generally be at least one rotation in length (Buongiorno and Gilless 2003).
To simplify the model, we assume all forest management activities (harvests) are going to take place at the same time during a given planning period. The harvests occur at the midpoint of the period. For example, if the harvest decision is for period 1, all

The data description of the target forest
To develop harvest planning for the target forest, we need forest resource data, economic data and data of current rotation age. Below is a summary of the data.
The target areas: All national forests in South Korea. Forest type: Red pine forests Rotation ages: 50 years, 60 years, 70 years respectively In Korea, the government sets the rotation age for national forests. The rotation age of the red pine forests are 60 years. However, this study assumes three different rotation ages, 50 years, 60 years and 70 years, respectively, to calculate various harvest strategies under different rotation ages. From this assumption, we can compare the changes in harvest strategies and forest CO 2 dynamics based on different scenarios. Rotation ages will influence harvest decisions because it affects the harvest volume per unit period.
Forest data: The initial age-class distribution of the target forests (by area) is shown in the following Table 2. Current data shows that major age classes of the target forest are 30-50 years (70%). This, however, seems to remain a favorable condition for carbon uptake since fast growing young growth would be more effective at capturing carbon. However, the current unbalanced distribution cannot guarantee sustainable carbon sequestration since most forest stands will turn into old growth after a few decades without appropriate treatment. The effects of an aging forest also include: (1) an increasing sensitivity to forest mortality from insects and disease outbreaks; (2) a decreasing timber productivity; (3) an increase in fuel loads, possibly resulting in a higher potential of wild fires (K€ ohl et al. 2017). In combining the production of timber and carbon sink, a more balanced age class distribution is necessary in order to sustain the high carbon sequestration capacity of young trees and high storage capacity of mature trees over time (Garcia-Gonzalo et al. 2007;Routa et al. 2012).
To simplify, we assume that the forests only produce wood products. However, more factors can be included in a future study and this is one of the strengths of linear programming compared to other techniques (McDill 1999) (Table 3).
Economic data: The economic data are used for calculating the profits (cost and return) for each management alternative (McDill 1999) and Table 4 shows the economic data for the target forests. The regeneration and harvest costs include labor and other management costs and the detailed information of the costs is shown in the Table in appendix.
3.3. The objective function for the model 3.3.1. The profit maximizing objective function The objective function maximizes the present value of the net revenue, from N years forest management planning horizon. The mathematical form of the objective function is the following Equation (1) where X ap ¼ the number of areas(ha) cut from initial age-class a (where a¼ age-class, 1,2,3 … M) in period p (where p ¼ 1,2,3,.N and p ¼ 0 means no harvest during the planning horizon). For example, X 31 implies the number of assigned areas (ha) to cut from initial ageclass 3 in period 1. c p ap ¼ Objective function coefficient that is the present value of the net revenue of assigning one area (ha) to the variable X ap : The objective function coefficients imply the discounted net profit per hectare for each variable that could be calculated using the following Equation (2). The expression 10 Ã p À 5 in Equation (2) implies the midpoint of the period P because we assume harvests occur in the midpoint of the period P where p w ¼ the stumpage price v ap ¼ the harvest volume per area for hectare assigned to the variable X ap e ¼ the regeneration cost per ha h ¼ the harvest cost per ha r ¼ interest rate

Loss of carbon from timber removal
Forest resource assessments should include expanded analyses of environmental issues such as CO 2 storage in forests (Richard 1992). The rate of carbon sequestration depends on the tree species, basic density of wood, biomass expansion index and carbon fraction (IPCC 2003). The total carbon losses due to harvests can be calculated by following equation (IPCC 2003) where V ¼ Volume of tree removal (m 3 ); WD ¼ Wood basic density; BEF ¼ Biomass expansion factor; R ¼ Root Ratio (CO 2 in roots); CF ¼ Carbon Fraction: Biomass)Carbon (IPCC ¼0.5); 44/12 ¼ CO 2 Fraction: Carbon(C))CO 2 ; V implies volume of tree removal in forests. The wood basic density is the ratio between the dry weight of wood and the green volume of the same wood which indicates the amount of actual wood substance present in a unit volume of wood (Zobel and Jett 1995). The Biomass Expansion Index (BEF) quantifies carbon stock in forests, which is calculated from the ratio of aboveground biomass and minimum DBH (Sanquetta et al. 2011). R implies ratio of the below-ground biomass to above-ground biomass, which is 0.26 for red pine. R can be set to zero if no changes of below-ground biomass allocation patterns are assumed. Carbon Fraction factor (CF) is used to convert biomass to carbon by multiplying it. The coefficients for red pine are summarized in Table 4.
The Korea Forest Service provides a standard forest carbon storage table, which estimates carbon storage in domestic forests by forest type based on Equation (3). Table 5 shows the yearly carbon storage of the red pine forest by unit area. As seen in the Table, the forest carbon sequestration by unit area is maximized in forest age 30, and then continually declines.
3.3.3. The profit maximizing objective function considering carbon loss by harvest Equation (4) implies the final form of the objective function considering the carbon sequestration value in forests. The objective function was generated by incorporating Equation (1); economic profit from harvesting and Equation (3); carbon loss due to tree removal and additional carbon gain from reforestation.
(4) c c ap is the objective function coefficient which is the present value of the net revenue of assigning one unit area (ha) to the variable X ap considering carbon value. The term ðD Â BEF Â CF Â 44=12Þ Á v ap implies total carbon loss due to timber removal. p c is carbon cost in market. The term D Â BEF Â CF Â 44=12 in the left  Koo et al. (2015).  side of the equation is obtained from Equation (3) but set with R ¼ 0 since we assume no changes of belowground biomass allocation patterns. We assume the p c is 10000 KOW/tc based on the literature review (Lee at al. 2010). Equation (4) also considers additional carbon gain from reforestation. When a harvest occurs, the area is reforested and the amount of new forest area is equal to the harvesting area. Therefore, ½D Â BEF Â CF Â 44=12 Á v 1p Á p c implies additional carbon sequestration from forest age-class one, which is newly generated. The term p c is necessary to convert carbon value into money term and v 1p implies the volume of new trees that are age class 1. Equation (4) simply indicates the economic benefit from cutting trees minus the loss from carbon release due to harvesting in monetary terms. This objective function is used to evaluate the optimal amount of harvest area that satisfies the maximizing profit from timber, considering carbon loss due to tree removal.

Area constraints
The area constraints simply imply that we cannot manage more areas (ha) than we have. The restriction for the model is specified by this set of constraints. The N þ 1 possible prescriptions for each targeted area are: cut in period 1, cut in period 2, cut in period 3, cut in period 4 … cut in period N and do not cut the trees in the areas during the planning periods. Therefore, the sum of the areas allocated from the analysis area to each potential prescription must be no more than the total area that we plan to manage (McDill 1999). The area constraints for the linear programming model follow Equation (5).
where A a ¼ the total number of hectare in initial age a 3.4.2. Harvest fluctuation constraints Minimum Harvest Constraints. We have set the minimum harvest constraints to meet the government's forest plan. The governmental plan projects that the volume to be harvested from national forests will increase to 1,500,000 m 3 by the year 2020 (Korea Forest Service 2008). Thus, the minimum harvest constraints could follow Equation (6). Equation (6) implies that the harvest level in each period should be more than 105,000 m 3 . The number 105,000 is generated by the equation, 1,500,000 Â 7%. Based on data from Table 1, we assume that the red pine forests cover 7% of the total national forests.
Harvest Fluctuation Constraints. The harvest fluctuation constraints are required to widely protect the fluctuating harvest level from one period to the next; furthermore, following constraints will limit harvest level that will be allowed to fluctuate from one period to the next (McDill 1999). We assume that the harvest level does not fluctuate from one period to the next by more than 15%. For example, we want to ensure that the harvest level in period 2 is not less than 15% below or more than 15% above the harvest level in period 1. Likewise, the harvest level in period 3 should not be less than 15% below or more than 15% above the harvest level in period 2. This can be explained in the following Equation (7).
H 2 ! 0:85H 1 ; H 2 1:15H 1 H 3 ! 0:85H 2 ; H 3 1:15H 2 :::::: The first line in the equation implies that the harvest level in period 2 is at least 85% and at most 115% of the harvest level in period 1. In the second line, we can see that the harvest level in period 3 is at least 85% and at most 115% of the harvest level in period 2 and so on. Here, H p represents variables rather than parameters. In order to use harvest fluctuation constraints like Equation (7), we need to introduce specific harvest accounting constraints. A harvest accounting constraint sums up the harvest level for a period and expresses this sum as a variable including H 1 , H 2 , H 3 … and so on. The harvest accounting constraints for this model can be expressed as the following Equation (8) and the constraint requires that the total harvest is greater than or equal to a minimum harvest target for the period p (McDill 1999) The constraints from Equation (8) can be expressed in terms of the variable H p , which is the total harvest volume at period p.

Average ending age constraint
We need to design the linear program to leave a specific age-class distribution at the end of the planning horizon. Our purpose is to achieve a "normal forest" through forest management action at the end of the planning horizon. For this, we have an evenly distributed age-class distribution in mind, but this approach is too restrictive. To enhance the potential of the model to achieve other goals, we introduce a set of the ending age constraints for the target forests. To calculate the average age of a forest, we use the following Equation (9) Age ¼ where Age ¼ the average age of the forest Area i ¼ the area in the ith unit of the forest, and Age i ¼ the age of the ith unit of the forest.
To formulate a constraint for average age of the forest at the end of the planning horizon, the term Area i in Equation (9) can be replaced with the variable X ap , which represents the areas in different blocks of the forest at the end of the planning horizon (McDill 1999). The following equation represents the average age of the forest at the end of the planning horizon.
where Age N ¼ the target minimum average age of the forest in N (end of the planning horizon) TotalArea ¼ the total area of the forest Age ap N ¼ the age in year N of areas in initial age-class a, which are planned to be cut in period p.
If we rearrange Equation (10), we can generate Equation (11), which is the general form of the ending average age constraint for our linear program model (Table 6).
The parameters Age ap N are determined by the following rules. If we consider the age of areas assigned to P ¼ 0 (do-not-cut prescription), it will be 60 years older after the end of the planning horizon under a 60-year planning horizon. Thus, the average age of the area (ha) is 5 years old if p ¼ 0 and a ¼ 1 at the beginning of the planning horizon, and their average age will be 65 years old at the end under the 60-year planning period. Likewise, if p ¼ 0 and a ¼ 2, the average age of the area is 15 years old at the beginning of the planning horizon, and they will be 75 years old at the end of the planning horizon. Given the areas assigned to cut in period 1 (p ¼ 1), they will be 55 years old at the end, regardless of their initial age. Table 7 summarizes the values of the Age ap 60 parameters. The parameter Age 60 (the target minimum average age of the forest in year 60) is related to the ideal forests (normal forests) that we would like to achieve through the LP solution. If the rotation age is 60 years old, the average age of the normal forests will be approximately a 30-year old (half of the rotation) since all age classes are represented in the same quantities in normal forest (Oldeman 2012). Actually,the basic rule of thumb to calculate Age N is (rotation age þ 1)/2 (McDill 1999).

Non-Negative constraints
The non-negative constraints are necessary in the model because the harvested area cannot have a negative value.

The complete linear program model for Profit-Maximization
To combine Equation (1) 6   1  65  55  45  35  25  15  5  2  75  55  45  35  25  15  5  3  85  55  45  35  25  15  5  4  95  55  45  35  25  15  5  5  105  55  45  35  25  15  5  6  115  55  45  35  25  15  5  7  125  55  45  35  25 15 5 Where X ap ¼ Variable, the number of areas to be harvested from initial age-class a (a ¼ 1,2,.6) in planning period p (p ¼ 0,1,2, … 6 and p ¼ 0 implies no harvest during the planning horizon) C ap c ¼ Coefficient, the present value of the profit from one-unit area (ha) to the variable X ap A a ¼ Coefficient, the total number of areas in initial age-class a v ap ¼ Coefficient, the harvest volume for each area assigned to the variable X ap H p ¼ Variable, the minimum harvest target for the planning period p Age N ¼ Coefficient, the average age of the forest in year N TotalArea ¼ Coefficient, the total area of the forest Age ap N ¼ Coefficient, the age in year N of areas in forest

Forest regulation scenarios
The linear programming (LP) harvest scheduling will be affected by various factors such as rotation age and planning horizon. The rotation age is important for managing desired forest structure when seeking optimized commodities and production goals (Bettinger et al. 2010). The rotation age also affects the amount of carbon stock in forests to capture these changes. We assume four different scenarios with single rotation management by rotation age and planning periods: 1. Baseline Scenario 2. Rotation age 50, planning horizon 50 3. Rotation age 60, planning horizon 60 4. Rotation age 70, planning horizon 70 All scenarios assumed single-rotation management and assumed forests are immediately re-established after harvesting. The harvest prescription for the baseline scenario is that forest stands are cut at harvest age (60 years by current law) and immediately re-established by planting, but the harvest prescription would follow the same harvest constraints from the LP model (minimum harvest constraint and harvest fluctuation constraint). Scenarios 2, 3 and 4 are established based on the LP the model, but each scenario reflects rotation age related differences. For example, the rotation age affects the ending age constraints because the basic rule of thumb to calculate the target minimum average age is (rotation age þ 1)/2.

LP solution
The following Tables and Figures indicate the results of the optimal solution from LP. Figure 2 shows the projected age-class distribution at the end of the planning horizon in various scenarios. The first graph in Figure 2 shows the age-class distribution at the beginning of the planning periods and the other graphs show the final age-class distribution following prescriptions from LP solutions. The age-class distribution at the end of the planning horizon is more balanced than the age-class distribution at the initial stage. The final age-class distribution does not follow completely uniform distribution (normal forest) under all scenarios. However, the results are satisfied with the minimum average ending age requirement. Generally, the LP model prescribes retaining additional areas of young growth rather than achieving normal forest that is composed of an equal area of forestland in each age-class. Figure 3 shows the changes in projected age-class distribution at the end of each planning period, under scenario 2, 3, and 4. The age-class distribution of the target forest changes over time, according to harvest prescription, and meets the average ending age requirement at the end of the planning horizon.
Tables 8 summarizes the results for the optimal harvest prescription by different scenarios. Tables 8 would let the forest owners know how many areas (ha) of what age-class will be harvested at a given period. Table 8 informs the harvest schedule in an intuitive way. For example, the forest owner harvests 54140.68ha from age-class 5, 15374.7ha from age-class 6, and 7687.35ha from age-class 7 at the first planning period under the scenario2 . We can interpret scen-ario3 and scenario4 in the same way. Table 8 also tells us that we need to harvest approximately 9332 ha of pine forest each year to meet our management purpose under scenario 2. Under scenario 3, approximately 8013 ha will be harvested each year. We harvest approximately 6868 ha each year under scenario 4 Tables 9 shows costs and revenues from the forest prescription for each period, by scenario. The costs for each period are from the harvest and from replanting the harvested areas, which is a function of the area planted and the harvested volume. Table 10 compares overall net revenues for each scenario. Scenario 2 provides the biggest net revenue, while the baseline scenario provides the least revenue among all scenarios. The shorter rotation age and planning periods generate greater net revenue because the total amount of the harvest each year tends to increase, as the rotation period would be shorter. Scenario 2 provides approximately twice as much revenue as the baseline scenario. Based on the LP solution, when comparing scenarios, scenario 2 generates 2.77% more revenue than scenario 3 and 5.28% more revenue than scenario 4. This is because the shorter the rotation periods, the more trees the LP solutions derive to harvest by the unit period to reach the management purpose (balanced age-class distribution). Also, another reason for generating less revenue in the longer rotation is because the objective functions consider the increasing discount rate according to the flow of time. Figure 4 shows the changes in yearly CO 2 sequestration in forests under different scenarios. "No treatment" means we keep the forest in its natural state, without any intervention. We can estimate the amount of CO 2 sequestration using Equation (3), then divide by the trees' age to get a yearly sequestration rate (Shodor Education Foundation, Inc. 1999). To estimate the yearly CO 2 sequestration rate in trees, we use a standardized yearly forest carbon sequestration table, developed by Korea Forest Service. In terms of carbon sequestration, harvesting and replanting scenarios are much superior to no treatment, and forest management scenarios (S2, S3, S4) from LP show a better carbon sequestration performance than the unsystematic forest plan (Baseline). Under no treatment and baseline scenarios, the net carbon sink in the forest would decrease as time progresses. Otherwise, even if the net carbon sink tends to decrease during the short term, it would rebound in the long term under scenarios from LP solutions. Table 10 shows the summary of the carbon sequestration performance under each scenario. The projected yearly carbon sequestrations in scenario 2 is 95% greater than baseline, and 261% greater than no treatment. The projected yearly carbon sequestrations in scenario 3 are 75% greater than baseline, and 257% greater than no treatment. The projected yearly carbon sequestrations in scenario 4 are 45% greater than baseline, and 215% greater than no treatment. Thus, progress toward a balanced age-class distribution in national forests could enhance forest carbon sequestration by 45% to 95% compared to baseline.

Discussion
We developed an adequate forest regulation with a single cut cycle in Korean national forests according to: (1) economic benefit from timber and (2) changes in net carbon sequestration regarding age class distribution of target forests. All national red pine forests are examined using four different regulation scenarios: (1) baseline, (2) 50 years rotation age and planning horizon, (3) 60 years rotation age and planning horizon and, (4) 70 years rotation age and planning horizon.
The harvest prescriptions that optimized the purpose of management are calculated under four different scenarios. Additionally, changes in yearly carbon sequestration from LP solutions are compared with baseline and no-treatment scenarios. To achieve our forest management goal of sustainable carbon storage and timber production, we introduced the concept of "normal forest." The simple definition of normal forest is a forest with an equal number of areas in each age class (McDill 1999). The normal forest provides sustainability to guarantee an even flow of timber products in perpetuity. Through harvest prescription from LP, the forest age structure at the end of the planning period is more balanced compared to the baseline scenarios. However, the solutions from LP did not achieve normal forests with perfectly even aged distribution. Instead it produced a left-skewed age-class  distribution curve due to the cost management ruling out the achievement of a normal forest as an optimal solution. This means that achieving a normal forest cannot produce the optimum solution to maximize profit. The results from our LP model also confirm that the forest management activities will enhance yearly carbon sequestration in forests for all scenarios compared to baseline and without treatment. The yearly carbon sequestration and economic benefits are maximized under the shortest rotation age (50 years), primarily due to a shorter rotation and planning horizon, where the LP prescription tends to cut more volume of trees per unit period compared to the longer rotation. With the forest management, the forests sequestrate additional 1.8, 1.5 and 0.9 million tons of CO 2 , under scenarios 2, 3, and 4, respectively, compared to baseline. However, it is hard to ensure that 50 years rotationis the optimal rotation age for the target forests, since we do not consider the biodiversity conservation benefit. A previous research showed that incorporating the biodiversity components such as the minimum viable population for birds in the optimization model led to a longer rotation age compared to the carbon rotation age (Nghiem 2014). Another study from Koskela et al. (2007) found that promoting biodiversity preservation prolonged rotation age using the simulation model. Also, a longer rotation could improve the soil condition (Wu et al. 2015) and the resilience of the forest to disturbance, disease and insect outbreaks (Moriarty 2008). Our LP model provides forest managers and policy makers a tool for establishing sustainable forest management plans considering both economic and carbon sustainability; however, the model has several limitations. First, the objective function does not allow for the possibility of interaction between species since we only considered the one forest type, red pine in particular, as the target species. Second, the objective function needs to reflect more realistic forest practices such as thinning. A thinning procedure will generate additional costs and affect trees' growth and volume. New results will be derived in the event that thinning is included in the model. In addition to carbon sequestration of trees, other factors such as biodiversity should be considered for a more sophisticated model. The model assumes that the  carbon price is fixed during the planning horizon, yet, carbon prices could change over time. Thus, volatility of the carbon price should be taken into account during studies. All in all, there is room to improve or develop better models: (1) expanding target forest species, (2) expanding the target area to private forests, and (3) including other values such as biodiversity in objective function.

Disclosure Statement
No potential conflict of interest was reported by the authors.