A two-stage stochastic programming model for scheduling replacements in sow farms

This paper presents the formulation and resolution of a two-stage stochastic linear programming model with recourse for sow farms producing piglets. The proposed model considers a medium-term planning horizon and specifically allows optimal replacement and schedule of purchases to be obtained for the first stage. This model takes into account sow herd dynamics, housing facilities, reproduction management, herd size with initial and final inventory of sows and uncertain parameters such as litter size, mortality and fertility rates. These last parameters are explicitly incorporated via a finite set of scenarios. The proposed model is solved by using the algebraic modelling software OPL Studio from ILOG, in combination with the solver CPLEX to solve the linear models resulting from different instances considered. The article also presents results obtained with previous deterministic models assessing the suitability of the stochastic approach. Finally, the conclusions drawn from the study including an outlook are presented.


Introduction
Traditional pig production in Spain was based on small familiar farrowing -to finish farms, but this is undergoing a rapid change as in the rest of the European Union (EU).Nowadays, production is being concentrated in bigger and specialised pig production units.In Spain, commercial pig production tends to be divided into three different phases according to the final product and the activities involved.The first phase relates to farms producing piglets, the second one to those producing feeder pigs and the third one to those producing fattened pigs.This division provokes a specialisation in farming activities by phase and gives additional efficiency gains as Rowland et al. (1998) already pointed out.The first phase represented by sow farms producing piglets is the most important because of the complexity of the reproduction process and the caring that piglets need before being sold or transferred to a different unit in the next phase.
The sow herd structure is central to maintain a steady production over time and replacement is the most crucial decision, not only due to consequences on piglet production, but also for being the main production cost (Jalvingh et al., 1992;Huirne et al. 1993).Furthermore, EU regulations concerning pig welfare have reduced profitability margins and increased competition.As consequence the interest in the optimisation or in improvements of management strategies by implementing suitable decisions in sow farms have been increasing constantly.
Improvements in modelling the decision-making process on sow farms to represent fairly the system have been done and thus, advances have been obtained in solving or circumventing methodological problems related to complex models (e.g.Kristensen, 1988Kristensen, , 1993)).Most of them are related to Markov chain and simulation models (see Plà, 2007) but none with stochastic programming in a finite time horizon.Researchers have the benefit of advances in computing, database and solving software which enables farming systems to be described in greater detail and with greater ease (Kingwell, 1996).Taking advantage of this situation the aim of the model presented in this paper is to provide a complimentary analytical tool to specialists and advisers serving to better decisions around the culling and replacement of sows.This is done scheduling purchases of gilts and culling less productive sows in the way that herd structure is preserved and productivity levels maintained.In the present paper, this scheduling is the result of solving an optimization model that can also represent explicitly the uncertainty present in the system.
After the seminal papers of Dantzig (1955) and Beale (1955), most of the references concerning optimization models in the presence of uncertainty come under the name of Stochastic Programming that allows to explicitly incorporate the uncertainty of the parameters in the model formulation, see Birge and Louveaux (1997), Ruszczynski and Shapiro (2003) and Wallace and Ziemba (2005), among others.In particular, a twostage stochastic program with recourse is an important class of models in stochastic programming and widely used in multiperiod planning problems.In such models, two kinds of decision variables there exist.The first-stage decisions represent proactive decisions whose values are not conditioned by any particular realization of the uncertain parameter.In this paper, these decisions are simply related to the purchase of gilts and replacement of sows in the more immediate planning period.On the other hand, the second-stage or recourse variables represent reactive decisions made in recourse or response to compensate for the decision made in the first stage after materializing the uncertainty and correspond to the remaining decision variables in the proposed model.
The objective of this paper is precisely to formulate and solve a two-stage stochastic programming model for scheduling replacements in sow farms.The model includes the productive and reproductive behaviour of a group of breeding sows over time where piglets are the commercial product.Hence, the herd model is mainly focused on reproduction and replacement management of sows.This model maximizes profits carrying out an efficient occupancy of farrowing facilities because these are the most expensive and as a consequence have a high impact in (or bound) the production.
The remainder of this paper is organized as follows.The description of the problems related to operations in sow farms is presented in section 2. This is followed by a deterministic model in section 3 that will serve as a base for the formulation of a more elaborated model incorporating parameters under uncertainty using a two-stage stochastic programming model in section 4. Results of both models and discussion are presented in section 5. Conclusions and future work are drawn in section 6.

Scheduling replacement and medium time decisions
Normally Spanish sow farms have three different facilities: breeding-control, pregnancy and farrowing facility.The breeding facility is where sows are inseminated and controlled in order to confirm the pregnancy (more or less three weeks after the insemination).Once the pregnancy is positively confirmed they are moved to the pregnancy facility.Otherwise, it is considered that conception has failed and they remain in the breeding facility for subsequent re-inseminations, according to a number of attempts that must be determined optimally.The farrowing facility is where farrowing and weaning operations are done.So, before farrowing (normally one week), pregnant sows are moved to the farrowing facility and sows remain there until weaning (normally 3 weeks after farrowing).
Different operations are performed on a sow farm.In general a weekly basis period is adopted to rationalise the daily work.Thus, for instance inseminations are scheduled Tuesdays and Thursdays, weanings on Wednesdays and so on.Also, purchases of gilts and culling of animals are not carried out daily, they are also grouped one day a week to send and receive breeding animals.Replacement is especially important because it determines future productivity of the herd.This is so because along the age-structure of the herd gilts and old sows are less productive than young or medium age sows.In addition, the age-structure of the herd is also related with sanitary aspects as the sensibility to diseases outbreaks or passive immunization by contact between young and mature sows.
Scheduling makes no sense in the long term because the population dynamics in a commercial sow herd is high (regular annual replacement rates are easily around 50%).Also, seasonal variations in reproductive performance and changes in pig meat demand are observed and may affect net revenues.Then scheduling should be modified depending on actual state of the herd and expected market conditions.When scheduling, farmers take into account future possible variations, and scheduling consider a medium term planning horizon.
Different replacement models concerning sow farms or other livestock species can be found in the literature (Jalvingh et al., 1992;Kristensen, 1993;Plà, 2007).However, most of them are developed for research purposes and then the details with practical relevance for daily operations on farms are left aside.

The deterministic model
In this section, a deterministic model used to determine the optimal purchase and replacement policy for a given planning horizon is briefly described, a detailed version can be found in Rodríguez et al. (2008).The deterministic model assumes that we know with certainty all the parameters present in the model.This assumption is very restrictive but will serve as a base for a more complex extension, in the next section, that replaces some of the uncertain parameters through a finite set of scenarios in the model formulation.In the model, the knowledge of population dynamic is crucial because the model represents the future state of the sow farm in terms of a present state.
More precisely, the life span of the sow is divided into different states in a way that all possible transitions of sows evolving from one state to another can be represented.
Assuming a given scenario for all the parameters in the life span of the sow, the deterministic model considers the following notation.The proposed deterministic model maximizes the total profit of the production plan and the non-negative feasible solutions must satisfy a set of constraints that mainly concern the population dynamic behavior and capacity constraints, given by the following optimization problem: The objective function in (1) represents the maximization of the total profit, that is, the addition of profits for each period of time and reproductive cycle per decision variable.
Profit is the difference between incomes obtained by sales and the different costs of production incurred.Incomes consider the sales of piglets and sales to the slaughterhouse of replaced sows while, the costs of production include feeding costs of sows and piglets, labour, insemination and veterinary expenses.According to incomes and costs affecting each decision variable coefficients of the objective function are set.
For instance, the coefficient rzt,1 involves feeding, insemination and labour but also takes into account the purchase cost of a new gilt.Constraints (2)-( 5) describes the initial stock at t=1.They account for the initial herd distribution of sows over different states and can be adapted to any particular situation either a starting farm or an existing one.
Constraints ( 6)-( 13) represent the herd dynamics, i.e., the flow of sows throughout the different states over time.More specifically, constraint (6) refers to the number of sows at the first insemination state.Constraints ( 7)-( 9) refers to states of sows being under control for pregnancy or a new insemination, the so-called repetition.Constraints (10) and ( 11) refer to the flow among gestation states.Constraints ( 12) and ( 13) refer to the flow among lactation states.
Constraints ( 14) -( 16) are the facilities capacity constraints.Three different facilities are considered depending on the state of the sow with respect to the reproductive cycle, i.e. breeding-control, gestation and farrowing facilities.
Constraints ( 17) and ( 18) refer to how purchases are limited.Then, a smooth variation in purchases from week to week was allowed (17) and a minimum and maximum amount of new gilts purchased ( 18) is fixed per week.
Constraints ( 19)-( 22) determine the imposed inventory of animals at the end of the planning horizon, representing the continuity of the farm beyond the end of the finite time horizon considered and imposed values that tend to the stationary or long term optimal decisions level, according to some complementary studies and methodologies used by the authors (Pla et al., 2008) Constraints ( 23) and ( 24) refer to culled sows.Constraint (25) refers to gestating sows suffering an abortion and culled from the herd.

Stochastic programming model
The proposed optimization model of the previous section includes, in practice, some parameters often known with uncertainty.The classical approach to face the uncertainty consists in replacing the stochastic parameters by theirs expected values in a deterministic model as the previous one.However, the optimal solution achieved in this way might not be sufficiently representative of the reality and it does not take into account the variability of these parameters with respect to their expected value.Thus, among the different methodologies for production planning under uncertainty (Mula et al., 2006), the previous deterministic model is reformulated as a stochastic optimization program, which maximize the expected value of the farm's profit.In what follows, a two-stage stochastic optimization model is stated to the corresponding production planning problem, where we look for a first-stage optimal decisions mainly related to the purchase of gilts and replacement of sows in the first T1 weeks of the planning horizon and a second-stage decisions related to optimal policy in the rest of the planning horizon as a recourse policy.The proposed model is based on previous papers by Escudero et al. (1993), Albornoz and Contesse (1999), Gupta and Maranas (2003), Alonso-Ayuso et al. (2005) and Albornoz and Canales (2006).In spite of the fact that all the decision variables are defined by scenarios, the decisions related to the more immediate planning period T1 satisfy an additional set of constraints, known as the non-anticipativity constraints.These constraints impose to concerned decision variables a value that does not depend on any particular scenario realization.
These decisions are called here-and-now decision variables and guarantee identical first-stage decisions (in particular sow replacement and purchase of gilts the T1 period) for all the scenarios considered.Therefore the decision variables for the rest of the The objective function in (26) represents the maximization of the expected profits.As in the deterministic model constraints ( 27)-( 30) describe the initial stock at t=1.The initial stock is the same for each scenario since only express the initial stocks of animals in farm considered.Constraints ( 31)-( 38) represent herd's dynamics.Constraints ( 39)-( 41) are the facilities' capacity constraints.Constraints ( 42) and ( 43) refer to how purchases are limited.Constraints ( 44)-(47) determine the imposed inventory of animals at the end of the planning horizon, where each scenario tends to its own long term optimal decisions level.Constraints ( 48) and ( 49) refer to culled sows.Constraint (50) refers to gestating sows suffering an abortion and culled from the herd.Constraints ( 51)-( 56) refer to the non-anticipativity constraints, therefore only decision variables in the first stage, i.e. with t  T1, are involved.

Computational results
In order to illustrate the suitability of the proposed deterministic model ( 1)-( 25) and the corresponding stochastic extension ( 26)-( 56), a case study is presented.Basic parameters of the study were taken from standard values under Spanish conditions and recorded in the BD-Porc databank (national record keeping system hosted at http://www.irta.es/bdporc/,accessed 14 May 2008), and do not correspond to a specific farm.An initial herd size of 2330 sows was considered a regular size.A maximum lifespan of 8 parities was allowed.The effective capacity of lactation facility was of 500 crates.The maximum number of allowed insemination was three, beyond that number it was considered as a culling reason, just like an abortion.The time horizon was 52 weeks (approximately one year).Models were implemented and solved using the algebraic modelling language ILOG OPL 6.1 and the solver CPLEX 11.2 respectively on a laptop computer (Intel Centrino Duo T5600 at 1.83 GHz and 1Gb RAM).

Basic example. Deterministic model.
Strictly speaking, decision variables of either models representing the number of sows ought to be integer and non-negative variables.However, given the computational time consumed for calculations in preliminary tests when all decision variables related to sows were considered integer, as well as the pertinent changes in the model's constraints to be included, make the pure integer model inappropriate for practical purposes (Rodriguez et al. 2008).As a consequence, only those decision variables corresponding to the first four periods where declared as integer and the rest as real variables.Beyond this limitation, these four periods represent the roller horizon where decisions must be implemented before new environmental changes could be appreciated or taken into account.Specific parameters of the linear programming model ( 1)-( 25) are detailed in Appendix A. Figures corresponding to the actual instance that was solved are presented in Table 1.The optimum average reward of the farm was 699 thousands of euros.The optimal replacement policy indicates to keep a sow until the end of the 6th cycle.Results provided the scheduling of the purchase and optimal replacement policy week by week.
The occupancy rate of the farrowing facilities was more than 0.95 over the time horizon of planning (Figure 1) and reached the full occupancy before the half time horizon was reached.This shows the rational behaviour of taking the maximum profit of lactation facilities given that piglets represent the main source of income for the farm.Occupancy Rate

Herd Size Rate
The initial (2)-( 5) and final ( 19)-( 22) herd distribution of the 2330 sows along different reproductive states was fixed.These distributions were important because no value resulted in a feasible program.Final herd distribution was selected from a nearer distribution to the ideal steady-state distribution of the herd (Plà et al., 2008) while initial herd distribution was arbitrary selected from those that made solvable the problem.These distributions may affect the herd size over time as shown in Figure 1 and the pattern of the purchase scheduling.

Sensitivity analysis
It is known that variation over time in the dynamic parameters (αt,c,g βt,c,r and γt,c) of a realistic biological system like a pig farm can be very high.Therefore, to prepare the extension of the model into a stochastic linear programming model and to value the impact of the uncertainty of dynamic parameters on model performances two additional cases were considered.The optimistic case, where dynamic parameters were increased a 5%, and the pessimistic case where dynamic parameters were reduced by 5%.
Results concerning the purchase scheduling are shown in Figure 2. The case showing a high productivity requires lower purchases against the case with low productivity and more purchases of gilts needed to maintain lactation facilities near full occupation.The occupancy rate of lactation facilities in the three cases was more than 0.95 along the horizon planning.However, the productivity of this facility varies among cases due to herd structure and composition of this occupancy.For instance, it was observed that the herd with low productivity has a higher rate of sows in gestation state occuping the lactation facility than the others cases.This is a logical result since more pregnant sows and less lactating sows lead to a low production.Comparing the optimal solution achieved in the deterministic model ( 1)-( 25) for the three cases considered, it can be concluded that the uncertainty inherent to the model it is not corresponded with that observed in real systems.As is shown in Table 2 changes of 5% in the dynamic parameters provoke changes of more that the 25% in the maximum revenue.The overall profit ranges from 699 to 888 thousands of euros.
Therefore, it would seem appropriate to extend the model into a new one dealing properly with this uncertainty.

Basic example. Stochastic model.
Stochastic model formulation requires the generation of a set of scenarios S. To illustrate and assess the suitability of the stochastic approach three scenarios were defined in this example.Then, the uncertain parameters (αt,c,g,s βt,c,r,s and γt,c,s) were considered to be modelled by scenario.Therefore, the optimistic, normal and pessimistic scenario were defined in correspondence with the values of high, average and low productivity respectively.Time horizon was of 52 weeks as with the deterministic example and T1={1,2,3,4}.
The resolution of this formulation ( 26)-( 56) give an optimal profit (RP) of 664 thousands of €/year.The results shown lactation facilities occupancy was maximized (see figure 4) as can be observed in Figure 1 for the deterministic model.However in that case the different behaviour for each scenario is observed and reveals the extra effort in the pessimistic scenario to take the maximum profit of lactation facilities.The optimistic scenario instead reaches the maximum occupancy of the lactation facility sooner.Concerning the herd size behaviour (see figure 5) shows how scenarios with high productivity need to maintain a lower size than the rest showing the regulation role of lactation facilities.Furthermore, depending on the initial and final inventory, the optimistic scenario shows a greater capability to reach a steady state sooner.With respect to the scheduling of purchases, again it is shown how the scenario affects the need for a supply of more gilts to the farm in the worst scenarios.Scenarios affecting negatively production require a higher replacement rate of sows which is translated in more gilts being purchased (Figure 6).Furthermore, if purchase scheduling of optimal solutions is compared for the first four periods, the stochastic model shows a better behaviour under practical point of view with lesser variations than the deterministic model (see Figures 2 and 6).This is a direct consequence of the inclusion of uncertainty by scenarios at the second stage of the stochastic model.In addition, just to analyse the importance of time horizon on outcomes for the first 52 weeks different instances for T=78, 104, 130 and 156 were solved (Table 3).It is observed that the time horizon has an very little influence on the first 52 weeks because in all instances the objective function never reports differences greater than a 0,08%.
Even less is the impact on the expected profit for the first stage period (0,02% as maximum).Inspecting the solution of the deterministic models with respect to the stochastic one, with the first 4 weeks as the first-stage, we compute the expected value of perfect information (EVPI), defined through the following expression: being Φ s the optimal value of the deterministic model ( 1)-( 25) when it was solved (separately) for each scenario s in Ω and RP the optimal value of the stochastic model ( 26)-( 56).For our study the EVPI=694-664 = 30 thousands of euros.EVPI measures the value of knowing the future with certainty.This is how much the farmer would be ready to pay this year to obtain perfect information about the population dynamic behaviour.
Additionally, the Value of the Stochastic Solution (VSS) was computed.Roughly speaking, it measures how good or bad results to use the optimal solution of the stochastic model instead of the deterministic one.Then, the Value of the Stochastic Solution is defined as VSS=RP-EEV, where EEV is the expected value assuming expected yields and expected parameters fixing the optimal values at the first stage..In our case, the VSS= 664-647=17 thousands of euros, this is the cost of ignoring uncertainty in choosing a decision.

Conclusions
In this study we formulated a two-stage stochastic programming model with recourse, for planning a piglet production farm.The approach is shown to be suitable to deal with the uncertainty present in the system and provides better solutions than any other of the deterministic models presented for managing the replacement of sows.The stochastic model incorporates the uncertainty associated with the biological system and provides an optimum replacement decision for the immediate time horizon (first stage) which does not depend on each particular scenario considered in the problem.Incorporating integer variables only at the first stage a better representation of the system was achieved and computational problems derived from complex pure integer programming models were avoided.Future uncertainties are taken into account through a finite set of scenarios and considered in the second stage.
It is demonstrated the different response of the herd not only when specific farm parameters are different but also when scenarios are uncertain.This shows the need for analytical tools for consultants and pig specialists to give better advice to their costumers.Thus, the use of stochastic programming in this paper constitutes a valuable tool to support decision-making in the field of pig production.
Finally, we intend to extend this methodology developing a refined set of scenarios and considering more than two stages.Such extensions would require the development of suitable numerical strategies of models with integer recourse variables, which are, in general, more complex than the ones examined in this paper.

Apendix A
Table 4 Daily feed intake of sows according to their physiological state (Kristensen and Sllested, 2004).

State
Feed (FEs) Comments Set and indexes C ={c} Set of number of cycles, Sg = {g} Set of number of gestation week, Sl = {l} Set of number of lactation week, T = {t} Set of periods, Nr ={r} Set of repetitions of insemination, Sr = {k} Set of insemination waiting week, Parameters rpt = price (€/head) offered per piglet at period t,, rut = price (€/head) offered by the slaughter for a replaced sow at period t, rzt,c = cost (€/week) to keep a sow in the insemination state at period t and cycle c, rzrt,c,k = cost (€/week) to keep a sow in the breeding-control state at period t, cycle c and waiting week k, rxt,c,g = cost (€/week) to keep a sow in the gestation state at period t, cycle c and gestation week g, rlt,c,l = cost (€/week) to keep a sow in the lactation state at period t, cycle c and lactation week l, αt,c,g = survival rate of gestation at period t, cycle c and gestation week g, βt,c,r = survival rate of insemination at period t, cycle c and waiting for the insemination attempt r, γt,c = number of piglets at period t and cycle c, cb = number of boxes in breeding facility, cp = number of boxes in pregnancy facility, cf = number of boxes in farrowing facility, z0c = initial stock of animals in insemination state at cycle c, zr0c,r,k = initial stock of animals in control state at cycle c, waiting for the insemination attempt r and week k, x0c,g = initial stock of animals in gestation state at cycle c and gestation week g , y0c,l = initial stock of animals in lactation state at cycle c and lactation week l, zf = number of animals in insemination at the end of the planning horizon, zrfr,k = number of animals in control at the end of the planning horizon waiting for the insemination attempt r and week k, xfg = number of animals in gestation at the end of the planning horizon at cycle c, yfl = number of animals in lactation at the end of the planning horizon at lactation week l, lz = lower bound of gilt purchase, uz = upper bound of gilt purchase, dz = maximum variation between gilt purchase of two consecutive periods, Decision Variables Yt,c,l = number of sows in lactation state at period t, cycle c, lactation week l, Xt,c,g = number of sows in gestation state at period t, cycle c, gestation week g, Zt,c = number of sows in insemination state at period t, cycle c, ZRt,c,r,k = number of sows in control state at period t, cycle c, waiting for the insemination attempt r, at the week k, ULt,c = number of replaced sows at the end of lactation state at period t, cycle c, UZt,c,r = number of replaced sows at the end of the insemination state at period t, cycle c, waiting for the insemination attempt r, ABt,c = number of sows with abortion at period t, cycle c,, The model includes the uncertainty in the dynamic behaviour parameters, future price and cost parameters modeled by means of a given finite set of scenarios.Let consider additional sets and indexes to those given in the previous section: S ={s} finite set of scenarios, T1T subset of T corresponding to the periods of the first stage, and the following notation for those parameters now defined by scenarios: ps = probability for scenario s.rpt,s = price (€/head) offered per piglet at period t and scenario s, rut,s = price (€/head) offered by the slaughter of a replaced sow at period t and scenario s, rzt,c,s = cost (€/week) to keep a sow in the insemination state at period t, cycle c and scenario s, rzrt,c,k,s = cost (€/week) to keep a sow in the breeding-control state at period t, cycle c, waiting week k and scenario s, rxt,c,g,s = cost (€/week) to keep a sow in the gestation state at period t, cycle c, gestation week g and scenario s, rlt,c,l,s = cost (€/week) to keep a sow in the lactation state at period t, cycle c, lactation week l and scenario s, αt,c,g,s = survival rate of gestation at period t, cycle c, gestation week g and scenario s, βt,c,r,s = survival rate of insemination at period t, cycle c, waiting for the insemination attempt r and scenario s, γt,c,s = number of piglets at period t cycle c and scenario s.Associated with each scenario sS, there is a given weight or probability ps.Once the scenarios are settled, all the non-negative decision variables in the model are defined by scenarios according to the following notation: Yt,c,l,s = number of sows in lactation state at period t, cycle c, lactation week l and scenario s, Xt,c,g,s = number of sows in gestation state at period t, cycle c, gestation week g and scenario s, Zt,c,s = number of sows in insemination state at period t, cycle c and scenario s, ZRt,c,r,k,s = number of sows in control state at period t, cycle c, waiting for the insemination attempt r, at the week k and scenario s, ULt,c,s = number of replaced sows at the end of lactation state at period t, cycle c and scenario s, UZt,c,r,s = number of replaced sows at the end of the insemination state at period t, cycle c, waiting for the insemination attempt r and scenario s, ABt,c,s = number of sows with abortion at period t, cycle c and scenario s, planning period T-T1, are called wait-and-see or recourse variables, whose values depend on the corresponding scenario realization that provide the flexibility needed to deal with uncertainty according to the number of periods to be included in the first stage.The resulting model is actually a two-stage stochastic linear program with recourse, whose extended deterministic equivalent program is given in the extensive form by the following optimization problem:

Figure 1
Figure 1 Representation of the behaviour occupancy rate of Farrowing facilities and herd size over time.Deterministic Model (Base).

Figure 4 .
Figure 4. Representation of the behaviour of the occupancy rate of Lactation Facilities with 3 scenarios.

Figure 5 .
Figure 5. Representation of the Herd Size Behaviour regarding three scenarios.

Figure 6 .
Figure 6.Representation of the Scheduling Purchase Behaviour regarding three scenarios.
has been transferred from kg to FEs (1 FE=NE of 1kg barley 13MJ DE).

Table 2 .
Deterministic model: expected profit of the 3 cases in thousands of €.

Table 3
Report of the size of the Stochastic Linear Programming Model Figure 7. Purchase of gilts at first stage.

Table 5
Unitary prices of different concepts related with costs and incomes in sow

Table 6
Main transition probabilities considered in the basic model, (1-(t,c,10)) represents the probability of abortion.Tabla 7 Expected number of piglets weaned per sow per cycle,  t,c,