THE COSTATE VARIABLE IN A STOCHASTIC RENEWABLE RESOURCE MODEL

NATURAL RESOURCE MODELING Volume 19, Number 1, Spring 2006 THE COSTATE VARIABLE IN A STOCHASTIC RENEWABLE RESOURCE MODEL KENNETH S. LYON Economics De...
Author: Tracey Stewart
1 downloads 1 Views 716KB Size
NATURAL RESOURCE MODELING Volume 19, Number 1, Spring 2006

THE COSTATE VARIABLE IN A STOCHASTIC RENEWABLE RESOURCE MODEL KENNETH S. LYON Economics Department Utah State University 3530 Old Main Hill Logan, UT 84322-3530 E-mail: [email protected] SAKET PANDE Center for World Food Studies (SOW-VU) Vrije Universiteit, Amsterdam The Netherlands E-mail: [email protected]

ABSTRACT. In this paper we discuss the costate variable in a stochastic optimal control model of a renewable natural resource, which we call a fishery. The role of the costate variable in deterministic control models has been discussed extensively in the literature. See for example Lyon [1999], Clark [1990, pp. 102 107] and Arrow and Kurz [1970, pp. 35 37]; however, there is little discussion of this variable for stochastic models, even though the costate variable has similar roles in the two models. In both models the costate variable is a shadow value of the associated state variable, and as such has the role of rationing the use of the state variable. In addition, as has been shown in Lyon [1999], in natural resource problems the costate variable can be partitioned into a scarcity effect and a cost effect. We show that this same partitioning can be done in the stochastic renewable resource problem. We discuss and contrast the similarities and differences in these concepts for deterministic and stochastic models. In addition, we present a numerical example to help solidify the results. KEY WORDS: Costate variable, stochastic dynamic programming, deterministic dynamic programming, renewable resource model, simulation.

Introduction. In this paper we discuss the costate variable in a stochastic optimal control model of a renewable natural resource, which we call a fishery. The role of the costate variable in deterministic control models has been discussed extensively in the literature. See for example Lyon [1999], Clark [1990, pp. 102 107] and Arrow and Kurz Received by the editors on April 19, 2004, and in revised form on Nov. 22, 2004. c Copyright 2006 Rocky Mountain Mathematics Consortium

45

46

K.S. LYON AND S. PANDE

[1970, pp. 35 37]; however, there is little discussion of this variable for stochastic models, even though the costate variable has similar roles in the two models. In both models the costate variable is a shadow value of the associated state variable, and as such has the role of rationing the use of the state variable. In addition, as has been shown in Lyon (1999), in natural resource problems the costate variable can be partitioned into a scarcity effect and a cost effect. We show that this same partitioning can be done in the stochastic renewable resource problem. We discuss and contrast the similarities and differences in these concepts for deterministic and stochastic models. In addition, we present a numerical example to help solidify the results. Below we first identify discrete time deterministic and stochastic optimal control models and manipulate them to identify the desired concepts. This includes the identification of the costate variable. We use discrete time rather than continuous time, because the stochastic elements require less overhead in discrete time. We then discuss the costate variable. After this we present a numerical stochastic fishery model to illustrate the concepts. The final section contains the summary and conclusions. The discrete time models. The models are of a renewable resource, which we will refer to as a fishery. In the deterministic model the objective is to maximize the present value of the resource and in the stochastic model it is to maximize the expected value of the resource. In both models we assume that current values of the variables are known with certainty; however, in the stochastic model the law of motion for the fish stock will have a random component. The real interest rate, r, is assumed to be constant to simplify the exposition, and we define β = 1/(1+r), the discount factor. To achieve the objective we maximize the present value of the net surplus stream, s(ht ), and the expected value of the net surplus stream, E[s(ht )], in the two models respectively, where ht is the harvest on the time horizon t = 1, 2, . . . , T . We let xt give the time path, and E[xt ] give the expected time path of the resource stock in the two models respectively on the same time horizon. In both the initial condition for x is x0 = x0 , and in both models the harvest cost function is written c(ht , xt ) with the following characteristics: marginal harvest costs are positive and increase with the rate of harvest, ch (h, x) > 0,

STOCHASTIC RENEWABLE RESOURCE MODEL

47

chh (h, x) ≥ 0. In addition, harvest costs decrease as the resource stock increases, cx (h, x) ≤ 0, which is the usual relationship. This is the result of greater abundance lowering harvest costs. The demand function in inverse form is written as D(h) with the usual negative slope, D (h) < 0. Net surplus can be written  s(ht , xt ) =

0

ht

D(v) dv − c(ht , xt ).

The deterministic model. For the deterministic model the objective functional can be written (1)

W =

T 

β t s(ht , xt ) + β T +1 S(xT +1 )

t=0

which is maximized subject to (2)

(3)

(4)

xt+1 = xt + g(xt ) − ht x 0 = x0

given

ht , xt ≥ 0 for t = 0, 1, . . . , T

with T given,

where g(x) is the growth, recruitment, or reproduction function for the resource stock, fish, S is the terminal value function and T is the end of the current time horizon. The graph of the growth function, g(x), is assumed to have an inverted “U” shape such as that related to the logistic curve and is assumed to be differentiable. Subtracting from this growth the harvest, h, yields the net growth, xt+1 − xt . The terminal value function, S, can be thought of as a bequest function for the fish stock that is left to some distant generation. The end of the current time horizon, T , is assumed to be sufficiently distant that the system will have evolved to the stationary state, the state where xt+1 − xt = 0. We will concentrate upon internal solutions, i.e., we will assume that the nonnegativity constraints are satisfied. We will, therefore, not be studying all of the possible solutions to models of this type, since it is

48

K.S. LYON AND S. PANDE

possible for the solution harvest to be equal to zero over an interval of time. We do this to keep the presentation from becoming too complex. The necessary conditions can be found using Bellman’s principle of optimality or using the Lagrange multiplier method. We present both to show the relationship between the two. Define Vt (xt ) = Max {ht }

T 



β t−t s(ht , xt ) + β T +1 S(xT +1 )

t=t

subject to equations (2) and (4) with xt given. Using the Bellman equation (principle of optimality) we can write (5)

Vt (xt ) = Max {s(ht , xt ) + βVt+1 (xt+1 )} ht

subject to equations (2) and (4) with xt given. Substituting equation (2) into equation (5) yields: (5 )

Vt (xt ) = Max {s(ht , xt ) + βVt+1 (xt + g(xt ) − ht )}. ht

The necessary conditions for an internal solution, i.e., where equation (4) is satisfied, are: (5a)

(5b)

(5c)

∂Vt+1 (x∗t+1 ) ∂s(h∗t , x∗t ) −β = 0, ∂ht ∂xt+1

t = 0, 1, . . . , T − 1,

∂s(h∗T , x∗T ) − βS  (x∗T +1 ) = 0, ∂hT x∗t+1 = x∗t + g(x∗t ) − h∗t ,

x∗0 = x0 ,

t = 0, 1, . . . , T,

where the asterisk indicates optimum value. envelope theorem, we have (5d)

In addition, by the

∂Vt+1 (x∗t+1 ) ∂Vt (x∗t ) ∂s(h∗t , x∗t ) = +β [1 + g  (x∗t )]. ∂xt ∂xt ∂xt+1

STOCHASTIC RENEWABLE RESOURCE MODEL

49

Equations (5a) (5d) will be put into a correspondence with equations from the Lagrange multiplier method, which we now present. In the Lagrange multiplier method (Chow [1997]) we define: (6) L =

T 

β t s(ht , xt )+β t+1 λt+1 (xt+g(xt )−ht−xt+1 )+β T +1 S(xT +1 )

t=0

where the λs are the Lagrangian multipliers. We ignore the nonnegativity constraints and examine an internal solution where the Lagrangian multipliers associated with these constraints are all zero. The necessary conditions for an internal solution of the problem in equation (6) are: ∂s(h∗t , x∗t ) − βλ∗t+1 = 0, ∂ht

(6a)

(6b)

(6c) (6d)

x∗t+1 = x∗t + g(x∗t ) − h∗t , λ∗t =

t = 0, 1, . . . , T

x∗0 = x0 ,

∂s(h∗t , x∗t ) + βλ∗t+1 [1 + g  (x∗t )], ∂xt

t = 0, 1, . . . , T

t = 0, 1, . . . , T

λ∗T +1 = S  (x∗T +1 ).

A comparison of equations (5) and (6) yields the conclusion that λ∗t+1 =

∂Vt+1 (x∗t+1 ) ∂xt+1

because they are necessary conditions for the same maximization problem. The λs are costate variables and are as usual the shadow values of the resource stock. We proceed now with the identification of the cost and scarcity effects in the costate variable. To achieve this we use λ∗0 to simplify the writing of the subscripts; however, we could use any t from 0 to T − 1. We use equation (6c) recursively as follows: ∂s(h∗0 , x∗0 ) + βλ∗1 [1 + g  (x∗0 )], ∂x0   ∂s(h∗1 , x∗1 ) ∂s(h∗0 , x∗0 ) λ∗0 = +β + βλ∗2 [1 + g  (x∗1 )] [1 + g  (x∗0 )], ∂x0 ∂x1 λ∗0 =

50

λ∗0

K.S. LYON AND S. PANDE   ∂s(h∗ , x∗ ) ∂s(h∗1 , x∗1 ) ∂s(h∗0 , x∗0 ) 2 2 = +β +β ∂x0 ∂x1 ∂x2   ∗  ∗  ∗ + βλ3 [1 + g (x2 )] [1 + g (x1 )] [1 + g  (x∗0 )].

We now shorthand ∂s(h∗t , x∗t )/∂xt = ∂s∗t /∂xt and 1 + g  (x∗t ) = 1 + gt∗ , and continue ∂s∗0 ∂s∗ + β 1 (1 + g0∗ ) ∂x0 ∂x1 ∗ ∂s + β 2 2 (1 + g0∗ )(1 + g1∗ ) ∂x2 ∂s∗ + β 3 3 (1 + g0∗ )(1 + g1∗ )(1 + g2∗ ) + · · · ∂x3 T +1 +β (1 + g0∗ )(1 + g1∗ ) · · · (1 + gT∗+1 ) λ∗T +1   t−1 T ∂s∗0  t ∂s∗t  ∗ ∗ λ0 = + β (1 + gk ) ∂x0 t=1 ∂xt k=0   T T +1 ∗ (1 + gk ) λ∗T +1 . +β λ∗0 =

(7)

k=0

Note that ∂s(h∗t , x∗t )/∂xt = −cx (h∗t , x∗t ) where the sub-x indicates t−1 partial derivative. Also note that k=0 (1+gk∗ ) can be viewed as giving the time path of an additional unit of x originating at time zero, t = 0, which we call ξt,0 . Note that g  (x) has been called the own biological rate of interest by Quirk and Smith [1970]. It is the discrete annual growth rate or physical return to one unit of the resource. For the more general case where the new unit originates in time τ we have

ξt,τ =

t−1 

(1 + gk∗ ).

k=τ

Using these definitions we can rewrite equation (7) as

(8)

λ∗0 = −

T  t=0

β t ξt,0 cx (h∗t , x∗t ) + β T +1 ξT +1,0 λ∗T +1 .

STOCHASTIC RENEWABLE RESOURCE MODEL

51

For the more general case we have (8 )

λ∗τ = −

T  t=τ

β t ξt,τ cx (h∗t , x∗t ) + β T +1 ξT +1,τ λ∗T +1 .

In equation (8 ), β T +1 ξT +1,τ λ∗T +1

(9a)

is the scarcity effect given that cx ≡ 0, and (9b)



T  t=τ

β t ξt,τ cx (h∗t , x∗t )

is the cost effect given that T is sufficiently large that β T +1 is nil. These results parallel those given in Lyon [1999] for the continuous time case, and these effects have exactly the same interpretation as in the continuous time case. When cx (ht , xt ) = 0 for all x(cx ≡ 0) it is easy to see that the cost effect disappears; thus, the value of the additional unit of the resource stock at time τ is the present value of ξT +1,τ times the value of one unit at T +1. The additional unit of resource at τ grows or shrinks to ξT +1,τ units in time T + 1. When the cost effect is nonnull, we can select T to be arbitrarily large so that equation (9a) is nil. Since cx is negative, the cost effect is positive and identifies the cost saving associated with the additional resource stock. The corresponding terms will be identified below for the stochastic costate variable. The stochastic model. We now develop the parallel results for the stochastic model. We consider the class of “admissible” control laws (policies) that consist of a finite sequence of functions π = {μ0 , μ1 , . . . , μT }, where μ : St → Ct and such that μt (xt ) ∈ Ut (xt ) for all xt ∈ St . Here Ct , St are spaces of elements ht and xt respectively. The control ht is constrained to take values from a given nonempty subset Ut (xt ) of Ct , which depends on the current value of the state variable. It is assumed that for some admissible policy vector, the optimal expected value could be achieved. The problem of maximizing the expected value of the net surplus stream then has the objective functional

 T E β t s(μt (xt ), xt ) + β T +1 S(xT +1 ) , (10) E0 [W ] = εt :t=0...T

t=0

52

K.S. LYON AND S. PANDE

which is to be maximized subject to xt+1 = xt + g(xt ) − μt (xt ) + εt

(11)

x 0 = x0

(12)

given

μt (xt ) ∈ Ut (xt ) ⊆ Ct , xt ∈ St ≥ 0 for t = 0, 1, . . . , T with T given,

(13)

where εt ∈ Dt (Dt is the space of elements εt ) is a random variable that is independently and identically distributed through time with mean zero and density function Pt (ε). We assume that at time t there is certainty about all variables except εt , which becomes known by time t + 1. The term g(xt ) + εt is the stochastic growth function over the time period t to t + 1. The term g(xt ) gives the mean growth since the mean value of εt is zero. However, the above formulation is equivalent to (Bertsekas [1987, p. 14]) (10 )

E0 [W ] =

E

εt :t=0...T

 T

β t s(ht , xt ) + β T +1 S(xT +1 ) ,

t=0

which is to be maximized subject to (12) and (11 )

xt+1 = xt + g(xt ) − ht + εt

(13 )

ht ∈ Ut (xt ) ⊆ Ct , xt ∈ St ≥ 0 for t = 0, 1, . . . , T

with T given. Since εt ∼ Pt (ε), it follows from (11) εt ∼ Pt (xt+1 − xt − g(xt ) + μt (xt )) ≡ Pt (xt+1 | xt ). If the expectation operator is defined as in (10), then   E [·] = ·Pt (ε) dε = ·P (xt+1 | xt ) dxt+1 = Et [· | xt ]. εt

STOCHASTIC RENEWABLE RESOURCE MODEL

53

As in the deterministic problem, the necessary conditions can be found using the Bellman equation or using the Lagrange multiplier method. As above we present both to show the relationship between the two methods. Define for any t with 0 ≤ t ≤ T vt (xt ) = Max Et

 T

{ht ∈Ut (xt )} t=t

β

t−t

s(ht , xt ) + β

T +1

S(xT +1 ) xt

subject to equations (11) and (13) with xt given. Using the Bellman equation for t ∈ [0, T ] we can write (14)

vt (xt ) = Max Et [{s(ht , xt ) + βvt+1 (xt+1 )} | xt ] ht

subject to equations (11) and (13) with xt given. We define vT +1 (xT +1 ) ≡ S(xT +1 ). Substituting equation (11) into equation (14) yields: vt (xt ) = Max Et [{s(ht , xt ) + βvt+1 (xt + g(xt ) − ht + εt )} | xt ] ht

or (14 ) vt (xt ) = Max {s(ht , xt ) + βEt [vt+1 (xt + g(xt ) − ht + εt ) | xt ]}. ht

The necessary conditions for an internal solution are:

∂vt+1 (xt+1 ) ∂s(h∗t , xt ) (14a) − βEt xt = 0, ∂ht ∂xt+1 (14b)

(14c)

t = 0, 1, . . . , T − 1

∂s(h∗T , xT ) − βET [S  (xT +1 ) xT ] = 0 ∂hT xt+1 = xt + g(xt ) − h∗t ,

x0 = x0 ,

t = 0, 1, . . . , T.

Applying the envelope theorem to equation (14 ) yields (14d)

∂vt+1 (xt+1 ) ∂vt (xt ) ∂s(h∗t , xt ) = + βEt xt [1 + g  (xt )]. ∂xt ∂xt ∂xt+1

54

K.S. LYON AND S. PANDE

We now identify the costate variables by comparing these conclusions to those from the Lagrange multiplier method. We define: (15) L = E0

 T

t

β s(ht , xt ) + β

t+1

 Λt+1 (xt + g(xt ) − ht + εt − xt+1 )

t=0



T +1

S(xT +1 ) | x0

where the Λ’s are Lagrangian multipliers. The necessary conditions for an internal solution of the problem in equation (15) are: (15a)

(15b)

(15c)

(15d)

∂s(h∗t , xt ) − βEt [Λ∗t+1 (xt+1 ) | xt ] = 0, ∂ht xt+1 = xt + g(xt ) − h∗t , Λ∗t (xt ) =

x0 = x0 ,

t = 0, 1, . . . , T

t = 0, 1, . . . , T

∂s(h∗t , xt ) + βEt [Λ∗t+1 (xt+1 ) | xt ] [1 + g  (xt )], ∂xt t = 0, 1, . . . , T,

Λ∗T +1 (xt+1 ) = ET [S  (xT +1 ) | xT ].

Equations (14a) (14d) and (15a) (15d) are for the same maximization problem; hence, we can conclude from equations (14d) and (15c) that ∂vt (xt ) . Λ∗t (xt ) = ∂xt The Λ’s are costate variables and are shadow values of the resource stock. The stochastic nature of the resource stock, however, does add some new information. We proceed now with the identification of the cost and scarcity effects in the costate variable. To achieve this we use Λ∗0 (x0 ) to simplify the writing of the subscripts; however, we could use any t from 0 to T − 1. We use equation (15c) recursively as follows: Λ∗0 (x0 ) =

∂s(h∗0 , x0 ) + βE0 [Λ∗1 (x1 ) | x0 ] [1 + g  (x0 )] ∂x0

STOCHASTIC RENEWABLE RESOURCE MODEL

Λ∗0 (x0 ) =

Λ∗0 (x0 ) =

55

∂s(h∗0 , x0 ) ∂x0 

 ∂s(h∗1 , x1 ) + βE0 + βE1 [Λ∗2 (x2 ) | x1 ] [1 + g  (x1 )] x0 ∂x1 × [1 + g  (x0 )]

∂s(h∗0 , x0 ) ∂x0



 ∂s(h∗1 , x1 ) ∂s(h∗2 , x2 ) + βE0 + βE1 + βE2 [Λ∗3 (x3 ) | x2 ] ∂x1 ∂x2  × [1 + g  (x2 )] x1  × [1 + g  (x1 )] x0 [1 + g  (x0 )]

(16) Λ∗0 (x0 ) =

∂s(h∗1 , x1 ) ∂s(h∗0 , x0 ) + β[1 + g  (x0 )]E0 x0 ∂x0 ∂x1



∂s(h∗2 , x2 ) + β 2 [1 + g  (x0 )]E0 [1 + g  (x1 )]E1 x1 x0 ∂x2



+ β 3 [1 + g  (x0 )]E0 [1 + g  (x1 )]E1 [1 + g  (x2 )]

∂s(h∗3 , x3 ) × E2 x2 x1 x0 + · · · . ∂x3

The growth of one unit of additional stock originating at the beginning of time period t over the period would be [1 + g  (xt )]. However, if the system has only reached the beginning of time period t − 1 and not t, one’s expected growth of one unit of additional stock injected at the beginning of time period t − 1 would be (by the end of time period t) [1+g  (xt−1 )]Et−1 [[1+g  (xt )] | xt−1 ]. Thus, by backward induction, the expected growth (by the end of time period t) of one unit of additional stock introduced at the beginning of time period zero would be E0 [ξt,0 | x0 ]     = [1+g  (x0 )]E0 [1+g  (x1 )]E1 · · · [1+g  (xt−1 )]Et−1 [1+g  (xt )] | xt−1   · · · | x1 | x0 t ∈ [1, T ].

56

K.S. LYON AND S. PANDE

Thus terms involving (1 + g  (·)) give the expected time path of an additional unit of the resource stock originating at time zero, E0 [ξt,0 | x0 ]. We have ∂s(h∗t , xt )/∂xt = −cx (h∗t , xt ) where the sub-x indicates partial derivative. Similar to the treatment given above, one additional unit of the stock injected at the beginning of time period t would yield a return of [1 + g  (xt )]Et [−cx (h∗t+1 , xt+1 ) | xt ] from expected cost savings in the next time period. By iteration, we therefore define

(17a)

   Et,0 = [1 + g  (x0 )]E0 [1 + g  (x1 )]E1 · · · Et−1 [1 + g  (xt )]    × Et [cx (h∗t+1 , xt+1 ) | xt ] | xt−1 · · · | x1 | x0 t ∈ [1, T − 1]

and

(17b)

   ET,0 = [1 + g  (x0 )]E0 [1 + g  (x1 )]E1 · · · ET −1 [1 + g  (xt )]    × ET [Λ∗T +1 (xT +1 ) | xT ] | xT −1 · · · | x1 | x0 .

Et,0 defines the expected return on a unit of additional stock injected at the beginning of time period 0 due to cost saving in time period t + 1. For the more general case where the new unit originates in time τ , we further define   Et,τ = [1 + g  (xτ )]Eτ [1 + g  (xτ +1 )]Eτ +1 · · · [1 + g  (xt )]   × Eτ [cx (h∗t+1 , xt+1 ) | xt+1 ] · · · | xτ +1 | xτ t ∈ [τ + 1, T − 1] and   ET,τ = [1 + g  (xτ )]Eτ [1 + g  (xτ +1 )]Eτ +1 · · · [1 + g  (xT )]   × ET [Λ∗T +1 (xT +1 ) | xT ] · · · | xτ +1 | xτ . By the law of iterated expectation, see Appendix A, equations (17a),

STOCHASTIC RENEWABLE RESOURCE MODEL

(17b) could be rewritten as

(18a)

 Et,0 = [1 + g  (x0 )]E{x1 ,x2 ,... ,xt+1 } [1 + g  (x1 )] [1 + g  (x2 )]  · · · [1 + g  (xt )]cx (h∗t+1 , xt+1 ) | x0

 t = E{x1 ,x2 ,... ,xt+1 } (1 + gk )cx (h∗t+1 , xt+1 ) x0 k=0

= E{x1 ,x2 ,... ,xt+1 } [ξt+1,0 cx (h∗t+1 , xt+1 ) | x0 ] t ∈ [1, T − 1] and (18b)

 ET,0 = [1 + g  (x0 )]E{x1 ,x2 ,... ,xT +1 } [1 + g  (x1 )] [1 + g  (x2 )]  · · · [1 + g  (xT )]Λ∗T +1 | x0 = E{x1 ,x2 ,... ,xT +1 } [ξT +1,0 Λ∗T +1 (xT +1 ) | x0 ].

For the more general case Et,τ = E{xτ +1 ,... ,xt+1 }

 t

(1 + gk )cx (h∗t+1 , xt+1 ) xτ

k=τ

= E{xτ +1 ,... ,xt+1 } [ξt+1,τ cx (h∗t+1 , xt+1 ) | xτ ] t ∈ [τ + 1, T − 1], and ET,τ = E{xτ +1 ,... ,xT +1 } [ξT +1,τ Λ∗T +1 (xT +1 ) | xτ ]. Using these definitions we can rewrite equation (16) as

(19)

Λ∗0 (x0 ) = −

T 

β t E{x1 ,... ,xt+1 } [ξt,0 cx (h∗t , xt ) | x0 ]

t=0

+ β T +1 E{x1 ,... ,xT +1 } [ξT +1,0 Λ∗T +1 (xT +1 ) | x0 ]. For the more general case we have

(20)

Λ∗τ (xτ ) =

T  t=τ

β t E{xτ +1 ,... ,xt+1 } [ξt,0 cx (h∗t , xt ) | xτ ]

+ E{xτ +1 ,... ,xT +1 } [β T −τ +1 ξT +1,0 Λ∗T +1 (xT +1 ) | xτ ].

57

58

K.S. LYON AND S. PANDE

In equation (20) (21a)

E{xτ +1 ,... ,xT +1 } [β T −τ ξT +1,0 Λ∗T +1 (xT +1 ) | xτ ]

is the expected value of the scarcity effect given that cx ≡ 0, and (21b)

T  t=τ

β t E{xτ +1 ,... ,xt+1 } [ξt,0 cx (h∗t , xt ) | xτ ]

is the expected value of the cost effect. As in the deterministic model the scarcity effect exists only when the cost effect is null. The difference is that the term is expected value. The time path of the additional unit is given by a Markov process and, since xT +1 is a random variable, so is the terminal costate variable. In the cost effect, the expected value of the summation is likewise considered. These have similar explanations. To emphasize the differences between the deterministic and stochastic cost effects we present them for time 0 in an expanded form:

(9b )

− cx (h∗0 , x∗0 ) − βcx (h∗1 , x∗1 )(1 + g  (x∗0 ))

− β 2 cx (h∗2 , x∗2 )(1 + g  (x∗0 ))(1 + g  (x∗1 )) − β 3 cx (h∗3 , x∗3 )(1 + g  (x∗0 ))(1 + g  (x∗1 ))(1 + g  (x∗2 )) + · · ·

(21b ) − cx (h∗0 , x0 ) − E{x1 } [βcx (h∗1 , x1 )[1 + g  (x0 )] | x0 ] − E{x1 ,x2 } [β 2 cx (h∗2 , x2 )[1+g  (x0 )][1+g  (x1 )] | x0 ] − E{x1 ,x2 ,x3 } [β 3 cx (h∗3 , x3 )[1+g  (x0 )][1+g  (x1 )][1+g  (x2 )] | x0 ] + · · · (21b )

−cx (h∗0 , x∗0 ) − βE1,0 − β 2 E2,0 − β 3 E3,0 + · · · .

This form emphasizes the fact that starting with t = 0 the resource stock is a stochastic variable, and that all terms that include xt for t > 0 are given as expected values, by law of iterated expectations. We next present a fishery numerical example to illustrate the stochastic nature of these concepts.

STOCHASTIC RENEWABLE RESOURCE MODEL

59

Fishery numerical example. We now present a numerical stochastic fishery example and its solution. We use logistic population growth with a stochastic element and maximize expected net surplus. This is a standard pedagogical fishery model with stochastic growth. This stochastic dynamic programming, SDP, problem does not have a closed form solution, which is typical for this type of problem. We solve this problem point-wise for a numerical policy function using value function iterations on a grid. We divide the state space, fish stock, into a 6000 point even spaced grid and numerically solve the Bellman equation contraction at each of these points. The result of this is a discrete numerical policy function, xt+1 = μ(xt ). The procedure we use is developed in detail in Bertsekas [1987, Chapter 5, especially pp. 188, 189]. This solution technique has the advantage that, as the step size in the grid approaches zero, the solution on the grid approaches the true solution. The problem is to find a time independent numerical decision rule (policy function), μ(xt ), that solves equation (10) subject to equations (11) (13), where D(h) = a0 − b0 h,

a0 , b0 > 0

c(h, x) = a1 h + 0.5b1 h2 − a2 x + 0.5b2 x2 ,

a1 , b1 , a2 , b2 > 0

so (22)

s(h, x) = (a0 − a1 )h − 0.5(b0 + b1 )h2 + a2 x − 0.5b2 x2 ,

and (23)

g(x) = γx − ηx2 ,

γ, η > 0.

The details of how we executed our value function iterations are discussed in Appendix A, and the resulting policy function is displayed in Figure 1. In Figure 1 the vertical axis is expected net investment. This is growth minus the harvest, g(xt ) − ht . We use this policy function and the normal density function with an expected value of ε of 0, E[ε] = 0, and a variance of ε of 1/5, E[ε2 ] = 1/5, to generate the simulation results. We simulate from x0 using (24)

xt+1 = μ(xt ) + εt

60

K.S. LYON AND S. PANDE

Expected Net Investment

0.15

0.1

0.05

0

14

16

18

20 22 Current Fish Stock

24

26

FIGURE 1. Policy function for the fish stock.

(25)

ht = g(xt ) + xt − μ(xt ) = γxt − ηx2t + xt − μ(xt )

(26)

E[Λt+1 | xt ] = β −1 (a0 − a1 − (b0 + b1 )ht )

Figure 2 shows the time path of the fish stock for 100 simulations and the time path with εt = 0, which is given by the white line down the middle of the simulation paths. The white line can also be thought of as the time path of the expected value of the stock. We start at a stock level that is less than the long run optimal level; hence, the fish stock grows through time. This is achieved, of course, through optimal harvests that allow the growth. The simulation time paths are clustered around this expected time path. Figure 3 shows the harvest time paths for the simulations and the path with εt = 0. Figure 4 shows the costate variable time path for the same computer runs. The three figures show a lot of similarity, which follows from the fact that the optimal harvest and the optimal costate variable are both well-behaved functions of the stock, equations (25) and (26). This diagram illustrates that the costate variable is stochastic. An alternative way to depict this variability is to view the variance of the

STOCHASTIC RENEWABLE RESOURCE MODEL

30

25

Fish Stock

20

15

10

5

0

0

50

100

150

200 Time

250

300

350

400

FIGURE 2. Optimal time paths of the fish stock.

2

1.8

1.6

1.4

Harvest

1.2

1

0.8

0.6

0.4

0.2

0

0

50

100

150

200 Time

250

300

350

FIGURE 3. Optimal time paths of the harvest.

400

61

62

K.S. LYON AND S. PANDE

14

12

Costate Variable

10

8

6

4

2

0

0

50

100

150

200 Time

250

300

350

400

FIGURE 4. Optimal time paths of the costate variable.

simulated values at each t. This concept is developed in Appendix A and is presented in Figure 5. Note that during the first few years the variance increases and then stabilizes around a long run value. This same feature is also noted in Figure 4 where the band of simulated values gradually increases and then stabilizes. Summary. In this paper we have discussed the costate variable in a stochastic optimal control model of a renewable natural resource, which we call a fishery. The role of the costate variable in deterministic control models has been discussed extensively in the literature; however, there is little discussion of this variable for stochastic models, even though the costate variable has similar roles in the two models. In both models the costate variable is a shadow value of the associated state variable, and as such has the role of rationing the use of the state variable. In addition, as is shown in Lyon [1999], in natural resource problems the costate variable can be partitioned into a scarcity effect and a cost effect. We have shown that this same partitioning can be done in the stochastic renewable resource problem. In the deterministic problem the costate variable is the value of an additional unit of the resource stock at time τ . It is the present value of

STOCHASTIC RENEWABLE RESOURCE MODEL

63

1

Variance

0.8

0.6

0.4

0.2

0

0

50

100

150

200 Time

250

300

350

400

FIGURE 5. Time path of the variance of the costate variable.

the return stream of an additional unit of the resource. The additional unit at time τ either grows or shrinks through time and will have a stream of cost savings associated with it if costs decrease as the resource stock increases. In the stochastic model the same characterization holds, except in an expected sense. The additional unit of the resource stock at time τ will have a stochastic time path; hence, the costate variable will be an expected value, and the cost effect will be the present value of a stochastic stream of cost savings. The numerical example illustrates the stochastic nature of the variables in the model, and the figures show that the optimal paths change from simulation to simulation. In this example the costate variable has only a cost effect, and, as can be seen, the time path of an additional unit or the marginal unit will be stochastic. The cost effect is, therefore, the expected value of a stochastic stream of cost savings. Appendix Policy function. In this Appendix we identify the procedure we use to generate the point-wise policy function, μ(xt ), for the numerical example. The procedure we use is developed in detail in Bertsekas

64

K.S. LYON AND S. PANDE

[1987, Chapter 5, especially pp. 188, 189] or Christiano [1990]. This solution technique has the advantage that, as the step size in the grid approaches zero, the solution on the grid approaches the true solution. We treat T is large without limit, giving us an infinite horizon problem, and write the Bellman equation as (A1)

v(xt ) = Max Et [{s(μ(xt ), xt ) + βv(xt+1 )} | xt ]. μ∈Ut

The state space is discretized into n points and denoted as S = {x1 , x2 , . . . , xn } = {1, 2, . . . , n}. For μ an admissible control, we define the transition probability, pij (μ) as (A2)

pij (μ) = P (xt+1 = j | xt = i, μt = μ),

i, j ∈ S, μ ∈ U.

Thus, pij (μ) is the probability that state xj will occur given that the current state is xi , and the admissible control μ is applied. The Bellman contraction mapping, T , for the grid point xi can be written (A3)

n  T (v)(i) = max s(μ, i) + β pij (μ)v(j) . j=1

Iterations on equation (A3) are performed until norm (Xk+1 − Xk ) < 0.001, where the vector Xk = [x1 , x2 , . . . , xn ] is the solution at iteration k. The result is the numerical policy function, μ(xt ), reported in the body of the paper. The parameter values for equations (22) and (23) are: a0 = 35, b0 = 10, a1 = 10, b1 = 2.5, a2 = −0.1, b2 = 0.003125, γ = 0.1, η = 0.00125 and β = 0.95123. Variances. We generate 100 simulations using the above policy function and report these in the body of the paper. From a cross section of these at each t we can calculate the variance of the generated series. For example, if λst,i is the calculated value of the costate variable

STOCHASTIC RENEWABLE RESOURCE MODEL

65

for time t in simulation i, then the variance of these simulated values for time t is given by the usual formula (A4)

var (λst ) =

100 

(λst,i − λst )2 /100

i=1

where the superscript s is for simulated value and the super-bar is for mean value. Law of iterated expectations (Ramanathan [1993, p. 88, Theorem 5.5]). From the law of motion of the stock, xi+2 is conditionally independent of {xi , xi−1 , . . . , xτ } given xi+1 , by the Markov property, (A5)

P (xi+2 | xi+1 ) = P (xi+2 | xi+1 , xi , . . . , xτ ).

By the product law of probability, P (xi+2 , xi+1 , . . . , xτ +1 | xτ ) = P (xi+2 | xi+1 , xi , . . . , xτ )P (xi+1 | xi , xi−1 , . . . , xτ ) · · · P (xτ +1 | xτ ). From (A5) we therefore get, (A6) P (xi+2 , xi+1 , . . . , xτ +1 | xτ ) = P (xi+2 | xi+1 )P (xi+1 | xi ) · · · P (xτ +1 | xτ ). The right-hand side of equation (17) could generally be represented as,   (A7) · · · f (xτ +1 , xτ +2 , . . . , xt+1 )P (xt+1 | xt )P (xt | xt−1 ) · · · P (xτ +1 | xτ ) dxt+1 dxt · · · dxτ +1 . Using (A6) we transform (A7) to   · · · f (xτ +1 , xτ +2 , . . . , xt+1 ) × P (xt+1 , xt , xt−1 , . . . , xτ +1 | xτ ) dxt+1 dxt · · · dxτ +1 = E{xτ +1 ,xt+1 } [f (xτ +1 , xτ +2 , . . . , xt+1 ) | xτ ].

66

K.S. LYON AND S. PANDE

The right-hand side of the above equation is what is used from equation (18) onwards. REFERENCES K.J. Arrow and M Kurz [1970], Public Investment, the Interest Rate, and Optimal Fiscal Policy, Johns Hopkins Press, Baltimore, MD. D.P. Bertsekas [1987], Dynamic Programming: Deterministic and Stochastic Models, Prentice Hall, New York. G. Chow [1997], Dynamic economics: Optimization by the Lagrange Method, Oxford Univ. Press, Oxford. L.J. Christano [1990], Linear-Quadratic Approximation and Value-Function Iteration: A Comparison, J. Bus. Econ. Statist. 8, 99 113. A.N. Kolmogorov and S.V. Fomin [1970], Introductory Real Analysis (Silverman, translator and ed.), R.A. Dover Publ. Inc., New York. K.S. Lyon [1999], The Costate Variable in Natural Resource Optimal Control Problems, Natur. Resource Modeling 12, 413 426. J.P. Quirk and V.L. Smith [1970], Dynamic Economic Models of Fishing, in Economics of Fisheries Management: A Symposium (A.D. Scott, ed.), Institute of Animal Resources Ecology, University of British Colombia, Vancouver. R. Ramanathan [1993], Statistical Methods in Econometrics, Academic Press, New York.

Suggest Documents