## Abstract

This paper investigates the early viral dynamics of foot-and-mouth disease (FMD) within infected pigs. Using an existing within-host model, we investigate whether individual variation can be explained by the effect of the initial dose of FMD virus. To do this, we consider the experimental data on the concentration of FMD virus genomes in the blood (viral load). In this experiment, 12 pigs were inoculated with one of three different doses of FMD virus: low; medium; or high. Measurements of the viral load were recorded over a time course of approximately 11 days for every 8 hours. The model is a set of deterministic differential equations with the following variables: viral load; virus in the interstitial space; and the proportion of epithelial cells available for infection, infected and uninfected. The model was fitted to the data for each animal individually and also simultaneously over all animals varying only the initial dose. We show that the general trend in the data can be explained by varying only the initial dose. The higher the initial dose the earlier the development of a detectable viral load.

## 1. Introduction

Foot-and-mouth disease (FMD) is an infectious viral disease (genus: *Aphthovirus*, family: Picornaviridae) that affects cloven-hoofed animals, such as cattle, pigs, sheep and goats. It has an almost worldwide distribution with an endemic and occasionally epidemic occurrence. If an outbreak occurs, it can have a significant economic impact on the areas affected. For example, in the UK, the 2001 FMD outbreak was estimated to have caused a loss of £3.1 billion to the agricultural sector and food chain, with a similar amount lost by the tourism industry (Thompson *et al.* 2002).

The spread of FMD has already been modelled at different levels, in particular farm-to-farm transmissions and within-farm transmissions. A general discussion of these types of models can be found in Woolhouse (2004). One area of FMD dynamics that has received little attention is a mathematical description of the viral dynamics within a host. An understanding of the within-host viral dynamics can be used to aid the study of FMD spreading within a population, providing more accurate information, such as incubation and infectious periods. Quantitative descriptions of FMD viraemia have already been reported (Cottral & Bachrach 1968; Sutmoller & McVicar 1976; Alexandersen *et al.* 2001, 2003*b*; Quan *et al.* 2004), as well as the clearance of virus from the circulation (Sutmoller & McVicar 1976).

Within-host modelling of viral dynamics has previously been directed at human diseases, in particular human immunodeficiency virus (Perelson *et al.* 1996, 1997; Corfec *et al.* 2000; Nowak & May 2000), hepatitis B virus (Herz *et al.* 1996; Nowak *et al.* 1997; Nowak & May 2000; Lewin *et al.* 2001; Perelson 2002; Ciupe *et al.* 2007), hepatitis C virus (Perelson *et al.* 2005), influenza A virus (Baccam *et al.* 2006) and human T-cell leukaemia virus (Wodarz *et al.* 1999). These models have mostly focused on the chronic nature of the diseases, as well as chemotherapeutic interventions. In the models, the life cycle of the viruses is restricted to the circulation or closely associated organs, and there are no divisions between the circulation and organs. There have been no reports on the interaction of virus in the circulation and virus replication within organs. By contrast, pharmacokinetic models of the movement of drugs between different compartments within the body are well established (Riviere 1999).

In this paper, we examine the relationship between the initial dose of FMD virus (FMDV) injected into the blood and the subsequent within-host dynamics. In particular, we question whether individual variation is mostly due to the effect of the initial dose of FMDV. The basis of this study is data from an experiment where 12 pigs were inoculated with one of three different doses of FMD virus: low; medium; or high. Readings of the viral load were recorded over a time course of up to 11 days for every 8 hours. Our model then adds to the simple inspection of, and the statistical analysis of, the raw data in three important ways. First, the model suggests a possible mechanistic explanation for the observed behaviour. Second, it suggests that the wide range of dynamics seen in the experimental data can be largely explained by the differences solely in the initial dose. Third, the model makes testable predictions that can be investigated in subsequent experiments.

The model which we have chosen to use for this study was developed earlier (Quan 2005) to reflect the current knowledge of FMDV biology. We do not consider here alternative and more speculative models of the within-host dynamics, but show that this simple model can be used to explain the patterns shown in the experimental data by varying only the initial dose while also providing estimates of within-host parameters.

## 2. Material and methods

A more detailed account of the experiment is documented in Quan *et al.* (2004) and Quan (2005).

### 2.1. Experiment overview

Twelve female white Landrace pigs, weighing 20–30 kg, were split into three groups of four animals. The first group was housed in two separate boxes, each box containing two animals that were physically separated from each other by double, solid, watertight, wooden partitions (30 cm apart), so that direct contact between the animals was not possible (Alexandersen *et al.* 2002*a*,*b*). The group was inoculated with a low dose of FMDV (10^{6.1} FMDV genomes).

The second group was treated in the same way as the first except for receiving a medium inoculation (10^{7.1} FMDV genomes).

The third group received a high inoculation (10^{8.1} FMDV genomes), and the pigs were all housed in the same box without any separation between animals. This was due to the limitation of space in the isolation unit. Each animal had serum collected every 8 hours for up to 11 days.

### 2.2. Isolation unit

The animals were housed in a secure isolation unit located at the Institute for Animal Health, Pirbright Laboratory. The unit consisted of five boxes in a row connected by a ‘clean’ corridor on one side and a ‘dirty’ corridor on the other. The personnel entrance to each box was from the clean corridor through the double doors with a shower unit between the doors. Before leaving the box, the personnel underwent disinfection with FAM (Evans Vanodine International, Preston, UK), followed by showering. Each box was self-contained with its own high-efficiency particulate air filtration, so that cross contamination was not possible.

### 2.3. Inoculation of animals

The animals were inoculated intravenously in the anterior vena cava/brachiocephalic vein with 0.5 ml diluted inoculum, either 1 in 10, 1 in 100 or 1 in 1000. This was equal to 10^{8.1}, 10^{7.1} or 10^{6.1} FMDV genomes, respectively. The inoculum used was isolated from a pig with FMDV strain O UKG 34/01 at Cheale Meats Abattoir, Brentwood, Essex, on 20 February 2001 (Alexandersen & Donaldson 2002; Alexandersen *et al.* 2002*a*,*b*, 2003*a*,*c*). The inoculum contained 10^{7.2} TCID_{50} in the primary bovine thyroid (BTY) cells or 10^{9.4} FMDV genomes ml^{−1}, and was diluted in HEPES tissue culture medium just before use. After inoculation, the dilutions were tested for genome quantity and infectivity using fluorogenic reverse transcription polymerase chain reaction (RT-PCR) and viral titration in the BTY cells.

### 2.4. Collection of samples

Serum was collected by venipuncture (anterior vena cava/brachiocephalic vein) into plain serum tubes with no additive (Vacutainer; Becton Dickinson, Plymouth, UK). Samples were taken before (to act as controls) and then approximately 5 min after inoculation to estimate the true nature of the initial dose of FMDV. Thereafter, the samples were taken every 8 hours: between 7.00 and 08.00, 15.00 and 16.00, and 23.00 and 24.00.

The virus was taken to be well mixed within the blood after 5 min. Taking the blood volume to be 1/13 of the pig weight (Chang *et al.* 2002) gives a volume of 1.923 l for a 25 kg pig. The mean cardiac output of a typical pig is 4.5 l min^{−1} (Åberg *et al.* 2004). This gives an estimate of (1.923/4.5)×60=25.6 s for the blood to circulate.

### 2.5. Processing of samples

A detailed description of RNA isolation, reverse transcription, PCR and quantitation of the serum samples from the experiment can be found in the previous work (Quan *et al.* 2004).

### 2.6. Model assumptions

The main assumptions in the model are as follows.

— There is a homogeneous mixing within each compartment of the model.

— The rate of exit of virus from each compartment can be approximated as a constant exit rate. This includes all immune-mediated clearance of virus in the blood and virus lost to the epithelial cells. Although the exit rate of virus will vary in practice due to a varying immune response, here we are only modelling the early dynamics of the virus spread.

— In some replicates, there is a clear peak in the viral load. Here, we are more concerned with the dynamics up to this peak than explaining the peak itself. In the model, a constraint on virus replication rates is provided only by reduced availability of cells to infect and not, for example, due to an increasing immune response.

— There is a reduced rate of infection of epithelial cells at low concentrations of virus in the interstitial space,

*Y*. This is consistent with the experimental data from Ohlmann*et al.*(1995), showing a reduced ability of FMDV to produce viral proteins at low concentrations. To achieve this, we use a term that converges to 0 as*Y*tends to 0, and converges to 1 as*Y*increases. We used a standard function in biological modelling, the sigmoidal function*Y*^{2}/(*m*^{2}+*Y*^{2}), where*m*is the value at which the infection rate is at half efficiency. The exponent of this function was set to a constant as it did not contribute to fitting the model. The surplus parameter was avoided to obtain the best fit for this simple model against the small dataset available. The value 2 was chosen as this was the most commonly used value.

### 2.7. The model

The model is a set of differential equations that describe the interactions between the virus concentration in the blood, *X*, and interstitial space, *Y*, and the proportion of available epithelial cells for infection that are either infected, *F*, or uninfected, *U*. The proportion of cells varies from 0 to 1; to obtain the number of cells either infected or not, the variables *U* and *F* must be multiplied by the total number of epithelial cells available for infection. A schematic of the model is shown in figure 1 and the variables *X*, *Y*, *U* and *F* are defined in table 1. The model differential equations arewhere the initial condition *X*(0) is a fitted parameter, *Y*(0)=0, *U*(0)=1 and *F*(0)=0.

#### 2.7.1. Constant values.

The ratio of volume of the blood to interstitial space, parameter *r*, was calculated from the volume of plasma (6.91% body weight) to interstitial space (13.09% body weight) (Pond & Hought 1978). Therefore, parameter *r* was fixed at a value of 0.53 in all our fits.

The lifespan of an infected cell, *l*_{f}^{−1}, was calculated as the approximate time for cell attachment and intracellular infection time of FMDV between 2.5 and 3.5 hours (Cartwright *et al.* 1957). For our fits, we chose to fix *l*_{f} at 0.33, so that *l*_{f}^{−1} is approximately 3.03 hours.

#### 2.7.2. Model dynamics.

Initially, the proportion of infected cells, *F*(0), is zero. The initial concentration of the FMDV genomes in the interstitial space, *Y*(0), is also zero. The initial concentration of the FMDV genomes in the blood, *X*(0), is a non-zero value given by intravenous inoculation, and is treated as a parameter due to the variable nature of the procedure. The virus then spreads as follows.

(i) There is a natural flow in both directions between the concentration of FMDV in the blood,

*X*, and the interstitial space,*Y*, governed by the parameters*k*_{xy},*k*_{yx}and the volume ratio,*r*, between the two compartments.(ii) Once there is a concentration of FMDV genomes in the interstitial space,

*Y*, the proportion of uninfected cells,*U*, changes to the proportion of infected cells,*F*, at a rate given by the parameter*h*_{y}. (The term*Y*^{2}/(*m*^{2}+*Y*^{2}) acts as a delay mechanism, so that the cell infection rate for low values of*Y*is very low. The higher the value of*m*the greater the delay.)(iii) FMDV is lost from the interstitial space,

*Y*, at a constant rate given by the fitted parameter*l*_{y}. This is the only exit point from the model of FMDV in the blood or interstitial space and therefore models*all*removal of viruses from the system. This includes loss into the epithelial cells and any kind of immune-mediated clearance.(iv) The proportion of infected cells,

*F*, bursts at a constant rate*l*_{f}. When the cells burst, FMDV genomes are released back into the interstitial space, increasing the concentration of FMDV by*b*per proportion of total cells. Furthermore, some of the FMDV genomes that are released infect the surrounding uninfected cells at a rate given by the parameter*h*_{f}.

The modelled dynamics are illustrated in figure 2. Sensitivity analysis of each parameter was investigated by scaling the parameter with different factors and plotting the resultant graphs.

### 2.8. Solution of differential equations

The system of differential equations was solved using a fourth-order Runge–Kutta method. An adaptive step-size control was used to ensure the accuracy of the solution, which was necessary due to the rapid changes in the derivative of the solution (Press *et al.* 1986).

### 2.9. Model fitting methods

Parameters were estimated using maximum-likelihood methods, assuming lognormal errors in the data and left-censoring observations below a minimum detection threshold of 10^{2} genomes ml^{−1}. The likelihood for the data, *D*, and parameter set, *P*, is given bywhere *D*_{i} is the observed concentration; *x*_{i} is the expected concentration given by the model; and *φ* and *Φ* are the probability and cumulative density functions for a normal distribution, respectively. The number of data points is given by *n*. The observations are left censored as indicated by *c*_{i}, where *c*_{i} is 0 if *D*_{i}>10^{2} and 1 otherwise. The standard deviation of the measurement error, *σ*, was unfortunately not available, and so was estimated by considering it as a parameter that must be fitted (Stewart *et al.* 1998).

The likelihood surface was sampled using Markov chain Monte Carlo with 200 000 iterations after a 200 000-iteration burn-in period (see Gilks *et al.* 1995). Proposal distributions were given by normal distributions, where the standard deviations were adjusted before the burn-in period to achieve good acceptance rates using a Metropolis–Hastings acceptance rule. The sample with the highest likelihood was taken as the estimate of the maximum likelihood. The Nelder–Mead downhill simplex method (Press *et al.* 1986) was used to search the parameter space to find suitable regions in the parameter space to initialize the chain. The 95% CI of a parameter was taken as the range between the 2.5 and 97.5 percentiles of the estimated cumulative distribution function of that parameter. The model parameter analysis was calculated on a Unix platform using GNU-complied C++ code.

#### 2.9.1. Model fitting.

The data from each animal were fitted to the model individually. The model was also fitted *simultaneously* to all the 12 animals, where the parameters are constant across all animals, except for the initial dose of FMDV inoculated, *X*(0), which was fitted to each animal individually.

The likelihood for the common fit, ℒ_{com}, is therefore given bywhere *D*_{i} and *X*_{i} are the data and the initial dose, *X*_{i}(0), for animal *i*, respectively, and *P*′ denotes the model parameters without the initial dose.

## 3. Results

In our study, we wished to explore the early within-host dynamics of FMDV in pigs using a within-host model. We can make some observations of the data, as shown in table 2, before considering the model.

The reading recorded at 0.1 hours gives an indication of the *effective dose* that each animal received, since this was taken shortly after inoculation. First, consider the low-dose animals, UQ10–UQ13. Animal UQ11 had a smaller reading at time 0.1 hours, which subsequently has a late rise in the viral load. Second, in the medium-dose animals, UQ14–UQ17, animal UQ17 has a lower reading at time 0.1 hours than any other animals, which again shows a late rise in the viral load. Interestingly, this initial value is still higher than all the readings from the low-dose animals. The high-dose animals, UQ18–UQ21, all have high readings at time 0.1 hours.

The early, smaller peak in the data for the animal UQ15 indicates an early round of viral replication and is frequently seen in naturally infected animals. In this case, FMDV replicates locally before disseminating systematically. In this animal, it is possible that epithelial cells were infected directly when the needle was pushed and pulled out of the skin. However, we do not attempt to model this feature with our simple model, which aims to still capture the basic dynamics.

The viral load readings could not be reliably detected for the values below a threshold of 10^{2} FMDV genomes ml^{−1}. This creates some uncertainty in the data, in particular with animals UQ11 and UQ17. The viral load readings for these animals are mostly low; therefore, the dynamics of the virus reproduction are mostly unobservable. The best-fit parameters in terms of the maximum likelihood for UQ11 are perhaps unrealistic, having a large delay in effective virus production and an ineffective clearance of virus after the peak. This fit is due to two factors in the data: first, the large number of values below the threshold of 10^{2}, including the initial value; second, the high later values with no clear peak. The model fit of UQ17 shows a delayed peak, which is likely to be the true behaviour of UQ11, which is supported by the model fitted simultaneously to all datasets.

### 3.1. Model fits

The viral load measurements for each animal (UQ10–UQ21) over the time course of the experiment are reported in table 2, and are also shown graphically in figures 3 and 4.

Unfortunately, there are no data over the full time course of the experiment for every animal. Many animals were killed or died before the end of the experiment. Some animals clearly reached peak viral load: UQ10; UQ14; and UQ15. The data for other animals suggest that the peak viral load may have been reached: UQ13; UQ16; UQ20; and UQ21. The results for the remaining animals are more uncertain.

#### 3.1.1. Model fitted individually.

The best-fit parameter values for each animal are shown in table 3 and the corresponding 95% CIs in figures 5 and 6, calculated as described in §2.9. The resultant viral load curves of the best-fit parameters for each animal are illustrated in figure 3.

#### 3.1.2. Common fit model.

The model was fitted to all the 12 animals simultaneously, where the initial dose of FMDV was fitted to each animal, as described in §2.9. The common maximum-likelihood parameter estimates found are listed in table 3.

The fitted initial values of *X* for each animal are shown in figure 5*b*. The 95% CIs are illustrated with the individually fitted parameters in figure 6 (rightmost value). The viral load curves for the common fit are shown in figure 4. Note that the general shape of the curves is identical only to the timing of the peak differing.

Although there will be differences between the animals (such as weight), the model fitted over all animals simultaneously shows that the overall trend can be represented quite accurately by only varying the initial dose of FMDV (figure 4). Some compromise must be made between the individual model fits in order to obtain a good overall fit to all datasets. Therefore, the common fit of the model for each animal will not be as accurate as the model fitted individually. Nevertheless, the common fit captures the basic dynamics of each animal correctly, which is discussed further within §3.2.

The two model-fitting approaches were compared using Akaike's ‘An Information Criterion’ with the second-order information correction (AICc) (Sugiura 1978; Hurvich & Chih-Ling 1989, 1995), showing that the common fit model is the better approach, because it is using only 20 parameters compared with 96 parameters of the individual fits. This test accounts for the large difference in the number of parameters, as well as considering the small amount of data. (The individual fits have an AICc value of 825.7 and the common fit a value of 640.1.)

### 3.2. Model parameter analysis

In this section, we consider each model parameter in turn: the role they have in the model dynamics and how they differ when fitted to different animals. The manner in which the model dynamics alter when a parameter is scaled is shown in figure 7.

#### 3.2.1. Flow between interstitial and blood compartments, *k*_{xy}, *k*_{yx} (h^{−1}) and *r*.

If the flow of virus is decreased from the blood (*k*_{xy}) to the interstitial space then few viruses are able to infect the cells creating a delay in the virus production. As more viruses are retained in the blood, an increase in the virus in the blood is also observed (figure 7*a*). Increasing the virus flow in the other direction (*k*_{yx}) increases the virus in the blood since few viruses are cleared via the interstitial space. This then increases the infection rate of the cells as the virus re-enters the interstitial space (figure 7*b*).

The model may be simplified since the variables for the blood, *X*, and the interstitial space, *Y*, are in equilibrium with one another, related by the equation , although the results presented here refer to the full model (see appendix A). The sensitivity of changing the ratio between the blood and the interstitial space, given by constant *r*, is shown in figure 7*c*. Parameters *k*_{xy} and *k*_{yx}, which are fitted, can account for any variability in *r*. Therefore, we are not too concerned with the individual values of *k*_{xy}, *k*_{yx} and *r*.

#### 3.2.2. Exit rate of FMDV genomes, *l*_{y} (h^{−1}).

The parameter *l*_{y} controls the exit of virus from the interstitial space and therefore also controls the exit of virus from the blood. This parameter has two main roles: first to create the initial dip in the viral load, and second to create a decrease in the viral load after the peak.

As the value of *l*_{y} is increased, the virus in the interstitial space is cleared quicker delaying the production of virus within the proportion of cells. The magnitude of the peak of the viral load also decreases since the rate at which virus is produced within the proportion of cells is met sooner by the clearance rate (figure 7*d*). For all of the individual fits, except UQ11 and UQ19, the upper bound of the 95% CI is unbounded. This is due to the exit rate to the blood, *k*_{yx}, also being a variable and can be set as high as required, so that all of the viruses are not lost through clearance. Unfortunately, this means that only a lower bound is found for *l*_{y} for these animals, showing that an almost instant clearance of virus from the interstitial space cannot be ruled out. The exit rate for the model fit of UQ11 is virtually zero, shown by the flat peak. The model fit of UQ19 relies less on virus production via the interstitial space and is therefore limited in magnitude.

#### 3.2.3. Interstitial space infection rate, *h*_{y} (genomes^{−1} ml h^{−1}).

The parameter *h*_{y} controls the rate at which the proportion of cells becomes infected by the free virion population in the interstitial space. As the value of *h*_{y} is decreased the production of virus is delayed (figure 7*e*), and if it were set to zero the virus would quickly disappear. Since the virus produced from the proportion of cells is returned to the interstitial space, this creates a positive feedback loop and the virus increases exponentially in the proportion of cells and the interstitial space. Lower values of *h*_{y} create a larger initial dip and delay the onset of virus production (as mentioned above for the animal UQ14).

We can calculate typical values for the interstitial space infection rate of the proportion of available epithelial cells for infection, given by the expression *h*_{y}(*Y*^{2}/(*m*^{2} + *Y*^{2}))*YU*, using the common fit model 24 hours before the peak in the viral load, at the peak and 24 hours after the peak, as shown in table 3 (animal UQ10 was used with a peak time of 138.64 hours). Since these values are for the proportion of available epithelial cells for infection, these values must be multiplied by the total number of available epithelial cells for infection to obtain values of the number of individual cells infected per hour.

#### 3.2.4. Cell-to-cell infection, *h*_{f} (cell^{−1}).

The parameter *h*_{f} controls the rate at which the proportion of infected cells infects the proportion of available uninfected cells. Increasing *h*_{f} causes the proportion of available uninfected cells to become infected quicker and so exhausting the supply of cells to infect sooner (figure 3*f*). If this parameter is set to zero, then cell-to-cell infection cannot occur. Therefore, the rate of infection of the proportion of cells by free virions in the interstitial space (*h*_{y}) must be sufficiently high to ensure that enough of the cells become infected to increase virus in the system. If *h*_{y} is too low, then the virus exits the system via the interstitial space quicker than it can be produced and so decreases to zero. The proportion of infected cells, *F*, increases exponentially provided that cells first become infected from the interstitial space, and, second, the proportion of available cells for infection does not decrease too low.

As above, we calculate typical values for the cell-to-cell infection rate of the proportion of available epithelial cells for infection, given by the expression, *h*_{f} *l*_{f} *FU*, as shown in table 3. The overall infection rate given by the first two terms of d*F*/d*t*, (*h*_{y}(*Y*^{2}/(*m*^{2} + *Y*^{2}))*YU* + *h*_{f}*l*_{f}*FU*), is given in table 4. This shows that infection from the interstitial space contributes 50–60 per cent.

#### 3.2.5. Reduced infectious ability, *m* (genomes ml^{−1}).

The parameter *m* represents a virus level below which there is a much reduced infectious ability, delaying the production of virus in the cell system (figure 7*g*). The maximum value of the 95% CIs of *m* for each of the high-dose animals is similar to one another. This is due to all of these datasets having few or no values below the detection threshold of 10^{2} genomes ml^{−1} (figure 6). Since virus levels in the high-dose animals are high, the reduced infectious term, *Y*^{2}/(*m*^{2} + *Y*^{2}), has less of an effect than it does for lower dose animals with lower virus levels. This is why higher doses have a larger confidence interval for *m*.

The common fit model shows that *m* is a necessary parameter, and hence justifies the inclusion of the term *Y*^{2}/(*m*^{2} + *Y*^{2}), having a lower 95% CI greater than zero. The value of *m* for the common fit model enables the correct dynamics to be displayed for all animals by only varying the initial condition. The parameter, *m* (genomes ml^{−1}), is estimated by the 95% CI for the common fit model: (2.11×10^{−1}, 1.68×10^{2}).

#### 3.2.6. Rate infected cells die, *l*_{f} (h^{−1}).

The sensitivity to the constant *l*_{f} apparent in figure 7*h* is considerable, but these graphs correspond to cell attachment and intracellular infection times of 24.24, 12.12, 6.06, 1.52, 0.76 and 0.38 hours for (1/8)*l*_{f}, (1/4)*l*_{f}, (1/2)*l*_{f}, 2*l*_{f}, 4*l*_{f} and 8*l*_{f}, respectively. These times are outside the calculated approximate time of between 2.5 and 3.5 hours.

#### 3.2.7. Released FMDV genomes, *b* (genomes ml^{−1}cell^{−1}).

When the amount of virus released per proportion of cells is increased, the infection rate of the proportion of cells also increases, consequently, exhausting the available cells for infection more quickly (figure 7*i*).

#### 3.2.8. Initial dose, *X*(0) (genomes ml^{−1}).

Decreasing the initial dose delays the production of virus, causing a viral load peak of the same magnitude later (figure 7*j*).

## 4. Discussion

The parameter *X*(0) controls the initial quantity of virus in the system. A higher initial concentration of FMDV in the blood results in an earlier peak of virus in both the data and in the model. The model peak always has the same dynamics for any value of *X*(0), which can be seen clearly from the common fit graphs in figures 4 and 7*d*. From the 95% CI plots of *X*(0) from the individual fits, as shown in figure 5*a*, we can see that this fitted parameter is in good agreement with the initial dose given to each animal. In particular, it reflects that animals UQ10–UQ13 were given a low dose, animals UQ14–UQ17 a medium dose and animals UQ18–UQ21 a high dose. This is the only parameter that seems to reflect the fact that different animals were given different doses. From the 95% CI plots of *X*(0) from the common fit model, as shown in figure 5*b*, we also see the relationship with the initial dose. Although there is no clear difference between the low- and medium-dose animals, the timing of the peaks is similar.

From the graphs in figure 4, it appears that the model predicted peak values of *X* are identical for each animal. When the peak value is calculated for each animal, we see that all the peak values are equal at 1.79×10^{8} FMDV genomes ml^{−1}. Unfortunately, it is not possible to write down an analytical solution for the peak value of *X* in terms of the given parameter set, since we must solve the system of differential equations numerically. However, in all the analyses performed in producing the best-fit common parameters, the peak value has always been the same across all animals.

In conclusion, although there is some variation in the dynamics of virus spread between animals, there is a general pattern that can be explained in terms of the initial dose of FMDV in the blood. The model predicts that, when the initial dose is varied, only the timing of the virus spreading is varied, with the dynamics remaining the same: *the higher the initial dose the earlier the development of detectable viral load*. This was shown by fitting the model to all 12 animal datasets only allowing the initial dose to vary between animals, as shown in figure 4.

The use of a model was important in this exercise providing estimates for missing data and providing a mechanistic explanation for the observed data. Furthermore, estimates for the parameters governing the within-host dynamics of FMDV were obtained. Some of the authors are currently involved in the ongoing experimental work to estimate the infectious and incubation periods of animals with FMDV, together with the implications for vaccines, an important component of this is the study of within-host dynamics.

## Acknowledgements

Animal experimentation was approved by the Institute for Animal Health (IAH) ethical review board under the authority of a Home Office project licence in accordance to the Home Office Guidance on the Operation of the Animals (Scientific Procedures) Act of 1986 and associated guidelines.We would like to thank the BBRSC for funding this research.

## Appendix A. Central compartment and interstitial space

In the model, the blood, *X*, and interstitial space, *Y*, are connected by a constant flow from *X* to *Y* with rate *k*_{xy}. Similarly, there is a constant flow from *Y* to *X* with rate *k*_{yx}, where *r* acts as a conversion factor for the concentration of FMDV between the different-sized spaces. Considering only the interactions between *X* and *Y* we have the following:These equations have a set of steady states (i.e. *X* and *Y* remain constant when d*x*/d*t* and d*y*/d*t* equals zero) given by the line . It is easily shown that this set of steady states is *stable*. The equations have the following solution:Therefore, as *t*→∞, and . Furthermore, since , we know that *X* and Y converge exponentially fast towards the stationary point (which is on the line *k*_{yx}*Y*=*rk*_{xy}*X*) along the line *Y*=*r*(*X*(0)−*X*)+*Y*(0).

When we consider *X* and *Y* in the model, we observe that *Y* has extra interactions that increase and decrease the value of *Y*. However, the interactions of *Y* with *X* remain the same; therefore, when *Y* changes, *X* converges exponentially fast to the corresponding value on the line . This can be observed in figure 2, where the values of *X* and *Y* track one another.

From these observations of *X* and *Y*, we can conclude that it would be sufficient to model *X* and *Y* with only one variable. For example, we can drop the variable *X* from the model and achieve the same behaviour. The initial values of *Y* are then given by the original model initial values of *X* as .

## Footnotes

- Received October 3, 2008.
- Accepted October 24, 2008.

- © 2008 The Royal Society