Agroclimatic Impact on Gastrointestinal Infection: A Mathematical Model

Identifying the correct dosage and time are key factors to successful implementation of anthelmintic. Comparing differential evolution of infection between anthelmintic treated animals against untreated ones, we present a mathematical model that first calibrates data collected and analyzed over an extended period of 10 years (2011-2019), and then predicts the dynamical evolution of gastrointestinal parasites in livestock, focusing specifically on Stronglye & Coccidia oocysts, the two prime negative contributors to cattle health, measured using the standard Eggs-Per-Gram (EPG) index. The model incorporates information about all three critical regimes of infection low infection regime (< 50 EPG), medium infection regime (50− 100 EPG) and high infection regime (> 100 EPG), including fatally large doses of infection (> 500), and probabilistically estimates the variation in animal weight due to infection propagation. A key success of our model is its ability to accurately predict the appropriate anthelmintic treatment times for cattle from a numerical solution of the model presented. The generic model can be applied to other agroclimatic conditions and can serve as a major diagnostic tool for anthelmintic strategizing.


Introduction
Infection by gastrointestinal nematodes (GIN) is ubiquitous in nature among ruminants, causing major economic losses due to constrained growth of ruminants leading to reduction in weight, hence meat volumes, and reduced milk production [1][2][3]. Grazing ruminants are infected by different species of GIN with varying pathogenicity attributes and geographical distribution [4]. The control of GIN in ruminants is largely based on the use of anthelmintics, alongside appropriate grazing management strategizing. However, indiscriminate use of anthelmintic treatments in ruminants resulted in anthelmintic (antimicrobial) resistance which is now a burning topic of discussion [5][6][7]. The increased prevalence of anthelmintic resistance has led to the search for alternative control strategies. Advanced statistical methodology is required to measure such resistance, especially in conditions where the worm burden is low. Faecal egg counts (FEC) based on the number of eggs per gram (EPG) of faeces, have been considered as gold standard test to assess the level of infection, that is directly related to both the worm burden and the fecundity of female adults in the host [8][9][10][11]. FEC have been used to measure genetic resistance to GIN, although in natural infections they can be quite variable both within and between populations [1].
Since most GIN parasites have free-living stages, their development and survival (transmission dynamics) are inextricably linked to environment [12]. As a result, environmental perturbations caused by climate and associated anthropogenic and environmental change is believed to have profound impact on parasite phenology, host-parasite dynamics and host population dynamics, considerably changing the quantitative approach based on mathematical modeling [13][14][15][16][17][18]. Empirical models, such as species distribution models, can be useful in identifying potential drivers of change [16,18], particularly where detailed data and knowledge of the system are unavailable. Moreover, the response of host and parasite to change in ambience, e.g. environment, is often non-linear or infection threshold dependent [13]. For example, there are species specific optimal temperature and moisture requirements for the development and survival of the free-living stages of common GINs of ruminants. Above and below these optima, the parasite growth rate slackens [19].
Alongside emerging infection kinetic preemptive models [20], mathematical modeling [21][22][23][24] is increasingly preferred as a non-invasive alternative to experimentation. It has the potential to compare existing strategies and predict future dynamics. In order to be able to assess the effectiveness of such control strategies, a stochastic (i.e. probabilistic, population-based) model allowing for individual response fluctuations is required. Stochastic model [22] can help to evaluate strategic uncertainties by simulating identical scenarios enabling a direct comparison of treatment effectiveness against anthelmintic regime, and to identify potential interactions, thereby aiding in the feasibility analysis of novel control strategies. Epidemiological modeling of infection propagation have started using nonlinear incidence [25] as also periodicity [26,32].
This article presents a mathematical model that analyses and quantifies the impact of agroclimatic variations on animal health, and hence can be integrated to mainline parasite control strategies. The model itself relates ideologically to an earlier model developed by Berk et al. [13] but is independent in spirit otherwise. The allusion to Berk's model is to account for the interactions between gastrointestinal parasites on animals to predict parasite infra-populations i.e. populations within individual hosts. By introducing variations in growth and resistant traits amongst animals, along with an epidemiological transmission layer, the model can quantitatively outline variations in stocking rate and larval pasture contamination. The model composition is entirely guided by data modeling unlike more phenomenological modeling structures [27]. The model is later validated against the prevailing management control strategies, such as reduced stocking rate, dosage and move strategy, together with strategic anthelmintic drenching, using field data collected by the team from anthelmintic treated farms.

Materials and Methods
The study uses complementary data versus mathematical modeling to analyze the diurnal proliferation of strongyle and coccidia parasitic infection and their impact on animal health. Only faecal samples of animals were used in full compliance with the Institute of the Animal Ethics Committee(IAEC) of the Indian Council of Agricultural Research (ICAR) norms and regulations, as outlined in project proposal ERS/IAEC/2019-20/INSTPROJ-001. No anesthesia or euthanasia was involved at any stage.
Unlike previous studies that relied more on empirical qualitative modeling [3,4,10,12,13], our analytical model is intrinsically data driven with data collected over extensive time periods, and hence is capable of depicting the real spatiotemporal characteristics of the EPG epidemic outbreak. For the purpose of developing the continuum model, field data are statistically analyzed, specifically focusing on the size growth rate of goats as a function of anthelmintic treatment. Towards this, our key descriptor is the growth rate of infection due to strongyle and coccidia, information on which we obtained through detailed data collection, follow-up sampling and then statistical analyzes. First, samples need to be collected from chosen sites that would allow for maximum variation of agroclimatic conditions so as to study the reaction of the animals to anthelmintic treatment over a wide spectrum of climatic possibilities. Second, statistical analysis, by definition, requires a large body of data sets, that for practical reasons, are often limited to a finite number of observations. The field data collection was constrained by numerous factors but we were still able to collect a statistically viable dataset to minimize sampling errors. A question that has hitherto remained unanswered is the contribution of Coccidia in Strongyle propagation and vice versa. Infection kinetics of multiple noninteracting pathogens has been a much studied issue of late [33] that could also be of interest for helminthic propagation.
The primary target of our model is to first analyze and then quantify the impact of rainfall fluctuations with diurnal variation, specifically targeting quantification of the proliferation of EPG (Strongyle & Coccidia) infection levels, to analyze implications on animal health, measured through variations in animal weight. To ensure a wide geographical coverage over varying agroclimatic conditions, we collected data separately from model farms in the Haripal-Hooghly region of West Bengal in India. The zone is traditionally renowned for high humidity and heavy rainfall with large levels of EPG infection during the pre-monsoon and monsoon period. Using standard tools from statistical physics and non-linear dynamics, the model analyzes specific features of the EPG infection dynamics and its climate dependence, ranging from utterly low infectious regimes (50-100) to high infection levels (> 250). The model can also be used to make probabilistic predictions of fatally large doses of infection (> 500), including predictions for peak infection times and agroclimatic parameters leading to such values.
A model based approach to EPG control strategizing has two major advantages over traditional field data based analyzes. First, such a probabilistic modeling based approach can go beyond collected data which are limited to the times and zones of study and predict what might happen in alien conditions where such farming has not yet been attempted but which yet might be cost productive. Second, modeling allows for a (practically) costless prediction of life expectancy of fodder in these regions and allows for optimization of the sustenance cost of these animals without having to further invest in a planning layout. The mathematical structure relies on Bayesian probabilistic mathematics that can potentially go beyond non-Markovian data modeling regimes.

Field Data Collection Methods
Goat farming in rural Bengal, particularly in the Indo-Gangetic planes, is the worst victim of gastrointestinal (GI) parasitism. GI parasitic infestation is responsible for severe economic loss, especially to economically challenged farmers due to loss of production, poor growth, loss of body condition and expensive treatment regime and prophylactic measures. Again emergence of anthelmintic resistance among the parasite population against common anthelmintic strategy has limited our therapeutic armaments. In rural Bengal, goats are mostly reared under free ranging system where sharing of common grazing field and common feeding areas or water bodies provide a convenient opportunity for cross-transmission of the parasite from other goats (homologous) and even from other species like sheep (heterologous transmission). The infective larvae are resistant to environmental factors and can survive for a long period [28,29]. Moreover, there is constant decrease in the availability of grazing pasture for the animals with rapid urbanization of rural bengal. Thus the grazing management scheme has been left with little opportunity to tackle the problem of GI parasitism in goats. Therefore, we require a simple and inexpensive tool to combat the ever-increasing problem of goat health deterioration and ill-productivity in goats especially due to parasitic infestation.
GI parasite load in animals can be conveniently predicted by the quantification of faecal egg amount using Mc Master method [30]. Monitoring changes in faecal egg volume can be used as a quantifier of the level of parasite infestation in goats. Faecal discharge of eggs [31] and its maturation to infected larvae affect several atmospheric factors and related climatic conditions. Our previous study [32] indicated that GI parasite infestation in animals throughout their productive year may be enumerated from seasonal variation in faecal egg burden. As mentiond in [32], faecal egg burden can vary depending on the relationship between microclimate environment and infective larva maturation. Likewise, the present study (under project proposal ERS/IAEC/2019-20/INSTPROJ-001) was undertaken to predict or forecast the disease occurence based on the meteorological data in agroclimatic zones of indo-gangetic belt of West Bengal. Some of the meteorological parameters like rainfall, humidity, temperature etc. were taken into account to compare the infection load with the corresponding meteorological parameters.

Choice of Sites
Using goats as samples, a wide range of data have been collected between 2011-2015 around the Haripal-Hooghly farming belt in West Bengal, India. Three villages were chosen for this study within the latitude-longitude specifications named Jejur (Latitude 22 • 88 , Longitude 88 • 12 ), Kasipur (Latitude 22 • 88 , Longitude 88 • 10 ) and Bhagaldighi (Latitude 22 • 81 , Longitude 88 • 10 ). which are carefully chosen so as to maximize the width of agroclimatic conditions, thereby ensuring that the study encompassed the full range of ambient variations. While none of the villages could singularly represent the whole array of climatic variations (although they are reasonably close by), the set of all 3 villages together, is ideally placed to represent an ensemble averaged version of all possible climate factors. Of the 3 villages, Jejur and Bhagaldighi were chosen for strategic anthelmintic treatment along with mineral mixture supplementation while Kasipur was left for "controlled" study, that is, the goats were left to nourish on their own without any intervention. Such an approach is tailored to ensure maximum statistical sampling of data while minimizing errors arising due to statistical fluctuations. Other affecting factors include easy accessibility and approachability factors, along with zones showing extremities in the infection levels (both on the high as well as on the low sides) subject to agroclimatic variations.

Criteria for sample collection
All data were collected from 3 model farms where the animals were treated with the correct anthelmintic dosage and their growth recorded against animals from non-protected farms without any specific treatment regimen. From each of these 3 villages, 20 farmers were randomly selected such that each of them had access to at least 4 goats. In other words, 80 samples were collected from each village, every month, over a span of approximately 4 years. The animals were evaluated each month to check their weight gain (or loss) for each condition, treated and controlled. The animals were consistently treated between months May-August each year between 2011-2015 and their weight statistics analyzed. This will explain the presence of 9 points corresponding to the respective months as the time identifier (x-axes of Figures 1 & 2).

Experimental design
This is the method for determining Eggs-per-Gram of faeces (EPG) in the faecal sample [7? ]. The test was performed by mixing 0.25 gm of faecal sample in saturated sugar solution (454 gm of table sugar in 355 ml of very hot water) in a test tube. The test tube was filled up to the top with the same solution and then placed a cover slip onto the meniscus and wait for 5 minutes. The cover slip was then placed in a slide and examine under microscope for counting the total number of eggs present in the entire cover slip area. The number of eggs counted, multiplied by 4 was the total number of eggs per gram of faeces (EPG).

Discussions
Our analysis is a complementary study of data and mathematical modeling. The uniqueness of this approach lies in the fact that empirical conclusions drawn from data modeling naturally leads to the mathematical model, that is then suitably reparameterized against real data to ensure acceptable tolerance level (uncertainty) of the conclusions drawn therefrom.

Data modeling
In this section, we tabulate the major conclusions derived from the collected data and use these qualitative results to structure a minimalist mathematical model that is capable of predicting EPG level with time variation (seasonal fluctuation) and as a function of agroclimatic parameters like rainfall and temperature.
In the following, we study the qualitative variations of Strongyle and Coccidia infection levels (expressed in percentage changes) as a function of time (seasonal pattern), subject to changes in rainfall and temperature. The dependence on soil salinity, pH, etc. will also affect the outcome. For the purpose of the present article, our focus will solely be on the rainfall dependent infection variation though.
Strongyle parasitic infection, in percentage, shows oscillatory decay with time, where time is represented by months of a year, starting with highest infection time in July, identied as "0", followed by August ("1"), then September ("2") going up to March ("8"). The maximum intensity level has been recorded during the monsoon period (July-August) with the minimum during winter. The transitions in between show a continuous (and decaying) pattern (Figure 1a).
The Coccidia parasitic infection is decidedly dissimilar in nature to its Strongyle counterpart. Interestingly, although the infection level for this parasite is an oscillating function of time (alike the Strongyle dominated phase), yet the mean profile is found to be largely unperturbed with time variations. This aside, the Coccidia infection shows two peaks, as shown in Figure 1b, indicative of augmented infection levels during monsoon and winter. It remains an open questions though why such double peaked behaviour is not observed in the Strongyle phase, an answer that we hope to generate from our mathematical model.  The time variation of the weight gain (or loss) when treated with anthelmintic (or not) shows a slow oscillatory pattern with a longer time period (ca 4 months) compared to relatively fast diurnal oscillations due to Strongyle (stronger) and Coccidia (weaker), both with circa 2 months' periods. For all treated farms in the three villages, in October, the animals were seen to gain weight when the Strongyle infection was high while the Coccidia infection was low. Another example is in December when the reverse effect was noted: gain in weight with high levels of Coccidia infection but low levels of Strongyle infection. In the phases when neither of these parasites show a high-low combination, the animal weight keeps oscillating with time.
This clearly suggests a pattern, the animals are found to gain weight when the high infection level of one parasite is mitigated by the low infection level of the other.
As to the order of magnitude, given the fact that for goats, Strongyle is a potentially stronger infection than Coccidia, the weight gain maximizes when Strongyle is the lowest with Coccidia on a relative high (e.g. December) and vice versa (e.g. October). Another key feature of note is the percentage change in animal weights compared to the relative infection levels. Weight gain is naturally a slower process than a fast rising infection level. Again, this shows how vital timely anthelmetic treatment is in deciding the quality of the animal's life cycle.
The feature explained above is independent of the region chosen (which is why we can take the liberty of using this as a general feature in our mathematical model) as can be seen from Fig 2. All three villages, Bhagaldighi, Kashipur and Jejur, albeit representing widely varying agroclimatic conditions, show similar qualitative temporal variations of infection levels and associated weight gain (or loss). Figure 2b compares the relative fluctuations between low infection level, defined as EPG between 50-100, against high infection levels (EPG greater than 100), and then against the overall weight gain  (or loss). The infection patterns clearly show robust oscillatory fluctuations and spans the entire range of infection levels (<50, 50-100, etc), as shown in Figure 2 for Bhagaldighi, Figure 3a for Kashipur and Figure 3b for Jejur. For all three model farms, the weight gain/loss is relatively optimal compared to the dynamics of infection.
The decaying oscillatory pattern of the percentage infection level S(t) versus time t for Strongyle in Figure 1a can be mathematically extrapolated to the temporal dynamics S(t) = 6 + 18 exp(−0.214t) cos(0.1t), whereas that for Coccidia, the same dynamics reads C(t) = 18 + 3 sin(2t).

Mathematical Model
Our target here is to structure a model that will combine information from the anthelmintically treated model farms, e.g. Bhagaldighi, Jejur, Kashipur, where the impact of the GI parasites on animal weight for fully treated animals (shown in Figures 4a & 4b), with those from partially treated animal data (show in Figures 1, 2, 3), and completely untreated ones (shown in Figures 5a & 5b), will be strategically used to design a predictive mathematical model that can quantify the impact of treatment and anticipate correct doses. To do this, we use generic empirical inputs from the model farms that depict the time rate of change of the average weight (W) of goat (our choice of cattle) due to the combined effects of Strongyle, Coccidia infection. In doing so, we assume that the time varying weight W(t) of an animal will be the multivariate sum of impacts from Strongyle and Coccidia (assumed to be "mutually exclusive"). The corresponding changes in infection level and their impact on the animal weight can then lead to 3 possible situations: a) no (or very low: <50) infection level, represented by W 0 (t) = W 0S (t) + W 0C (t), where W 0S and W 0C respectively represent the Strongyle and Coccidia dominated contributions in the weight loss; b) the Strongyle dominated phase when Strongyle infection is high and the Coccidia level is low, represented by W 1 (t) and c) the Coccidia dominated phase when the Coccidia infection is high and the Strongyle level is low, represented by W 2 (t). The above segregation is based on our analysis in the previous section that confirms that the average animal weights increase under (anthelmintic treated) conditions (b) and (c) while (a) is an omnipotent condition that is always there (or is always possible, for example through aggravated health care).
Based on the above hypothesis, if W(t; S, C) represents the average weight of the goat at time t with a time dependent Strongyle infection level given by S(t) and Coccidia infection level given by C(t), the rate of change of animal weight due to both parasitic infection can be modelled as the following dynamical system For Bhagaldighi data, the specific values of these parameters are given by A 1 s = 6, A 2 s = 18, A 3 s = 0.214, A 4 s = 0.1 and A 1 c = 18, A 2 c = 3, A 3 c = 2; for other model farms, the qualitative pattern remains unchanged, only the absolute values change. Similarly, W 0 (t), W 1 (t) and W 2 (t) are both qualitatively where a 1 and a 2 are empirical fit parameters to be extracted from real data. For our data sets, we get a 1 = 13 and a 2 = 4. The above extrapolations lead to the following The dynamics relating parasite affected weight (loss), as detailed in equations (4), with anthelmintically treated regime (equation (1a) leads to Since Strongyle dominates Coccidia in strength, in general | A s A c | > 1. As a generic formulation, our analysis leads to the following epidemic model We can solve Equation (5) to estimate the animal weight W(t) at time t as affected by GI infection:

Gestation timeline
The model presented in equation (7) is a generic parasitic infection model that analyzes the impact of GI parasites on animal weight. With knowledge of the nature and virulence of the (Strongyle or Coccidia) infection, W 1 and W 2 can be empirically evaluated (equation (6)), together with the individual infection dynamics of Strongyle and Coccidia. The latter will remain qualitatively unchanged as in equation (5) but will be parametrically modified depending on location and other agroclimatic parameters like rainfall, salinity, etc.
Solution of the time t c from the model in equation (5), i.e. dW dt = 0 gives us probabilistically accurate time periods for administering appropriate anthelmintic dosage. For Bhagaldighi, using parameters (A 1 , we can numerically estimate the best Strongyle anthelmintic treatment times, reflecting infection gestation times, as t c = 23, 86, 142, all in units of days, indicating that anthelmentric treatment needs to start latest before 23 days from the start of July and ideally 142 days before July. The latter number (142 days) corresponds to a time of approximately 4.5 months prior to July, suggesting March as the ideal time for starting the treatment unlike April, that is generally used. On the other hand, for Coccidia treatment, the model predicts critical time period values as follows t c = 5, 43, 128, indicating that Coccidia treatment can start relatively later but the ideal start is again in March, about 128 days prior to July. Using data for other regions, we can similarly calculate the correct treatment times for those regions.

Conclusions
We present a mathematical model that is capable of explaining the effects of agroclimatic fluctuations on EPG (Strongyle & Coccidia) parasitic infection levels. The model is motivated by real data collected from the Haripal-Hooghly belt of West Bengal India, zones that are traditionally renowned for high humidity and heavy rainfall with large levels of EPG infection during the summer and monsoon period. The model explains a broad range of features of the gastrointestinal (GI) infection dynamics and its climate dependence (measured through parametric variations), ranging from utterly low infectious regimes (50-100) to high infection levels (> 100). The major incentive of such a model is (practically) costless prediction of life expectancy of fodder in these regions and how to optimize the sustenance cost of these animals without having to invest in a planning strategies, that for developing economies, is often a practical challenge.
A major strength of this mathematical approach is the ability to go beyond collected data which are starkly limited by the times and zones of study whereas the mathematical model can probabilistically predict what might happen in alien conditions where such farming has not yet been attempted but which yet might be cost productive. Equation (6) represents such a preemptive model that will apply to other agroclimatic conditions where the only change will be in parameter values.
A highlight of our model is its ability to accurately predict the best treatment times required for successful anthelmintic treatment of GI infected animals. Our model first identifies the monsoon period, starting with July, as the key high infection regime (identified as peaks in Figures 1 & 2) and then predicts the correct anthelmintic intervention timings. We show that treatment for Strongyle infection should start earlier than Coccidia infection but both need to start ideally by March (unlike April, as is commonly practised), and no later than 3 weeks prior to July for Strongyle while for Coccidia, this can wait until the last week of June. These dates can change for other regions but that can be easily calculated from this model. Such precise quantification allows simple time lining of best treatment practice.
The model can be further extended to make probabilistic predictions of fatally large doses of infection (> 500), including predictions for peak infection times and aggro-climatic parameters leading to such values. Work is presently underway to ascertain such critical thresholds as also to integrate other agroclimatic factors, like rainfall, humidity, and salinity in what will be a multivariate approach through machine learning utilities.