Sensitivity analysis of a stochastic simulation model for foot and mouth disease

A spatial and temporal Monte-Carlo simulation model was developed to analyse the epidemiology and control of foot and mouth disease (FMD). Animal, people and vehicle contacts as well as airborne and local spread represented the FMD virus transmission between farms housing cattle, pigs or sheep. Contacts were explicitly modelled by routes, airborne transmission by the Gaussian Dispersion model and local spread by distance dependent transmission probabilities. Control measures were implemented according to the EU Directive (2003/85/EG). A sensitivity analysis with a two-level fractional factorial design was used to examine the robustness of the simulation model to extreme input values. The influence of eleven input parameters and interactions between them were estimated: ability of airborne spread, duration of the incubation period, time from infection until infectivity, time from onset of clinical signs until diagnosis, farm density, type of index case, number of farms visited per route, visiting interval, type of the animal sales, control strategy, and delay until start of control strategies. The considered parameters as well as certain two-factor interactions between them showed a significant impact on the epidemic duration and the number of infected and culled farms. Particularly, the parameter airborne spread, farm density, number of farms visited per route and control strategy influenced the course of the epidemic. The consideration of airborne spread as well as the implementation of contacts between farms with routes allowed a detailed analysis of these transmission paths.


Introduction
Foot and mouth disease (FMD) is a highly contagious OIE (Office International des Epizooties) listed disease that can spread rapidly and widely among all cloven-hoofed animals (DONALDSON et al. 2001, NISSEN et al. 2003, KITCHING et al. 2005).Simulation models are a feasible mean to gain better insight into the processes of the epidemiology of FMD virus and to evaluate control strategies against FMD, particularly for countries that are free from the disease (GREEN and MEDLEY 2002, GARNER et al. 2007, KARSTEN et al. 2007).Epidemiological models for FMD have been developed for different countries (FERGUSON et al. 2001, BATES et al. 2003a, PEREZ et al. 2004, MENACH et al. 2005, CARPENTER et al. 2007, GARNER et al. 2007, HARVEY et al. 2007, KOBAYASHI et al. 2007a, WARD et al. 2009).However, it is difficult or even wrong to extrapolate experience of one country to another (SCHOENBAUM andDISNEY 2003, DUBÉ et al. 2007).Farming practices, farm densities, farm sizes, and contact patterns are different.
In the present study, a stochastic simulation model with special focus on the simultaneous estimation of the airborne spread and the modelling of contacts between farms is presented.People and vehicles moving from one farm to another were modelled as a routing problem.A sensitivity analysis was performed to test the models' robustness to extreme values of the input parameter and to identify risk factors that are to be included in further applications of the model.

Material and methods
The concept of the stochastic simulation model was based on a spatial and temporal Monte-Carlo approach.The modular structure of the program was implemented in the object oriented language C++.Routines from the NAG C library (NAG 2002) were used to generate random numbers.

Description of the simulation model
An individual-based approach was chosen to model the spread of FMD among farms housing pig, cattle or sheep.A farm consisted of one or several units that represent different livestock (dairy cows, beef cattle, suckling cows, farrowing herd, finishing herd, sheep) and were treated separately to include livestock-dependent susceptibility to infection (DONALDSON 1979), virus emission (SØRENSEN et al. 2000), and contact structures (DICKEY et al. 2008, THORNLEY andFRANCE 2009).Farms were located using x-and y-coordinates to take the influence of the farms' spatial arrangement within a region into account that affects the character and extent of disease spread as well as the effectiveness of spatially-targeted control strategies (CHOWELL et al. 2006, BECKETT andGARNER 2007).The farm data in this analysis was derived from the Farm Structure Survey 2003 (FORSCHUNGSDATENZENTRUM 2005).The region was located in Northern Germany and consisted of 729 farms: 462 dairy cow, 412 beef, 165 suckling cow, 49 farrowing, 59 finishing and 79 sheep units.
The FMD epidemic started with the infection of all animals at the index-case farm which was set by the user of the simulation program.The course of the disease was modelled for each livestock unit at farm level by subsequently moving through a predefined order of infection states, one at a time: ›susceptible‹, ›infected‹ (carries virus but unable to transmit virus), ›infectious‹ (carries virus and able to transmit virus), ›clinical signs‹ (first animals on a farm showing signs), ›diagnosed‹ (confirmed presence of FMD virus on a farm) and ›culled‹.The flowchart (Figure 1) illustrates the development of the disease from infection until culling and shows in which infection states farms could be infected or could be a source of infection for other farms, respectively.
According to other FMD simulation models, the transmission paths of FMD were classified as airborne spread, direct and indirect contacts as well as local spread.Consequently, the spread of the FMD virus was estimated simultaneously.
Airborne spread was modelled using a Gaussian Dispersion model, which calculated the virus concentration on farms located downwind up to a distance of 10 km (SØRENSEN et al. 2000).As input parameters for the Gaussian Dispersion model, livestock-specific virus emission was considered.Infection of animals on the farms located downwind was assessed using the concept of minimum infectious doses (FRENCH et al. 2002).Exposure to doses below that quantity will not, whereas doses above this will, with increasingly probability, result in infection and disease (SØRENSEN et al. 2000, FRENCH et al. 2002).For the parameterization of livestock-specific virus emission rates and minimum infectious doses quantities published and applied by SØRENSEN et al. (2000) and ( 2001) were used.

 
The modelling of direct and indirect contacts between farm units was represented by animal, personnel and vehicle contacts.As direct contacts, three types of animal distribution were modelled: integrated production system (dairy or farrowing unit on same farm as the accordant beef or finishing unit), direct contract sale (distribution of animals to constant purchasers) and trade (random allocation by an animal hauler).The three most important indirect contact types, the bulk milk tanker, the veterinarian and the animal feed vehicle, were implemented.These people or vehicles visited farms in a certain order (called route), consequently only subsequent farms of an infectious farm were at risk of being infected.As a reasonable approximation of these routes, the farms within a route were arranged so that the travelling distance was minimised (routing problem).This was the most common optimisation criterion under economical and ecological aspects.Routing problems were NP-hard (NP: Non-Deterministic Polynomial) and solving algorithms were requested to give optimal solutions in a reasonable amount of time.Exact solving algorithms were usable only for a small numbers of farms (n) because the shortest route had to be extracted from the immense number of n! possible solutions.Therefore, the heuristic of the simple and fast Nearest-Neighbour Algorithm was used to give an approximate solution (GRÜNERT and IRNICH 2005).The procedure of the algorithm was illustrated in an example (Figure 2).
The route started and ended at the depot (D).Starting the algorithm, the distance between the depot and all units was calculated and the unit with the shortest distance was chosen to be the first stop (Figure 2a).Then, all the distances from this unit to all others were calculated and again the closest one was added to the route (Figure 2b).This procedure continued until the last unit had been added (Figure 2c) and the route was completed with a direct return to the depot (Figure 2d).Local spread summarized the spread of disease when no clear linkage other than geographical proximity was found (YOON et al. 2006).It was estimated by distance-depending transmission probabilities and restricted to a maximum distance of 1 km.Different control measures were implemented based on the EU Directive (2003/85/EG) (ANONYMOUS 2003) and German Regulations (ANONYMOUS 2004): the immediate culling of all animals of an infected farm, the establishment of protection and surveillance zones on all farms surrounding an infected farm (0-3 and >3-10 km, respectively), contact tracing (forward and backward during last 3 weeks), preventive culling and emergency vaccination (both in variable circles).All control strategies could be established after a user-defined time of delay, modelling the time needed for arrangements before the control strategies were fully implemented.

Sensitivity analysis of the simulation model
To test the simulation model for its robustness model to extreme values of the input parameters KLEIJNEN (1995a) and(1995b) recommended a sensitivity analysis.The present model has a large number of input parameters.KLEIJNEN (2009) described that simulation experiments are well-suited to the application of sequential designs instead of inferior 'one shot' designs (change one parameter at a time).(Fractional) factorial designs resulted in smaller variances (higher accuracy) and furthermore interactions between the input parameters were to be estimated (KLEIJNEN 1995b).Qualitative expert knowledge existed to appraise the outcome of the simulation model for the sensitivity analysis (KLEIJNEN 1995b).Following the study of KARSTEN et al. (2005), a sensitivity analysis with a two-level fractional factorial design was applied (BOX and HUNTER 1961a, b).The results were independently evaluated by five German epidemiologists whose expertise was the research on animal diseases in theory and field.
The sensitivity analysis was divided into three two-level fractional factorial designs (Model I, II and III) to accomplish the effort of a full factorial design, but to include the estimation of interactions.In each model, the basic parameters airborne spread, incubation period, time from infection until infectiousness and time from onset of clinical signs until diagnosis were included representing the variations in the epidemiology of the FMD virus.Moreover, each design includes two or three additional parameters to describe specific environmental circumstances.In Model I farm density and type of index-case farm were considered, in Model II farms per route, visiting interval, and proportion of trade to contract sale and in Model III control strategy and delay until its start.Consequently, using the three different fractional factorial designs the number of required scenarios was reduced from 2 048 (2 11 ) to 96 (2 6-1 +2 7-2 +2 6-1 ) and each scenario was repeated 250 times.
Table 1 shows the values for the two levels (lower and upper) for each parameter.The upper level of each parameter was defined to be the higher value.Airborne spread was included or not (yes, no).The assumed farm densities described mean densities in Northern Germany (0.67/1.05 farms/km 2 ).As index farms, a beef and a dairy farm were chosen as beef farms had fewer and different contacts of animals, vehicles and people than dairy farms.

Statistical analysis
To analyse the dynamical development of the FMD outbreak, epidemic curves of six exemplarily chosen scenarios are shown.The total size of the simulated FMD epidemics was described by three response variables: the mean epidemic duration in days, the mean number of infected farms and the mean number of culled farms.A linear model with a negative binomial distribution and a log-link function provided the best fit to the data.Estimation of the main factors as well as the two-factor interactions is carried out using the GLIMMIX procedure of the SAS program package (SAS 2006).

Results
In Figure 3 six epidemic curves are presented which show the typical positively skewed distribution of epidemic curves.The much smaller outbreak curves A and B are generated when no airborne transmission is assumed.In scenarios A, B and C the control strategies start with only one day delay, in the others with three days.The one day delay curves have smaller peaks and tails, so that the epidemics are eradicated within a shorter period of time.The peak of the curve E is smaller but later in time than that of curve F. Additional preventive culling in scenario E decreases the number of infected farms and slows down the development of the epidemic compared to F.

Basic parameters
The results for the three basic parameters airborne spread, incubation period and time from the onset of clinical signs until diagnosis are identical for all three models: they show a highly significant influence on the response variables (Table 2).The contrast coefficients (difference of upper level and lower level) are positive (Table 3) indicating that the upper level of each of these parameters increases the mean epidemic duration as well as the mean number of infected and culled farms.The contrast coefficients for airborne spread are much higher than for the other parameters.The basic parameter 'time from infection until infectivity' shows negative contrast coefficients.The longer the time period until infectivity, the fewer farms are infected or culled, respectively.Although contrast coefficients for the epidemic duration are only significant for Model I, tendencies are identical in all three models.<.0001 <.0001 <.0001 <.0001 <.0001 <.0001 <.0001  <.0001 <.0001 <.0001 <.0001 <.0001 0.0003 0.0014

Model I
Both, the farm density and the type of index-case, show a highly significant influence on the response variables.If an epidemic starts on a dairy farm or in a region with a high farm density, more farms are infected and the epidemic lasts longer.In Table 3, the coefficients for the number of infected farms (1.09 and 3.7) are much higher than for the epidemic duration (0.08 and 0.31).So the epidemic duration has less variation then the number of infected farms.A significant interaction between the airborne spread and the farm density is observed.For both levels of airborne spread, more farms are infected if a high farm density is present.But the difference in the number of infected farms between the lower and the higher farm density is smaller if the virus can not spread via air than if airborne spread is possible.In a region with a high farm density, airborne FMD virus transmission has more severe consequences than in an area with a low density.Figure 4a shows that the farm density significantly interacts with the index type as well.In low density regions an epidemic results in more infected farms if it starts on a beef farm and not on a dairy farm.In a region with a high farm density the reverse is the case.

Model II
In Model II the number of visited farms per route shows a significant influence on the mean number of infected farms.the route is longer, more farms are infected (+0.3 of the LSMs on log scale, Table 3).The parameters ›visiting interval‹ and ›proportion of contract sales to trade‹ show the same tendencies.The more frequent the farms are visited or the more animals are sold via trade, the higher are the epidemic duration and the mean number of infected farms.Moreover, the number of farms per route shows a significant interaction with the airborne spread (Figure 4b).
If the routes comprise of many farms, more farms are infected.This influence is more important if airborne spread is present than if no airborne spread occurs.

Model III
The delay of the establishment of control measures and the control strategy are highly significant factors (Table 2).The longer the time period until the control measures are established, the higher are the mean epidemic duration and the mean number of infected and culled farms.The control strategy with additional preventive culling and contact tracing results in shorter epidemics and fewer infected farms, but a higher number of culled farms (Table 3).The control strategy shows a significant interaction with airborne spread.Airborne transmission decreases the difference in the response variables between the two control strategies (not shown).Furthermore, the delay of establishment interacts significantly with the control strategy for all response variables (Table 2).More farms are infected if only a protection and a surveillance zone are implemented.A three-day delay has a negative influence on the success of both control strategies.Contact tracing and preventive culling causes a reduction in the number of infected (Figure 4c) and culled farms (Figure 4d).

Discussion
The present experimental design comprises three statistical models what implicates that interactions between the additional parameters in each model can not be estimated.However, to evaluate the models' robustness to changing input parameters, the design is appropriate.The definition of two levels (lower and upper) for each input parameter allows the estimation of the models reaction to extreme values of the parameters.An additional level, e.g.medium, obtains the opportunity of a more detailed analysis, particularly for the interactions.A factional factorial design for three parameters is feasible but often very complex and includes many confounded parameters.

Basic parameters
In general, the three models show conform results regarding the basic parameters.Conform to other recent simulation studies (GARNER andBECKETT 2005, HESS et al. 2008), airborne spread is the most important influence.However, due to the fact the airborne spread is uncontrollable and the definition of the two levels of airborne spread (possible or not) a high difference in the results was expected.This is realistic as also in real epidemics the FMD virus types differ in their ability of airborne spread.While the strain of the epidemic 1967/68 in the UK could easily be spread via air, strain from 2001 in the UK could not.Moreover, the airborne spread is only possible under certain weather conditions so both scenarios have to be considered.Furthermore, the significant interaction between airborne spread and the control strategy indicates that the optimal control strategy depends on the ability of the FMD virus to be transmitted by air.If the weather conditions allow the transmission of FMD virus by air, the implementation of a control strategy including preventive culling and contact tracing is to be applied to induce a rapid and long range elimination of susceptible animals.

Model I
The results of the present analysis confirm findings of recent studies that farm density has a significant influence on the epidemic duration as well as on the number of infected and culled farms (MORRIS et al. 2001, GERBIER et al. 2002, RIVAS et al. 2003).In areas with a high farm density the distance-dependent virus transmission such as local and in particular airborne spread becomes more important as more farms are within the grasp of an infected farm.As the ability of FMD virus to be transmitted via air varies between different FMD virus types, a fast identification of the virus type in a real epidemic is essential for good predictions of further virus spread and decisions which control strategies are to be established.The results for the parameter ›type of the index-case farm‹ are conform to other simulation results and analysis from recent outbreaks (BATES et al. 2003b, WARD and PEREZ 2004, KOBAYASHI et al. 2007a, b, DICKEY et al. 2008).Dairy herds as index cases result in larger epidemics than beef farms.Dairy farms have more contacts concerning the sale of animals or the daily visit of the bulk milk tanker.KOBAYASHI et al. (2007a, b) concluded that dairy herds are to be given preferential treatment in allocating control measures.Especially in densely populated areas the allocation of control measures to dairy herds is even more important than in sparsely populated regions.

Model II
As expected, the more farms are visited per route, the higher is the number of farms infected.On longer routes an infectious farm often has more subsequent farms susceptible to infection.This influence reduces if the FMD virus can be transmitted via air.Many farms already become infected by airborne transmission so that the number of farms per route is less important.The visiting interval and the proportion of animal sales via trade or contract selling have no significant influence.A possible explanation is that assuming visiting intervals of 14 or 21 days an infected farm sells animals only once during the infectious period.Consequently, it has no influence if these animals are sold randomly by trade or directly by contract selling.On the contrary, farms that always receive their animals from one farm have only this single farm as infection source, while farms with changing suppliers have multiple chances of infection via animal purchasing.For a better understanding, further analysis should include smaller visiting intervals of e.g.seven days and larger differences for the ratio trade to contract selling.In shorter visiting intervals, multiple visits during the infectious period are less likely.

Model III
As in studies by BATES et al. (2003b) andHARVEY et al. (2007), the control strategy has a significant impact on course of the disease.The application of preventive culling and contact tracing reveals a faster eradication of the epidemic but at the expense of a high number of culled farms.With every day of delay until control measures are put in place, more farms are infected and preventive culling becomes more effective in reducing the epidemic size, but again at the expense of the high number of culled farms.Studies by YOON et al. (2006) on the 2002 FMD outbreak in Korea are in line with the present findings.The additional preventive culling and contact tracing reduce the number of susceptible animals within the region and allow a faster and more efficient identification of infected farms.In this analysis, the effect of preventive culling and contact tracing is not separated, but it can be suggested that both have a significant influence (KISS et al. 2005, KOBAYASHI et al. 2007a, b).However, it is shown that additional preventive culling and contact tracing are most effective if they are established immediately.Furthermore, if airborne spread is possible, preventive culling is favourable.Without airborne spread, movement restrictions and contact tracing can be sufficient taking the economical and ethical consequences into account.
In conclusion, the models outcome is conform to expectations.In further studies, the simultaneous estimation of FMD virus transmission via airborne spread is to be taken into account.It is the most important basic parameter and shows significant interactions with the parameters farm density, farms per route and control strategy.The modelling of contacts between farms by routes is an alternative to distance-dependent estimation and allows a more detailed analysis of this transmission path.

Figure 2
Figure 2 Example of using the Nearest-Neighbour Algorithm to find the shortest route from one depot (D) through 9 farms presented in Figures a-d Beispielanwendung des Nearest-Neighbour Algorithmus zur Bestimmung der kürzesten Route mit neun Betrieben ausgehend von einem Depot (D)

Figure 4
Figure 4 Selected two-factor interactions (LS-Means on log scale) for the number of infected farms, (a) index × density in Model I, (b) air × route in Model II, (c) delay × control in Model III, and for the number of culled farms (d) delay × control in Model III Ausgewählte zweifach Interaktionen aus Modell I, II und III

Table 1
Parameters and their levels in the fractional factorial designs Parameter mit den jeweiligen Faktorstufen des fraktioniert faktoriellen Designs PZ protection zone (0-3 km), SZ surveillance zone (>3-10 km), PC preventive culling (<1 km), CT contact tracing, a ANONYMOUS 2003),b BATES et al. 2003a, c CARPENTER et al. 2004, d FERGUSON et al. 2001, e GLOSTER  et al. 2003, f HUGH-JONES and WRIGHT 1970, g KEELING et al. 2001, h KOBAYASHI et al. 2007a, i NIELEN et al. 1996,  j PEREZ et al. 2004, k SCHOENBAUM and DISNEY 2003, l THORNLEY and FRANCE 2009, m YOON et al. 2006 a C cattle, P pig, S sheep, V veterinarian, A animal feed vehicle, M bulk milk tanker, D dairy farm, F farrowing farm,

Table 2 P
-values of the main factors and two-factor interactions (Type III likelihood ratio statistics) for the mean epidemic duration in days (Dur) and mean number of infected (Infec) and culled (Cull) farms P-Werte der Hauptfaktoren und signifikanter zweifach Interaktionen für die mittlere Epidemiedauer und die mittlere Anzahl infizierter und gekeulter Betriebe Sensitivity analysis of a stochastic simulation model for foot and mouth disease