A mathematical model describing the localization and spread of influenza A virus infection within the human respiratory tract
Christian Quirouette, Nada P. Younis, Micaela B. Reddy, Catherine A.A. Beauchemin
AA mathematical model describing the localization and spread ofinfluenza A virus infection within the human respiratory tract
Christian Quirouette , Nada P. Younis , Micaela B. Reddy , Catherine A.A. Beauchemin , , ∗ Department of Physics, Ryerson University, Toronto, Canada Array BioPharma Inc., Boulder, USA Interdisciplinary Theoretical and Mathematical Sciences (iTHEMS), RIKEN, Wako, Japan ∗ Corresponding author: [email protected]
August 23, 2019
Abstract
Within the human respiratory tract (HRT), viruses diffuse through the periciliary fluid (PCF) bathingthe epithelium, and travel upwards via advection towards the nose and mouth, as the mucus escalatorentrains the PCF. While many mathematical models (MMs) to date have described the course of influenzaA virus (IAV) infections in vivo, none have considered the impact of both diffusion and advection on thekinetics and localization of the infection. The MM herein represents the HRT as a one-dimensional trackextending from the nose down to a depth of 30 cm, wherein stationary cells interact with the concentrationof IAV which move along within the PCF. When IAV advection and diffusion are both considered, theformer is found to dominate infection kinetics, and a 10-fold increase in the virus production rate isrequired to counter its effects. The MM predicts that advection prevents infection from disseminatingbelow the depth at which virus first deposits. Because virus is entrained upwards, the upper HRT sees themost virus, whereas the lower HRT sees far less. As such, infection peaks and resolves faster in the upperthan in the lower HRT, making it appear as though infection progresses from the upper towards the lowerHRT. When the spatial MM is expanded to include cellular regeneration and an immune response, it cancapture the time course of infection with a seasonal and an avian IAV strain by shifting parameters in amanner consistent with what is expected to differ between these two types of infection. The impact ofantiviral therapy with neuraminidase inhibitors was also investigated. This new MM offers a convenientand unique platform from which to study the localization and spread of respiratory viral infections withinthe HRT.
Mathematical models (MMs) of influenza A virus (IAV) infection kinetics are mainly based on ordinarydifferential equations (ODEs) that describe the time evolution of infection (virus as a function of time) withthe implicit assumption that all cells see all virus and vice-versa, as discussed in several reviews [10, 23, 29].A different approach from ODE MMs is to use agent-based MMs which treat the dynamics of each cell, andeven each virus, individually and track local cell-virus interactions [5, 7, 71]. Such MMs are computationallyintensive, restricting the number of cells, and therefore the area, which can be represented, and often lacksupport or validation in the form of experimental localization data [4, 7, 50]. Infection localization can alsobe modelled by dividing the human respiratory tract (HRT) into compartments corresponding to differentregions of the HRT, where each compartment has different parameters based upon physiological differences[62].IAV infection spread faces a tremendous physical, spatial barrier in the form of the thick mucus layerwhich lines the HRT, trapping and carrying free virions rapidly upwards towards the nose and throat andout of the HRT [26]. This physical clearing process is typically taken somewhat into account in non-spatial1 a r X i v : . [ q - b i o . CB ] A ug igure 1: Representation of the human respiratory tract by the mathematical model.
The MMEqn. (1) represents the HRT as a one-dimensional track, as illustrated. The MM considers only the virusthat is located in, and diffusing within the PCF, which also moves along with the PCF upwards at a fixedadvection speed. What happens to the virus when it reaches the bottom (reflective) or top (absorptive)of the MM-represented HRT is indicated. MM (1) describes the interactions between the moving virusconcentration in the PCF, and the stationary cells that carpet the HRT.ODE MMs via a term for the exponential clearance of virions. This simplification, however, does not accountfor the fact that not all cells will have equal exposure to the virus: cells that are downstream of the mucusflow (higher in the HRT) will have greater exposure to virus than those upstream, located deeper withinthe HRT. To our knowledge, the effect of virus advection due to entrainment by the mucus layer on IAVinfection localization and spread within the HRT has never been evaluated.In this work, a spatiotemporal MM for the spread of IAV infection within the HRT using partial differentialequations (PDEs) is constructed. Through its one-dimensional representation of the HRT, the MM is usedto study the effect of viral transport modes on the course of an IAV infection in vivo. Via the addition ofonly two additional spatial parameters whose value is relatively well-established, namely the rate of diffusionand advection of virions, the MM able to produce a richer range of IAV kinetics. The effect of cellularregeneration and a simplified immune response are also considered. The PDE MM is also used to exploredifferences in the kinetics of IAV infection, in the absence and presence of antiviral therapy, for IAV infectionwith either a seasonal or an avian-adapted strain, with the latter being typically more severe [56].
In the proposed spatial MM, the HRT is represented as a one-dimensional tract running along the x -axisindicating the depth within the HRT, with x = 0 cm located at the top of the HRT (nose), and x = 30 cmterminating somewhere within the bronchi [60, 70], as illustrated in Figure 1. It is an extension of the‘standard’ MM for IAV in vitro [55, 59, 68] which adds: (1) the diffusion of virions through the periciliaryfluid (PCF) which lies between the cells’ apical surface and the thick mucus blanket which lines the airways;(2) their advection due to upwards entrainment of the PCF by the mucus layer; and (3) their effect on theone-dimensional, depth-dependent fraction of non-motile (stationary) cells in various stages of infection. Thespatial MM is formulated as 2 T ( x, t ) ∂t = − βT ( x, t ) V ( x, t ) ∂E ( x, t ) ∂t = βT ( x, t ) V ( x, t ) − n E τ E E ( x, t ) ∂E i ( x, t ) ∂t = n E τ E E i − ( x, t ) − n E τ E E i ( x, t ) i = 2 , , ..., n E ∂I ( x, t ) ∂t = n E τ E E n E ( x, t ) − n I τ I I ( x, t ) (1) ∂I j ( x, t ) ∂t = n I τ I I j − ( x, t ) − n I τ I I j ( x, t ) j = 2 , , ..., n I ∂V ( x, t ) ∂t = p n I (cid:88) j =1 I j ( x, t ) − cV ( x, t ) + D PCF ∂ V ( x, t ) ∂x + v a ∂V ( x, t ) ∂x The fraction of uninfected, target cells, T ( x, t ), located at a depth x at time t are infected at rate βV ( x, t ),proportional to the concentration of virions, V ( x, t ), at that depth and time. Newly infected cells are inthe eclipse phase, E i ( x, t ), i.e. they are infected but are not yet producing virus. After an average time τ E ,infected cells leave the eclipse phase to become infectious cells, I j ( x, t ), producing virions at constant rate p for an average time τ I until they cease virus production and undergo apoptosis. As in [55, 68], the eclipseand infectious phases are each divided into n E = n I = 60 age classes so that the time spent by cells in eachphase follows a normal-like distribution, consistent with biological observations [36].The MM assumes virions are released from stationary infected cells into the PCF, so V ( x, t ) is theconcentration of virus in the PCF. The produced virions diffuse through the PCF at rate D PCF , and aretransported upwards (towards the nose) within the PCF at speed v a , as the PCF is entrained by the mucusblanket pushed along by the beating cilia of the ciliated cells that line the HRT [47, 69]. When the virusin the PCF moves beyond the bottom or top of the HRT, it either ‘bounces back’ (reflective) or is lost(absorptive), respectively. Absorption of virions into the mucus blanket (out of the PCF), their loss of viralinfectivity over time (in the PCF), and other modes of non-specific virion clearance are all taken into accountvia a single exponential viral clearance rate term, cV ( x, t ). Virions contributing to the infection in this MMare those remaining in the PCF, in direct contact with infectable cells.In the spatial MM simulator, infection is initiated by a spatially localized, spray-like virus inoculumdeposited at depth x d within the HRT. The spray-like inoculum is represented by a Gaussian centred at thesite of deposition, with a standard deviation of 0 . × the size of a large cough droplet [80]. Atsites far from x = x d , V ( x, t = 0) ≈
0. The baseline values of the spatial MM’s parameters were all set tovalues estimated in Baccam et al. [1], obtained by fitting a non-spatial MM to patient data from experimentalprimary infections with the influenza A/Hong Kong/123/77 (H1N1) virus. Table 1 lists the initial conditionsand parameters used. A complete description of the spatial MM is provided in the Methods.
Infection of cells by the IAV, and the subsequent release of virions, occurs almost exclusively at the apicalsurface of the epithelium [15,51]. When advective motion of the PCF is neglected, virions can be assumed toundergo Brownian motion and the virus concentration distribution in the PCF will be governed by diffusion.Figure 2(a–c) shows IAV infection kinetics in the presence of diffusion alone, where quantities are shownaveraged over all spatial sites as a function of time (see Methods). The infection takes place at a much slowerpace in the spatial (diffusion only) MM than in the non-spatial MM which implicitly assumes an infinitediffusion coefficient. Figure 2(d–f) shows the spatial extent of the infection dissemination (with diffusiononly) through the HRT at certain times post-infection. Virus becomes available to target cells gradually asit diffuses away from the site of initial infection, x d = 15 cm, and the infection wavefront moves outwards3able 1: Default initial conditions and parameter valuesSymbol Description Value ∗ [Source]— Spatial simulator parameters — D PCF
IAV diffusion coefficient in PCF at 37 ◦ C 10 − m / s [37] v a advection speed of PCF 40 µ m / s [47]∆ x size of one spatial grid site 100 µ m ( N x = 3000 sites) [see Methods]∆ t duration of one time step 2 . t = ∆ x/v a ) [see Methods]— Initial conditions — x d deposition depth of virus inoculum 15 cm (cid:104) V ( x, t = 0) (cid:105) initial virus inoculum averaged over x . × − TCID / mL [1, *] T ( x, t = 0) initial fraction of uninfected target cells 1 . ∀ xE i , I j ( x, t = 0) initial fraction of infected cells 0 . ∀ x — Infection parameters — τ E duration of eclipse phase 8 h ( n E = 60) [*] τ I productively infected cell lifespan 20 h ( n I = 60) [*] c virus clearance rate 0 .
22 h − [1] β infection rate of cells by virus 1 . × − (TCID / mL) − · h − [1] p virus production rate 8 . × (TCID / mL) · h − (11 × p Baccam ) [*]* These values are used in all simulations unless otherwise stated. The spatially averaged initial inoculum, (cid:104) V ( x, t = 0) (cid:105) , equals the value in [1] (see Methods). Values for τ E and τ I were chosen so as to lie nearthe middle of the ranges of [6 ,
10] h and [10 ,
40] h, respectively, obtained from MMs of IAV infections invitro [55, 68]. The value for p is discussed in Section 2.2, and that for p Baccam is explained in Methods.from that site symmetrically towards the two ends of the HRT at x = 0 and x max . S1 Video depicts thespatiotemporal course of the infection, in the presence of diffusion only, in the form of a video.Along the length of the HRT, a layer of mucus about 0 . µ m–5 µ m thick covers the PCF [26]. Thecollective motion of the underlying epithelial cells’ beating cilia, dubbed the mucociliary escalator, drivesthis mucus layer upwards. It also leads to an upward advection of the PCF, at a speed similar to that of themucus layer [47], entraining any virus in the PCF upwards at that speed. Given the advection speed of themucus and PCF ( v a ≈ µ m / s), any newly produced or deposited virion would be cleared from the HRTin less than ∼
12 min (30 cm /v a ). Thus, it is not surprising that adding advection to the spatial MM, whilestill using the same parameter values, results in a subdued, low viral titer, slow-growing infection which stillhas not peaked by 7 dpi. In contrast, viral titer in patients infected with influenza A/Hong Kong/123/77(H1N1) virus in Baccam et al. [1] peaks 2–3 dpi, which the non-spatial MM with these same parametersreproduces well.Figure 3(a–c) shows the IAV infection kinetics in the presence of both diffusion and upward advection asthe virus production rate, p , is increased from the base value estimated by Baccam et al. [1] for their non-spatial MM, p Baccam (see Methods). Increasing the virus production rate to 11 × p Baccam yields an infectionthat peaks at ∼ p = 11 × p Baccam , see Table 1) isused in the remainder of this work so that the spatial MM qualitatively reproduces the approximate timingof viral titer peak in an IAV infection.Using the new value for the virus production rate, Figure 3(d–f) shows the spatial extent of the infectiondissemination in the presence of both diffusion and advection in the spatial MM. It illustrates the protectiveeffect of advection: preventing the infection from travelling much beyond its initial deposition depth ( x d =15 cm), as seen from the target cell depletion shown in Figure 3(d). Figure 3(g–i) explores the effect ofvarying this depth of deposition of the initial virus inoculum, x d . When the inoculum deposits lower in theHRT (as x d increases), the fraction of HRT consumed by the infection increases, resulting in higher viral titerpeak and total virus yield. However, the virus titer curve is largely insensitive to the depth of deposition for4 . . . . . . F r a c t . t a r g e t ce ll s ( T ) (a) . . . . . . F r a c t .i n f ec t i o u s ce ll s ( I ) (b) V i r u s ( TC I D / m L ) (c) . . . . . . . . . . . F r a c . t a r g e t ce ll s ( T ) (d) . . . . . . . . . . . F r a c .i n f ec t i o u s ce ll s ( I ) (e) . . . . . V i r u s ( TC I D / m L ) (f) Figure 2:
IAV infection kinetics in the presence of diffusion alone (no advection). (Top) Timecourse (averaged over space) of the infection for the fraction of cells in the (a) target/uninfected or (b)infectious state, and (c) the infectious virus concentration, obtained using the non-spatial ODE MM (dashed)or the spatial (diffusion only) MM (solid). (Bottom) Localized fraction of cells in the (d) target or (e)infectious state and (f) the infectious virus concentration as a function of depth within the HRT shown atspecific days post-infection (dpi). 5 . . . . . . F r a c . t a r g e t ce ll s ( T ) (a) . . . . . F r a c .i n f ec t i o u s ce ll s ( I ) (b) p Baccam p Baccam p Baccam p Baccam p Baccam p Baccam V i r u s ( TC I D / m L ) (c) . . . . . . . F r a c . t a r g e t ce ll s ( T ) (d) . . . . . . F r a c .i n f ec t i o u s ce ll s ( I ) (e) V i r u s ( TC I D / m L ) (f) . . . . . . F r a c . t a r g e t ce ll s ( T ) (g) . . . . . . . . F r a c .i n f ec t i o u s ce ll s ( I ) (h) x d = 1 cm x d = 5 cm x d = 10 cm x d = 15 cm x d = 20 cm x d = 29 cm V i r u s ( TC I D / m L ) (i) . . . . . . F r a c . t a r g e t ce ll s ( T ) (j) . . . . . . F r a c .i n f ec t i o u s ce ll s ( I ) (k) x = 1 cm x = 5 cm x = 10 cm x = 12 cm x = 13 cm x d = 14 cm V i r u s ( TC I D / m L ) (l) Figure 3:
IAV infection kinetics in the presence of diffusion and advection. (a,b,c) Time course(averaged over space) of the infection for the fraction of cells in the (a) target/uninfected or (b) infectiousstate, and (c) the infectious virus concentration, obtained using the ODE MM (dashed) or the spatial MM(solid), as the rate of virus production, p , is varied. (d,e,f) Localized fraction of cells in the (d) target or (e)infectious state and (f) the infectious virus concentration as a function of depth within the HRT shown atspecific days post-infection (dpi). (g,h,i) Same as (a,b,c) but varying the depth of deposition of the initialvirus inoculum ( x d ). (j,k,l) Similar to (a,b,c) but rather than being averaged over space, the infection timecourse is shown at specific, spatially localized depths ( x ). Unless otherwise noted, p = 11 × p Baccam and x d = 15 cm. 6 . . . . . F r a c .i n f ec t i o u s ce ll s ( I ) (a) − − (b) AmplitudeSpeed (cm/d)(14 . − µ p )(cm)
12 13 14 15Depth (cm) − − − F r a c . I ce ll s ( × − ) (c) (diff+adv)-(adv only) Figure 4:
A closer look at infection resolution: the impact of deposition depth and diffusion .(a) A zoomed-in continuation of Figure 3(e) showing the localized fraction of infectious cells as a function ofdepth within the HRT, at times beyond 7 dpi. (b) The (nearly exponential) decay over time of the amplitudeand speed of the Gaussian pulse of fraction of infected cells shown in (a), and of its centre location (mean, µ p ) relative to its final, asymptotic location of 14 . D PCF = 0.depths greater than ∼
15 cm.In the presence of both diffusion and advection in the spatial MM, the initial virus inoculum travelsquickly from its deposition site at a fixed speed up the HRT, leaving in its wake a small, equal fraction ofinfected cells at all sites above the deposition point ( ∀ x < x d ). As the newly infected cells begin to releasevirus which is entrained up the HRT, uninfected cells are only exposed to virus produced by cells lower than(upstream of) their own depth (greater x value), such that cells higher in the HRT (smaller x value) areexposed to the most virus. Consequently, infection in the upper HRT proceeds much faster than that inthe lower HRT. Figure 3(j–l) shows the localized fraction of target and infectious cells, and infectious virusconcentration over time at specific sites or depths along the HRT. Despite infection kinetic parameters beingthe same at all sites, infection progresses faster at sites located higher in the HRT, downstream of the virusflow.Interestingly, the more rapid progression of the infection higher in the HRT makes it appear as thoughthe infection is moving downwards, from the top (nose, x = 0) towards the bottom ( x = 30 cm) of the HRT.Both the target cell depletion wavefront in Figure 3(d) and the infected cell “pulse” in Figure 3(e) appearto be moving to the right (towards the lower HRT) as time advances. This can be seen more clearly in thevideo (S2 Video) which depicts the spatiotemporal course of the infection, in the presence of both diffusionand advection. Cells in the upper HRT are consumed first, after which the pulse of infected cells appears toslowly move backwards, down the HRT, back to the initial deposition site, where it then slowly dissipates.Figure 4(a) provides a more “zoomed-in” view of what happens to the pulse in Figure 3(e) at times beyond7 dpi. The pulse is shrinking and slowing down as it approaches the initial inoculum deposition depth of x d = 15 cm. Figure 4(b) quantitatively displays the (nearly exponential) decay of both the amplitude andspeed of the pulse, as its peak’s position ( µ p ) approaches its final, asymptotic position of ∼ . in the fraction of cells infected) in light of typical, in vivo experimental variability ( ∼ During an uncomplicated IAV infection, viral loads typically fall below detectable levels between 6 dpi–8 dpi [2,13]. An experimental study in which mice were sacrificed at various times over the course of an IAVinfection [61] found that by 3 dpi, ciliated cells were rarely observed on the tracheal surface, and insteadthe latter was mostly composed of a layer of undifferentiated cells. By 5 dpi, signs of repair were apparent,by 10 dpi the tracheal surface was covered with ciliated and nonciliated cells, and by 14 dpi the surface wasindistinguishable from the uninfected surface. Another experimental study in hamsters suggests that 6 d–7 dafter mechanical injury, most of the epithelium is comparable to that of control hamsters [42]. A similarstudy in guinea-pigs suggests that 5 d after mechanical injury, the epithelium was ciliated and differentiated,and by 15 d it was indistinguishable from that of control guinea-pigs [25]. Together, these studies suggestthat the process of cellular regeneration following damage caused by an IAV infection would at least begin,if not be well underway, before the IAV infection has fully resolved.Following airway injury, numerous factors are thought to trigger the start of the repair process [18]. Oncetriggered, it is believed neighbouring epithelial cells will try to stretch to cover the denuded area, and divide soas to regain their normal shape while maintaining coverage. As damage becomes more significant, progenitorcells, mainly basal cells exposed to the PCF as the epithelial cells above them detach or are removed byapoptosis or damage, also undergo proliferation and differentiation until the epithelium is restored to itspre-injury state [64]. With these processes in mind, cellular regeneration was implemented in the spatialMM by replacing the equation for uninfected, target cells in MM Eqn. (1) with ∂T ( x, t ) ∂t = − βT ( x, t ) V ( x, t ) + r D T ( x, t ) D ( x, t − τ D ) , (2)where D ( x, t ) = 1 − T ( x, t ) − n E (cid:88) i =1 E i ( x, t ) − n I (cid:88) j =1 I j ( x, t ) is the fraction of dead cells, i.e. the extent of the damage or injury. As such, cellular regeneration in Eqn.(2) proceeds at a rate that is greater in the presence of greater damage (higher D ) and of greater fractionof cells available to regenerate (higher T ). Parameter r D sets the scale of the regeneration rate (which alsodepends on T and D ), and τ D is the regeneration delay such that the current regeneration rate depends onthe fraction of target cells ( T ) currently available to repopulate the area, and on the amount of damage ( D )that was perceived some time τ D ago. The delay between damage and regeneration accounts for the timerequired for the damage to activate appropriate signalling pathways, and for both division and differentiationto take place so that newly regenerated cells are susceptible to infection. This same equation has been usedby others to represent target cell regeneration during an IAV infection [12, 28], but the authors therein didnot include a delay ( τ D = 0).An experimental study of cell regeneration in hamsters [42] suggests that in many cases, between 18 h–24 hfollowing a mechanical injury, cells have already stretched or migrated to cover the complete denuded areaand cell division has begun. A similar study in guinea-pigs [25] also suggests that at 15 h after mechanicalinjury, cells have migrated to cover the denuded area and cell proliferation is underway. This suggests thatthe delay between damage and the start of cellular division is around 15 h–24 h. Herein, this delay —namely the time elapsed between damage and the start of both cellular division and differentiation to anextent that newly regenerated cells are susceptible to infection — was chosen to be τ D = 1 d, based on8 . . . . . . F r a c . t a r g e t ce ll s ( T ) (a) r D = 0 . − r D = 0 .
75 d − r D = 1 d − r D = 1 . − r D = 2 d − . . . . . . F r a c . t a r g e t ce ll s ( T ) (b) τ D = 0 d τ D = 1 d τ D = 3 d τ D = 5 d τ D = 7 d Figure 5:
Target cell regeneration following a MM-simulated mechanical injury.
Regenerationfrom mechanical injury was simulated using Eqn. (2) in the absence of infection, i.e. V ( x, t ) = E i ( x, t ) = I j ( x, t ) = 0 and D ( x, t ) = 1 − T ( x, t ). Initially, T ( x, t = 0) = 0 . D ( x, t = 0) = 0 .
99, and D ( x, t <
0) = 0,representing an injury inflicted at t = 0 which removed 99% of all target cells. The effect of varying (a)the regeneration rate ( r D ) or (b) the regeneration delay ( τ D ) are shown. Unless varied, τ D = 1 d and r D = 0 .
75 d − .these experimental studies. To select an appropriate value for the regeneration rate ( r D ), HRT epithelial cellregeneration following mechanical injury was simulated using Eqn. (2). Results are shown in Figure 5 fordifferent regeneration rates and delays: the regeneration rate determines the steepness of the regeneration,while the regeneration delay sets its timing. These parameters should be easily identifiable if experimentaldata was available. An intermediate value of r D = 0 .
75 d − was chosen to ensure that with a delay of 1 d,regeneration is well underway by 5 d–8 d, and completely resolved by 12 d–14 d [25, 42, 61].Figure 6 shows the infection kinetics in the presence of cellular regeneration. At this stage, and overthe range of values considered here, the spatial MM is insensitive to the choice of regeneration rate or itsdelay. This is largely due to the fact that, in the absence of an immune response, the HRT above theinoculum deposition point is decimated by the infection at a rate much faster than the density-dependentregeneration can counter. As damage increases, the number of target cells available to replenish the lost cellsalso decreases dramatically, preventing any significant regeneration. To get a more realistic view of the roleand impact of cellular regeneration on infection kinetics, the protective role of the immune response mustbe considered. Past MM for IAV infection have incorporated one or more immune response components with varying degreesof success [11,12,23]. Difficulties in incorporating an immune response when modelling IAV infections includethe complex nature of the interactions (large networks of cells and signals), and the lack of appropriate data(quantity and quality) to inform the MMs [29]. The simplified immune response considered herein comprisesan innate response based on interferon (IFN), a humoral response represented by antibodies (Abs), anda cellular response embodied by cytotoxic T lymphocytes (CTLs). Since our aim is to display the MM’srange of kinetics — rather than identify parameters, analyze data, or challenge hypotheses — a parsimoniousapproach was preferred wherein key immune components are represented using empirical curves that broadlyreproduce the scale and timing of these experimentally measured quantities (see Methods). Figure 7 showsthe empirical MM curves against their corresponding experimental time courses, for IFN, Abs, and CTLs,during experimental, in vivo IAV infections.Although IFN is known to have many effects [40], herein its effect is to reduce the virus production rate p [1, 30], via a resistance parameter, f , analogous to the IC used to describe antiviral resistance. Theresistance is defined such that if f = 0 .
8, the virus production rate is halved when IFN concentration is9 . . . . . . F r a c . t a r g e t ce ll s ( T ) (a) . . . . F r a c .i n f ec t i o u s ce ll s ( I ) (b) r D = 0 . − r D = 0 .
75 d − r D = 1 d − r D = 1 . − r D = 2 d − V i r u s ( TC I D / m L ) (c) . . . . . . F r a c . t a r g e t ce ll s ( T ) (d) . . . . F r a c .i n f ec t i o u s ce ll s ( I ) (e) τ D = 0 d τ D = 1 d τ D = 3 d τ D = 5 d τ D = 7 d V i r u s ( TC I D / m L ) (f) Figure 6:
The effect of cellular regeneration on the MM prediction of IAV infection kinetics.
Results are shown with a regeneration rate r D = 0 .
75 d − and a delay τ D = 1 d, unless otherwise specified. − − − I F N ( s c a l e d ) (a) HoshinoIwasakiFritz (human)Saenz (horse) − − − A b s ( s c a l e d ) (b) Iwasaki IgMIwasaki IgAIwasaki IgGMiao IgGMiao IgM − − − CT L s ( s c a l e d ) (c) MiaoMcLaren (ferret)
Figure 7:
Time courses of key immune response components during IAV infections in vivo.
Empirical MM curves (solid, black) are shown against experimental measurements (dashed, coloured) for (a)interferon (IFN) [27,38,41,65], (b) antibodies (Abs) [41,49]; and (c) cytotoxic T lymphocytes (CTLs) [48,49],taken over the course of in vivo IAV infections in mice, unless otherwise specified.10 . . . . . . F r a c . t a r g e t ce ll s ( T ) (a) . . . . . . . F r a c .i n f ec t i o u s ce ll s ( I ) (b) no IFN f = 1 f = 0 . f = 0 . f = 0 . f = 0 . V i r u s ( TC I D / m L ) (c) . . . . . . F r a c . t a r g e t ce ll s ( T ) (d) . . . . . F r a c .i n f ec t i o u s ce ll s ( I ) (e) no Ab k A = 50 h − k A = 100 h − k A = 200 h − k A = 500 h − k A = 10 h − V i r u s ( TC I D / m L ) (f) . . . . . . F r a c . t a r g e t ce ll s ( T ) (g) . . . . . . . F r a c .i n f ec t i o u s ce ll s ( I ) (h) no CTL k C = 0 . − k C = 1 h − k C = 10 h − k C = 50 h − k C = 200 h − V i r u s ( TC I D / m L ) (i) Figure 8:
Spatial MM-predicted IAV infection in the presence of key immune response com-ponents. (a–c) The effect of interferon (IFN) as IFN resistance, f , is decreased (sensitivity is increased).(d–f) The combined effect of IFN and antibodies (Abs) as the rate of infectious virus neutralization by Abs, k A , is increased. (g–i) The combined effect of IFN+Abs and cytotoxic T lymphocytes (CTLs) as the rate ofinfected cell killing by CTLs, k C , is increased. Unless otherwise noted, f = 0 . k A = 500 h − .11 V i r u s ( [ V ] ) (a) Seo Hoshino V i r u s ( [ V ] ) (b) Iwasaki V i r u s ( [ V ] ) (c) Yap Wells Kris V i r u s ( TC I D / m L ) (d) V i r u s ( TC I D / m L ) (e) V i r u s ( TC I D / m L ) (f) Figure 9:
MM predictions vs experimental data of immune response knockout experiments.
Experimental (top row) or MM-simulated (bottom row) viral titer time course for IAV infections with a fullimmune response (solid lines) or with one immune response component experimentally or mathematicallydisabled or knocked-out (dashed lines), respectively. Either (a,d) the IFN response [38, 66], (b,e) the Abresponse [41] or (c,f) the CTL response [43, 77, 81] were disabled. In these experimental studies, the animalmodel is mice, except Seo et al. [66] which was conducted in pigs. Experimental viral load is measured inTCID / mL for [66], in TCID / mouse (of lung homogenate) for [38], in EID / mL for [81], [43] and for [77]and in pfu / mL for [41]. The MM results are shown with f = 0 . k A = 500 h − and k C = 50 h − , ordisabled with F ( t ) = 0, A ( t ) = 0 or C ( t ) = 0 in (d,e,f), respectively.80% of its peak value (see Methods). Figure 8(a–c) shows the effect of IFN in the MM as resistance toIFN, f , is lowered (increased sensitivity). The initial viral titer peak is reduced due to IFN presence, withlower resistance (smaller f ) having a greater impact. However, once IFN starts to decay, its effect rapidlydissipates, leading to a rise, or even a rebound, in the viral load. On its own, IFN does not lead to infectionresolution in this MM.The effect of Abs herein, like elsewhere [9,28,30,44,49], is to enhance infectious virus clearance such that c becomes c + k A A ( t ), where k A represents the neutralization rate of infectious virus by Abs, A ( t ). Figure 8(d–f) shows the effect of IFN+Abs in the MM as the neutralization rate, k A , is increased. Low binding affinityAbs ( k A ≤
200 h − ) cannot clear the infection. As k A increases, viral titer decay rates increase, leading toinfection resolution within 11 dpi–16 dpi. However, the time to infection resolution remains longer than the6 dpi–8 dpi typically observed for IAV infections in humans [2, 13].The effect of CTLs herein, like elsewhere [9, 28, 44, 49], is to increase the rate of loss of infected cellsexpressing IAV peptides on their MHC-1. Infected cells begin expressing IAV peptides ∼ k C C ( t ), where k C representsthe killing rate of infected cells by CTLs, C ( t ). Figure 8(g–i) shows the effect of IFN+Abs+CTLs in theMM as the killing rate, k C , is increased. Higher CTL killing rates ( k C ) cause the infection to resolve earlier.12able 2: Parameter values used to reproduce the IAV infection kinetics in Figure 10.IAV strain β ([ V ] − · d − ) p ([ V ] · d − ) τ E (h) τ I (h) IFN, f Abs ∗ , A CTLs, k C (h − )seasonal 1 . × − . × . . × − . × − . × . . × − —avian (late) 1 . × − . ×
12 100 5 . . × − —* The rate of infectious virus neutralization by Abs was fixed to k A = 500 h − for all 3 infections.Figure 9 shows the MM predictions and experimental data of immune response knockout experiments,wherein one component of the immune response is suppressed or neutralized. The MM prediction of an IFNknockout experiment is in good agreement with the experimental data of IFN knockout experiments [38, 66]at early time as both suggest IFN acts early and helps to reduce the initial viral titer peak. The experimentalAb knockout experiment [41] suggests Abs help reduce viral load at a later time in the infection but doesnot affect the initial viral titer peak. The MM prediction also suggests Abs help reduce the viral load ata later time but also affects the initial viral titer peak. This could be better reproduced by the MM byhaving a lower initial amount of Abs ( A ) so that Abs appear later and thus act later. Finally, the MMprediction of the CTL knockout experiment is in good agreement with the experimental data of CTL knockoutexperiments [43, 77, 81] as both suggest that CTLs act only at a late stage in the infection, helping reduceinfection duration. Generally, the experimental knockout appear to show a greater impact for the knock-outs than predicted by the MM. This is likely due to the intricate interaction network between the differentimmune response components in vivo which makes it difficult to experimentally disrupt one componentwithout affecting another. The difference in infection severity between patients naturally infected with seasonal IAV strain or avian-origin H5N1 strains is thought to be the result of a number of different possible factors, including more rapidinfection and replication kinetics, cell tropism, lack of pre-existing immunity, and/or an aberrant immuneresponse [14, 16, 19, 39, 66–68, 72, 74, 75]. Figure 10(a) shows viral load measurements from pharyngeal swabsof patients (each point is a single measure from a single patient) naturally infected with either a seasonal(H3N2 or H1N1) IAV strain or a strain of the H5N1 subtype. This study wherein these measures weretaken over the same time period and following the same methodology, such that viral loads for infectionwith seasonal and H5N1 strains can be readily compared, is the only one of its kind. Since patients wererecruited into the study days after illness onset, and began receiving antiviral therapy immediately uponhospital admission, untreated infection kinetics leading up to and after admission are not available. Themeasured pharyngeal viral load (determined via qRT-PCR) is ∼ or PFU) which was not measured as part of this study, although the ratio of infectious to total IAVis known to change over the course of an infection, and inconsistently so between experiments [58]. Anotherfeature is that H5N1-infected patients reported to the hospital 3 d–4 d later after illness onset than thoseinfected with seasonal IAV strains. It is unclear whether this is due to a slower progression of infection withan H5N1 strain, or because patients infected with H5N1 strains mostly came from remote provinces andtook longer to reach the city hospital than those infected with seasonal strains which came from the cityitself or neighbouring provinces.From the data in Figure 10(a), two hypothetical infection time courses or portraits were developed forinfection with IAV H5N1 strains, and one for seasonal strains, shown in Figure 10(b). The portrait forinfection with a seasonal strain peaks at ∼ –10 TCID / mL at ∼ ∼ V i r u s ( c D NA c o p i e s / m L ) (a) seasonal strainavian strain V i r u s ( TC I D / m L ) (b) seasonal strainavian strain - early peakavian strain - late peak V i r u s ( TC I D / m L ) (c) k E = 50 h − (early peak) k E = 50 h − (late peak) no CTL (early peak)no CTL (late peak) Figure 10:
Kinetics of infection in patients naturally infected with a seasonal or avian IAVstrain. (a) Total viral load measurements (cDNA/mL via qRT-PCR) from throat (pharyngeal) swabs ofpatients which naturally contracted infections with either a seasonal (H3N2 or H1N1) or avian (H5N1) IAVstrain (data taken from [19]). The shaded polygons trace out the extent of the data points. (b) Proposedtime courses for in vivo IAV infections with either a seasonal (grey) or H5N1 (black) IAV strain. The solidand dashed black lines represent two plausible time courses for infection with an IAV H5N1 strain: a morerapid rise to higher viral titers than for infections with seasonal strains (solid); or a more moderate rise,similar to that seen for infections with seasonal strains, that grows to higher titer due to lack of rapid andeffective immune control (dashed). MM parameters for these curves are listed in Table 2. (c) The possibleeffect of a delayed CTL response, i.e. one which peaks at 12 dpi with a killing efficacy of k C = 50 h − (grey)in both the early (solid) and late (dashed) infection time courses which are otherwise shown without a CTLresponse (black). 14ith a seasonal strain, but that growth is sustained over a longer period and thus peaks later, reaching highertiters than infection with a seasonal strain. For infection with an H5N1 strain, viral titer peak was chosen tobe 100-fold higher, as observed for total viral loads, and the reported 3 dpi–4 dpi delay in hospital admissionwas captured either as a longer, sustained infection in the early time course, or as a longer delay to reachpeak titer in the late time course. The IAV H5N1 portraits were constructed from the MM of the seasonalportrait by shifting parameters controlling cell-virus interactions and immunity in keeping with differencesbelieved to be responsible for that shift. These parameters are presented in Table 2. The shaded areas,which depict the extent of the data points in Figure 10(a) where time is based on days since illness onset, areshown time-shifted in Figure 10(b) where time is days post-infection. For a seasonal infection, the (pale grey)area was shifted to one day later since there is generally a delay of one day between infection and symptomonsets [13]. For avian infections, the (dark grey) area was shifted to 4 days later since reports suggest ∼ p ), increased virusinfectivity ( β ), and a shorter eclipse phase ( τ E ), consistent with shifts in these parameters estimated fromin vitro infections of A549 cells with seasonal (H1N1) versus H5N1 and H7N9 IAV strains [68]. In contrast,the late time course relies on lower virus infectivity and a longer eclipse phase, along with a higher virusproduction rate. The early and late portraits also depend on 8- and 10-fold longer periods of virus production(infectious cell lifespan, τ I ), respectively, compared to that for cells infected with a seasonal strain. Incapturing differences in immunity, both the early and late portraits of infection with IAV H5N1 strains relyon 10-fold increased resistance to the effect of IFN ( f ), decrease in pre-existing immunity in the form ofa 200-fold decrease in the neutralizing effect of Abs, and the absence of a CTL response. The resistance ofH5N1 strains to the effect of IFN- α and - γ has been reported in vitro in porcine cell cultures and in vivoin pigs [66]. The decrease in the neutralizing efficacy of Abs is captured in the simplistic Ab MM utilizedherein as a 200-fold decrease in their initial quantity ( A ) which means a longer time ( ∼ k A )of Abs. A recruitment study of patients naturally infected with avian H7N9 IAV strains reported thatthe patient group with the earliest recovery, had a prominent CTL response by 10 d post admission and aprominent Ab response, 2–3 d later [76]. In contrast, the patient group with the most delayed recovery hadan Ab and CTL response that remained low even after 20 d post admission. This suggests that an infectionwith an H5N1 strain, characterized by a prolonged viral shedding period, would have delayed Ab and CTLresponses. Because the portraits in Figure 10(b) only go up to 12 dpi, the effect of these delayed CTLs wouldnot be apparent. Figure 10(c) illustrates what the time course of infection with an avian strain could looklike in the presence of a delayed CTL response, but one which would be as effective (same k C ) as that forinfection with a seasonal strain. It shows an important role for CTLs in controlling the persistent sheddingin these infections, consistent with the findings that a delayed CTL response correlated with slower recoveryfrom infection with H7N9 strains [76].Figure 11 explores the impact of antiviral therapy with neuraminidase inhibitors (NAIs), administeredat various times post-infection, on the time courses for infection with seasonal or avian IAV strains. SinceNAIs block the release of newly produced virions, their effect is implemented in the MM, as elsewhere[8,11,20–22,31,45,54], as a reduction in the virus production rate, namely (1 − ε NAI ) p , from the time treatmentis administered, where ε NAI ∈ [0 ,
1] is the drug efficacy. The efficacy of NAIs was set to ε NAI = 0 .
98 to studythe case of treatment with a high efficacy. Table 3 quantitatively compares the impact of NAI treatment forvarious endpoints: reduction in resolution time, peak viral titer, and area under the viral titer curve (AUC).Human volunteer studies in patients experimentally infected with seasonal (H1N1) IAV strains and treatedwith NAIs report that early treatment (24–32 h post infection) is effective in reducing viral load [3, 32, 33],and one such study reports that delayed treatment (50 h post infection) is less effective [33], in line withthe MM results. A study of H5N1-infected patients recruited days after illness onset, and beginning NAItreatment immediately after admission, reports that late treatment (4–8 d after illness onset) is still effectivein some patients, reducing viral load and thus possibly also reducing time to infection resolution. This wouldbe consistent with some avian strain-infected patients experiencing an infection time course similar to theearly portrait with moderate to no benefit from delayed NAI treatment, and others experiencing infection15 V i r u s ( TC I D / m L ) seasonal (a) No treatment Treatment at 2 dpi Treatment at 3 dpi Treatment at 6 dpi V i r u s ( TC I D / m L ) avian early (b) V i r u s ( TC I D / m L ) avian late (c) Figure 11:
MM-predicted NAI antiviral therapy efficacy in patients infected with a seasonalor avian IAV strain.
The viral titer time course for IAV infections under antiviral therapy initiated atvarious times post infection (see legend) with NAIs, captured as decreasing virus production rate. The effectof treatment is shown in the context of infection with (a) a seasonal IAV strain; or with an avian IAV strainbased on either the (b) early or (c) late infection portraits. The drug efficacy, ε , was chosen to be 0 .
98. Theother parameters used in the MM are those listed in Table 2.Table 3: Measures of the impact of NAI treatment under various conditions. ∗ t adm (dpi) ∆[Resolution time (d)] ∆[AUC (d · TCID / mL)] ∆[Viral titer peak(TCID / mL)]—Seasonal —2 7 . − . . . / . = 10 . . / . = 10 . . − . . . / . = 10 . . / . = 10 . . − . . . / . = 10 . . / . = 10 . — Avian (early) —2 ND −
17 = ND 10 . / . = 10 . . / . = 10 . −
13 = ND 10 . / . = 10 . . / . = 10 . −
12 = ND 10 . / . = 10 . . / . = 10 . — Avian (late) —2 ND −
15 = ND 10 . / . = 10 . . / . = 10 . −
17 = ND 10 . / . = 10 . . / . = 10 . −
16 = ND 10 . / . = 10 . . / . = 10 . ∗ For NAI therapy administered at various times post infection ( t adm ), the change (∆) in 3 commonlyreported endpoints are shown for the viral titer curves in Figure 11. The change in resolution time iscomputed as a difference (∆ = without − with), whereas that for the area under the viral titer curve (AUC)and viral titer peak corresponds to the fold-change (∆ = without / with). ND stands for not determined incases where infection resolution ( V ( t ) < / mL) was not achieved without treatment.16ore consistent with the late portrait and benefiting from NAI therapy even when treatment initiation isdelayed. Mathematical models (MMs) for the course of an influenza A virus (IAV) infection in vivo typically assumethe infection is spatially homogeneous, i.e. that all cells see all virus, and vice-versa, instantly over all space.This simplification is reasonable in vitro in cell culture infections whose spatial extent is often no more thanone or two square centimetres [37], but it is unclear whether it remains appropriate at the scale of the entireHRT. With simplicity and parsimony in mind, the spatial MM introduced herein represents the HRT as aone-dimensional track that extends from the nose down to a depth of 30 cm. It implements two modes ofviral transport: advection of virus upwards towards the nose, and diffusion of the virus within the periciliaryfluid.When diffusion alone is considered, infection in the spatial MM proceeds at a slower pace than in thenon-spatial MM, but virus production is sustained for longer. The diffusion initially acts to spatially restrictthe number of cells available to the infection, and then releases these cells progressively as the diffusinginfection wave reaches ever further in the HRT. When advection is added to diffusion, the former dominatesthe infection kinetics, requiring a 10-fold increase in the virus production rate to restore the timing andlevel of viral titer peak to that in the non-spatial MM. This is noteworthy, firstly, because it shows thatadvection constitutes an effective physiological mechanism to suppress infection. Secondly, it shows thatuse of a non-spatial MM to analyze infection data possibly underestimates the virus production rate, andconsequently also the total amount of virus produced over the course of an infection. Since the more virusare produced, the more mutations accumulate [24, 31, 57], underestimating the virus production rate meansunderestimating the rate and likelihood of emergence of drug resistance.Interestingly, the MM suggests that the depth at which the initial virus inoculum deposits also plays amajor role in determining the extent of HRT involvement in the infection. Specifically, the MM predicts notarget cell consumption below the depth at which the initial virus inoculum deposits. While deposition atdepths lower than ∼
15 cm did not affect the viral titer time course, initial depositions above that point bothdelayed and reduced peak titer. As such, higher virus production rates would be required to compensatefor shallower deposition depths to obtain equivalent viral titer time courses. This is one more way a non-spatial MM might underestimate the rate of virus production, and thus the likelihood of antiviral resistanceemergence.A more unexpected finding was that although kinetically in the spatial MM the infection deposits atsome depth and is dragged upwards by advection, the infection appears to be moving from the upper tothe lower HRT. This is because the upper HRT is downstream of virus advection and sees most of thevirus produced, while the lower HRT is upstream and sees little or none. Therefore, infection progresses,peaks and resolves faster in the upper HRT and slower in the lower HRT. It is this difference in time scalewhich makes it appear as though the infection is moving from the upper to the lower HRT. Experimentalstudies of mice infected with IAVs engineered to be fluorescent or bioluminescent observe infection spreadfrom the upper towards the lower regions of the lung [46, 73], consistent with the MM’s predictions. In theabsence of sufficient immune control, this could serve to maintain the infection with, wherein the slower,lower HRT infection could re-ignite more rounds of infection in the upper HRT, depending on the delicatetiming between infection time course in the lower HRT and cell regeneration in the upper HRT.The spatial MM was expanded to include density-dependent cell regeneration, and a simplistic immuneresponse comprising IFN acting to down-regulate the rate of virus production by infected cells, Abs neu-tralizing infectious virions, and CTLs killing infected cells. In the presence of this immune response, theMM-simulated IAV infection was well controlled and resolved fully by 8 dpi. The spatial MM could alsolargely reproduce the key features of in vivo immune response knockout experiments [23]. In the MM-simulated infection (see Figure 8g), 10% of the HRT was involved in the infection by 3 dpi, around 6 dpithe fraction of cells involved in the infection peaked at 30%, and by 10 dpi 10% of the epithelium had yet toregenerate, although viral titers had fallen below the detection limit (see Figure 8i). These numbers align17ell with a 1989 report in Russian cited by Bocharov and Romanyukha in their 1994 IAV MM paper [9]of 10% damage at symptom onset, 30–50% of the upper airway destroyed at the peak of the disease, andresolution of disease at a time when as much as 10% of the normal epithelium is still damaged.The full spatial MM was applied to simulate the kinetics of infection with either a mild, seasonal IAVstrain or a severe infection with an avian strain such as H5N1. Since complete, untreated, infection kinetictime courses for infection with avian IAV strains are not available, two different portraits were constructedto represent plausible time courses: one peaking early with sustained, high titers over several days, andanother rising slowly, similarly to seasonal infections, but continuing on to reach high titers 3–4 days later.Differences in the time courses for infection with a seasonal vs avian IAV strain could be captured by the MMby shifting the parameters in ways that are consistent with known or expected differences between infectionswith seasonal and avian strains. The portraits were used to evaluate the impact of antiviral therapy withneuraminidase inhibitors. Treatment was most effective for the seasonal and late peak avian IAV strainportraits when administered early during the infection at 2 dpi. Later treatment was of limited effectivenessfor the seasonal IAV strain portrait, but was still effective in the late peaking portrait of infection with anavian IAV strain.The spatial MM presented herein, despite its simplicity, offers a novel and interesting opportunity to studythe spatial localization and dissemination of IAV infections within the HRT. In the future, the addition ofan explicit, dynamical immune response — to replace the empirical equations used herein as stand-ins forimmune response components — would enable the use of this MM to reproduce the kinetics of re-infection,similar to the work done by others [12, 79].
S1 Video.
Spatiotemporal course of IAV infection in the presence of diffusion only (no advection).
S2 Video.
Spatiotemporal course of IAV infection in the presence of diffusion and advection.
This work was supported in part by Discovery Grant 355837-2013 (to CAAB) from the Natural Sci-ences and Engineering Research Council of Canada ( ), Early Researcher AwardER13-09-040 (to CAAB) from the Ministry of Research and Innovation of the Government of Ontario( ), and by Interdisciplinary Theoretical and Mathemat-ical Sciences (iTHEMS, ithems.riken.jp ) at RIKEN (CAAB).
CAAB received financial support in the form of a research contract and a consultancy fee from F. Hoffmann-La Roche Ltd. in the early stages of this project. MBR was employed by F. Hoffmann-La Roche Ltd. whenfirst participating in this work, and is currently employed by Array BioPharma Inc.. Neither company hadany role in the study design, data collection and analysis, or decision to publish.18
Methods
The Euler method was used to numerically solve the cell population equations, namely T n +1 m = T nm − ∆ t βT nm V nm E n +11 ,m = E n ,m + ∆ t (cid:20) βT nm V nm − n E τ E E n ,m (cid:21) E n +1 i,m = E i,m + ∆ t (cid:20) n E τ E E ni − ,m − n E τ E E ni,m (cid:21) i = 2 , , ..., n E I n +11 ,m = I n ,m + ∆ t (cid:20) n E τ E E nn E ,m − n I τ I I n ,m (cid:21) (3) I n +1 j,m = I nj,m + ∆ t (cid:20) n I τ I I nj − ,m − n I τ I I nj,m (cid:21) j = 2 , , ..., n I where T nm = T ( x m , t n ), x m = m ∆ x , t n = n ∆ t , and ∆ x and ∆ t are the chosen spatial and temporal stepsizes, respectively. We define N x = x max / (∆ x ) = 3000, the number of spatial boxes or sites that makeup the simulated HRT such that ∆ x = 0 . / (3000 sites) = 100 µ m / site, the diameter of ∼ N x was chosen by verifying that choosing a larger number of sites (smaller ∆ x ) did not affectthe solution. When presenting the fraction of target cells or infectious cells, or the virus concentration, as afunction of time only, averaged over space, those are computed as T ( t ) = T n = 1 N x N x (cid:88) m =1 T nm I ( t ) = I n = 1 N x N x (cid:88) m =1 n I (cid:88) j =1 I nj,m V ( t ) = V n = 1 N x N x (cid:88) m =1 V nm For the virus concentration, for simplicity, the diffusion and advection terms are each treated separately,and are applied before the production and clearance terms are considered. To solve the diffusion term, theCrank-Nicolson method was used, namely ∂V ( x, t ) ∂t = D PCF ∂ V ( x, t ) ∂x V n +1 m − V nm ∆ t = D PCF (cid:34) ( V n +1 m +1 − V n +1 m + V n +1 m − ) + ( V nm +1 − V nm + V nm − )(∆ x ) (cid:35) − αV n +1 m − + (2 + 2 α ) V n +1 m − αV n +1 m +1 = αV nm − + (2 − α ) V nm + αV nm +1 where α = D PCF (∆ x ) / (∆ t ). The rate of viral diffusion in the PCF was estimated based on the Stokes-Einstein equation for IAV diffusing in plasma at body temperature as D PCF ≈ − m / s [6, 37]. Anabsorbing boundary condition was used at the top of the HRT, V (0 , t ) = 0 (or V n = 0), to capture virusflow out through the mouth and nose. A reflective boundary condition, V ( x max + ∆ x, t ) = V ( x max , t ) (or V nN x +1 = V nN x ), was used at the bottom of the HRT. The bottom boundary condition becomes irrelevantonce advection is introduced as the flow at x = x max becomes negligible. With these boundary conditionsin place, the diffusion of virus over a step size ∆ t can be expressed as19 V n +11 V n +12 ... V n +1 N x − V n +1 N x = α − α . . . − α α − α . . . ...0 . . . . . . . . . 0... . . . − α α − α . . . − α α − − α α . . . α − α α . . . ...0 . . . . . . . . . 0... . . . α − α α . . . α − α V n V n ... V nN x − V nN x The virus advection term can be written as follows ∂V ( x, t ) ∂t = v a ∂V ( x, t ) ∂xV n +1 m − V nm ∆ t = v a V nm +1 − V nm ∆ x (4) V n +1 m = (cid:18) − v a ∆ t ∆ x (cid:19) V nm + (cid:18) v a ∆ t ∆ x (cid:19) V nm +1 . This simplistic numerical scheme is known to lead to a dispersion (diffusion) of the solution [52]. However,the spatial MM simulator imposes v a ∆ t/ (∆ x ) = 1, such that Eqn. (4) simplifies to V n +1 m = V nm +1 , atrivial translation of the solution, which ensures the latter will remain stable and will not disperse. Thespeed of the PCF is believed to be approximately the same as that of the mucus blanket and is estimatedto be v a ≈ µ m / s [47]. This requires the time step for the spatial MM simulator to be ∆ t = ∆ x/v a =(100 µ m) / (40 µ m / s) = 2 . V n +11 = V n ), and the boundary condition at the bottom is chosen suchthat V nN x = 0, i.e. there is no virus beyond the end of the HRT. Once diffusion and advection have beenapplied, the Euler method is used to solve the remaining terms of the virus equation, namely V n +1 m = V nm + ∆ t (cid:20) p n I (cid:88) j =1 I nj,m − cV nm (cid:21) (5)The droplet-like, initial IAV distribution V ( x, t = 0) is represented by a Gaussian centred at the site ofdeposition, V m = V ( x m , t = 0) = V ∗ √ πσ exp (cid:18) − ( x m − x d ) σ (cid:19) , (6)where x d is the site of deposition, σ = 0 . V ∗ is chosen so that (cid:104) V ( x, t = 0) (cid:105) = N x (cid:80) N x m =1 V m =7 . × − TCID / mL, i.e. the spatial average of the initial virus inoculum concentration in the spatial MMherein is equal to V ( t = 0) in the non-spatial ODE MM by Baccam et al. [1].While most parameters were taken directly from [1], some were adapted. Whereas in Baccam et al. [1] T + E + I = 4 × corresponds to the number of cells, in the spatial MM herein, fraction of cells in eachstate are considered instead such that T + E + I = 1. As such, what we consider the Baccam et al. [1] virusproduction rate is converted from that reported in their paper as p Baccam = (cid:18) .
046 (TCID / mL)d · cell (cid:19) × (cid:0) × cells (cid:1) × (cid:18) (cid:19) = 7 . × (TCID / mL) · h − . Once advection is introduced (see Results), a value of 11 × p Baccam is used so that the virus titer will peakat roughly 2–3 dpi, consistent with that observed in Baccam et al. [1].20 .2 Empirical, mathematical representation of the immune response
The time course for interferon (IFN) is represented as F ( t ) = 2e − λ g ( t − t p ) + e λ d ( t − t p ) , (7)where λ g and λ d are the growth and decay rates of IFN respectively, and t p is the time of peak. Since F ( t ) ∈ [0 , F ( t ) represents the fractional amount of IFN, relative to peak IFN, at time t . In the MM, λ g = 2 d − , λ d = 1 d − and t p = 3 dpi to match experimental data. Its effect in the MM is to attenuatevirus production, p , captured as (1 − ε IFN ) p = (cid:18) − FF + f (cid:19) p, (8)where f is the amount of IFN required to reduce the virus production rate to one half its normal value,i.e. p → p when F = f . In the MM, F ( t ) is normalized so as to have a maximum value of 1 ( F max = 1).As such, if f = 0 .
1, then p → p when F = 0 . F max = 0 . A ( t ) = 11 + (cid:16) A − (cid:17) e − αt , (9)where A is the initial amount of Abs, and α is the growth rate of Abs. Since A ( t ) ∈ [0 , A ( t ) representsthe fractional amount of Abs, relative to peak Abs, at time t . In the MM, α = 0 .
75 d − and A = 2 × − to match experimental data. The effect of Abs in the MM is to enhance virus clearance, c , captured as c + k A A ( t ) (10)where k A (h − ) represents the binding affinity of Abs.The time course of cytotoxic T lymphocytes (CTLs) is represented as C ( t ) = 2e − λ (cid:48) g ( t − t (cid:48) p ) + e λ (cid:48) d ( t − t (cid:48) p ) , (11)where λ (cid:48) g and λ (cid:48) d are the growth and decay rates of CTLs, respectively, and t (cid:48) p is the time of peak. Since C ( t ) ∈ [0 , C ( t ) represents the fractional amount of CTLs, relative to peak CTLs, at time t . In the MM, λ (cid:48) g = 2 d − , λ (cid:48) d = 0 . − and t (cid:48) p = 8 dpi to match experimental data. The effect of CTLs in the MM is toincrease the rate of loss of IAV-infected cells which CTLs recognize as infected. Newly infected cells begin topresent IAV peptides on their MHC-1 for CTL recognition ∼ τ E / E i =[1 ,n E / wouldnot be presenting peptides, while those in E i = n E / ,n E would be recognizable and thus could be killed byCTLs. As such, killing of recognizably-infected cells by CTLs in the MM is captured as ∂E i ( x, t ) ∂t = · · · − i = 1 , , ..., n E ∂E i ( x, t ) ∂t = · · · − k C C ( t ) E i ( x, t ) for i = n E , ..., n E (12) ∂I j ( x, t ) ∂t = · · · − k C C ( t ) I j ( x, t ) for j = 1 , , ..., n I where k C (h − ) provides the scale for the rate of infected cell killing by CTLs.21 eferences [1] P. Baccam, C. Beauchemin, C. A. Macken, F. G. Hayden, and A. S. Perelson. Kinetics of influenza Avirus infection in humans. J. Virol. , 80(15):7590–7599, 2 August 2006. doi:10.1128/JVI.01623-05 .[2] S. Balasingam and A. Wilder-Smith. Randomized controlled trials for influenza drugs and vaccines:A review of controlled human infection studies.
Int. J. Infect. Dis. , 49:18–29, August 2016. doi:10.1016/j.ijid.2016.05.013 .[3] L. Barroso, J. Treanor, L. Gubareva, and F. G. Hayden. Efficacy and tolerability of the oral neu-raminidase inhibitor peramivir in experimental human influenza: Randomized, controlled trials forprophylaxis and treatment.
Antivir. Ther. , 10(8):901–910, 2005.[4] A. L. Bauer, C. A. A. Beauchemin, and A. S. Perelson. Agent-based modeling of host-pathogen systems:The successes and challenges.
Inform. Sciences , 179(10):1379–1389, 29 April 2009. doi:10.1016/j.ins.2008.11.012 .[5] C. Beauchemin. Probing the effects of the well-mixed assumption on viral infection dynamics.
J. Theor.Biol. , 242(2):464–477, 21 September 2006. doi:10.1016/j.jtbi.2006.03.014 .[6] C. Beauchemin, S. Forrest, and F. T. Koster. Modeling influenza viral dynamics in tissue. In H. Bersiniand J. Carneiro, editors,
Proceedings of the th International Conference on Artificial Immune Systems(ICARIS 06) , number 4163 in Lecture Notes in Computer Science, pages 23–36. Springer-Verlag BerlinHeidelberg, 2006. doi:10.1007/11823940_3 .[7] C. Beauchemin, J. Samuel, and J. Tuszynski. A simple cellular automaton model for influenza A viralinfections.
J. Theor. Biol. , 232(2):223–234, 21 January 2005. doi:10.1016/j.jtbi.2004.08.001 .[8] N. F. Beggs and H. M. Dobrovolny. Determining drug efficacy parameters for mathematical models ofinfluenza.
J. Biol. Dyn. , 9(sup1):332–346, 9 June 2015. doi:10.1080/17513758.2015.1052764 .[9] G. A. Bocharov and A. A. Romanyukha. Mathematical model of antiviral immune response III. InfluenzaA virus infection.
J. Theor. Biol. , 167(4):323–360, 21 April 1994. doi:10.1006/jtbi.1994.1074 .[10] A. Boianelli, V. K. Nguyen, T. Ebensen, K. Schulze, E. Wilk, N. Sharma, S. Stegemann-Koniszewski,D. Bruder, F. R. Toapanta, C. A. Guzmn, M. Meyer-Hermann, and E. A. Hernandez-Vargas. Modelinginfluenza virus infection: A roadmap for influenza research.
Viruses , 7(10):5274–5304, 12 October 2015. doi:10.3390/v7102875 .[11] P. Cao and J. M. McCaw. The mechanisms for within-host influenza virus control affect model-basedassessment and prediction of antiviral treatment.
Viruses , 9(8):E197, 26 July 2017. doi:10.3390/v9080197 .[12] P. Cao, A. W. Yan, J. M. Heffernan, S. Petrie, R. G. Moss, L. A. Carolan, T. A. Guarnaccia, A. Kelso,I. G. Barr, J. McVernon, K. L. Laurie, and J. M. McCaw. Innate immunity and the inter-exposure inter-val determine the dynamics of secondary influenza virus infection and explain observed viral hierarchies.
PLoS Comput. Biol. , 11(8):e1004334, 18 August 2015. doi:10.1371/journal.pcbi.1004334 .[13] F. Carrat, E. Vergu, N. M. Ferguson, M. Lemaitre, S. Cauchemez, S. Leach, and A.-J. Valleron. Timelines of infection and disease in human influenza: A review of volunteer challenge studies.
Am. J.Epidemiol. , 167(7):775–785, 1 April 2008. doi:10.1093/aje/kwm375 .[14] M. Chan, C. Cheung, W. Chui, S. Tsao, J. Nicholls, Y. Chan, R. Chan, H. Long, L. Poon, Y. Guan,and J. Peiris. Proinflammatory cytokine responses induced by influenza A (H5N1) viruses in primaryhuman alveolar and bronchial epithelial cells.
Respir. Res. , 6(1):135, 11 November 2005. doi:10.1186/1465-9921-6-135 . 2215] M. C. Chan, R. W. Chan, W. C. Yu, C. C. Ho, W. Chiu, C. Lo, K. M. Yuen, Y. Guan, J. M. Nicholls,and J. M. Peiris. Influenza H5N1 virus infection of polarized human alveolar epithelial cells and lungmicrovascular endothelial cells.
Respir. Res. , 10(1), 30 October 2009. doi:10.1186/1465-9921-10-102 .[16] C. Y. Cheung, L. L. M. Poon, A. S. Lau, W. Luk, Y. L. Lau, K. F. Shortridge, S. Gordon, Y. Guan,and J. S. M. Peiris. Induction of proinflammatory cytokines in human macrophages by influenza A(H5N1) viruses: A mechanism for the unusual severity of human disease?
Lancet , 360(9348):1831–1837,7 December 2002. doi:10.1016/s0140-6736(02)11772-7 .[17] T. Chotpitayasunondh, K. Ungchusak, W. Hanshaoworakul, S. Chunsuthiwat, P. Sawanpanyalert, R. Ki-jphati, S. Lochindarat, P. Srisan, P. Suwan, Y. Osotthanakorn, T. Anantasetagoon, S. Kanjanawasri,S. Tanupattarachai, J. Weerakul, R. Chaiwirattana, M. Maneerattanaporn, R. Poolsavatkitikool,K. Chokephaibulkit, A. Apisaranthanarak, and S. F. Dowell. Human disease from influenza A (H5N1),Thailand, 2004.
Emerg. infect. Dis. , 11(2):201–209, February 2005. doi:10.3201/eid1102.041061 .[18] L. M. Crosby and C. M. Waters. Epithelial repair mechanisms in the lung.
Am. J. Physiol. Lung CellMol. Physiol. , 298(6):L715–L731, June 2010. doi:10.1152/ajplung.00361.2009 .[19] M. D. de Jong, C. P. Simmons, T. T. Thanh, V. M. Hien, G. J. D. Smith, T. N. B. Chau, D. M. Hoang,N. V. V. Chau, T. H. Khanh, V. C. Dong, P. T. Qui, B. V. Cam, D. Q. Ha, Y. Guan, J. S. M. Peiris, N. T.Chinh, T. T. Hien, and J. Farrar. Fatal outcome of human influenza A (H5N1) is associated with highviral load and hypercytokinemia.
Nat. Med. , 12(10):1203–1207, October 2006. doi:10.1038/nm1477 .[20] H. M. Dobrovolny, M. J. Baron, R. Gieschke, B. E. Davies, N. L. Jumbe, and C. A. A. Beauchemin.Exploring cell tropism as a possible contributor to influenza infection severity.
PLoS One , 5(11):e13811,23 November 2010. doi:10.1371/journal.pone.0013811 .[21] H. M. Dobrovolny and C. A. A. Beauchemin. Modelling the emergence of influenza drug resistance: Theroles of surface proteins, the immune response and antiviral mechanisms.
PLoS One , 12(7):e0180582,10 July 2017. doi:10.1371/journal.pone.0180582 .[22] H. M. Dobrovolny, R. Gieschke, B. E. Davies, N. L. Jumbe, and C. A. A. Beauchemin. Neuraminidaseinhibitors for treatment of human and avian strain influenza: A comparative study.
J. Theor. Biol. ,269(1):234–244, 21 January 2011. doi:10.1016/j.jtbi.2010.10.017 .[23] H. M. Dobrovolny, M. B. Reddy, M. A. Kamel, C. R. Rayner, and C. A. A. Beauchemin. Assessing math-ematical models of influenza infections using features of the immune response.
PLoS One , 8(2):e57088,28 February 2013. doi:10.1371/journal.pone.0057088 .[24] J. W. Drake. Rates of spontaneous mutation among RNA viruses.
Proc. Natl. Acad. Sci. U.S.A. ,90(9):4171–4175, 1 May 1993. doi:10.1073/pnas.90.9.4171 .[25] J. S. Erjef¨alt, I. Erjef¨alt, F. Sundler, and C. G. A. Persson. In vivo restitution of airway epithelium.
Cell Tissue Res. , 281(2):305–316, August 1995. doi:10.1007/bf00583399 .[26] J. V. Fahy and B. F. Dickey. Airway mucus function and dysfunction.
N. Engl. J. Med. , 363(23):2233–2247, 2 December 2010. doi:10.1056/NEJMra0910061 .[27] R. S. Fritz, F. G. Hayden, D. P. Calfee, L. M. R. Cass, A. W. Peng, W. G. Alvord, W. Strober, and S. E.Straus. Nasal cytokine and chemokine response in experimental influenza A virus infection: Results of aplacebo-controlled trial of intravenous zanamivir treatment.
J. Infect. Dis. , 180(3):586–593, September1999. doi:10.1086/314938 .[28] B. Hancioglu, D. Swigon, and G. Clermont. A dynamical model of human immune response to influenzaA virus infection.
J. Theor. Biol. , 246(1):70–86, 7 May 2010. doi:10.1016/j.jtbi.2006.12.015 .2329] A. Handel, L. E. Liao, and C. A. A. Beauchemin. Progress and trends in mathematical modelling ofinfluenza A virus infections.
Curr. Opin. Syst. Biol. , 12:30–36, December 2018. doi:10.1016/j.coisb.2018.08.009 .[30] A. Handel, I. M. Longini, and R. Antia. Towards a quantitative understanding of the within-hostdynamics of influenza A infections.
J R Soc Interface. , 7(42):35–47, 6 January 2010. doi:10.1098/rsif.2009.0067 .[31] A. Handel, I. M. Longini, Jr., and R. Antia. Neuraminidase inhibitor resistance in influenza: Assessingthe danger of its generation and spread.
PLoS Comput. Biol. , 3(12):e240, December 2007. doi:10.1371/journal.pcbi.0030240 .[32] F. G. Hayden, J. J. Treanor, R. F. Betts, M. Lobo, J. D. Esinhart, and E. K. Hussey. Safety andefficacy of the neuraminidase inhibitor GG167 in experimental human influenza.
JAMA. , 275(4):295–299, January 1996.[33] F. G. Hayden, J. J. Treanor, R. S. Fritz, M. Lobo, R. F. Betts, M. Miller, N. Kinnersley, R. G. Mills,P. Ward, and S. E. Straus. Use of the oral neuraminidase inhibitor oseltamivir in experimental humaninfluenza: Randomized controlled trials for prevention and treatment.
JAMA. , 282(13):1240–1246, 6October 1999. doi:10.1001/jama.282.13.1240 .[34] F. G. Hayden, A. R. Tunkel, J. J. Treanor, R. F. Betts, S. Allerheiligen, and J. Harris. Oral LY217896for prevention of experimental influenza A virus infection and illness in humans.
Antimicrob. AgentsChemother. , 38(5):1178–1181, May 1994. doi:10.1128/aac.38.5.1178 .[35] T. T. Hien, N. T. Liem, N. T. Dung, L. T. San, P. P. Mai, N. van Vinh Chau, P. T. Suu, V. C. Dong,L. T. Q. Mai, N. T. Thi, D. B. Khoa, L. P. Phat, N. T. Truong, H. T. Long, C. V. Tung, L. T. Giang,N. D. Tho, L. H. Nga, N. T. K. Tien, L. H. San, L. V. Tuan, C. Dolecek, T. T. Thanh, M. de Jong,C. Schultsz, P. Cheng, W. Lim, and P. Horby. Avian influenza A (H5N1) in 10 patients in Vietnam.
N.Eng. J. Med. , 350(12):1179–1188, 18 March 2004. doi:10.1056/NEJMoa040419 .[36] B. P. Holder and C. A. A. Beauchemin. Exploring the effect of biological delays in kinetic modelsof influenza within a host or cell culture.
BMC Public Health , 11(Suppl 1):S10, 25 February 2011. doi:10.1186/1471-2458-11-S1-S10 .[37] B. P. Holder, L. E. Liao, P. Simon, G. Boivin, and C. A. A. Beauchemin. Design considerations inbuilding in silico equivalents of common experimental influenza virus assays.
Autoimmunity , 44(4),June 2011. doi:10.3109/08916934.2011.523267 .[38] A. Hoshino, H. Takenaka, O. Mizukoshi, J. Imanishi, T. Kishida, and M. G. Tovey. Effect of anti-interferon serum of influenza virus infection in mice.
Antiviral Res. , 3(1):59–65, March 1983.[39] S.-M. Hsieh and S.-C. Chang. Insufficient perforin expression in CD8+ T cells in response tohemagglutinin from avian influenza (H5N1) virus.
J. Immunol. , 176(8):4530–4533, 15 April 2006. doi:10.4049/jimmunol.176.8.4530 .[40] A. Iwasaki and P. S. Pillai. Innate immunity to influenza virus infection.
Nat. Rev. Immunol. , 14(5):315–328, May 2014. doi:10.1038/nri3665 .[41] T. Iwasaki and T. Nozima. Defense mechanisms against primary influenza virus infection in mice. I.The roles of interferon and neutralizing antibodies and thymus dependence of interferon and antibodyproduction.
J. Immunol , 118(1):256–263, January 1977.[42] K. P. Keenan, T. S. Wilson, and E. M. McDowell. Regeneration of hamster tracheal epithelium aftermechanical injury.
Virchows Arch. B Cell Pathol. Incl. Mol. Pathol. , 43(3):213–240, 1983.2443] R. M. Kris, R. A. Yetter, R. Cogliano, R. Ramphal, and P. A. Small. Passive serum antibody causes tem-porary recovery from influenza virus infection of the nose, trachea and lung of nude mice.
Immunology ,63(3):349–353, March 1988.[44] H. Y. Lee, D. J. Topham, S. Y. Park, J. Hollenbaugh, J. Treanor, T. R. Mosmann, X. Jin, B. M.Ward, H. Miao, J. Holden-Wiltse, A. S. Perelson, M. Zand, and H. Wu. Simulation and prediction ofthe adaptive immune response to influenza A virus infection.
J. Virol. , 83(14):7151–7165, July 2009. doi:10.1128/JVI.00098-09 .[45] L. E. Liao, S. Kowal, D. A. Cardenas, and C. A. A. Beauchemin. Exploring virus release as a bottleneckfor the spread of influenza A virus infection in vitro and the implications for antiviral therapy withneuraminidase inhibitors.
PLoS One , 12(8):e0183621, 24 August 2017. doi:10.1371/journal.pone.0183621 .[46] B. Manicassamy, S. Manicassamy, A. Belicha-Villanueva, G. Pisanelli, B. Pulendran, and A. Garc´ıa-Sastre. Analysis of in vivo dynamics of influenza virus infection in mice using a GFP reporter virus.
Proc. Natl. Acad. Sci. U.S.A. , 107(25):11531–11536, 22 June 2010. doi:10.1073/pnas.0914994107 .[47] H. Matsui, S. H. Randell, S. W. Peretti, C. W. Davis, and R. C. Boucher. Coordinated clearance ofpericiliary liquid and mucus from airway surfaces.
J. Clin. Invest. , 102(6):11251131, 15 September 1998. doi:10.1172/JCI2687 .[48] C. McLaren and G. M. Butchko. Regional T- and B-cell responses in influenza-infected ferrets.
Infect.Immun. , 22(1):189–194, October 1978.[49] H. Miao, J. A. Hollenbaugh, M. S. Zand, J. Holden-Wiltse, T. R. Mosmann, A. S. Perelson, H. Wu, andD. J. Topham. Quantifying the early immune response and adaptive immune response kinetics in miceinfected with influenza A virus.
J. Virol. , 84(13):6687–6698, July 2010. doi:10.1128/JVI.00266-10 .[50] H. Mitchell, D. Levin, S. Forrest, C. A. A. Beauchemin, J. Tipper, J. Knight, N. Donart, R. C. Layton,J. Pyles, P. Gao, K. S. Harrod, A. S. Perelson, and F. Koster. Higher level of replication efficiencyof 2009 (H1N1) pandemic influenza virus than those of seasonal and avian strains: Kinetics fromepithelial cell culture and computational modeling.
J. Virol. , 85(2):1125–1135, January 2011. doi:10.1128/JVI.01722-10 .[51] R. Mora, E. Rodriguez-Boulan, P. Palese, and A. Garc´ıa-Sastre. Apical budding of a recombinantinfluenza A virus expressing a hemagglutinin protein with a basolateral localization signal.
J. Virol. ,76(7):3544–3553, April 2002. doi:10.1128/jvi.76.7.3544-3553.2002 .[52] K. W. Morton and D. F. Mayers.
Numerical solution of partial differential equations: An introduction .Cambridge University Press, Cambridge, UK, 2nd edition, 11 April 2005.[53] A. F. Oner, A. Bay, S. Arslan, H. Akdeniz, H. A. Sahin, Y. Cesur, S. Epcacan, N. Yilmaz, I. Deger,B. Kizilyildiz, H. Karsen, and M. Ceyhan. Avian influenza A (H5N1) infection in eastern Turkey in2006.
N. Engl. J. Med. , 355(21):2179–2185, 23 November 2006. doi:10.1056/NEJMoa060601 .[54] J. Palmer, H. M. Dobrovolny, and C. A. A. Beauchemin. The in vivo efficacy of neuraminidase inhibitorscannot be determined from the decay rates of influenza viral titers observed in treated patients.
Sci.Rep. , 7:40210, 9 January 2017. doi:10.1038/srep40210 .[55] E. G. Paradis, L. T. Pinilla, B. P. Holder, Y. Abed, G. Boivin, and C. A. A. Beauchemin. Impactof the H275Y and I223V mutations in the neuraminidase of the 2009 pandemic influenza virus invitro and evaluating experimental reproducibility.
PLoS One , 10(5):e0126115, 20 May 2015. doi:10.1371/journal.pone.0126115 .[56] J. S. M. Peiris, M. D. de Jong, and Y. Guan. Avian influenza virus (H5N1): A threat to human health.
Clin. Microbiol. Rev. , 20(2):243–267, April 2007. doi:10.1128/CMR.00037-06 .2557] A. S. Perelson, L. Rong, and F. G. Hayden. Combination antiviral therapy for influenza: Predictionsfrom modeling of human infections.
J Infect Dis. , 205(11):1642–1645, June 2012. doi:10.1093/infdis/jis265 .[58] S. M. Petrie, T. Guarnaccia, K. L. Laurie, A. C. Hurt, J. McVernon, and J. M. McCaw. Reducinguncertainty in within-host parameter estimates of influenza infection by measuring both infectious andtotal viral load.
PLoS One , 8(5):e64098, 15 May 2013. doi:10.1371/journal.pone.0064098 .[59] L. T. Pinilla, B. P. Holder, Y. Abed, G. Boivin, and C. A. A. Beauchemin. The H275Y neuraminidasemutation of the pandemic A/H1N1 virus lengthens the eclipse phase and reduces viral output of infectedcells, potentially compromising fitness in ferrets.
J. Virol. , 86(19):10651–10660, October 2012. doi:10.1128/JVI.07244-11 .[60] O. G. Raabe, H. C. Yeh, G. M. Schum, and R. F. Phalen. Tracheobronchial geometry: Human, dog,rat, hamster. Report, Lovelace Foundation, 1976.[61] R. Ramphal, W. Fischlschweiger, J. W. Shands, Jr., and P. A. Small, Jr. Murine influenzal tracheitis: Amodel for the study of influenza and tracheal epithelial repair.
Am. Rev. Respir. Dis. , 120(6):1313–1324,December 1979. doi:10.1164/arrd.1979.120.6.1313 .[62] L. A. Reperant, T. Kuiken, B. T. Grenfell, A. D. M. E. Osterhaus, and A. P. Dobson. Linking influenzavirus tissue tropism to population-level reproductive fitness.
PLoS One , 7(8):e43115, 28 August 2012. doi:10.1371/journal.pone.0043115 .[63] P. D. Reuman, D. I. Bernstein, M. C. Keefer, E. C. Young, J. R. Sherwood, and G. M. Schiff. Efficacy andsafety of low dosage amantadine hydrochloride as prophylaxis for influenza A.
Antivir. Res. , 11(1):27–40,February 1989.[64] J. R. Rock, M. W. Onaitis, E. L. Rawlins, Y. Lu, C. P. Clark, Y. Xue, S. H. Randell, and B. L. M.Hogan. Basal cells as stem cells of the mouse trachea and human airway epithelium.
Proc. Natl. Acad.Sci. U.S.A. , 106(31):12771–12775, 4 August 2009. doi:10.1073/pnas.0906850106 .[65] R. A. Saenz, M. Quinlivan, D. Elton, S. MacRae, A. S. Blunden, J. A. Mumford, J. M. Daly, P. Digard,A. Cullinane, B. T. Grenfell, J. W. McCauley, J. L. N. Wood, and J. R. Gog. Dynamics of influenzavirus infection and pathology.
J. Virol. , 84(8):3974–3983, April 2010. doi:10.1128/JVI.02078-09 .[66] S. H. Seo, E. Hoffmann, and R. G. Webster. Lethal H5N1 influenza viruses escape host anti-viralcytokine responses.
Nat, Med. , 8(9):950–954, September 2002. doi:10.1038/nm757 .[67] K. Shinya, M. Ebina, S. Yamada, M. Ono, N. Kasai, and Y. Kawaoka. Avian flu: Influenza virusreceptors in the human airway.
Nature , 440(7083):435–436, 23 March 2006. doi:10.1038/440435a .[68] P. F. Simon, M.-A. de La Vega, E. Paradis, E. Mendoza, K. M. Coombs, D. Kobasa, and C. A. A.Beauchemin. Avian influenza viruses that cause highly virulent infections in humans exhibit distinctreplicative properties in contrast to human H1N1 viruses.
Sci. Rep. , 6:24154, 15 April 2016. doi:10.1038/srep24154 .[69] D. J. Smith, E. A. Gaffney, and J. R. Blake. Modelling mucociliary clearance.
Respir. Physiol. Neurobiol. ,163(1–3):178–188, 30 November 2008. doi:10.1016/j.resp.2008.03.006 .[70] M. D. Stoneham. The nasopharyngeal airway. Assessment of position by fibreoptic laryngoscopy.
Anaes-thesia , 48(7):575–580, July 1993. doi:10.1111/j.1365-2044.1993.tb07119.x .[71] M. C. Strain, D. D. Richman, J. K. Wong, and H. Levine. Spatiotemporal dynamics of HIV propagation.
J. Theor. Biol. , 218(1):85–96, 7 September 2002. doi:10.1006/jtbi.2002.3055 .2672] C. I. Thompson, W. S. Barclay, M. C. Zambon, and R. J. Pickles. Infection of human airway epitheliumby human and avian strains of influenza A virus.
J. Virol. , 80(16):8060–8068, August 2006. doi:10.1128/JVI.00384-06 .[73] V. Tran, L. A. Moser, D. S. Poole, and A. Mehle. Highly sensitive real-time in vivo imaging of aninfluenza reporter virus reveals dynamics of replication and spread.
J. Virol. , 87(24):13321–13329,December 2013. doi:10.1128/JVI.02381-13 .[74] D. van Riel, V. J. Munster, E. de Wit, G. F. Rimmelzwaan, R. A. M. Fouchier, A. D. M. E. Osterhaus,and T. Kuiken. H5N1 virus attachment to lower respiratory tract.
Science , 312(5772):399, 21 April2006.[75] D. van Riel, V. J. Munster, E. de Wit, G. F. Rimmelzwaan, R. A. M. Fouchier, A. D. M. E. Osterhaus,and T. Kuiken. Human and avian influenza viruses target different cells in the lower respiratory tract ofhumans and other mammals.
Am. J. Pathol. , 171(4):1215–1223, October 2007. doi:10.2353/ajpath.2007.070248 .[76] Z. Wang, Y. Wan, C. Qiu, S. Quinones-Parra, Z. Zhu, L. Loh, D. Tian, Y. Ren, Y. Hu, X. Zhang, P. G.Thomas, M. Inouye, P. C. Doherty, K. Kedzierska, and J. Xu. Recovery from severe H7N9 disease isassociated with diverse response mechanisms dominated by CD8+ T cells.
Nat. Commun. , 6:6833, 13May 2015. doi:10.1038/ncomms7833 .[77] M. A. Wells, P. Albrecht, and F. A. Ennis. Recovery from a viral respiratory infection: 1. Influenzapneumonia in normal and T-deficient mice.
J. Immunol. , 126(3):1036–1041, March 1981.[78] T. Wu, J. Guan, A. Handel, D. C. Tscharke, J. Sidney, A. Sette, L. M. Wakim, X. Y. X. Sng, P. G.Thomas, N. P. Croft, A. W. Purcell, and N. L. La Gruta. Quantification of epitope abundance revealsthe effect of direct and cross-presentation on influenza CTL responses.
Nat. Commun. , 10(1):2846, 28June 2019. doi:10.1038/s41467-019-10661-8 .[79] A. W. C. Yan, S. G. Zaloumis, J. A. Simpson, and J. M. McCaw. Sequential infection experi-ments for quantifying innate and adaptive immunity during influenza infection.
PLoS Comput. Biol. ,15(1):e1006568, 17 January 2019. doi:10.1371/journal.pcbi.1006568 .[80] S. Yang, G. W. Lee, C.-M. Chen, C.-C. Wu, and K.-P. Yu. The size and concentration of dropletsgenerated by coughing in human subjects.
J. Aerosol. Med. , 20(4):484–494, 2007. doi:10.1089/jam.2007.0610 .[81] K. L. Yap and G. L. Ada. Cytotoxic T cells in the lungs of mice infected with an influenza A virus.
Scand. J. Immunol. , 7(1):73–80, 1978. doi:10.1111/j.1365-3083.1978.tb00428.xdoi:10.1111/j.1365-3083.1978.tb00428.x