SENSITIVITY AND UNCERTAINTY ANALYSES FOR A SARS MODEL WITH TIME-VARYING INPUTS AND OUTPUTS

This paper presents a statistical study of a deterministic model for the transmission dynamics and control of severe acute respiratory syndrome (SARS). The effect of the model parameters on the dynamics of the disease is analyzed using sensitivity and uncertainty analyses. The response (or output) of interest is the control reproduction number, which is an epidemiological threshold governing the persistence or elimination of SARS in a given population. The compartmental model includes parameters associated with control measures such as quarantine and isolation of asymptomatic and symptomatic individuals. One feature of our analysis is the incorporation of time-dependent functions into the model to reflect the progressive refinement of these SARS control measures over time. Consequently, the model contains continuous time-varying inputs and outputs. In this setting, sensitivity and uncertainty analytical techniques are used in order to analyze the impact of the uncertainty in the parameter estimates on the results obtained and to determine which parameters have the largest impact on driving the disease dynamics.


Introduction.
A novel disease known as severe acute respiratory syndrome (SARS) was reported in March, 2003 by the World Health Organization (WHO) [31].This respiratory disease, caused by a previously unknown coronavirus, SARS-CoV [6,14,18,22,32], spread to 32 countries and regions, causing 774 deaths and 8,096 infections globally [30].The high mortality and morbidity associated with the virus, coupled with its rapid global spread, prompted a coordinated global effort spearheaded by WHO and aimed at determining effective control strategies for 528 R. MCLEOD, J. BREWSTER, A. GUMEL, AND D. SLONOWSKY combatting the spread of SARS.Thankfully, these efforts resulted in the elimination of all cases of SARS by August 2003 [7,30].
Numerous mathematical models have been designed and used to evaluate control strategies against SARS in some of the SARS-stricken regions (see, for instance, [5,7,9,11,12,13,21,28,29,33]).These models, generally in the form of systems of deterministic or stochastic differential and difference equations, contain many parameters.The study of the dynamics of such a new emerging disease presents important challenges not only because of the uncertainty in the estimates of the model parameters, but also because control measures are gradually refined as more data (and further knowledge about the disease) become available.
This paper, which focuses on the mathematical model presented in [7], aims to investigate the effect of the uncertainty associated with the parameter estimates used in [7].In this work we (i) extend the model in [7] by incorporating continuous time-dependent functions into the model to reflect the gradual refinement of control measures over time and (ii) apply sensitivity and uncertainty analysis (SUA) techniques to the resulting model containing continuous inputs and outputs.
It is known that the behaviour of nonlinear and multidimensional mathematical models can be explored using SUA (see [23] and the references therein).While an uncertainty analysis quantifies the variability in the outcome of the model attributable to the uncertainty in the values of the associated input parameters, a sensitivity analysis enables the determination of the most influential parameters driving a model.
In the context of disease epidemiology, SUA has been used (notably by Blower and co-workers [2,3,4,24]) to gain insights into the transmission and control dynamics of many human diseases such as HIV and tuberculosis.For example, in [2], SUA was used to investigate a model for the evolution of drug-resistant HIV in San Francisco.This model featured two continuous time-varying outcome variables (prevalence and transmission of drug resistance).Our paper incorporates continuous time-varying input functions (Section 3) into the underlying mathematical model in [7], and we use SUA techniques similar to those in [2,3,4,24] to assess the epidemiological significance of the input parameters, including those associated with the new functions that we introduce to model the gradual refinement of anti-SARS control measures.
In this paper, SUA techniques are used to identify the key epidemiological parameters affecting the dynamics (persistence or elimination) of SARS.Owing to the uncertainty in parameter estimation associated with a new emerging disease such as SARS (where some of the key biological and epidemiological features are not precisely known at the beginning of the epidemic), a detailed SUA is imperative to ascertain the effect of such uncertainties on the results obtained from associated mathematical models.
Unlike most epidemiological models, the modified model considered here involves time-varying, continuous, control parameters in addition to a time-dependent response.The time-varying control parameters are used to model the progressive refinement of quarantine and isolation measures for SARS-infected individuals.Furthermore, the model includes a time-varying component representing the gradual implementation of hygienic precautions by health workers and others in close contact with infected individuals in isolation.The existence of a time-dependent response implies that functional data is obtained as output.Functional data arise in many different fields, with the characteristic feature that the output consists of (typically smooth) curves [19,20].In the SUA literature, however, little use is made of the functional analytic concepts described in [19] and [20].Although the analytical techniques used in this paper do not fully exploit the functional nature of the data, they do provide a starting point and the motivation for future research on the development of SUA techniques for epidemiological models with functional output.
The paper is organized as follows.The underlying SARS transmission model is briefly described in Section 2. Section 3 presents time-dependent control strategies used to further refine the model of Section 2. SUA techniques employed in the analysis of the modified model are described in Section 4.
2. Mathematical model.The model proposed in [7], which subdivides a given SARS-affected region into six mutually exclusive subpopulations (susceptible (S), asymptomatically infected (E), quarantined (Q), symptomatic (I), isolated (J) and recovered (R) individuals), consists of the following differential equations: In ( 1)-( 6), the total population at time t, denoted by N (t), is given by The associated model parameters and variables are described in Tables 1 and 2. Since the above model monitors human populations, all parameters and state variables are assumed nonnegative for all t.
For further details about the model and its associated parameters, the reader may refer to [7].
The model has a disease-free equilibrium (DFE) given by Following [7], let with By linearizing the model around the DFE with p = 0, the following result can be proven.
Lemma 1.The DFE, (7), of the model ( 1)- (6), is locally asymptotically stable if The threshold quantity R c is called the control reproduction number [1], which represents the average number of new infections generated by a single infective when SARS control measures are in place.Lemma 1 shows that introducing a small number of infectives into the population will not lead to a major epidemic if R c < 1 (i.e., in this case, SARS will be eliminated if R c < 1).Similarly, SARS will persist in the population if R c > 1.For a review of reproduction numbers, see [8].
2.1.Parameters for control measures.In the absence of a rapid diagnostic test, therapy or vaccine, SARS control measures were based on quarantine (of exposed individuals), isolation (of symptomatically infected individuals) and stringent hygienic precautions during isolation (prescribed for health care personnel).In the model ( 1)-( 6), parameters corresponding to the health control measures (HCMs)quarantine, isolation and hygienic precautions-are denoted by γ 1 , γ 2 and J , respectively (see Table 2).It should be noted that, in the model, individuals in quarantine are assumed to be asymptomatically infected.
In [7], γ 1 and γ 2 were set to 0 on the day (February 23, 2003) the index case was reported in the Greater Toronto Area (GTA) and switched to γ 1 = 0.1 and γ 2 = 0.5 a few weeks later ( March 30, 2003), signifying the full implementation of the quarantine and isolation programs.
One key feature of the 2003 SARS outbreak was that, subsequent to the introduction of quarantine and isolation measures, a large number of nosocomial infections were reported.These infections were attributed to inadequate hygienic precautions taken by health care workers who came in contact with isolated infectives.The modification factor, J , was introduced into the model in [7] to reflect the level of hygienic precautions undertaken during isolation.Perfect hygienic precautions would mean J = 0, implying that further SARS infections were not arising through contact with isolated individuals.The transmission coefficient, J β, associated with SARS transmission during isolation was larger before stringent hygienic measures were put in place.In [7], J was set to 0.36 prior to April 30, 2003, after which J was switched to 0.
Although γ 1 , γ 2 and J are modelled as discontinuous step functions in [7], as described above, such functions are somewhat simplistic in reflecting the progressive refinement of HCMs over time.Rather, it seems more plausible to model γ 1 , γ 2 and J as continuous, time-varying functions.In this manuscript we construct such input functions, and analyze the resulting continuous, time-varying, output functions.

Development of time-varying control functions.
As mentioned in Section 2, it is assumed that the aforementioned SARS control measures were ineffective at the beginning of the epidemic, but were progressively refined over time.In contrast with modelling γ 1 , γ 2 and J via simple step-functions, as in [7], this paper introduces continuous, time-varying functions for better representing the gradual implementation of HCMs.Notationally, γ 1 (t), γ 2 (t) and J (t) are now time-dependent, and their functional forms are described in detail below.
As a starting point for constructing the time-varying control functions, consider the following cumulative distribution function, defined on : where e ≥ 0. The function H e is a standardized version of a control function, where the parameter e determines how quickly the anti-symmetric function H e moves from its initial value of 0, at z = −1, to its final value of 1, at z = 1.To allow for the possibility of initial and final values other than 0 and 1 and for alternative ranges for the support of the distribution, the following continuous time-dependent function The functional forms (with inflection point at t = c) for γ 1 (t), γ 2 (t) and J (t) are now defined by setting for i = 1, 2 and and Although the functions above may appear somewhat complex at first glance, they are actually quite simple and easy to interpret.They represent a natural way of moving from an initial value of a to a final value of b in a continuous antisymmetric manner between times t = r and t = s.The parameter e reflects the speed at which implementation of the relevant control measure takes place between times r and s.Other functional forms could have been used instead; in fact, any suitably parameterized cumulative distribution function for a continuous random variable on the interval (−1, 1) could have been used in place of H e .The function need not be antisymmetric.Moreover, the standardized function could be defined on an infinite range, such as (0, ∞), with appropriate modifications being made to the definition of h.This would permit time-varying functions analogous to exponential decay to be used.However, it is thought that the functions used in this paper provide a suitable representation of the implementation of SARS control measures, while being flexible enough to permit a SUA to be performed on the model parameters.
3.1.Parameters for the time-dependent functions.Prior to implementation of quarantine, γ 1 (t) = a 1 = 0.That is, no one is quarantined for t < r 1 , where the parameter r 1 represents the time (number of days since onset of epidemic) at which health care professionals first introduce quarantine measures for asymptomaticallyinfected individuals.The parameter s 1 represents the number of days (since onset of epidemic) until quarantine is implemented throughout the region of interest.At t = s 1 , the function γ 1 (t) = b 1 , the limiting (final) value for γ 1 (t).It is assumed that the maximum attainable quarantine rate is reached at time t = s 1 .The parameter e 1 determines the "shape" of γ 1 (t) between days r 1 and s 1 .For example, if e 1 = 0, then γ 1 (t) will be linear between days r 1 and s 1 .Furthermore, e 1 > 0 implies that implementation of HCMs are slower near days r 1 and s 1 , and faster near the inflection point, (r 2 + s 2 )/2.Shape parameter for J (t) between r J and s J From a health policy standpoint, it is also useful to consider the reparametrization d 1 = s 1 − r 1 , the number of days required to fully adhere to quarantine HCMs.An analysis (see Section 4) investigating the time required for the health-care system to respond completely and effectively to the SARS epidemic will shed light on how delays in implementation could significantly increase the number of SARS cases.
Table 4 provides ranges (and nominal values) of the parameters associated with the function γ 1 (t).While the parameter ranges were selected so that the SUA analysis could investigate a wide variety of plausible disease scenarios, the nominal values correspond to reasonable baseline values for the parameters based on the GTA estimates contained in [7].The functional form for γ 1 (t) is depicted in Figure 1 using the nominal values of a 1 , b 1 , r 1 , s 1 and e 1 (given in Table 4).
The parameters associated with the function γ 2 (t) are similarly defined (see Table 3), although several distinctions are noted below.First, before isolation measures are introduced (t < r 2 ), it is assumed that a small number of individuals may voluntarily "isolate" themselves at hospitals because of rapidly failing health.
This possibility is incorporated by letting a 2 vary between 0 and 0.2 (Table 4) so that a 2 ≥ a 1 .In line with [7], it is assumed that the isolation rate is greater than the quarantine rate so that b 2 > b 1 .From the preceding discussion, it is reasonable to expect that control measures will be introduced and subsequently implemented sooner for symptomatic individuals than for asymptomatic individuals (so that r 2 ≤ r 1 and s 2 ≤ s 1 ).
The infectiousness and contact rates between a susceptible and a SARS-infected individual in the asymptomatic, quarantined and isolated classes are denoted by E , Q and J , respectively.As in [7], we set E = Q = 0 since small values for these modification parameters yield similar simulation results.However, the infectiousness and contact rate between a susceptible and a SARS-infected individual in the isolated class is modelled as a function of time using J (t).The modification We assume that J (t) is initially positive by setting J (t) = a J until t = r J , the time at which proper hygienic precautions (gloves, face masks, etc.) are first introduced.When hygienic precautions are fully implemented (at time t = d J ), then J (t) = b J = 0. 4. SUA of the refined SARS model.This section employs analytical tools from the field of SUA to identify key parameters in the refined SARS model, given by ( 1)-( 6), with γ 1 (t), γ 2 (t) and J (t) as defined in Section 3. We begin the section by describing techniques for sampling the input parameter values from their proposed ranges (Tables 4 and 5).This is followed by application of uncertainty analysis techniques to the refined SARS model.An uncertainty analysis is conducted to quantify how the uncertainty in the choice of input parameter values produces variability in the response(s).The uncertainty analysis is followed by a sensitivity analysis to ascertain those parameters in the model which are most influential (this enables the ranking of input parameters in terms of their effect on the response).
It is worth stating that replacing the control parameters γ 1 , γ 2 and J in the autonomous model ( 1)-( 6) with their respective time-dependent formulations from Section 3 makes the model nonautonomous.Nonetheless, since the control reproduction number for the autonomous case (R c ) governs the persistence or effective control of SARS near the DFE, it is instructive to analyze the dynamics of this threshold when time-dependent functions are used.In other words, this study monitors, as output, the behavior of R c as a function of time (R c (t)).
Before assessing the effect of parameter uncertainty on the response, a methodology for selecting (sampling) input parameter values has to be chosen.Many sampling techniques are available to practitioners, but each must be judged in the light of model assumptions, resource availability and other practical considerations.
For the SARS model under consideration, a given sampling approach must ensure reasonable coverage of the high-dimensional parameter space while remaining parsimonious in terms of total sample (run) size.To accomplish this, we employ a sampling technique known as Latin hypercube sampling (LHS) [2,3,15,23,24,25].In terms of its space-filling qualities, LHS is typically superior to simple random sampling (SRS) in cases where many parameters are present in a model and limitations exist on the number of runs that are computationally practical.The superiority of LHS is particularly realized when the response of interest is a monotonic function of each of the input parameters.
LHS has been extended to include other criteria in the sample selection process to produce samples with additional desirable space-filling properties [16,17,25,26,27].These refinements to LHS are not required here, because of the ease with which relatively large samples can be selected in this problem.
Parameters are treated as random variables in LHS.Thus, probability distributions must be assigned to each parameter prior to the selection of the parameter values.Our ensuing analysis assumes that the parameters are either independently uniformly distributed, pairwise uniform on constrained (triangular) regions, or conditionally specified.In the context of the refined SARS model, ranges for the parameters are provided in Tables 4 and 5.It should be noted that the parameters a 1 , b J , e J , µ, E and Q are fixed.Using the preceding methodology, samples of size n = 1, 000 and n = 2, 000 were generated.Although the SUA in the subsequent sections will focus only on the results obtained for n = 1, 000, the plots, analyses and conclusions were similar for the case with n = 2, 000.As noted above, it was relatively easy to obtain a large sample in this problem, so it was possible to compare the results from different sample sizes.The sample size of n = 1, 000 was chosen because the curves were smooth and the results were stable by that point.In general, the required sample size will depend on the number of parameters and the characteristics of the problem.However, it is of interest to note that [2,3] also use n = 1, 000 in their simulations.4.1.Uncertainty analysis.For each of the 1,000 runs of the LHS, a "curve", R c (t), is generated as the response.A question that arises is how to conduct a SUA when the output has such functional characteristics.In our uncertainty analysis, we adopt a cross-sectional approach, except when analyzing the response t * below.Although the cross-sectional approach does not fully exploit the functional analytic nature of the output, it will be shown that this approach leads to useful conclusions in this problem.This is largely because the curves of interest are monotone decreasing for the parameter combinations under consideration.The cross-sectional approach would not work nearly as well if the curves were more complex.The approach used here provides a starting point for conducting a SUA for a model with functional output and serves as a springboard for future research.Fig. 2 displays the boxplots of the response, y = R c (t), for days t = 0, 10, . . ., 100.At a given day, each boxplot displays the 25 th (Q 1 ) and 75 th (Q 3 ) percentiles of R c (t).The 25 th and 75 th percentiles are denoted by the lower and upper horizontal lines on a box, respectively.The horizontal line within a box denotes the median value (50 th percentile) of R c (t).The "whiskers" protruding from each box extend to the most extreme values for R c (t), which are no more than 1.5(Q 3 − Q 1 ) away from the box.Any value for R c (t) plotted beyond the whiskers is classified as an outlier.Finally, the solid line depicts the line R c (t) = 1.Although no theoretical analyses of the nonautonomous model are reported here, it is reasonable to infer that values of R c (t) less than one correspond to cases where initial SARS outbreaks are effectively controlled, and values greater than one correspond to cases where SARS will persist (in line with the theoretical results for the autonomous case given by Lemma 1).
The boxplots of R c (t) are monotonically decreasing over t = 0, 10, . . ., 100.Most R c (t) curves become less than one at around 40 to 60 days after the onset of the epidemic, although some R c (t) values are less than one at t = 0. Figure 2 further demonstrates that there are no changes in the R c (t) curves between t = 0 and t = 10.Also, there is very little movement in these curves between t = 10 and t = 20.Recall that health control measures are not implemented until t = 15 (Table 4); therefore, changes in values for the time-dependent functions (γ 1 (t), γ 2 (t) and J (t)) comprising R c (t) do not occur until after this time.Lastly, the presence of outliers in Figure 2 indicates that certain parameter configurations result in values of R c (t) that are markedly higher than the "typical" values inside of the boxes.Under these parameter configurations, SARS would persist for a longer period of time.
We now consider the distribution of the response, y = t * , the first day at which R c (t) = 1. Figure 3 presents a boxplot and histogram of t * .The median number of days (since the onset of epidemic) until R c (t) = 1 is just over 40.Aside from the values for which t * = 0, the histogram exhibits a roughly symmetric distribution.4.2.Sensitivity analysis.Sensitivity analysis employs quantitative methods for determining which parameters have the largest impact on a model's response variables.From the previous section, it should be recalled that the response variables are R c (t) at fixed time points (days) and t * , the first day at which R c (t) = 1.
Many sensitivity analytical techniques are available [23]; some of the more popular techniques include various graphical methods (scatterplots, cobweb plots), correlation measures (partial correlation, rank correlation), regression methods (stepwise regression, rank regression, nonparametric regression) and more advanced approaches such as variance-based methods (sensitivity indices, total effect indices).4.2.1.Partial rank correlation coefficients.In line with [2,3,4,10,23,24], our sensitivity analysis begins with the calculation of partial rank correlation coefficients (PRCCs).This means that, in our case, the PRCCs must be calculated between each of the k = 19 nonconstant parameters and the responses R c (t) at each time t and t * .Calculating PRCCs between inputs and a response is a useful exploratory technique for ascertaining the importance of individual parameters, or main effects.The attraction of PRCCs also lies in their computational simplicity and in their ability to parallel conclusions obtained by more advanced sensitivity procedures.
To calculate PRCCs, we first convert the input parameter values to their respective rank (1 to n) with a response.The reasons for using ranks in place of the actual parameter and response values are two-fold.First, unlike the actual measurements, the ranks for each parameter and output automatically have the same scale and range.More importantly, whereas the output (R c (t) or t * ) is monotone in each parameter, the relation between any particular parameter, x i , and the output may be highly nonlinear.In this scenario, the usual partial correlation coefficient, being suited for assessing a linear relationship between a parameter and a response, may not properly reflect the strength of the influence of an x i on the output.Once PRCCs have been calculated, parameters can be ranked in descending order of epidemiological importance, according to the magnitude of their PRCCs.Parameters having PRCCs approaching the bounds −1 or +1 signify a stronger impact on the output (R c (t)).
Note that, in this section, we consider d 1 = r 1 − s 1 and d 2 = r 2 − s 2 rather than s 1 and s 2 , explicitly.As mentioned in Section 3.1, this reparametrization will assist in the interpretability of the model and the development of useful health policy strategies.
PRCC curves between each of the 19 nonconstant parameters (Tables 4 and 5) and R c (t) are plotted continuously as a function of time, t (Figure 4).(It is worth emphasizing that this approach is similar to that used in [2], in which PRCCs were calculated for different years, when modelling the transmission of drug-resistant HIV.)To reduce "clutter," we use separate graphs when plotting the PRCC curves for the time-independent parameters (Figure 4(d)) and the parameters corresponding with the time-dependent functions, γ 1 (t), γ 2 (t) and J (t) (Figures 4(a)-(c)).
Figure 4(d) shows that for all t, β has the largest positive PRCC of the timeindependent parameters.This large PRCC for β implies that as the infectiousness and contact rate between a susceptible and symptomatic individual increases, so too will the length of time required before SARS is eliminated.Justification for such a strong relationship between β and R c (t) may also be obtained directly from the equation for R c (t) in (8), from which it is evident that β is (positively) linearly related to R c (t) (since β could be factored out from each numerator in ( 8)).Of  the remaining time-independent parameters, the most influential are κ 1 (rate of development of clinical symptoms in the asymptomatically infected population) and σ 2 (rate of recovery in the isolated class).The parameter σ 2 has PRCCs of about −0.4 until (approximately) t = 25 days, at which point its correlation with R c (t) drifts toward zero.The parameter κ 1 has PRCCs of approximately 0.4 when t > 60 days.For t > 40 days, the PRCCs of all remaining time-independent parameters tend toward zero (suggesting their marginal significance in driving the disease dynamics).
For the 12 nonconstant parameters associated with γ 1 (t), γ 2 (t) and J (t), Figures 4(a)-(c) suggest three distinct stages in the transmission and control dynamics of SARS.In the early stages of the outbreak (approximately 0 -30 days), the two parameters (other than β) that most influence R c (t) are a 2 , the rate of self-reporting (self-isolation) and a J , the infectiousness of SARS prior to the implementation of hygienic precautions during isolation.In the intermediate stage (approximately 30 -60 days), the influence of a 2 and a J on R c (t) decreases.During this second stage, there is also a marked increase in the influence of d 1 (the number of days taken to fully implement quarantine HCMs) and d J (the number of days needed to fully adhere to proper hygienic precautions) on R c (t).In the third and final stage (approximately 60 days and beyond), the influences of d 1 and d J wane while the influence of b 1 (the asymptotic quarantine rate) and b 2 (the asymptotic isolation rate) become prominent.Figure 5 displays the PRCCs of the time-independent and time-dependent parameters with t * .The parameters most influential in determining the time at which R c (t) = 1 are β, d J and r J .Although β is intrinsic to SARS, the parameters d J and r J can be controlled by health care professionals, and as suggested by Figure 5, an early and rapid implementation of stringent hygienic precautions will have a significant impact on the timely elimination of SARS.4.2.2.Contour plots of two-parameter interactions.Although PRCCs are suitable for determining the effect of individual parameters on a response, they fail to provide insight regarding the effect of parameter interactions on the model output.
As a means for investigating interactions amongst parameters, we construct twodimensional contour plots of (x i , x j ) for R c (t) and t * .Clearly, an investigation of all 19 C 2 = 171 possible two-parameter interactions for each response would be formidable.Instead, rank regression models [23,25] can be used to identify which two-parameter interactions to investigate.In rank regression analysis one replaces the original data with their corresponding ranks to form first-and second-order regression models.Those interactions having the largest (statistically significant) impact on R c (t) and t * can be identified and plotted.
Figures 6(a) and (b) display contour plots of R c (40) for the two-parameter interactions βd 1 and βd 2 , both of which were found to have a significant impact on R c (40).Given that we have no control over β, these two-parameter interactions suggest that, in the presence of a large β, it is imperative that health care personnel implement quarantine and isolation as quickly as possible, once such measures are deemed necessary. .Conclusions.In this paper we modify a deterministic model for the transmission dynamics of SARS to include the progressive refinement of anti-SARS control measures.This is accomplished by introducing monotone, continuous, timedependent functions that mimic the introduction and subsequent implementation of quarantine, isolation and hygienic precautions.Consequently, the model contains time-varying inputs and functional output.SUA techniques are applied to determine those parameters having the largest effect on the disease dynamics, as represented by R c (t).Our analysis demonstrates the existence of three distinct stages in the evolution of SARS transmission dynamics where, other than the transmission rate (β), (i) in stage 1 (0 -30 days of the epidemic), the self-isolation rate (a 2 ) and the infectiousness and contact rate prior to the implementation of hygienic precautions (a J ) are the most critical parameters in determining the fate of the epidemic; (ii) in stage 2 (30 -60 days of the epidemic), the number of days taken to fully implement quarantine (d 1 ), isolation (d 2 ) and hygienic precautions (d J ) are the most influential; (iii) in stage 3 (after 60 days), the maximum attainable rates of quarantine (b 1 ) and isolation (b 2 ) are the most critical parameters.
It is hoped that the conclusions drawn from this study will enhance our understanding of SARS and that the manner in which we have modelled the implementation of control measures will prove useful in other disease scenarios.

Figure 1 .
Figure 1.Functional form of γ 1 (t), the quarantine rate of asymptomatic individuals, with a 1 , b 1 , r 1 , s 1 and e 1 at their nominal values

Figure 2 .
Figure 2. Control reproduction number, R c (t) vs. number of days since onset of epidemic

Figure 3 .
Figure 3.The distribution of t * , the first day at which R c (t) = 1

Figure 4 .
Figure 4. Partial rank correlation coefficients, over time, of the time-independent and time-dependent parameters with the control reproduction number, R c (t)

Figure 5 .
Figure 5. Partial rank correlation coefficients of the timeindependent and time-dependent parameters with t * , the first day at which R c (t) = 1

Table 1 .
Descriptions of State Variables and Time-Independent Parameters

Table 3 .
The Time-Dependent Functions Shape parameter for γ 2 (t) between r 2 and s 2

Table 4 .
Ranges and Nominal Values for the Time-Dependent Functions

Table 5 .
Ranges and Nominal Values for the Time-Independent Parameters