Aggregating Harvest Activities in Long Term Forest Planning by Minimizing Harvest Area Perimeters

The study presents a new approach for aggregating stands for harvest in strategic forest planning. In fragmented landscapes this could benefit nature conservation as well as reduce costs. The approach is built on the idea of minimizing the outside perimeter of contiguous harvest areas. The formulation allows for the use of exact solution methods such as mixed integer programming. The method was tested in a landscape consisting of 2821 stands. The application showed that large and compact harvest areas were created with limited sacrifice of financial value. The mixed integer programs were in most cases solved within a couple of hours. The method needs to be tested on different landscapes with different degrees of fragmentation. It is also necessary to evaluate the long term consequences of the large clear cuts that appear to be a consequence of this problem formulation.


Introduction
Forest planning encompasses three features that, taken together, make the planning exercise a formidable challenge to the forest manager.First, there are diverse and sometimes conflicting goals.Second, forest planning characteristically covers very long time frames, in order to reflect the nature of forestry problems.Third, forest planning is complicated by the complexity of forest systems.Because of many different interacting processes, large areas, incomplete data, etc. there is typically a high degree of uncertainty when predicting the outcome of different economic and ecological variables.To cope with this forest planning is usually divided into a hierarchical structure with strategical, tactical and operational levels (Davis et al. 2000, Weintraub andCholaky 1991).At the highest level of planning, the focus is on long term goals such as targets for harvest volumes and regeneration strategies.These results are then transferred to planning levels operating on shorter time horizons.The location of harvest sites has essentially been the prerogative of the tactical level of the planning hierarchy (Boston and Bettinger 2001, Cea and Jofré 2000, Church et al. 2000, Covington et al. 1988, Davis and Martell 1993, Gustafsson et al. 2000, Jamnick and Walters 1993, Kirby et al. 1986, Nelson et al. 1991, Richards and Gunn 2000, Sessions and Sessions 1992).However, not considering the spatial layout of activities at the strategic level could lead to misleading results or plans that are impossible to implement (Cea and Jofré 2000, Clements et al. 1990, Davis et al. 2000, Daust and Nelson 1993, Nelson et al. 1991, Weintraub and Davis 1996).One reason for considering harvest site location, such as clustering of harvest sites, is the potential savings due to a reduction of road building, road maintenance and harvest machine transportation activities (Baskent andJordan 1991, Gustafsson 2000).The clustering of harvests could also be motivated by nature conservation since the practise of dispersing the cuttings is a major contributor to the reduction of forest interior (Barret et al. 1998, Franklin and Forman 1987, Wallin et al. 1994).If, instead, the harvests are clustered in time and space more interior habitat would be sustained on the landscape as a whole and the fragmentation of old forest would decrease (Gustafson 1996, Gustafson 1998).Thus the problem is here the opposite of the well known and thoroughly investigated problem with adjacency or green-up constraints (see e.g.Baskent and Keles 2005 for a review).
Unfortunately, when aggregation of harvests is included as a consideration in the strategic planning problem it will grow in complexity.One reason for this is that to represent aggregation in the models integer variables must be used since they allow one to specify if a unit is managed with a certain treatment.Another reason is that regulating aggregation would require that each stand is characterized not only by its state but of the state of the neighbouring stands as well.Further, it is not always obvious what measure should represent aggregation.Is it enough that harvest areas are within a certain distance from each other, or should some other requirement be satisfied?Because of this complexity, a number of different approaches have been suggested in the literature for expressing aggregation of harvests in the optimization.One approach is to include road accessibility and the costs associated with road constructing and maintenance in the long range planning model (e.g.Andersson and Eriksson 2007, McNaughton et al. 1998, Weintraub and Navon 1976).However, these models tend to be very large.Also, even though such models make harvests spatially coordinated they do not necessarily lead to the clustering of harvests motivated by nature conservation.Another approach builds on direct modelling of the aggregation of harvesting activities in pixels where aggregation is governed by a net present value criterion that rewards cost reductions due to harvest of adjacent pixels.Holmgen and Thuresson (1997) solved the problem of forming treatment units in an area with about 3000 pixels using a greedy heuristic.Lu and Eriksson (2000) used genetic algorithms for delineating harvest units on an area consisting of 10 000 pixels.Lind (2000) used simulated annealing for aggregating pixels for harvest in an area consisting of 6200 pixels.Yet another approach is to use criteria that are a direct measure of the degree of aggregation.Heinonen et al. (2007) suggested that the harvest could be clustered by maximizing the proportion of the boundaries between two adjacent pixels that were both cut during the same period and by minimizing the proportion of cut-uncut boundaries.In a case study they solved the resulting problem with threshold accepting for an area consisting of 4612 hexagons.Another approach building on the direct clustering of harvest activities was in 2003 suggested by Öhman and Lämås.They proposed a special criterion building on the harvested volume that can be included in the model formulation for temporal as well as spatial clustering of harvest activities.In a case study they solved the clustering problem with simulated annealing for a landscape consisting of 2600 stands.
In summary, it seems as there are relatively few studies dealing with the aggregation of harvest in forest planning.Those that do consider the phenomenon either based on the concept of road construction or use the concept of small pixels instead of stands or the problems are solved with a heuristic method.The drawback of including road networks into strategic planning is that the models tend to be large and require data of existing or potential road objects.The stands of a forest are normally represented by polygons; converting the polygons to a raster model with pixels would normally mean a major increase of the number of spatial elements in the problem.The application of metaheuristics has its own complications.There is a matter of choice of method, testing the parameter values of the algorithm, and, often, the need to rerun the problem to assess the quality of the solution (Glover and Kochenberger 2003).An alternative to heuristic techniques is to use an exact solution technique such as mixed integer programming (MIP) with the branch and bound algorithm.One limitation of the branch and bound has been the time required to locate the solution to a complex problem and, connected to that, the limited size of the problems that can be solved.However, recent developments in optimization software systems and computer hardware have increased the scope for solving large scale problems in a reasonable time (Atamtürk andSvavelsbergh 2005, Johnson et al. 2000).
The objective of this study is to study the applicability of an approach for aggregating stands for harvest in strategic forest planning.The approach is built on the idea of minimizing the outside perimeter of harvest areas.If two adjacent stands are harvested in the same period the total perimeter of the harvest areas would decrease compared to if two not adjacent stands with the same area and shape would be harvested.An advantage of using the perimeter as a criterion is that it makes it possible to formulate the problem as an MIP problem.The perimeter criteria has been used in earlier studies for mitigating the fragmentation of old forest (Fischer and Church 2003, Tòth and McDill 2008, Öhman and Wikström 2008).It would, thus, be of value to see what consequences on landscape patterns this concept has when applied to the aggregation of harvest activities.In contrast to the study by Fischer and Church (2003), whose mathematical formulation of the perimeter criteria is adopted, the problem dealt with here is dynamic, not a single period.The practical limitations of the approach for large scale long range planning are, thus, not known.The approach is tested in a case study in a landscape consisting of 2821 stands.

Model Formulation
The first part of the model is connected to traditional timber production whereas the second part will introduce the aggregation criterion.The first objective is to maximize the net present value (NPV) from forest management, subject to traditional long range forest planning constraints.The constraints include requirements on even harvest flow, minimum amounts of old forest, and a minimum ending inventory.This is accomplished by assigning a treatment schedule j, j ∈{1,…,J i }, to each stand i, i ∈{1,…,I}.A treatment schedule is a sequence of treatments for a stand from period 1 to period P, the end of the planning horizon.Consequently, the model formulation for the first part of the problem is an example of the standard model I formulation building on the treatment schedule concept (Johnson and Scheurman 1977).These treatments include regeneration, thinning and final felling.Associated with each treatment schedule j for each stand i is its NPV, D ij, the volume harvested in period p, H ijp , the area of old forest in period p, G ijp , and the amount of standing volume in period P, R ijP .D ij ,is here the sum over period 1 to P of discounted revenues and costs that include timber revenues and costs for harvesting and regeneration but do not include timber transportation or entry costs (the transportation of harvesters, road maintenance, etc.).Given that the decision variable x ij = 1 signifies that treatment schedule j of stand i is assigned to the stand, and otherwise not, the mathematical formulation of the first part of the problem then becomes Eq. 1 expresses the objective of maximizing NPV.
Eq. 2 requires that the flow in harvest volume between one period and the next deviates less than a fraction a from the average volume.To include consideration to biodiversity a demand of a certain amount of old forest area, G _ p , in each period p should be satisfied (Eq.3).Eq. 4 guarantees that the ending inventory is at least R _ .Eq. 5 keeps track of the area of treatment schedules assigned to each stand.To ensure that a stand is subject to one and only one treatment in each period, only integer solutions are valid (Eq.6).
The second part of the problem deals with the aggregation of harvests, i.e. thinning and clear cutting.To describe if two adjacent stands, i.e. stands that share a border, are treated with a schedule that has harvest in period p or not, a new indicator variable, z ilp , is introduced.This variable takes on the value 1 if stand i and the adjacent stand l are assigned a treatment schedule that has harvest in period p and the value of 0 otherwise.Let Y be the set of stand pairs i < l which share borders and M ip the set of treatment schedules for stand i that consists of a harvest in period p.In order for z ilp to be set to 1 the appropriate schedules for stand i and l need to be made active.This is ensured by Eqs.8-9 as follows: The aggregation of harvest can now be expressed by minimizing the total outside perimeter of harvests.This is expressed by minimizing the sum of the total perimeter of all stands that are harvested minus double the length of the common border between pairs of adjacent stands that are harvested in the same period, thus, where B i is the perimeter of stand i and S il is the length of border between stand i and stand l.
The formulation of perimeter in the model is, if terms are rearranged, in principle the same as in the model by Fischer and Church (2003) with the difference that the z-variables are here defined as continuous and not binary.
To solve the two-objective problem the model is converted to a single-objective problem by weighting the two objectives, Eqs. 1 and 10, together: with weights such that w 1 > 0 and w 2 > 0 and w 1 + w 2 = 1.

Case Study
The landscape contains 22 100 ha productive forest land divided in 2821 stands located in northern Sweden, Fig. 1.The forest data for the area was estimated by the k nearest neighbor method (kNN) (Reese et al. 2003) with stands delineated with the algorithm developed by Hagner (1990).The initial age class distribution can be seen in Fig. 1.The planning horizon was set to 100 years divided into 10-year periods, i.e.P was set to 10. NPVs were based on a 3 percent real discount rate.The harvest volumes were allowed to increase and decrease within 10 percent from the average harvest volume, i.e. a was set to 0.1.The old forest demand was set to be 5 percent of the productive forest land in each period, i.e.G _ p was set to be 1100 hectares for p = {1,…,10}.The demand on standing volume after harvest in the last period was set to be at least 2.4 million m 3 which is equal to the initial standing volume.Treatment schedules with NPV and future forest conditions were simulated using the GAIA stand simulation system (Eriksson 1983, Hoen andEid 1990).Each schedule consists of different timings for the allowed silvicultural activities, thinning and final felling with appropriate regeneration following the harvest.The minimum age for final felling for each stand was set according to the forestry act of Sweden regarding corresponding site productivity (Swedish National Board of Forestry 1994).Thinning treatments had an intensity of 30 percent removal of the basal area, except for deciduous trees where the intensity was 10 percent.Stand establishment activities after clear felling were given by a fixed program and did not vary between treatment schedules.These settings resulted in that 25 174 schedules were generated, corresponding to an average of 8.9 schedules for each stand.Growth functions according to Agestam (1985) were used.Revenues were computed with functions from Ollas (1980) and the year 2002 timber price list for the region.Silvicultural costs was assessed with data from Brunberg (2002), and machine costs were computed with functions from SLU (1989) and with unit prices adjusted such that the average cost level coincided with the data given in Brunberg (2002).
The stated problem was solved using a branch and bound algorithm, which is a standard algorithm for solving mixed integer problems.The model was formulated with the AIMMS® 3.5 optimization modeling system (Bisschop and Roelofs 1999), and solved from within AIMMS® with ILOG™ CPLEX 9.0.The software was run on a PC with a 2.8 GHz Pentium 4 processor and 1.5 GB of RAM.In the case study a convergence bound of 1 percent was used, meaning that the solution procedure stops if the solver can guar-   antee that the current best solution is within 1 percent of the global optimum.
To produce a good estimation of the effi cient frontier or trade-off curve the stated problem was solved for 5 different weighting combinations, see Table 1.

Results
The differences in the spatial layout of harvest activities between different weight combinations are clearly seen in the maps over the harvest activities in period 1 (see Fig. 2).This is the period for which the pattern would have the most immediate bearing to operational planning.With a fairly high weight on the perimeter measure the harvest activities become more aggregated and the number of small isolated harvest areas decrease.This is refl ected in the average as well as the maximum size of contiguous harvest areas (see Fig. 3a and b).The amount of perimeter over time for all different cases could be seen in Fig. 4. It follows a consistent pattern; the higher the weight on perimeter the lower the amount of perimeter in each period.
To evaluate whether the proposed model had any effect on the thinning activities the thinning proportion in each period for each case was recorded   5).It seems that with an increased weight on amount of perimeter the thinning proportion decreases.However, even if the thinning proportion is infl uenced the total harvest volumes seem to be unaffected except for case 5 (see Fig. 6).
The trade-off curve between NPV and the length of perimeter is shown in Fig. 7.It appears that it is possible to decrease the perimeter of the harvest areas with a moderate decrease in NPV.The reduction in amount of perimeter between weightings 1 and 4 was almost 52 percent while the corresponding reduction in NPV was only 2.6 percent (see Table 2).
Finally, the solution times are shown in Table 2. Except for case 1, the solution time increases rapidly from a matter of minutes in case 2 to almost six hours in case 5.The extremely long run time for case 1 is to us diffi cult to explain, but is of little practical consequence; in case the spatial pattern is of no interest, as in case 1, you would remove the latter part of the model and, most likely, relax the integer requirement on the x-variables.

Discussion
The major focus of this study was to study the applicability of the approach both from a model related, or technical, perspective and from the viewpoint of the resulting economic and landscape consequences.This was done by presenting and solving a long term forest planning problem aiming at both a high NPV and a low amount of total perimeter of harvest areas.The results from the case study indicate that the perimeter meas- ure makes it possible to achieve rather different harvest area aggregations with limited economic sacrifice, at least for the chosen landscape and with a modest weight on the perimeter criteria.This is consistent with results from studies where related measures such as shape index, have been used (Öhman and Lämås 2005).
The method gives the decision maker possibilities of generating a range of plans and selecting the plan with the most appropriate spatial pattern.The method appears capable of handling problems of the size you find in long range planning in Swedish forestry.This goes for the number of stands and the number of planning periods.Under the more realistic perimeter requirements the total time approaches a couple of hours.It may be regarded as long.Still, considering the very long range nature of this planning stage and the infrequency with which it appears, it is quite realistic for practical planning purposes.Stand was here used as a basic unit.A similar formulation with pixels would probably not be able to handle a problem of the size required at this stage of the planning process, at least not without the further development of the solution procedure.An obvious advantage of using the method is that it is based on MIP for which standard solution procedures exist.Thus, the deliberations associated with the heuristic methods used in many studies that have approached the aggregation problem are avoided.The exact solution method also means a good assessment of the quality of the solution in terms of the objective function is obtained.Another advantage is that global constraints, whether on harvest activities or on required states of the forest, are easily included.However, even though the method appears to be immediately available for practical use there are aspects that need to be discussed.
The procedure forms large aggregates of harvesting areas.In some cases they are very large, encompassing contiguous areas of several hundred hectares.This could be good or bad depending on what values are pursued.It is profitable in terms of economies of scale.It is also beneficial from the point of view of certain environmental values.The more clustered harvesting areas will over time create landscapes with more interior habitats (Gustafson 1996, Gustafson 1998).The creation of harvesting sites of a couple of or sev-eral hundred hectares may be less attractive from the point of view of other values such as esthetic, recreational and reindeer herding.In many countries there is also a legally binding upper limit on opening size, making the control of maximum size a decisive factor for practical use in those cases.One way of tempering the maximum size of harvest areas would be to introduce needed elements of an area restriction model (ARM; Murray 1999).The ARM would however considerably increase the complexity of the planning exercise both with preparation of cliques of the MIP model and in limiting the run time.
A clear management effect of increased aggregation is the reduction of thinning areas.More and more of the harvest is in the form of final felling as aggregation increases because larger volumes can be extracted from smaller areas.Under the current growth and yield and economic conditions this appears to be of no consequence.However, under other conditions it could.One option would be to leave out thinning from the harvest area.This would probably create the opposite problem, i.e. harvests are moved from final felling to thinning the more aggregation is emphasized.Also, from an economic point of view it is not advantageous to disentangle thinning from final felling.Another option that appears more promising could instead be to constrain the thinning quota of the total harvest.The requirement could for instance be derived from a run with the model without the aggregation objective.Yet another option is to have separate perimeter criteria for final felling and thinning, respectively.
The approach that is proposed here could be seen as an alternative to models that explicitly includes the road network, an aspect that traditionally belongs to tactical planning.The adaptation of harvest activities to the construction or use of roads does have obvious spatial implications (see e.g.Gustafsson et al. 2000).One drawback of not including the road network is that some of the economic logic is lost.If roads need to be constructed and if the whole harvesting area uses the same road segment is not assessed.The results of the model should therefore be more applicable in forest areas with already developed road networks, as in most of Sweden, than in areas with virgin forests.The most obvious advantage of not modeling the road network is that the model becomes smaller.Tactical problems often result in such large problems that they rely on a sequential procedure, i.e. first a long range problem is solved and after that the result of that model is fed into the tactical planning problem (Davis et al. 2000).Thus, this approach should be more amenable to the augmentation of issues related to long range planning.
As noted above, economic aspects related to roads and road use are not addressed.That means that the NPV does not capitalize on the economies of scale ensuing from aggregating harvest in the form of, for instance, reduced road construction and maintenance and transport of harvesting machines between harvesting sites.Thus, the cost in terms of NPV of creating aggregates could be even smaller than reported.
As with all long range planning, the result from planning performed with this model is of an indicative nature.The results of this planning model suggest, along with traditional non-spatial models, where the cutting blocks could be.It then becomes the object of tactical planning to delineate borders and sequence the harvests over the next few years such that undesired effects of large harvest areas can be mitigated and adjustments to road networks made.
Although promising, the results should be checked for other forest areas.It can then be ensured whether the limited cost of aggregation found here is an artifact of the particular area or if it has more general bearing.In Swedish forestry, final felling blocks are on average about 10 hectares; 50 hectares is considered large.The results of our model indicate that it could be optimal to create harvest areas that are much larger than this.Following this, it would be valuable to have a new assessment of the implications of large harvest blocks, asking questions like 'How is biodiversity affected in the long run?' Or 'Is it or is it not in agreement with recreational values to concentrate harvest activities in a few locations rather than spread them over the whole forest area?' Another issue that the results bring up is how to conduct tactical planning.If tactical planning is traditionally aiming at aggregating stands from an unwieldy list of stands into a limited number of harvesting tracts, the result of our model rather seems to require the manager transforming a limited number of aggregates into suitable harvest areas.It is to be hoped that the method presented here will, after further tests and developments, find its way into the forest manager's tool box for spatial planning.

Fig. 1 .
Fig. 1.The case study area with initial age class distribution.

Fig 2 .
Fig 2. The distribution of harvest area in (a) case 1, period 1 (b) case 2 period 1 (c) case 3, period 1 (d) case 4, period 1 (e) case 5, period 1.The red areas indicate stands treated with final harvest and the blue areas indicate stands treated with thinning.

Fig. 3 .
Fig. 3. a) Average size of harvest areas over the planning periods for each case.b) Maximum size of harvest areas over the planning periods for each case.

Fig. 5 .
Fig. 5.Thinning area in percent of total harvest area for cases 1-5 over the planning periods.

Fig. 4 .Fig. 7 .
Fig. 4. The trends in amount of perimeter for the different weight combinations.

Fig. 6 .
Fig. 6.The harvest volume over the periods for each weight combination.

Table 1 .
Weights for the two objectives.ε is a very small number.

Table 2 .
The solution times for all weight combinations and the reduction in NPV and in amount of perimeter between weight combination 1 and all other combinations.