Adaptive economic and ecological forest management under risk

Forest managers must deal with inherently stochastic ecological and economic processes. The future growth of trees is uncertain, and so is their value. The randomness of low-impact, high frequency or rare catastrophic shocks in forest growth has significant implications in shaping the mix of tree species and the forest landscape. In addition, the fluctuations of wood prices influence greatly forest revenues. Markov decision process models (MDPs) offer a rigorous and practical way of developing optimum management strategies, given these multiple sources of risk. Examples illustrate how such management guidelines are obtained with MDPs for combined ecological and economic objectives, including diversity of tree species and size, landscape diversity, old growth preservation, and carbon sequestration. The findings illustrate the power of the MDP approach to deal with risk in forest resource management. They recognize that the future is best viewed in terms of probabilities. Given these probabilities, MDPs tie optimum adaptive actions strictly to the state of the forest and timber prices at decision time. The methods are theoretically rigorous, numerically efficient, and practical for field implementation.


Background
Forest planning and decision making models are often deterministic. They assume that the future is known, or they reduce uncertain variables to their expected values, with the hope that the results will work, at least "on average". However, the future is highly uncertain a , and it can only be described in terms of probabilities, be they "objective probabilities" derived from factual data, or "subjective probabilities" reflecting a personal belief (de Finetti 1937).
There is a large literature concerning the role of risk in forest ecosystems (see e.g. Kant and Alavalapati 2014, p. 307-369). But there are few applications of Markov decision process models (MDPs) in forestry, a paucity also observed in the general operations research literature (White 1993). Nevertheless, MDPs are powerful and general methods that deserve more attention in forest management. This was recognized early on in Lembersky and Johnson (1975), and Lembersky (1976) who dealt with applications to the management of forest plantations. Their pioneering work was followed by several applications to the uneven-aged management of various forest types for financial objectives (Kaya and Buongiorno 1987), ecological objectives (Zhou and Buongiorno 2006), and combined financial and ecological objectives Rollin et al. 2005;Zhou et al. 2008a, b).
The purpose of this paper is to describe how stochastic simulation and MDPs can be used to deal with decision making under risk in forestry. The data and the specific applications described below deal with mixed loblolly pine (Pinus taeda L.)-hardwood stands in the southern United States. However, the methods are general and can be applied to a variety of forest ecosystems, ranging from evenaged monospecific plantations to mixed-species unevenaged forests, and from single stands to large forest areas. The objective throughout is to develop simple adaptive management guidelines that can be readily applied in the field. The paper is organized as follows. The next section presents the methods: stochastic forest growth models and their reduction to Markov chains which are then used for prediction and optimization. This is followed by the results which show first some of the consequences of (wrongly) predicting the evolution of forest ecosystems without taking into account the effects of random ecological, climatic, or economic shocks. Next are shown predictions of the long-term evolution of forest ecosystems with and without management and the implications for economic and ecological criteria. And, MDP results optimizing decisions adaptively for discounted or undiscounted economic and ecological objectives, and combinations thereof. The discussion and conclusion deal with limitations of the methods and potential improvements.

Methods
The guiding principle of the methods is to represent the evolution of a forest ecosystem subject to various impact variables, such as prices and catastrophic disturbances, with Markov chains: probability matrices describing the frequencies of transition of each variable between discrete states. After the transition probabilities between the system states have been determined, the Markov model is used to predict future system states subject to specific management policies. Furthermore, decision process models (MDPs), optimization techniques based on Markov chains, are used to determine the best management policy for a particular objective. A policy is a set of rules that prescribe a decision for each observed system state at decision time. Various adaptive strategies are obtained in this way depending on the objective function which may deal with discounted or undiscounted criteria and constraints that may limit the decision domain.

Markov forest model
The different possible states of a forest stand are described in a practical manner by the level of basal area in trees of different sizes and species. Basal area per hectare, the sum of the cross sectional area of the trees, typically measured at breast height, is one of the simplest measures of forest stand density and composition. To be practical, the number of tree size and species categories should be kept low b . For example Table 1 shows the definition of forest stand states used here, with the basal area in three tree size classes (small, medium, large), and two species categories (softwoods, hardwoods). With two levels of basal area, high or low, indicating below average or above average basal area in each species and size category, there are 2 6 = 64 possible stand states. The example stand state #13 in Table 1, denoted (001, 100) has high basal area in large pine trees and small hardwood trees, and low basal area in the other categories (Zhou and Buongiorno 2006).
In MDPs, the evolution of stands over time is described by T = [p(s'|s)], a matrix of yearly transition probabilities between stand states s and s'. These probabilities cannot be established directly by direct observation due to lack of sufficient data in each stand. Simulation methods are used instead. A typical stochastic simulation model of stand growth suitable for this purpose has the form: where y t = [y ijt ] is the vector of the number of trees per unit area of species i and size j at time t. G t is a matrix of parameters describing the probability that a tree of a particular species stays alive and grows into a higher size class from t to t + 1. The vector i t refers to recruitment, the number of trees of a particular species that enters the smallest diameter class from t to t + 1. Both G t and i t vary depending on the stand condition, y t , through the stand density that affects tree growth and recruitment. ε t is a vector of random disturbances. In the following applications the parameters of model (1) are estimated from observations on permanent sample plots in the mixed loblolly pine (Pinus taeda L.)hardwood forests of the Southern United States (Schulte et al. 1998) c . The differences between the deterministic predictions of model (1), and the observations on the plots give observations on the random shocks ε t due to ice storms, wind, insect outbreaks, abnormal weather, etc.… that have affected forest growth during the observation period (Zhou and Buongiorno 2004).
To get the transition probabilities between stand states p(s'|s) model (1) is used to predict the future stand state, s', of a random initial stand in state s, by bootstrapping a random shock ε t from the set of the observed shocks. This is repeated a sufficient number of times to obtain stable estimates of the transition probabilities p(s'|s) (Zhou 2005, p. 51).
Due to the short time typically covered by plot observations, the transition probability matrix T = [p(s'|s)] obtained in this manner reflects mostly small scale, high frequency disturbances. In addition, catastrophic, low The example state #13 has high basal area in large pine trees and small hardwoods, and low basal area in other categories. Source: Zhou and Buongiorno (2006).
frequency events, such as hurricanes must be added where needed, based on weather and geographic data. Table 2 shows the non-zero elements of the 64 × 64 transition probability matrix computed in this way and used in the results of this paper. The entry in row and column 000,000 means that if a stand has now low basal area in all species and size categories, the probability is 0.797 that it will be in the same state next year. The entry in row 000,000 and column 000,001 means that the same stand has a 0.027 probability of growing into a stand with high basal area in the largest size class. Table 2 also indicates that regardless of its current state, there is a 0.025 probability that the stand will be next year in state 000,000 due to a catastrophic event, a hurricane in this example (Zhou and Buongiorno 2006). The resulting T matrix is then used to predict future states: where Π t = [π st ] is the vector of probabilities of stand state s at time t, and the steady-state probabilities, Π = Π t = Π t + 1 are obtained by solving the system of equations: where Π = [π s ] is the vector of long-term steady state probabilities for each stand state s, and I is the identity matrix. As equation (3) indicates, the long-term probabilities are independent of the initial state. Each probability is also interpreted as the fraction of a large forest landscape that is in each state s in the long run when only natural disturbances occur. As indicated below, an equation equivalent to (3) can also be used to predict the steady state probabilities when a particular management regime is followed indefinitely.

Markov price model
The financial returns from forest management depend critically on the price of wood. Perfect prediction of future prices is not possible, but probable price levels and changes can be drawn from the behavior of past prices. In parallel with the procedure for forest stands, the method describes market fluctuations by transition probabilities between price levels. For example, based on past data prices may be classified as low, medium, or high. And the transition probabilities between price levels are obtained from a price model such as d : where Q t is the price in year t and α, β are parameters estimated by regression analysis of past data. The matrix of transition probabilities between price levels is M = [p(m'|m)] where p(m'|m) is the probability that next year's price level is m' given the current price level m.
These probabilities are obtained by simulations and repeated sampling from the distribution of the residual v t . An example of price transition probabilities used later in this paper is in Table 3 (Zhou 2005). Given a currently low price level (less than $84 · m −3 ), the probability that it will be low following year is 0.82. The probability that it will be medium (between $84 · m −3 and $94 · m −3 ) is 0.10, and the probability that it will be high (above $94 · m −3 ) is 0.08. Table 3 illustrates the high price autocorrelation (price stickiness): there is a high probability that the price stays at the same level from year to year.

System states and transition probabilities
Given transition probability matrices between stand states T = [p(s'|s)] and price levels M = [p(m'|m)], the transition probabilities of the forest-price system are: where i is the current system state (stand state s and price level m), and j is next year's system state (stand state s' and price level m'). For example, with the 64 stand states in Table 2 and the three price levels in Table 3, there are 192 system states and 192 × 192 transition probabilities.

Decisions and immediate rewards
Decisions consist in moving instantly from one stand state to another by harvesting some of the trees. It is assumed that market prices are exogenous, i.e. the decisions have no effect on prices. A decision produces a stand state with less basal area, less carbon sequestered in the stand, and a different composition of tree species and size. It also generates revenues that depend on the volume of timber harvested and on the price at decision time e . A policy consists of a set of decisions applied systematically to each stand-price state. If the policy is applied to an entire forest, decisions alter the distribution of stand states and thus the landscape diversity, and possibly the part of the forest that stays in an old-growth (late-seral) state. Ecological criteria include a variety of indices concerning the state of the stand and of the landscape induced by a policy. Although several diversity indices are possible (Magurran 1988), Shannon's index (Shannon 1948) was used here throughout. For example, applied to tree species diversity Shannon's index is: where H is the tree-species diversity of a stand state, f a is the fraction of trees of species a, and n is the number of species. The fraction f a is preferably based on the basal area of trees rather than their number, to give more importance to the large trees. With this index, species  diversity is between 0 and ln(n), and it is at the maximum when f a is the same for all species. The tree size diversity of a stand state is calculated in a similar way.
The forest landscape diversity is also expressed similarly with f s , the fraction of the forest landscape in state s. Another useful criterion is the fraction of the forested landscape in "old growth" state. For some forest types, old growth or late-seral stands have been defined in previous studies (Hummel and Calkin 2005). Here, old-growth is defined as the most frequent states that develop in the steady state predicted with equation (3) in the absence of harvest (Zhou and Buongiorno 2006). With the data in Table 2, the five old-growth states are 001,011; 001,111; 011,111; 101,111; and 111,111, which all have high basal area in the largest softwoods and hardwood trees.
The data on the expected timber volume and ecological criteria by stand state are obtained during the simulations with model (1) that give the transition probabilities. The immediate ecological reward of a decision is the characteristic of the state induced by the decision, for example, its tree diversity. The immediate financial reward is equal to the change in volume obtained by harvesting the stand from the current state to another, multiplied by the price level at decision time. The carbon sequestered was estimated from the amount of growing stock left after harvest, assuming 1.24 t · m −3 of CO 2 e for pulpwood and 1.57 t · m −3 for sawtimber (AFC 2014).
In sum, each immediate reward, monetary or ecological, is summarized by a vector V = [v ik ] where v ik is the immediate reward when the stand-price system is in state i and the decision is k (which changes the stand state or leaves it intact).

Consequences of management policies
A policy, c, is described by a particular decision matrix: where p c (k|i) is the probability of harvesting a stand to state k, given the stand-price state i under policy c. If the harvest policy is deterministic (the same action is taken Source: Zhou and Buongiorno (2006). 1 Basal area in pulpwood, small sawtimber, and large saw timber of pines (first three digits), or hardwoods (last three digits) 1 = higher than current average, 0 = lower than current average. all the time in each state) the probability of a particular decision is 0 or 1. To investigate the consequences of continuing the current management, independently of the price level, D c is based on the frequency at which harvests have changed the state of the sample plots during the observation period.
The transition probabilities between stand-price states given a particular policy are then obtained from: where p c (j|i) is the probability of ending in stand-price state j after one period from stand-price state i under policy c, and p(j|k) is the probability of ending in standprice state j after one period from stand state k. Substituting E c for T in equation (3) gives the probability vector where π c i is the steady state probability of standprice state i under policy c.

Predicting expected undiscounted rewards
The expected undiscounted reward, R, such as the expected tree species diversity, or the expected annual timber revenue, that results in the long run by following a specific policy is: where V 0 ¼ r c i Â Ã is a row vector of the immediate reward in stand-price state i under policy c. For financial returns r c i is the timber income when the stand-price state is i and the policy is c. The expected long-run diversity of tree species and size, the expected basal area and carbon sequestered, and the expected annual production under a particular policy are obtained in similar fashion. The expected long-run forest landscape diversity is based on equation (6) where p a is replaced by the steady state probability of each stand state, π c s , under the policy c. The expected interval between harvests (or cutting cycle), in years, is the inverse of the yearly probability of a harvest. Let h i c = 1 if the decision called for by policy c is a harvest when the stand-price state is i. Then, the annual probability of a harvest is: And the expected cutting cycle is: Predicting expected discounted rewards Some criteria are typically discounted, future outcomes receiving less weight than current ones. This is especially the case for timber revenues, but ecological criteria may also be discounted (Boscolo et al. 1997;Howarth 2009).
For a particular management policy, c, the present value of a criterion is obtained by solving the following system of equations (Hillier and Lieberman 2005, p. 785): where V c i is the present value of the criterion of interest over an infinite horizon for an initial stand-price state i, when the harvest policy c is followed, with an immediate reward r c i and a yearly interest rate ρ. In the results presented below, the interest rate was set at 3.08% per year, the real yield on AAA corporate bonds from 2001 to 2013, using the consumer price index as deflator (U.S. Government 2014). This choice implies a long-term, low-risk investment by a conservative investor.
Then, the average expected present value per ha over an entire forest with a particular initial distribution of stands and probability of price level is: where π 0 i is the initial probability of stand-price state i.

Optimizing management policies
In parallel with the predictions of the effects of specific management policies, the optimization of policies differs depending on whether the rewards are discounted or undiscounted.

Optimizing expected discounted rewards
The policy that maximizes an expected discounted reward over an infinite horizon, such as the net present value of the timber income, without any other constraint, is found by solving the following linear programming program, originally due to d'Epenoux (1963).
where y ik is the total discounted probability of being in state i and making decision k over an infinite horizon f . β j is the probability that the system starts in stand-price state j. This probability must be strictly positive, i.e. β j > 0 and X j β j ¼ 1. After finding the best solution, y Ã ik , the best decisions are obtained from: where D ik is the probability of making decision k when the stand-price state is i. The optimal policy is deterministic (D ik = 0 or 1), and it is independent of the initial condition, i.e. of the β j 's (Hillier and Lieberman 2005, p. 921), although the value of the objective function, the maximum discounted value of the reward, depends very much on the initial condition g .

Optimizing expected undiscounted rewards
Another linear program gives the policy that maximizes an unconstrained undiscounted reward (Manne 1960). Examples of undiscounted criteria are the expected species diversity observed yearly in steady state, or the expected annual wood harvest.
where z ik is the steady-state probability of stand-price state i and decision k.
Given the best solution, z Ã ik , the corresponding best decisions are: The best decisions obtained with equation (17) are still deterministic (D' ik = 0 or 1) and independent of the initial stand-price state. And, the maximum value of the objective function, such as the maximum expected undiscounted species diversity over an infinite time horizon, is also independent of the initial condition.

Multiple objectives
Policies that best meet multiple objectives simultaneously are obtained by modifying the objective function of models (14) and (16), or/and by adding constraints. In the results shown below the models were kept linear by expressing both the objective function and the constraints in undiscounted or discounted terms.
For example, in a discounted framework, a policy that maximizes the net present value of timber income, while keeping a specific fraction of the forest area in old growth state, is obtained with a model consisting of equations (14)-(15) and the additional constraint: where oik = 1 if state k is an old-growth state, oik = 0 otherwise. O* is the maximum, unconstrained, discounted value of the fraction of old growth, obtained with model (16) to (17), and α is the desired fraction of this maximum that must be maintained by the policy.
As an example of an undiscounted framework, the maximum expected annual harvest that keeps the sequestered carbon dioxide equivalent (CO 2 e) at a specific level is obtained by solving model (16)-(17) with r ik equal to the amount of timber harvest resulting from decision k in state i, and the additional constraint: where c ik is the amount of CO 2 e stored in the stand state that results from decision k in state i, and C* is the desired expected value of the stored CO 2 e.

Results
Effects of disturbances on predicted stand growth Figure 1 illustrates the importance of recognizing random shocks in predicting the evolution of forest ecosystems. The data result from long-term simulations of a forest stand of mixed hardwoods-loblolly pines in the south of the United States, growing without human intervention (Zhou and Buongiorno 2004). The figure shows long-term predictions of total stand basal area, and basal area of pine trees. The continuous line refers to predictions obtained with the deterministic part only of model (1), ignoring the random shocks, ε t . The jagged shows instead the predictions obtained with the random shocks reflecting the high-frequency, low impact disturbances derived from the plot data, without catastrophic events. While the total basal area is qualitatively similar for both models, the predictions of the basal area of pine trees are very different. Ignoring the random shocks, pines totally disappear in about 300 years. Instead, when the random disturbances are taken into account the basal area of the pine trees is never less than 20 m 2 · ha −1 , about a third of the total basal area. Forest scientists agree that stochastic disturbances "can be a major determinant of forest structure" (Oliver and Larson 1990). In the case of loblolly pine, recurrent fires are typical in the south of the United States. Lobllolly pines are more fire resistant than hardwoods. Fires also cause openings that favor pine regeneration. The simulations with stochastic shocks are thus more likely to predict correctly the long-term evolution of loblolly pine-hardwood stands than the deterministic version.

Effects of catastrophes on forest landscape
In addition to the high-frequency, low-impact disturbances, the forests are also subject to low-frequency catastrophes, such as hurricanes in the region investigated here. The long-term effects of catastrophic events on the forest landscape were predicted with the Markov model in Table 2, with and without the 0.025 probability of a catastrophe appearing in the first column. The fraction of the landscape that would be present in the long run (steady state) in each of the 64 possible states was predicted with equation (3). Figure 2a shows the prediction without catastrophic events, while Figure 2b shows the prediction with catastrophes. Both figures show the fraction of the landscape area occupied in the steady state by each of the five oldgrowth states, and the aggregate fraction occupied by the 59 other states. Without catastrophic disturbances, 81% of the landscape was in an old growth state, and 38% consisted of stands in state #16 (001,111) with high basal area in hardwoods in all three size categories, and high basal area in large pine trees. However, when catastrophic events were taken into account, only 29% of the landscape was in an old-growth state and the dominant old-growth state #16 occupied only 12% of the landscape area. In sum, the long-run expected diversity of the forest landscape measured with Shannon's index (6) was 80% higher due to catastrophic disturbances than it would have been without them.

Economic and ecological consequences of current management
The economic consequences of the current management, were expressed by the net present value (NPV) of timber income over an infinite horizon. The NPV was derived Figure 1 Simulated total basal area (a) and basal area of pine trees (b) in an uneven-aged loblolly pine stand with a deterministic model (smooth line) or a stochastic model with low-impact high-frequency shocks but no catastrophic event (jagged line). Source: Zhou and Buongiorno (2004).
with equation (12) with the initial distribution of stands reflected in the plot data and a low initial price level h . Other consequences, such as the annual wood production and the diversity indices were undiscounted expected values in steady state, obtained with equation (9), while equation (11) gave the expected length of the cutting cycle. The results in Table 4 compare the consequences of the current management with those arising from forest growth alone without management.
Continuing the current management indefinitely generated an NPV of $1,339 · ha −1 in timber revenues and an average expected production of 1.9 m 3 · ha −1 · y −1 over the entire forest. Due to this harvest the expected average basal area was 2 m 2 · ha −1 lower with the current management than it would be without human intervention.
Using the diversity indices without management as a basis of comparison, the current management increased tree species diversity slightly (1%). Like catastrophic events, Figure 2 Expected distribution of forest area in steady state, without catastrophic disturbances (a) and with catastrophic disturbances (b). management also increased landscape diversity, though to a lesser extent (3%). However, by harvesting large trees it decreased the expected tree size diversity by 2%.
The fraction of old-growth in the forest maintained by the current management was 11% lower than it would be in a natural forest, and the value of the average carbon sequestration was 42 t · ha −1 lower. All these outcomes stemmed from harvests in individual stands that occurred at average intervals of 7 years under the current management.

Optimizing financial objectives
The policy that maximized the expected NPV was obtained by solving model (14)-(15), with a uniform distribution of initial stand-price states (i.e. β j = 1/192 ∀ j). This best policy was then used to obtain with equations (12)-(13) the expected NPV given the stand states distribution in the region of interest, as reflected by the plot data, and the current price level ("low" in 2014 according to the definitions in Table 3).
The results in Table 5 show how the best decisions adapted to the stand state and price level observed at decision time. For example, for a stand in state #64 (111,111), at low price, the best decision was to cut the stand to state #52 (110,011) by reducing the basal area of the large pine trees and of the small hardwood trees. If instead the current price was average, the best decision was to cut the stand to state #3 (000,010) which meant reducing further the basal area of the small and medium pine trees. Last, if the price was high, the best decision was to reach stand state #1 (000,000) by reducing the basal area to a low level in all species and size categories. In contrast, under current management, the plot data revealed that a stand in this same state #64 (111,111) was left intact 85% of the time and cut to state #44 (101,011) only 15% of the time.
The expected consequence of this best financial policy for timber production and the various ecological indicators is in Table 6. The expected NPV and annual timber production were more than four times higher than with the current policy i . However, the expected basal area, and the diversity of tree species and size were all lower. The landscape diversity was less than a third of its level with the current policy, and no old growth was left in the entire landscape. The amount of CO 2 e sequestered was less than half, while the frequency of harvests was about double that with the current policy.

Optimizing ecological objectives
Two criteria were used as applications of policies that best meet ecological objectives: maximizing tree species diversity, or carbon sequestration. The results in Table 7 show that at most, tree species diversity could be 7% higher than with the current management. Maximizing tree species diversity would not hurt financial returns, as in fact the expected NPV was $713 · ha −1 higher than with the current management as it led to substantially higher harvest of sawtimber and pulpwood. This more intensive harvest, however, reduced the fraction of the forest in old-growth state by 5 percentage points. Accordingly, the landscape diversity was reduced by approximately 30%, and the amount of CO 2 e sequestered was 52 t · ha −1 less. The tree size diversity index was also reduced slightly. Figure 3a shows the expected fraction of the total forest area in four stand states, starting in the current condition and evolving over a century with the policy that maximized the expected tree species diversity. The states in Figure 3a dominated the steady-state landscape when state #1 (000,000), with low basal area in all species and size categories occupied 20% of the forest landscape, state #7 (000,110) with high basal area in small and medium size hardwoods occupied 14.5%, state #8 (000,111) with in addition high basal area in large hardwoods occupied 14%, and state #16 (001,111) with in addition high basal area in large pine trees occupied 13%. With the tree species maximization policy, state #1 quickly over reached its steady-state level, while the other states approached it more slowly. Figure 3b shows how the tree species diversity index of the entire forest evolved from its initial level to the steady state with the policy that maximized expected tree species diversity. The diversity index is relative to the steady-state diversity of the unmanaged forest. With this policy it took approximately 30 years for the species diversity index to reach 90% of its maximum steadystate value of 108%.

Maximizing carbon sequestration
The amount of CO 2 e sequestered in the trees living biomass in steady state could be increased by 42 t · ha −1   compared to the current management (last column of Table 7), but in doing so, wood production would be very low (0.2 m 3 · ha −1 · y −1 ). In fact, the level of CO 2 e sequestered in this way would be almost the same as that obtained by letting the forest grow naturally without any harvest (Table 4, last column). Indeed, the hands-off option maybe superior as it kept 29% of the landscape in old-growth state against only 15% when maximizing CO 2 e sequestration, at the small opportunity cost of $119 · ha −1 in foregone NPV from timber revenues j . The other criteria were similar to those obtained by natural forest growth. Figure 4 shows the results of a series of optimizations in which the objective was to maximize the NPV of timber revenues given its initial stands distribution and current low price, subject to increasing constraints on the fraction of old growth maintained on the forest k . Figure 4 documents the opportunity cost, in terms of the NPV of foregone revenues, of higher levels of old growth. While approximately $5,300 · ha −1 of NPV could be obtained without any old growth, this was reduced to about $1500 · ha −1 when 22% of the landscape was kept in old growth. The marginal cost of old growth increased at increasing levels of old growth, from a low of $368 · ha −1 from 0 to 5% of old growth to a high of $1145 · ha −1 for 21% to 22% of old growth.

Constrained financial and ecological objectives
As an example of application with multiple ecological objectives, Figure 5 shows the results of maximizing annual wood production (sawtimber and pulpwood), subject to constraints on the amount of CO 2 e sequestered in the growing stock. Both objectives were undiscounted and the constraint was as in equation (19). As indicated by Figure 5, there was no conflict between wood production and carbon sequestration up to approximately 100 t · ha −1  of CO 2 e. Beyond that point, wood harvest decreased rapidly at an increasing rate. Increasing carbon storage from 100 to 150 t · ha −1 decreased wood production by 1.2 m 3 · ha −1 · y −1 , and increasing it from 200 to 250 t · ha −1 decreased production by 3.3 m 3 · ha −1 · y −1 . As shown above, the maximum carbon storage achievable was about 245 t · ha −1 , and at that point there was hardly any wood harvest.

Discussion and conclusions
In a 1999 interview, G.B. Dantzig, the inventor of linear programming, remarked that "all problems that are solved under deterministic means have that fundamental weakness-they don't properly take uncertainty into account" (Dantzig 1999). Forestry planning is no different. Risk, uncertainty, and stochastic behavior are critical parts of how forest ecosystems work. As illustrated above, simulations of mixed loblolly pine-hardwood stands reveal how ignoring random disturbances due to biological or catastrophic shocks in forest growth models can lead to wrong ecological predictions, such as the disappearance of pines in this context. To these natural sources of risk must be added the high risk that stems from price fluctuations which complicates decision making with financial objectives. The methods suggested in this paper to handle risk in forest decision making follow the modeling approach outlined by Holling et al. (1986). The first step is to "bring the world to the laboratory" with possibly complex and non-linear stochastic models, and then simplify the models to allow for efficient optimization with Markov decision process models. MDPs are attractive for their simplicity while keeping the essence of planning problems under risk l . They recognize that the future is unknown and that it may be described only in terms of probabilities. And, they lead to best adaptive management policies whereby decisions depend entirely on the systems state, i.e. the state of knowledge of the decision maker at decision time.
The fact that with MDP models the future state of a forest ecosystem depends only on its current state should not be viewed as a shortcoming, but recognition of fact. Predictions can rely only on current knowledge. This current knowledge may include past behavior of the system (such as past tree growth) if necessary. Thus, Markov models are not necessarily "memory less" and apparent failure of Markov models to correctly represent forest growth (Binkley 1980;Roberts and Hruska 1986;Johnson et al. 1991) may be due to incomplete description of the current ecosystem state rather than to a shortcoming of the Markov model. Similarly, Markov models of price changes are very general, embracing random walk, rational expectations, autoregressive, and "any stochastic model in which the price is conditional on previous prices" (Taylor 1984). In sum, Markov chains allow the simplification of complex, multi-dimensional stochastic processes and make their optimization easier or at all possible (Holling et al. 1986;Insley and Rollins 2005).
There exists considerable potential and flexibility to enrich the models presented in this paper as they are relatively small compared with the size of the problems that can be solved with current linear programming software. Still, to be successful in forest ecosystem management the MDPs should adhere to the general principles of parsimony and simplicity. To this end, the number of system states must  be kept as small as possible, with a few state variables (tree species and size categories, prices) and a few levels of these variables.
In the same spirit of simplicity, multiple objectives such as optimizing financial returns subject to ecological constraints, of vice-versa, have been treated here with linear models. This requires expressing both the objective function and the constraints in either discounted or undiscounted criteria. While discounting financial returns from forestry has long been a standard procedure (Faustmann 1849), it is less so for ecological criteria. Nevertheless, it is plausible to give more weight to the present ecological characteristics of a forest ecosystem than to their future values. Indeed, it seems preferable to maintain desirable current states (such as old-growth stands) rather than lose them, even if they could be restored later (albeit after a long delay). It has also been argued on theoretical grounds that "the future benefits provided by a generic public good-environmental quality-should be discounted at a rate that is close to the market rate of return for risk-free financial assets" (Howarth 2009), a principle that has been applied for example in discounting future carbon sequestered in forests (Boscolo et al. 1997) The methods presented here can be expanded in many directions, in particular to deal with the new carbon markets and the highly stochastic prices for CO 2 e sequestered in forests. The treatment is parallel to the stochastic prices of timber. This would increase substantially the size of the problem, but still keep it well within the capabilities of current linear programming software.
Other issues may require deeper modifications, to deal with non-stationary processes. For example, "climate change can affect forests by altering the frequency, intensity, duration, and timing of fire, drought, introduced species, insect and pathogens outbreaks, hurricanes, windstorms, ice storms, or landslides" (Dale et al. 2001). This implies a change in the transition probabilities over time. Although non-stationary problems can in principle be converted to stationary ones by a reformulation (Bertsekas 1995, p. 167), efficient numerical methods are still elusive (Ghate and Smith 2013). Nevertheless, the consequences of changes in transition probabilities on the steady-state criteria can be readily explored with the methods described in this paper. In sum, given the conceptual generality and the well-developed theory of MDPs, coupled with the powerful solution techniques available, the MDP approach is well suited to deal with risk in forest ecosystem management, and to develop practical adaptive management policies with both economic and ecological objectives.
Endnotes a In this paper the terms "risk" or "uncertainty" are synonymous. Although Knight (1921) defines uncertainty as the absence of probabilities, it may be argued that probabilities, objective (data based) or subjective (opinion based) always exist to some degree. b While modern linear programming software places almost no limit (computationally) on the number of states, it is best for all practical purposes (identifying states of nature, doing the computation, and laying out the recommendations) to use the minimum number of states needed for sufficient accuracy and realism of applications. c For other applications, any stochastic model of forest growth can be used to compute the transition probabilities used in the MDP approach. d If the price series is not stationary, the price data are first de-trended, and the trend is added to the discount rate in present value calculations. e In the particular application described here, prices vary by tree size and species, and all the prices change in parallel with the probabilistic changes of the price index described in Table 3. f That is: y ik = z 0ik + dz 1ik + d 2 z 2ik + d 3 z 3ik + … where z tik is the probability of state i and decision k at time t and d = 1/(1 + ρ) is the discount factor (Hillier and Lieberman 2005, p. 921). g As applied to forestry, the MDP described by equations (10) and (11) is a generalization of (Faustmann 1849) formula recognizing that future stand states and prices are known only as probability distributions. The classical, deterministic Faustmann formula is a special case in which the transition probabilities are 0 or 1 (Buongiorno 2001). h According to Timber Mart South (2014), the average price of softwood sawtimber stumpage in the South of the United States in the second quarter of 2014 was $25 per short ton or approximately $19 · m −3 , placing it in the low range according to the definitions in Table 3.
i Previous results also show that the optimum adaptive policy derived here is superior in terms of NPV to an optimum fixed policy that converts stands to a chosen state at fixed intervals (Zhou et al. 2008b). j The reduction in old-growth fraction from 29% to 15% was due to the best decision calling for a harvest when the stand was in the old-growth state #32 (011,111). Although the harvest occurred only every 86 years on average, this was equivalent to a low frequency natural catastrophe reducing substantially the fraction of the landscape in oldgrowth state. k In Figure 4, the old growth fraction is the right-hand side of equation (18) with O* replaced by its annual constant perpetual equivalent, ρO */(1 + ρ).
l See Buongiorno and Gilless (2003, p. 337-371) for an introduction to Markov and MDP models in forestry. m Alternatively, models (10)-(11) and (12-13) can be extended to deal with discounted objective functions and undiscounted constraints, or the reverse, by introducing non-linear constraints (e.g. Rollin et al. 2005;Zhou 2005), but at the cost of the attendant numerical difficulties.