INTRODUCTION
After the identification of influenza A(H1N1)pdm09 virus in Mexico in April 2009, it spread rapidly worldwide, resulting in more than 16 900 laboratory-confirmed cases and 500–1000 deaths in over 67 countries, by mid-February 2011 [1]. The A(H1N1)pdm09 virus is presumed to spread in a spatio-temporal pattern similar to those of previous pandemics, but at an accelerated rate because of the frequent air travel in modern times [2]. The epidemiology of A(H1N1)pdm09 differs from both seasonal influenza epidemics and previous pandemics. Initially, most cases were clustered in households and schools; of these, more than half of the reported cases were school children aged between 5 and 18 years [Reference Yang3]. Public policy typically advocated in-home care for the ill, or school closures when there was high potential transmission risk at school. These policies were designed to diminish the frequency of close contacts [Reference Yang3]. However, human behaviour changes not only through policy recommendations, but also through fear of infectious disease; social contact avoidance is a potential mechanism by which disease transmission may be reduced.
Recently, Fenichel et al. [Reference Fenichel4] and Ferguson [Reference Ferguson5] used the concept of adaptive behaviour to model the transmission of infectious disease in an epidemiological system. People can be assumed to adapt their contact patterns, and these changes feed back to alter epidemic dynamics [Reference Epstein6]. Models of population dynamics have incorporated the influence of human behaviour on the spread of infectious diseases on behaviour-associated contact structures, model parameters, and individual disease-state-based behavioural changes [Reference Epstein6–Reference Perisic and Bauch11]. Thus, the contact behaviour for respiratory infectious agents plays an indispensable role in mathematical modelling.
Generally, social contact patterns can be quantified by questionnaire and survey responses to estimate daily contact number or a transmission rate matrix for specific populations who are at potential risk. High rates of influenza transmission have been detected for school children and teenagers in particular. Wallinga et al. [Reference Wallinga, Teunis and Kretzschmar12] found that school children and young adults experienced the highest incidence of infection and contributed the most to its further spread during the initial phase of an emerging respiratory-spread epidemic. Moreover, an age-specific contact matrix of daily contact number (C ij ) can be estimated through the same means to construct a more detailed transmission rate matrix for disease transmission modelling.
A previous study estimated the contact behaviour in non-pandemic periods [Reference Chen13]. We are not aware of studies which assess adaptive contact patterns during epidemic or pandemic periods. We, therefore, attempted to provide a reasonable adjusted factor of human contact behaviour changes which affected disease transmission.
Handel et al. [Reference Handel, Longini and Antia14] provided the information on the relationship point between the daily contact numbers and symptom levels of infectious individuals. They suggested that a sick person might reduce the frequency of their contacts with others, i.e. an increased symptom score might be associated with behavioural change. Thus, researchers can only capture the behaviour of a susceptible individual for healthy persons who have no intention of changing their behaviour. In order to extend the association between contact behaviour, symptom scores, and viral shedding, the published experimental human influenza infection data for A(H1N1)pdm09 [Reference Cowling15] could be used to understand the relationship between its contact properties, dose (viral load), and response (respiratory symptom scores).
The most well-known susceptible-infected-recovery (SIR) model is a basic and potentially powerful model in mathematical modelling of infectious diseases [Reference Nichol16]. The key epidemiological parameter of basic reproduction number (R 0) was also estimated for a comparison between adjusted and unadjusted behavioural changes. R 0 essentially determines the rate of spread of an epidemic and how intensive a policy will need to be control the epidemic. When R 0 > 1, it implies that the epidemic is spreading within a population, whereas R 0 < 1 means that the disease is dying out [Reference Anderson and May17]. We estimated R 0 in order to compare adjusted and unadjusted behavioural changes.
Taken together, the objective of the present study was to determine the changes of human contact behaviour that affect the indoor transmission of A(H1N1)pdm09 in school children. This work outlines a practical tool to determine how to implement control measures during epidemic periods by taking into account important contact behaviours and experimental human influenza infection data.
MATERIALS AND METHODS
The framework along with the computational algorithm employed in the present study is shown in Figure 1. Based on survey data [Reference Chen13] and an experimental human study [Reference Cowling15], we performed a (mathematical) modelling of the indoor transmission of A(H1N1)pdm09 in social contact structures.
Schematic showing research framework and flowchart used in this study.
Estimation of age-specific social contact behaviour
This study used data from our previous survey [Reference Anderson and May17], which focused on school children in junior high school (grades 7–9; children aged 13–15 years). Questionnaires were completed only after the participants and their parents (or legal guardian) supplied written informed consent. The results of the questionnaires and the study were fully anonymous. In brief, a total of 404 questionnaires (202 participants) were given to junior high-school students. The effective sample size was 274 questionnaires with a 67% response rate. The survey data collected was for gender, household size, health status, and contact duration and frequency for each individual. The period investigated by the survey was for 1 week, separated into to weekdays and weekend. The contact populations investigated were classified into three age groups (0–12, 13–39, ⩾40 years). Baseline information of the survey is provided elsewhere [Reference Chen13]. Table 1 provides essential information of sampling data along with contact numbers. The mean numbers of contact for grades 7–9 ranged from 9·44 ± 8·68 (mean ± s.d.) to 11·18 ± 7·98 per person day−1, with similar contact behaviour between school grades. Statistical analysis showed that contact numbers did not differ among the three grades (P < 0·05) (Wilcoxon's rank sum test).
Characteristics of sampling data (mean ± s.d.)*
n.a., Not available.
* Total sample size of questionnaire = 404, and effective sample size = 274. Thirty-nine questionnaires were not returned, 49 included incomplete data (lacking sampling date or basic individual information) and 42 questionnaires were only completed for one of the two days required.
† Health: cough, runny nose, headache, sneezing, fever.
We used a matrix to describe the contact frequency for school children interacting with different age groups. The simplified matrix of contact number (C ij ) can be arranged as
$${C_{{\rm ij}}} = \left( {\matrix{ {{C_{11}}} & {{C_{12}}} & {{C_{13}}} \cr {{C_{21}}} & {{C_{22}}} & {{C_{23}}} \cr {{C_{31}}} & {{C_{32}}} & {{C_{33}}} \cr}} \right),$$
where, i 1,2,3 = specific age in school grades 7–9 (1: 13 years; 2: 14 years; 3: 15 years); j 1,2,3 = specific age group contacted by school children (1: 0–12 years; 2: 13–39 years; 3: ⩾40 years).
Relationship of pandemic H1N1 viral titre to respiratory symptom score
A published study, assessing the comparative epidemiology of pandemic and seasonal influenza A [Reference Cowling15] provided information on viral titre data based on RT–PCR assay and culture throughout the course of illness for pandemic and seasonal influenza, especially for teenagers (<15 years). There were 24 index patients aged 0–15 years (48% of all index patients). The daily viral titres (tissue-culture infectious dose, log TCID50 ml−1) at day 0 (time of onset of acute respiratory illness) to day 10 were re-analysed.
Mean symptom scores were calculated from a composite of three groups of signs and symptoms of influenza – systemic, upper respiratory, and lower respiratory – and each ranged from 0 to 1, with higher scores indicating greater severity of symptoms [Reference Handel, Longini and Antia14]. The A(H1N1)pdm09 virus titre and symptom score dataset was analysed using Didger 4 software (Didger® v. 4·2, Golden Software Inc., USA). This study integrated the daily-based viral titre of children with daily symptom score relationship of A(H1N1)pdm09 to match the age groups in the questionnaire.
This study calculated total symptom score (TSS) as the summation of systemic, upper respiratory, and lower respiratory symptom scores, and expressed it as a function of volunteers' nasal influenza viral titre (V). Table Curve 2D software v. 5·01 (Systat Software Inc., USA) was used to perform curve fitting. Functional equations were fitted to determine the best-fitting correlation.
Estimation of normalized contact rate
Handel et al. [Reference Handel, Longini and Antia14] previously described a mapping technique that incorporates behaviour and viral load. They assumed that a sick person might reduce their frequency of contact with other persons; in general terms, an increasing symptom score changes behaviour. Handel et al. [Reference Handel, Longini and Antia14] also expressed the normalized contact rates (w) as a function of viral load (V) as. w = 1/(1+TSS(V). In the present study, we predicted the time-dependent normalized contact rate by a fitted virus dynamic model.
The nonlinear regression models were fitted to the experimental viral shedding data. The optimal fitted model was able to describe well the distributions of viral shedding dynamics. Based on these good fits, we employed response surface methodology to explore the relationship between viral titre, TSS, and normalized contact rate. Table Curve 3D software v. 4·0 (Systat Software Inc.) was used to perform the model-fitting techniques.
Transmission rate estimation
To assess the age-specific transmission rate in junior high-school students, we adopted the concept of an infectious contact rate, which is a function of social mixing patterns and transmission probabilities for a given social contact [Reference Cowling15]. The infectious contact rate equals the number of daily contacts multiplied by the transmission probability; the probabilities of an infected individual transmitting to a susceptible contact were estimated to be in the range of 0·025–0·087 [Reference Longini18].
To understand the behaviour of decreasing contact presumed to be exhibited by the infected children, we compared normalized contact rates vs. time for infected individuals in terms of the area under the curve (AUC). In this plot, the maximum normalized contact rate is equal to 1; this corresponds to normal contact behaviour. The quantification of normal contact behaviour can be seen as a rectangular area.
The AUC of the normalized contact rate was introduced to quantify the contact behaviour change under influenza infection. Next, we calculated the adjusted factor (p) by AUC over the rectangular area, which can modify the actual contact number for susceptible children after influenza infection.
We defined a parameter, p, which expresses the adjustment factor for contact behaviour change for the infected population. p can be estimated from AUC based on the normalized contact rate vs. time curve [Reference Hayden19, Reference Liao20]. p ranges from 0 to 1, corresponding to a low to high level of behavioural change in response to a pandemic threat.
In this study, we applied the concept from Nichol et al. [Reference Nichol16] to re-interpret the algorithm of the transmission rate (β) as:
$$\beta = C \times p \times q,$$
where C is the mean contact number for school children in each grade (contacts per day), p is the adjustment factor that can reflect a reduced contact ratio for an infected individual (−), and q is the transmission probability from an infected person given one contact (−) (Fig. 1e ).
We further constructed the transmission rate matrix, which can be seen as a ‘who acquires infection from whom’ (WAIFW) matrix [Reference Opatowski21]. The transmission rate matrix (β ij ) aligns with the integrated contact number matrix above, and can be written as
$${\beta _{{ij}}} = \left( {\matrix{ {{\beta _{11}}} & {{\beta _{12}}} & {{\beta _{13}}} \cr {{\beta _{21}}} & {{\beta _{22}}} & {{\beta _{23}}} \cr {{\beta _{31}}} & {{\beta _{32}}} & {{\beta _{33}}} \cr}} \right).$$
Model of infection dynamics
The SIR model is a simple and basic mathematical model of infectious disease [Reference Anderson and May17]. To explore the impact of contact behaviour on epidemiological processes during the pandemic period, this study used an SIR model embedding normalized contact rate (w) and transmission rate matrix (β ij ) to represent pandemic modelling
We defined a population size N in a given area, and divided N into three compartments: susceptible, S; infected and infectious, I; and recovered with immunity R [Reference Fenichel4]. We further integrated each element in the transmission rate matrix by calculating the dominant eigenvalue, which can represent the dominant transmission rate (β M) for the whole population under the social contact structure.
The SIR model can provide a basic description of the transmission dynamics by using a simple parameterized set of ordinary differential equations,
$$\displaystyle{{{\rm d}S} \over {{\rm d}t}} = - \displaystyle{{{\beta _{\rm M}}SI} \over N},$$
$$\displaystyle{{{\rm d}I} \over {{\rm d}t}} = \displaystyle{{{\beta _{\rm M}}SI} \over N} - \nu I,$$
$$\displaystyle{{{\rm d}R} \over {{\rm d}t}} = \nu I,$$
where β M is the transmission rate, and νis the rate at which an infectious individual recovers per unit time.
We assumed that the population size of school children was N = 34 individuals and the initial I(t = 0) = 1. R 0 took the classic form R 0 = β M/ν [Reference Fenichel4].
RESULTS
Age-specific social contacts
Only the covariates of household size (P < 0·05), survey date (P < 0·0001), and health status of participants (P < 0·05) showed significant difference within each group. The mean numbers of contacts were 4·27 ± 0·64 d−1 and 9·63 ± 1·21 d−1 during the weekend and weekday periods, respectively, and the highest contact frequency was in the 13–39 years age group (Fig. 2).
Contact frequency estimation across three estimated contact age groups (0–12, 13–39, ⩾40 years) by questionnaire among school children for (a) the weekend and (b) weekdays, respectively, during March 2010.
Estimated adjusted fraction by viral shedding and symptom scores
Figure 3a shows the optimal fitted model for viral shedding data (r 2 = 0·80). To mimic viral shedding dynamics, we integrated time-dependent symptoms after influenza illness onset. Our results showed that TSS peaked at ~2·1 at day 1 after illness onset (Fig. 3b ).
(a) Fitted time-dependent viral titres and (b) symptom scores including systemic, upper respiratory, lower respiratory and total symptom scores which were adopted from Cowling et al. [Reference Cowling15].
Table 2 summarizes the optimized fitted equation that best describe the dose–response relationship between symptom scores and viral shedding with fitted parameter values for systematic, upper and lower respiratory, and total symptom scores. The exponential function best describes the trends for viral shedding dynamics. Results showed that systemic score and TSS yielded the best predictability for viral shedding, with r 2 = 0·93 and 0·77, respectively (Fig. 4a, b ).
Best-fitted models describing the dose–response relationship between symptom scores and viral shedding for (a) systematic and (b) total symptom scores, respectively. (c) A response surface describing the relationship between viral titres, total symptom scores, and normalized contact rates.
Fitted exponential equation parameters for viral shedding dynamics and relationship of viral shedding and symptom scores
* Exponential function: V(t) = k 1exp[−k 2 t].
† Exponential function: SS(V) = k 1exp[k 2log(V)].
We further constructed a response surface to describe the relationship between viral titre, TSS, and normalized contact rate (Fig. 4c ). Virus shedding of 0–5 logTCID50 ml−1 allows for TSS within 0·6–2·0 and a reduction of normalized contact rate from 0·62 to 0·33 (Fig. 4c ). Based on the viral shedding dynamics, we can predict the contact behaviour change for an infected individual. Figure 5 shows the time-dependent normalized contact rate from 0 to 5 days after illness onset. In this period, the normalized contact rates were increasing from 0·41 to 0·63. The adjusted fraction (p) could then be calculated based on the mean ratio of AUC and rectangular area (0·519, 95% CI 0·455–0·659).
Time-dependent normalized contact rate with 95% confidence intervals from days 0 to 5 since illness onset.
SIR dynamic modelling
According to the matrix for contact numbers, we further estimated the transmission rates during the weekday and weekend periods. The transmission rates on weekdays were estimated as the matrix:
$${\beta _{{ij}}} = \left( {\matrix{ {0.087} & {0.415} & {0.075} \cr {0.042} & {0.388} & {0.050} \cr {0.035} & {0.323} & {0.091} \cr}} \right).$$
The transmission rates on the weekend were calculated as:
$${\beta _{{ij}}} = \left( {\matrix{ {0.083} & {0.166} & {0.148} \cr {0.042} & {0.142} & {0.119} \cr {0.029} & {0.192} & {0.107} \cr}} \right).$$
Results showed that the adjusted contact behaviour-based transmission rate (β M) was estimated to be 0·389 [95% confidence interval (CI) 0·261–0·534] on weekdays and 0·259 (95% CI 0·168–0·355) on the weekend, respectively. In addition, we calculated the recovery rate (ν) as 0·199 (95% CI 0·159–0·269). Table 3 summarizes the results of model parameterization.
Summary of estimated parameters with influenza infection-associated contact behavioural change for pandemic H1N1 2009
CI, Confidence interval; AUC, area under the curve.
* Estimated based on reported mean duration of viral shedding from Bhattarai et al. [Reference Bhattarai36] for pandemic influenza: 5·02 d (95% CI 3·73–6·29).
† Adopted from Nichol et al. [Reference Cowling15].
‡ Calculated whole rectangular area in Figure 5.
§ β M = calculated maximum eigenvalue from transmission rate matrix.
¶ R 0 = β M/ν.
To predict the population dynamics of influenza transmission in the classroom, we incorporated the estimated probability distributions of parameters for dominant transmission rates into the population dynamics model (Fig. 6). In 50-day simulations, the proportions of time-dependent infected vs. total number showed the greatest differences between adjusted and unadjusted behaviour groups. The dynamics of behaviour-adjusted populations were investigated by tier percentile of the dominant transmission rates of 0·25, 0·32, 0·36, 0·41, and 0·52 d−1 (2·5, 25, 50, 75, and 97·5 percentiles, respectively). Our results indicated that the peak infected population could be reduced by 38–69% by incorporating different levels of contact behaviour change (Fig. 6c ).
Population dynamics modelling of influenza trans- mission in school children with adjusted (a, c, e) and un- adjusted (b, d, f) behaviour change for susceptible, infected, and recovered individuals, respectively.
To investigate contact behaviour change-induced influenza transmission decreases in school children, we calculated the distribution of R 0 under different scenarios (Table 3). The results showed that contact behaviour change for an infected population can truly reduce R 0, with estimates for adjusted contact behaviour change of 1·854 (95% CI 1·146–2·781) and 0·860 (95% CI 0·617–1·174) on weekdays and the weekend, respectively, based on the probability distributions of dominant transmission and recovery rates. The results revealed that the social contact structure during the weekend may ease the disease-spreading potential in school children, i.e. the R 0 value may be <1. Yet, the disease will spread rapidly under an unadjusted contact scenario, for which the mean R 0 estimates were 3·3 and 2·2 on weekdays and the weekend, respectively.
DISCUSSION
Social contact among school children
Contact processes and age-specific transmission rates among populations for respiratory-spread infectious agents play an important role in the spread of disease, and thus necessitate the need for mathematical modelling. The transmission rate matrix is a classic method for expressing the transmission rate between age groups [Reference Wallinga, Teunis and Kretzschmar12, Reference Opatowski21]; however, early research could only assume contact patterns for modelling a priori. Thus, the strength of the present study lies in the fact that we provided real daily contact numbers for school children between and within each age group.
Based on the survey data of contact behaviour, our results found that the mean number of contacts for school-going young adults in the 13–19 years age group was higher than for the adult group. McCaw et al. [Reference McCaw22] indicated that for all encounters, the contact number varied for each age group, ranging between 1 d−1 (70–79 years age group) and 22 d−1 (40–49 age group). In addition, Mikolajczyk et al. [Reference Mikolajczyk23] conducted a questionnaire in a primary school in Germany and indicated that the mean number of contacts was 25·1 ± 16·5 d−1 (min–max: 0–78) for children and 7·5 ± 5·0 d−1 (min–max: 1–47) for adults. Therefore, the afore-mentioned results implicate that influenza control strategies should focus on school children within the same age group.
Contact behaviour response influenced by self-awareness of symptoms
Human social contact and disease-associated behavioural change play an important role in the spread of diseases, and understanding them can help to reinforce the necessity of other control strategy efforts. Funk et al. [Reference Funk, Salathé and Jansen7] indicated that disease states of individuals can change contact behaviour, and this could potentially modify model parameters and contact frequency. This conclusion led us to employ an adjusted contact ratio parameter for the infected population in the present study.
Indeed, contact rates may also be reduced by symptomology. Hayden et al. [Reference Hayden24] and Fritz et al. [Reference Fritz25] found that the time-course of TSS exhibited similar trends to virus dynamics after experimental influenza infection, Handel et al. [Reference Handel, Longini and Antia14] elaborated that exponential functions were capable of describing the relationship between viral titre and TSS. In light of this, we introduced exponential functions into our model, which best fit the relationship between viral titre and TSS and systemic symptom scores. A point of concern was that model reliability may be affected negatively by uncertainty surrounding the adjusted contact ratio. However, Handel et al. showed that normalized contact rate could be compared against symptom strength to reveal the reduction in contact.
Based on this conclusion, we decided to investigate the relationship between viral titre, TSS, and normalized contact rate. We were able to reduce the uncertainty of adjusted fractions by using virus dynamics-associated normalized contact rates. We used the AUC ratio to calculate the probability distribution of adjusted fractions during the infectious period. The adjusted fractions indicated that the infected population decreases by about 50%, paralleling diminished contact frequency with other people. This supports the notion that our study can effectively link behavioural dynamics to emerging epidemiological models.
Funk et al. [Reference Funk, Salathé and Jansen7] indicated that behavioural change could be associated with infectious disease dynamics. It can exert effects on (i) the disease state of the individual, (ii) the parameters of transmission rate and recovery rate, and (iii) the contact structure and frequency as they pertain to the spread of a disease. In our study, the adjusted fractions further affected the transmission rate and influenced the population dynamics. When the adjusted fractions were not taken into account, the disease spread rapidly. Indeed, under the worst cases, the mean R 0 was estimated as high as 3 for an infected population without behaviour change. Yet the inclusion of behaviour change intervention in our model allowed the infected individual to decrease the spread of disease. It should be noted that contact behaviour change of an infected population was unable to completely negate disease transmission on weekdays due to the distribution of R 0 being >1.
It is important to mention Brauer's [Reference Brauer26] suggestion that a susceptible population may decrease their contact frequency with an infected group. The major reason for this is assumed to be fear of disease infection [Reference Epstein6]. We only quantify the contact behaviour change for the infected population in our study due to existing difficulties in accurately modelling the behaviour change for the susceptible population; however, we intend to focus on constructing a complete framework to quantify the contact behaviour and structure during the period of emerging infectious disease in future studies.
Implications for control strategy
In consideration of contact behaviour change in the infected population, we can also determine the disease control efficacy for class suspension strategy by comparing the R 0 in weekdays and the weekend. We found that R 0 < 1 in the weekend, implying that school closures or class suspensions are probably an effective social distancing policy to control the transmission of pandemic influenza [Reference Gross, D'Lima and Blasius8, Reference Jefferson27]. School closure is also the best control strategy among non-pharmaceutical interventions, since this can reduce the contact frequency and disease spread among school children [28–Reference Cauchemez30].
In Taiwan, the so-called ‘3–2–5 intervention policy’ for class suspension was implemented to control disease transmission among all students aged <18 years during the A(H1N1)pdm09 pandemic period [Reference Hsueh31]. It stipulated that if three students in the same class were confirmed A(H1N1)pdm09 within 2 days, the class should be suspended for 5 days. Since the school-based vaccination programme started on 16 November 2009, and a 74·7% coverage rate for A(H1N1)pdm09 vaccination in students aged 7–18 years in Taipei City, the rule for class suspension has been expanded from 3–2–5 alone to include the ‘8–14 intervention policy’. This removes the obligation for class suspension if >80% of the students in a class have been vaccinated for more than 14 days. Along with the school-based immunization programme, the 3–2–5 intervention policy combined with the 8–14 policy was able to mitigate and suppress the rapid spread of the disease.
Although some benefits can be expected from the above-mentioned policies, there are still problems with respect to the best time that policies should be enforced or the plausible economic impacts of such decisions [Reference Inglesby32]. Borse et al. [Reference Borse33] indicated that school closure could cause different levels of economic impact on households in New York City. Furthermore, class suspensions or school closures alone may not inhibit the pandemic completely, and may need to be combined with other public health strategies to make disease control more efficient [Reference House34, Reference Wu35]. Due to a lack of recommendations regarding a class suspension threshold during the influenza epidemic, we could estimate a potential closure threshold by considering influenza characteristics. We can also model control efficiency of class suspension in conjunction with other public health interventions, and consider additional parameters such as the proportion of infectious individuals in class.
In conclusion, our study provides a better understanding for knowing how contact behaviour change can affect disease transmission dynamics in school children. In previous research, we conducted a questionnaire survey in order to quantify the daily contact numbers of various groups of school children. In this study, we furthered our knowledge of contact numbers and behaviour by successfully applying the questionnaire-derived transmission rate matrix to simulate contact behaviour changes during an influenza epidemic. Such a linkage facilitates the incorporation of human contact behaviour to mathematical transmission models for respiratory-spread infectious diseases that transmit from person to person via the environment.
DECLARATION OF INTEREST
None.
