Modeling effect of cover condition and soil type on rotavirus transport in surface flow

Nearly all animal housing facilities such as feedlots, barns, and pens are believed to be primary sources of microbial contamination of surface water resources. Microbial pathogens present in agricultural runoff impose a significant risk to health when acquired directly via the fecal-oral route or indirectly as a waterborne contaminant. As per Waterborne Disease and Outbreak Surveillance System (WBDOSS) of Center for Disease Control and Prevention (CDC), there have been 794 waterborne disease outbreaks causing more than 575,000 illnesses from 1971 to 2006 in the United States (Bhattarai 2011). Furthermore, rotavirus is a leading cause of death of children around the globe. Rotavirus infection kills more than 500,000 children under the age of five every year, most of them in developing countries (Feng et al. 2002). Even in a developed country like the United States, rotavirus was estimated to be the cause of about 60,000 hospitalizations and 37 deaths annually during 1993–2003 (Fischer et al. 2007). Apart from humans, rotavirus has a wide host range including many mammalian species like pigs (Ferrari & Gualandi 1986), lambs (Wray et al. 1981), foals (Conner & Darlington 1980), rabbits (Conner et al. 1988) and deer (Smith & Tzipori 1979).

The risk of pathogenic contamination is not only a human health hazard but also a potential threat to the continued existence of animal operations. To design and develop pathogen control practices, microbial transport processes in surface runoff need to be properly understood and quantified along with various factors that affect pathogen transmission to the environment. There have been several studies which looked into the transport of pathogens such as Cryptosporidium (Harter et al. 2000; Atwill et al. 2002; Tate et al. 2004; Trask et al. 2004; Koch et al. 2013; McLaughlin et al. 2013; Davidson et al. 2014) and bacteria (Smith et al. 1985; McMurry et al. 1998; Gagliardi & Karns 2000; Jamieson et al. 2002; Passmore et al. 2010). There are very limited studies which have investigated the fate and transport of viruses in soil and water systems (Powelson et al. 1990; Sim & Chrysikopoulos 2000; Chu et al. 2003; Ojha et al. 2011; Davidson et al. 2013, 2016; Gotkowitz et al. 2016).

It is very time and resource consuming to monitor the influence of multiple environmental factors on pathogen fate and transport. Therefore, modeling can provide a consistent quantitative approach to estimate pathogen loading under a wide range of environmental conditions. In terms of modeling pathogen transport phenomena, there are few studies in the literature that model the transport and fate of pathogens in water flow (Moore et al. 1989; Walker et al. 1990; Canale et al. 1993; Wilkinson et al. 1995; Fraser et al. 1998; Tian et al. 2002; Dorner et al. 2006; Parajuli et al. 2009; Bhattarai et al. 2011; Niazi et al. 2015; Liao et al. 2016). Pachepsky et al. (2006) reviewed the status and challenges in manure-borne pathogen transport modeling in the environment. Similarly, Jamieson et al. (2004) reviewed different approaches in microbial transport modeling in surface waters at the watershed scale. A few studies have estimated surface transport of bacteria using simple empirical equations which relate bacteria transport to runoff rate or depth (Moore et al. 1989; Fraser et al. 1998; Tian et al. 2002). However, the empirical transport equations used in these models have not been well tested and validated. The parameters do not have a physical basis and must be obtained through calibration. Surface transport of bacteria in SWAT (Sadeghi & Arnold 2002) and COLI (Walker et al. 1990) are computed based on the modified universal soil loss equation which is also an empirical model. As pathogen transport models move from a small laboratory setting to a large scale such as a watershed (Walker & Stedinger 1999; Dorner et al. 2006; Parajuli et al. 2009), the model prediction uncertainties are enhanced. Very few modeling studies have linked sediment transport to pathogen movement in water (Bhattarai et al. 2011; Pandey et al. 2016). The objective of this study was to model the fate and transport of infective rotavirus particles in surface flow using a process-based model. The results of this study can help in developing the guidelines for translating the model into practice for designing best management practices (BMPs) such as a vegetative filter strip (VFS). Previous studies have shown the promising potential of VFS in reducing sediment, nutrients, pesticides and pathogens from runoff (Atwill et al. 2002; Tate et al. 2004; Trask et al. 2004; Koch et al. 2013; Mclaughlin et al. 2013; Davidson et al. 2014).

The rotavirus overland transport experiments were completed using three soils (Catlin silt-loam, Alvin fine sandy-loam, and Darwin silty-clay) and three cover conditions (bare soil, Smooth Brome, and Fescue vegetation) at the University of Illinois. Three soil types used in this study span the range of sand and clay contents for soils typically found in Illinois. The methodology is described in detail in Davidson et al. (2016) with a condensed summary included below.

The rainfall simulator used in this study measured 1.07 × 0.66 m in size and was constructed using 12.7 mm PVC pipe. Since adequate droplets were not produced by the mister nozzles, a common fiberglass window screen with 1 mm square openings was placed 305 mm below the arrangement of nozzles. The nozzles sprayed water onto the screen and droplets were allowed to fall onto the surface below, more closely mimicking natural rainfall droplets. This small scale rainfall simulator can produce rainfall intensities from 23 to 90 mm/hr. In this study, 65 mm/hr rainfall intensity, which represents the 25-year event for east-central Illinois, was used for all soil and cover conditions.

Soil bed boxes (0.610 m long by 0.305 m wide by 0.015 m deep) used during experimental studies were made from 6.35 mm thick Plexiglas. A total of 11 holes with 4.76 mm diameter were drilled 76.2 mm above the bottom of the soil box (corresponding to soil surface location) to allow for the collection of surface runoff. The Plexiglas soil boxes were placed in a wooden frame which allowed the slope adjustment of the box frame. Two aluminum trays were placed on the down-slope end of the wood frame to collect surface flow samples.

Three different soils were collected and placed in the small-scale soil chambers. These soils represent the typical range of soil textures in Illinois. A moderately drained silt-loam (Catlin series: 24% sand, 50% silt, 26% clay) soil was obtained from an agricultural field in Champaign, IL. Similarly, a well-drained fine sandy loam (Alvin series: 60% sand, 25% silt, 15% clay) soil and poorly drained silty-clay (Darwin series: 5% sand, 50% silt, 45% clay) soil were collected from a farm near Newton, IL. Since the soil chambers were designed for 76 mm of soil and soil texture would remain somewhat constant in such a small depth range, the top 76 mm of soil was collected after removing all vegetation and organic material from the soil surface. The soil was oven-dried at 105 °C for 24 hours to reduce existing microbial activity that could interfere with sample analyses and then packed in three one-inch (25.4 mm) layers. Each layer was moistened and compacted with a plate. Rainfall was applied on a regular basis after placing three layers of soils in the chamber to simulate natural compaction. Two vegetation covers used in this study were Smooth Brome and Tall Fescue (purchased at Illini FS, Urbana, IL). The seeding rate for each grass was 16.8 kg/ha. It took approximately 4–6 months to achieve the desired maturity level of vegetation before use in transport studies. The surface runoff samples were collected and analyzed for rotavirus. Surface runoff samples were collected every minute after runoff began.

Group A Porcine rotavirus OSU strain (P9[7], G5) was obtained from the American Type Culture Collection catalog # VR-892 and passaged two additional times in cultured Ma104 cells Three colostrum-deprived newborn piglets were inoculated with 1 mL (3 × 10⁶ FFU) each, and the feces were pooled and collected for 2 days. Rotavirus was prepared by cultivation in Ma104 cells and purified as described by Rolsma et al. (1998), excluding the last CsCl gradient step.

Since the collected surface runoff samples contained a relatively high concentration of rotavirus particles, a dilution of the samples was necessary for analysis. A 1/5 dilution was prepared by mixing 200 μL of the sample with 800 μL of serum-free MEM and applying 500 μL directly to the plate wells. Similarly, soil core samples were prepared by combining 250 μL of the sample with 750 μL of serum-free MEM and applying 500 μL to the plate wells.

The focus forming unit (FFU) assay was performed following previously published methods (Rolsma et al. 1998). The rotavirus FFUs were enumerated by counting the number of the FFUs present in a given sample at a ×100 magnification using a Nikon TS 100 inverted microscope equipped with a computer-controlled electronic stage and Spot RT-slider CCD camera. The experimental results were used to parameterize and validate the model results.

The following possible states of virus fate and transport in overland flow were considered based on the stochastic framework proposed by Yeghiazarian et al. (2006) and modified from the process-based framework proposed by Bhattarai et al. (2011):

• 1.

  Attachment: viruses attach to soil particles with an attachment rate k₁₂

• 2.

  Movement: viruses attached to soil particles (soil + virus aggregate) are moved by surface flow velocity (u) with rate constant k₂₃

• 3.

  Detachment: viruses attached to soil particles become detached with the detachment rate k₂₁

• 4.

  Inactivation: viruses become inactivated with decaying rate k_d

• 5.

  Infiltration: viruses are removed by infiltration with the rate equal to f/D (f is saturated hydraulic conductivity, and D is flow depth).

The system of mass-balance equations was solved using the following initial and boundary conditions: (a) initial conditions: ⁠, where is the concentration of rotavirus applied to the system; and (b) boundary condition: These conditions reflected the assumptions that at time zero the rotavirus particles were unattached and freely moving in the water, and no continuous source of rotavirus exists at x = 0.

Since viruses were immobile in the C₂ stage, the flow velocity component in Equation (2) was considered zero. However, virus concentrations in the C₁ and C₃ stages were in mobile flow and moved downstream with surface flow. In this study, the surface velocity was obtained from the Water Erosion Prediction Project (WEPP) simulation results and assumed to be constant over the space domain for each iteration. The partial differential equations (PDEs) for C₁ and C₃ have the form of a wave equation (hyperbolic PDE). Since the flow velocity (v) was assumed constant, the grid for the finite difference method could be discretized using dx=v dt for each iteration. The system of Equations (1)–(3) was solved with the method of characteristics (Park et al. 2008; Bhattarai et al. 2011). In the presence of vegetation, an additional equation (similar to Equation (2)) was added to the model which could represent attachment to the vegetation with the rate k₁₄ and detachment from the vegetation with the rate k₄₁.

The modeling framework included the coupling of a soil erosion model (WEPP), with a pathogen transport model to capture the dynamics of sediment-bound pathogens. The WEPP is a physically based erosion prediction model applicable to erosion processes on single hillslopes and small watersheds (Nearing et al. 1989; Cochrane & Flanagan 1999). Outputs from WEPP simulations were used as input for the pathogen transport model (Bhattarai et al. 2011). The WEPP model is capable of predicting soil loss and sediment deposition from the overland flow on hillslope and the concentrated flow in small channels.

For the rotavirus transport study, a 20-minute single rainfall simulation on a laboratory tilting soil bed was performed. The input for the WEPP model included: climate data including the rainfall amount and duration; hillslope profile including area, shape, and inclination of the hillslope; soil types; and cover conditions, etc. The WEPP simulations provided: (a) average flow depth and average flow rate, which can be calculated from the Chezy coefficient and the peak runoff rate; (b) average soil loss per detachment area; and (c) the soil particle fraction exiting in the flow.

For modeling purposes, microorganisms and soil particles can be considered as two substances that interact with each other producing the third constituent: microorganism-soil particle aggregates. The strength of attachment of microorganisms to soil particles depends on the chemical characteristics of soils which are highly influenced by the particle size. The mathematical representation of inter-particle interactions can be obtained experimentally or from an erosion model (suspended sediment). The particle size distribution in the suspended sediment is predicted by the WEPP model.

Both of the vegetation parameters k₁₄ and k₄₁ were treated as calibration parameters in the absence of reliable parameter estimates in the literature.

The pathogen transport model developed by Bhattarai (2011) was combined with WEPP to predict the fate and transport kinetics of infective rotavirus particles. The model was validated using data provided by Davidson et al. (2016). The results are presented and discussed in the following sections.

It was observed that the volumes of surface runoff obtained from the bare-ground conditions were consistently higher than those from vegetated conditions. This was observed for all three soils and both surface cover types. These findings were expected since reduced runoff volumes from the vegetation cover were associated with more resistance to overland flow and longer contact times allowing a higher infiltration of water into the soil profile. Similarly, the rotavirus recovery also followed the same trend. More infective rotavirus particles were recovered from the bare surface condition compared to the vegetated surface conditions. In terms of bare soil beds, the highest amount of infective rotavirus particles was recovered from the bare Catlin soil and the lowest recovery was observed from the bare Alvin soil. Rotavirus recovery from the bare and vegetated Alvin soil beds was smaller compared to those from Catlin or Darwin soil beds (Davidson et al. 2016).

For the WEPP model, many input parameters were obtained from the experimental setup such as soil plot length of 0.61 m, rainfall intensity of 65 mm/hr applied for 20 minutes and three different soil types (Catlin, Darwin and Alvin). Hence, initial saturation level (antecedent moisture content) and soil hydraulic conductivity were used as the calibration parameters. The parameters were calibrated in such a way that Nash–Sutcliffe efficiency index or EI (Nash & Sutcliffe 1970) and coefficient of determination (R²) between observed and simulated runoff rate were maximized. Since the moisture holding capacities and hydraulic concavities of bare and vegetated surfaces are different, the parameters were calibrated separately for those two conditions.

The comparison of the experimental data with the WEPP simulation results for surface runoff from the bare ground and two vegetated conditions for Darwin also showed a good agreement (R² values of 0.78, 0.83 and 0.72 for values for bare, Brome and Fescue cover conditions, respectively). The observed surface runoff rate for the Darwin soil bed with Fescue vegetation was much lower than the rate for the Brome grass cover condition. The difference in factors such as root zone depth, surface cover and leaf area index between the Brome and Fescue grasses may be attributed to the decrease in surface runoff rate.

The statistical measure of goodness of fit between measured and WEPP simulated runoff rate for the three different cover and soil conditions is presented in Table 1. It was observed that both EI and R² values were high for the bare condition for all three soils. EI values for three vegetated cases were negative indicating the mean observed was a better estimate of the simulated results.

After calibrating the WEPP model, the simulated outputs were used to drive the pathogen transport model. The WEPP model simulations provided the input parameters: bed slope, rainfall intensity, flow velocity, flow depth, saturated hydraulic conductivity and soil loss for the rotavirus transport model. During the calibration process, the break-through curve from the rotavirus transport model was compared with the laboratory experiments to estimate different attachment factors k₁₂, k₂₃, k₂₁ for bare soil and additionally k₁₄ and k₁₄ for vegetated conditions.

For Alvin soil, the experimental data showed more than one peak in the rotavirus break-through curve for all three cover (bare and two vegetated) conditions. The rotavirus transport model could reproduce the first peak seen during the experiments; it could not replicate remaining peaks. Hence, the comparison of observation and model simulation resulted in low R² values (0.61 for bare and 0.26 for Brome cover conditions).

Although the relation provided by Ling et al. (2003) gave an initial parameter to drive the model for bare soil conditions, the one parameter approach for all soil and cover conditions did not produce satisfactory results in many cases. In the absence of an alternative reliable approach that can incorporate various soil and cover types in a partitioning coefficient, the parameters k₁₂ had to be calibrated for each case by comparing model results with the experimental data. The different values of attachment and detachment rate parameters (k₁₂, k₂₃, k₂₁ for bare soil and additionally k₁₄ and k₁₄ for vegetated condition) had to be used during simulation process to replicate the result obtained from the experiments.

Table 2 shows calibrated rotavirus partition coefficients for bare and vegetated surface conditions for Alvin, Catlin and Darwin soils. It was observed that the parameters k₁₂ and k₂₁ were more sensitive compared to the parameter k₂₃. The parameter k₁₂ affected the peak of the rotavirus recovery curve while the parameter k₂₁ dictated the shape of the curve. Parameter k₂₃ governed the formation of pathogen and soil particle aggregate which was not quantified in this study. In the case of vegetated surface, parameters k₁₂ and k₁₄ were two sensitive parameters, which governed how many pathogens got attached to soil surface and vegetation, respectively. Parameters k₂₁ and k₄₁ were another set of sensitive parameters which governed how many pathogens got released from soil surface and vegetation, respectively.

These simulated results indicated that one partitioning parameter based on soil clay content as used by Ling et al. (2003) might not be valid for rotavirus and similar microorganisms for varying environmental and topographic conditions. Future modeling attempts must include the dynamic nature of environmental and topographic condition for formulating the partitioning parameters.

A physically based approach for modeling rotavirus fate and transport in surface flow has been presented. Deterministic descriptions of the processes have been employed and a correspondence between them has been established and model parameters have been identified from the literature. Another feature of the model is the incorporation of an erosion model (WEPP) into the microorganism transport model, thereby establishing a link between the sediment and microbial transport, and highlighting the capability of the suspended sediment to transport microbes in the overland flow. Model results were compared with small-scale experimental data and the model performed well in most cases. Although the model replicated the kinetics of infective rotavirus recovery from surface runoff, the model under-predicted recovery when experimental data showed multiple peaks in the rotavirus break-through curve. The under-prediction in simulated results can be attributed to simplification of complex processes and model parameter uncertainty. The model results can be improved by training the model with more data which will help reduce model and parameter uncertainty.