Optimal Management of Pinus radiata Silvopastoral Systems Established on Abandoned Agricultural Land in Galicia ( North-Western Spain )

Timber production has been the main objective in forest production in Galicia for a long time. Nevertheless, factors such as fire risk and the need to obtain non-timber benefits make other production alternatives like silvopastoral systems worth of consideration. Integration of grazing in the production system not only diversifies products and benefits, but also decreases fire risk by enhancing fuel control. Nonetheless, few studies have examined the economic profitability of these systems. This article analyses the economics of silvopastoral systems established on abandoned agricultural soils afforested with Pinus radiata D. Don. Different tree planting densities, discounting rates, grass values and fire risk scenarios were analysed. The technique employed is based on the combination of an optimization algorithm and a simulator of stand growth and grass yield. The most profitable schedules were obtained with initial stand densities of 1500 trees per hectare. However, with high unit values of pasture production (high value of grass), schedules with an initial stand density of 500 trees per hectare were the most profitable. When the risk of fire was included in the analyses, silvopastoral systems were always more profitable than timber production systems. With an assumption that grazing reduces fire risk thinnings should be done earlier and heavier to reduce the expected losses due to fire and to promote grass production. This lengthens the pasture period. In general, rotation lengths of silvopastoral systems were shorter than in timber production.


Introduction
The main use of Galician forest is timber production.Many forests are managed to maximise biomass production for pulp and board industries.Steady income from these management regimes is hampered by both forest fires (Núñez Regueira et al. 2003) and low timber prices.Therefore, there is a need to search for alternative regimes that would make income-generation less risky.Silvopastoral systems could be an appropriate alternative.Grazing has sometimes been considered to promote fire damages because fire has been carelessly employed by shepherds to promote the growth of more palatable grass species.However, if prescribed fires are excluded, grazing reduces the risk of fire by diminishing the fuel loads in the forest (Rigueiro-Rodríguez et al. 2005).Another advantage is the multiplicity of products that make silvopastoral systems economically safer under market uncertainty than the traditional timber production oriented forestry (Anderson andSinclair 1993, Sharrow 1999).Furthermore, silvopastoral systems generate incomes much earlier (Sharrow 1999) than pure timber production systems.Finally, silvopastoral systems improve accessibility (Knowles 1991) and scenic value of the landscapes, enhancing their recreational use (Ruark et al. 2003, Alavalapati et al. 2004).
Despite the advantages showed by the implementation of silvopastoral systems, not much research has focused on the optimal design of these systems.One of the very few examples is the study of Muchiri et al. (2002b), which optimized the management of an agroforestry system composed of maize and Grevillea robusta.
This study is focused on silvopastoral systems established on abandoned agricultural lands.These lands are fertile (stand dominant heights up to 30 meters at 20 years for Pinus radiata D. Don).Accordingly, this study analyses only good sites on which economically viable agroforestry is possible (Hawke 1991).Pinus radiata has been the most common tree species employed in silvopastoral systems.The system is established by planting trees and sowing grass at the same time.Therefore, there is forage production already in the first year suitable for instance to sheep grazing (Rigueiro-Rodríguez et al. 2002).A study carried out at the Department of Crop Production at the University of Santiago de Compostela (Spain) found that silvopastoral systems of this type need an area of about 200-300 hectares to make their implementation profitable.On this scale the management costs of the silvopastoral system (veterinary costs, shepherd costs and other related costs) become affordable.It has been claimed (Adams et al. 2001) that the establishment of conifer plantations on this type of terrains may decrease soil fertility in the Spanish Atlantic region due to a pH reduction.This impoverishment of soil fertility leads to a change in the type of understorey vegetation, from herbaceous vegetation with low fuel loads towards more inflammable shrub communities (Rigueiro-Rodríguez et al. 2005, Mosquera-Losada et al. 2006).Therefore, it seems even more important to keep the herbaceous stratum at the understorey level by grazing in order to reduce the risk of fire as long as possible.
For the optimal management of silvopastoral systems, the influence of trees on pasture production must be known.The key factor of the success of the system is to achieve a compromise between the two sources of economic benefit.Grazing is possible when the tree canopy allows light to reach the understorey layer.Canopy cover is commonly used to set the limits for pasture production (Knowles et al. 1998).Literature suggests that canopy covers higher than 50% drastically decrease the pasture production (Rigueiro-Rodríguez et al. 1998).Other studies suggest a maximum canopy cover of 70% (Knowles et al. 1998).Also green crown lengths (Percival andKnowles 1983, 1988), horizontal projection of the crowns (Sibbald et al. 1994) or stand density (Pearson et al. 1995, Burner and Brauer 2003, Rozados-Lorenzo et al. 2007) are variables that have been used to predict pasture production.Canopy cover is difficult to measure in the field and predict in simulations.Green crown length and horizontal projection of the crowns are also problematic because they are not measured in normal inventories.Stand density (number of trees per hectare) is not a good predictor neither because, by itself, it does not give enough information about competition in the stand.Therefore, we decided to model the dependence of pasture production on stand basal area and site index, which are easily obtained from regular inventories.Using this model with a growth and yield model for tree stand dynamics, we were able to calculate the profitability of the system.The aim was to study how the economic profitability and optimal management of silvopastoral systems established on abandoned agricultural terrains depend on site quality, grass value and planting density of trees.The effect of fire risk on the profitability of the silvopastoral systems was also studied.Moreover, we analysed the effects of an assumption that grazing reduces fire risk by diminishing fuel loads and promoting less inflammable species.The analyses of the study are divided into three parts: effect of (I) stand density and grass price, (II) fire risk and (III) the influence of grazing on fire risk, on the economic profitability and optimal management of silvopastoral systems.

Simulation of the Tree Stand Dynamics
Silvopastoral systems have three components, tree stand, forage and livestock.We used the model of Castedo-Dorado et al. (2007) for even-aged P. radiata stands in Galicia to simulate stand development in different management schedules.In this model, the initial stand conditions are defined by three state variables: number of trees per hectare, stand basal area and dominant height.The model uses three transition functions to project each state variable for a given time period.It also includes a function for predicting the initial stand basal area when no inventory data are available.Once the state variables are known for a specific moment, a distribution function is used to estimate the number of trees in each diameter class by recovering the parameters of the Weibull function, using the moments of the first and second order of the distribution.By using a height-diameter function to estimate the height of the average tree in each diameter class, and a taper function, the total and merchantable stand volume are calculated.
The model for the dominant height development is as follows: with, 1 755 4 12 44 ) . ) (2) ) ) (3) where H 1 is the dominant height (m) at age T 1 (years), and H 2 is dominant height at age T 2 .
Reduction in the number of trees per hectare (natural mortality) is predicted with: where N 2 is the number of trees per hectare at age T 2 and N 1 is the number of trees per hectare at age T 1 .The following function was used for basal area initialization:

SI N T
where G is stand basal area (m 2 ha -1 ) at age T (years), N is the number of trees per hectare and SI is the site index (m), estimated using Eq. 1 at a reference age of 20 years.The function for basal area projection is: where G 2 is the stand basal area (m 2 ha -1 ) at a given projection age T 2 , and G 1 is stand basal area (m 2 ha -1 ) at age T 1 .The equation for predicting the arithmetic mean diameter, to be used to derive the diameter distribution with the parameter recovery approach, is: where d is the arithmetic mean diameter (cm) and d g the quadratic mean diameter (cm), The equation for predicting the height of a representative tree in each 1-cm diameter class is: / .
where h is the total tree height (m), d is diameter at breast height (cm), and D and H are, respectively, dominant diameter and dominant height of the stand.
Both uniform and low thinnings can be simulated as intermediate treatments.Uniform thinnings remove an equal percentage of trees from every diameter class.When a low thinning is simulated, the remaining number of trees in diameter class i (n i ) is calculated following the distribution independent approach proposed by Alder (1979): where N before is the total number of trees per hectare before low thinning, L is low-thinning intensity expressed as one minus the proportion of removed trees (1-N removed /N before ) and

Simulation of Grass Production
Pasture production is highly dependent on the stand development since pasture production is only possible when the canopy allows light to reach the understorey level.Site and stand characteristics were used as predictors to fit a model for the pasture production.We used data from an experiment in Castro de Riberas de Lea (Lugo) that consisted of the measurement of trees and pasture production in a silvopastoral system during seven years since tree planting (see Mosquera-Losada et al. (2006) and Fernández-Núñez et al. (2007)).The trial has plots of different site quality and two different planting densities (833 and 2500 trees ha -1 ).The pasture was a mixture of Lolium perenne L., Dactylis glomerata L., Trifolium repens and Trifolium pratense.The fitted regression model is as follows: where grass is the annual grass production (dry mass) (t ha -1 ), SI is site index of Pinus radiata stand (m) (dominant height at the reference age of 20 years) and G is the basal area (m 2 ha -1 ) of the tree stand.The R 2 of the equation is 0.425.This equation shows that the higher the basal area is, the smaller the grass production becomes (Fig. 1).
Better sites produce more grass, as expected.The pasture was considered to generate income only when the grass yield was higher than 0.3 t ha -1 .This is the minimum amount required for feeding one sheep per hectare per year (data provided by the Department of Crop Production of Santiago de Compostela).

Study Cases
The effect of stand density and unit value of pasture production was tested with two site indices for radiata pine, namely 29 and 25 meters at 20 years.These site indices are typical for silvopastoral systems.The studied planting densities ranged from the sparsest to the densest stockings used in forestry practise in the region: 500, 1500 and 2500 trees ha -1 .Two site indices with three planting densities resulted in six different initial stands.In every stand we tested different number of thinnings (0-2) and different unit value of pasture production.The revenues from pasture come from the animals fed by the grass.Grass production was converted into fed livestock (lamb and sheep) to calculate the income using data of the Department of Crop Production at the University of Santiago de Compostela.The data indicate that one ton of grass can feed three sheeps.Taking into account that each sheep delivers 1.6 lambs per year on an average, one ton of grass generates an annual income of about 200 € when silage making and all the related costs such as veterinary, shepherd and silage are considered.This is called as the unit value of pasture production.This value was varied to see the effect of market fluctuations; the used unit values of pasture production were 100, 200 and 400 € t -1 .Fire risk was assumed to have two components: probability of occurrence and damage.We tested four different probabilities of occurrence: 0, 1, 3 and 5%.When fire takes place we assumed that it ends the rotation prematurely and only a part of the growing stock volume can be harvested (salvaged).The proportion of salvaged timber describes the second component of risk: damage (salvage rate = 1 -damage rate).The proportion of timber that can be salvaged depended on the mean tree diameter (Pasalodos-Tato et al. 2009b) (Fig. 2): where d is the diameter at breast height measured in cm.In addition to loosing a part of timber in fire, the salvaged timber was depreciated by 25% (Arenas and Izquierdo 2007).This price reduction of salvaged timber was used in all optimisations.The optimizations were done for one site index (29 meters) and two different stand densities (500 and 1500 trees ha -1 ).
Grazing may reduce fire risk by reducing fuel loads and promoting less inflammable species (Rigueiro-Rodríguez et al. 2005, Rigueiro-Rodríguez et al. 2009).Even though the literature has mentioned this effect widely (Blackmore and Vitousek 2000, Elmore and Asner 2003, Casal et al. 2009) it is difficult to find quantitative information on it.Therefore, we used several reduction factors, namely 25, 50, 75 and 100% to reduce the probability of fire occurrence in every grazing year (grass yield ≥ 0.3 t ha -1 ).Two different stand densities (500 and 1500 trees ha -1 ), one site index (29 meters) and two different probabilities of fire occurrence (1 and 5%) were employed to analyse the effect of the reduction in fire risk due to grazing.

Objective Function
Soil expectation value (SEV) calculated with 3% discounting rate was used as the objective variable.The SEV was calculated as the net present value (NPV) of all future net incomes: where NPV is the net present value of one rotation, r is the discounting rate and R is the rotation length (years).The expression for the NPV is: Where I wt and I gt are the incomes and C wt and C gt are the costs derived from timber and pasture production in year t, respectively (see Tables 1  and 2).The incomes from timber production (I wt ) were calculated from: where s t is the proportion of salvage calculated from Eq. 13 (s t = 1 if there is no fire), J is the number of diameter classes, n j is the number of trees in diameter class j, P k is the unit price of timber assortment k and v kj is the volume of assortment k of a tree in diameter class j.The following road side timber prices were used: 90 € m -3 for grade I, 50 € m -3 for grade II and 18 € m -3 for grade III (see Pasalodos-Tato et al. 2009a, b).The unit price was reduced by 25% when fire ended the rotation.
The costs of the silvopastoral system (both C wt and C gt ) depended on site index (Table 1).Timber production costs were different when there was no grazing (Table 2).Regeneration cost (Tables 1 and 2) was assumed to be a linear function of the number of planted trees per hectare with the constant part representing the cost of site preparation and the variable part representing the planting cost per tree.In silvopastoral systems the regeneration cost is higher because of an additional cost of individual tree protectors to avoid the damages that sheeps can cause on the seedlings.
The tree harvesting cost was calculated from (based on Ambrosio et al. 2000): where HCost is harvesting cost (€ ha -1 ), ECost is entry cost (€ ha -1 ), V is the total harvested volume (m 3 ha -1 ), FCost is forwarding cost (€ m -3 ), S is slope (%), and ν is the mean volume of harvested trees (m 3 ).It was assumed that the entry cost of moving the machinery to the forest (ECost) is 200 € ha -1 .The forwarding cost was assumed to be 5 € m -3 and the slope was taken as 20%.

Integrating Fire Risk into Objective Function
In order to include fire risk in the calculation of SEV we used the approach developed by Bright and Price (2000).The method consists of the sum of all possible outcomes, weighted by their probabilities.The expression for the expected SEV was (see Pasalodos-Tato et al. 2009a, b): where p t is the probability that the stand burns in year t and survives the previous years, i.e., p t = (1p fire ) t p fire , where p fire is the annual probability of fire occurrence, and p R p R = (1p fire ) R  is the probability that there is no fire before the rotation age.NPV first is calculated from: where NPV t is the net present value if fire hits the stand at age t and ends the rotation prematurely, and NPV R is the net present value if there is no fire during the rotation (R).
A thinning intensity higher than 30% was assumed to make the stand sensitive to windthrow and snow breakage (Castedo-Dorado et al. 2009).Therefore, a penalty function was added to the SEV of the management schedule as a means to avoid too heavy thinnings.The objective function (OF) which was maximized in optimization was therefore where H% m is thinning intensity in percent of removed stand basal area in thinning m and M is the number of thinnings.The penalty function implies that the penalty of harvesting too much at a time increases from 0 to 10 000 € ha -1 when the thinning percentage increases from 30 to 100.

Decision Variables
Decision variables such as the number and intensity of thinnings, and rotation length define the management schedule.Optimizing a management schedule is equal to finding optimal values for decision variables.Due to the fact that the number of thinnings is not a continuous variable schedules that have a different number of thinnings must be treated as different optimization problems.In this study management schedules were optimized with 0, 1 and 2 thinnings, which are all feasible options for Pinus radiata silvopastoral system.
The simulated thinnings were combinations of uniform and low thinning.Therefore the management regime was defined by the number of thinnings and the following decision variables: The number of optimized decision variables was therefore 3 × M + 1 where M is the number of thinnings.

Optimisation Method
The optimisation algorithm used was the direct search method of Hooke and Jeeves (1961).This method uses a form of coordinate optimization and does not require explicit evaluation of any partial derivative of the objective function.The direct search method compares each new trial solution with the best obtained up to that time.The search has two components, the exploratory search and the pattern search.For a given base point, the exploratory search examines points around that base point in the direction of the coordinate axes (decision variables).The pattern search moves the base point in the direction defined by the given (current) base point and the best point found in exploratory search.

Profitability of Silvopastoral Systems
After running the optimisation for the six different initial stands (2 site indices with 3 planting densities) with three different thinning schedules (0, 1 and 2 thinnings), we chose that number of thinnings that gave the maximum SEV.The optimal schedules had one thinning with initial density 500 trees ha -1 and two thinnings with the other planting densities with both site indices (25 and 29 m).Silvopastoral system was always more profitable than mere timber production, planting density 1500 trees ha -1 being the most profitable (Fig. 3).The establishment of pasture improved profitability most with the lowest density, 500 trees ha -1 .The improvement was up to 50% with a unit value of pasture production of 200 € t -1 .SEV improved 15% with planting density 1500 trees per hectare and 4-7 % with 2500 trees per hectare.
The optimal rotation lengths without pasture were 40 and 42 years, respectively, for planting densities 500 and 1500 trees per hectare in site index 25 m, and 38 years for both densities in site index 29 m.In general, rotation lengths decreased with the inclusion of pasture.This decrease was more noticeable with lower planting densities (Fig. 4).

Effect of Grass Prices
When the value of grass was high (400 € t -1 ) schedules with lower planting densities became the most profitable, and optimal rotation lengths were shorter (Fig. 4).With grass value of 400 € t -1 and planting density of 500 trees ha -1 the most profitable alternative was to produce only forage since tree growing was no longer profitable.However, since silvopastoral systems were analyzed in this study, tree planting was forced in the solution although the landowner should not plant trees in this case.With a low unit value of grass production (100 € t -1 ) silvopastoralism was not the best alternative anymore with the highest planting density (2500 trees ha -1 ).

Effect of Fire Risk
Silvopastoral systems were always more profitable than timber production systems when the risk of fire was included in the analysis (Fig. 5).When the planting density was 500 trees per hectare, the improvement in SEV was 55% for 1% annual fire probability, 93% for 3% probability and 167% for 5% probability.The trend was the same with 1500 trees per hectare.When the annual probability of fire occurrence was 5% the profitability of the silvopastoral system was 40% higher than in timber production.The superiority was 24 and 14%, respectively, with annual fire probabilities of 3 and 1%.
The optimal rotation lengths were shorter with increasing fire risk.With the planting density of 500 trees per hectare the optimal rotation length decreased 5% and 11%, respectively, when the probability of fire was 1% and 5% (Fig. 6).With planting density of 1500 trees per hectare the reduction was from 3 to 6%.The higher is the fire risk the heavier and earlier the thinnings become.

Results when Grazing Reduces Fire Risk
The more grazing was assumed to reduce fire risk, the heavier and earlier the optimal thinnings became (Fig. 7).With risk reductions of 50% or more the optimal thinnings were so heavy and early that grazing could continue for most of the rotation (Fig. 7).However, this happened only with a high risk of fire (5% annual probability without fire reduction).With low planting density and high fire risk, increasing risk reduction due to grazing shortened optimal rotation lengths.When the annual fire probability was 5% and planting density was 500 trees per hectare, risk reductions of 25, 50, 75 and 100% resulted, respectively, in an increase of 6, 12, 20 and 30% in SEV (Fig. 8).The improvements were slightly smaller for planting density 1500 trees per hectare.When the annual fire probability was 1% the increase in SEV was not much, only 2% to 4%.Reduction in risk of fire Fig. 7. Development of stand basal area and annual grass yield in the optimal management schedule of silvopastoral systems when grazing is assumed to reduce the risk of fire.The planting density is 500 trees per hectare, annual fire probability is 5% and unit income from grass is 200 € t -1 .
Fig. 8. Soil expectation value of the optimal silvopastoral schedules for different planting densities (500 and 1500 trees per hectare) and probabilities of fire occurrence (1 and 5%) with unit income from grass of 200 € t -1 when grazing is assumed to reduce fire probability by 0, 25, 50, 75 and 100%.

Discussion
This study showed that with the current grass values, silvopastoralism is more profitable in Pinus radiata plantations than mere timber production on good sites in Galicia.The study presents a method that allows managers to optimise the design of silvopastoral systems established on abandoned agricultural land relating pasture production to stand basal area and site quality.
The study has some limitations, which are mainly related to the lack of information regarding the dynamics of silvopastoral systems.One shortcoming is the assumption that there is no difference between the competition effect of trees on understorey vegetation on improved and unimproved grass vegetation.Even though many studies show that grasses compete with trees in plantations, there are also studies (Ares et al. 2003) which suggest that the sowing of a pasture improves tree growth because grasses compete less than shrub species with trees.Therefore, when grass occupies the lower level of the forest, some modification may be required in growth and yield models, in order to account for the effect of grass on early tree growth.It would also be helpful to study the species dynamics of the understorey vegetation to better predict the influence of grazing on fire risk.
The results on the optimal number of thinning with different planting densities were logical: the denser the stands, the more thinnings.Some other results attract attention.The first one is that a schedule with an initial stand density of 2500 trees ha -1 is never the best, not even when pasture is not considered.This result has practical importance because dense stands present difficulties also from the practical point of view, e.g. when mechanising some operations such as thinnings and cleaning, since some machinery can not easily enter the stand.
Silvopastoral systems generate revenues much earlier than timber production systems.Silvopastoral systems have shorter optimal rotation lengths than timber production systems.This difference is greatest with high grass value and low planting density.In these cases grazing generates much incomes, most of which are obtained in the beginning of the rotations.Therefore, it is profit-able to shorten rotation lengths so as to have a productive pasture period soon again.
When the annual probability of fire was at least 1% silvopastoral systems were more profitable than timber production even when the unit value of pasture production was halved.The relative profitability of silvopastoral systems was enhanced with increasing probability of fire occurrence.This is because timber production needs many years to reach a reasonable mean annual income, and the probability that all production is lost due to fire is high for long production times.If the first fire occurs before the trees are merchantable, all the grass produced so far has generated income, but all benefits from tree planting are lost.
Planting density of 1500 trees per hectare is more profitable than 500 trees per hectare if grazing is considered.However, when the annual fire probability is 3% or 5% 500 trees per hectare becomes the best planting density.The reason is that with 1500 trees per hectare the proportion of timber revenues of the total income is high, which means that the potential losses due to fire are also high.With 500 trees per hectare the incomes from pasture dominate and the potential losses due to fire are therefore low.When fire occurs and the rotation ends prematurely, a new pasture is available already next year.
Further studies in the optimization of silvopastoral systems require the integration of fertilisation effect in the production models of both grass and trees because fertilization is commonly used to increase grass yields.Inclusion of amenity values, which are significant in silvopastoral systems, would also be an interesting research topic.The effect of grass on tree growth should also be studied.The effect of tree cover on grass yield also requires additional research, spatial modelling being an interesting option (Muchiri et al. 2002a, b).
the cumulative frequency distribution at diameter d i .The taper model proposed by Fang et al. (2000) fitted for P. radiata by Castedo et al. (2007) was used to calculate the stem volume of trees extracted in thinning operations or clear cuttings.The following top diameters were used: 35, 18 and 7 cm.The timber assortments therefore corresponded to the following over-bark stem diameters: (I) d ≥ 35 cm; (II) 35 cm > d ≥ 18 cm; and (III) 18 cm > d ≥ 7 cm.The following minimum piece lengths were assumed in this study: (I) 3.0 m; (II) 2.5 m; and (III) 1.0 m.If the piece was shorter, the volume was moved to the next (with a smaller minimum top diameter) timber assortment.

Fig. 1 .
Fig. 1.Dependence of grass production on stand basal area for different site indices.

Fig. 2 .
Fig. 2. Dependence of the salvage proportion of the breast height diameter.
For thinnings: -Stand age at the first thinning and number of years between the first and the second thinning.-Percentage of uniform thinning (% of number of trees) -Percentage of low thinning (% of trees removed after uniform thinning) For final felling -Number of years since the last thinning

Fig. 3 .
Fig. 3. Soil expectation value of the optimal silvopastoral schedules for different initial stand densities when different grass prices are considered.

Fig. 4 .
Fig. 4. Development of stand basal area and annual grass yield in the optimal management schedule for different silvopastoral systems for different stand densities and unit value of grass production when site index is 29 meters.

Fig. 5 .Fig. 6 .
Fig.5.Soil expectation value of the optimal timber production (black) and silvopastoral schedules (grey and white) for different planting densities (500 and 1500 trees per hectare) and annual probabilities of fire occurrence (1, 3 and 5%) when the unit value of grass production is 200 € t -1 .

Table 1 .
Years and costs of tending operations for silvopastoral systems.N is the number of planted trees per hectare.

Table 2 .
Years and costs of tending operations for a timber-production schedule.N is the number of planted trees per hectare.