Cross Sectional Study and Risk Factors Analysis of Francisella tularensis in Soil Samples in Punjab Province of Pakistan

Tularemia is an endemic zoonotic disease in many parts of the world including Asia. A cross-sectional study was conducted to determine genome-based prevalence of Francisella tularensis (Ft) in soil, assess an association between its occurrence in soil and likely predictors i.e., macro and micro-nutrients and several categorical variables, and determine seroconversion in small and large ruminants. The study included a total of 2,280 soil samples representing 456 villages in eight districts of the Punjab Province of Pakistan followed by an analysis of serum antibodies in 707 ruminants. The genome of Ft was detected in 3.25% (n = 74, 95% CI: 2.60–4.06) of soil samples. Soluble salts (OR: 1.276, 95% CI: 1.043–1.562, p = 0.015), Ni (OR: 2.910, 95%CI: 0.795–10.644, p = 0.106), Mn (OR:0.733, 95% CI:0.565–0.951, p = 0.019), Zn (OR: 4.922, 95% CI:0.929–26.064, p = 0.061) and nutrients clustered together as PC-1 (OR: 4.76, 95% CI: 2.37–9.54, p = 0.000) and PC-3 (OR: 0.357, 95% CI: 0.640, p = 0.001) were found to have a positive association for the presence of Ft in soil. The odds of occurrence of Ft DNA in soil were higher at locations close to a water source, including canals, streams or drains, [χ2 = 6.7, OR = 1.19, 95% CI:1.05–3.09, p = 0.004] as well as places where animals were present [χ2 = 4.09, OR = 2.06, 95% CI: 1.05–4.05, p = 0.02]. The seroconversion was detected in 6.22% (n = 44, 95% CI: 4.67–8.25) of domestic animals. An occurrence of Ft over a wide geographical region indicates its expansion to enzootic range, and demonstrates the need for further investigation among potential disease reservoirs and at-risk populations, such as farmers and veterinarians.


INTRODUCTION
Tularemia is caused by the bacterium Francisella tularensis (Ft), a category A classified select agent by the Center for Disease Control and Prevention (https://www.selectagents.gov/ SelectAgentsandToxinsList.html). Francisella tularensis (Ft) is a pleomorphic Gram-negative intracellular bacterium (Schulert and Allen, 2006) that has zoonotic implications across many parts of the globe (Oyston, 2008;Vogler et al., 2009). Four subspecies, namely, tularensis, mediasiatica, holarctica, and novicida, have been identified (Sjösted, 2005;Champion et al., 2009;Penn, 2015). The presence of Ft subspecies holarctica has been reported in Europe and Asia, whereas Ft subspecies tularensis has been reported in North America (Garaizar et al., 2006). Among the Asian countries, most of the outbreaks and cases have been reported in Turkey (607 cases in 2012), China (31 cases in 1986), and Iran (36 casesin 2013) (Esmaeili et al., 2014;Gürcan, 2014). The organism has a broad and complex host distribution that includes vertebrates, invertebrates, and environmental matrices such as soil, aerosols, and water (Kuske et al., 2006;Silvestri et al., 2016). Humans can acquire infection through inhalation, an arthropod bite, ingestion of contaminated food, or water, as well as through contact with infected tissues or fluid from animals (Oyston, 2008;Ulu-Kilic and Doganay, 2014;Silvestri et al., 2016). Individuals living in rural areas or those in close proximity to animals/disease reservoirs, especially farmers and veterinarians, are considered the most at-risk population for tularemia (Lévesque et al., 1995;Ulu-Kilic and Doganay, 2014). Though glandular, oculo-glandular, ulcero-glandular, typhoidal, and pneumonic symptoms are common in affected humans, the clinical signs and severity of disease depend on the entry route and infectivity dose (<10 CFU) (Helvaci et al., 2000;Bossi et al., 2004;Pechous et al., 2009). Just as in humans, clinical signs in animals are varied. Cats are more susceptible than dogs and remain mostly in non-clinical form however in some cases, symptoms may include fever, lymphadenopathy, anorexia, oral ulceration, hepatospleenomegaly, and dehydration (Gliatto et al., 1994;Woods et al., 1998).
Different subspecies of Ft have been reported worldwide with varying geographic distribution and disease potential. Occurrence of tularemia caused by the most virulent type (F. tularensis: biovar A) has been reported only from North America, while cases caused by the less virulent type (F. holarctica: biovar B) have been observed in Europe, North America, and Asia (Garaizar et al., 2006;Oyston, 2008;Esmaeili et al., 2014;Ulu-Kilic and Doganay, 2014). However, there is a paucity of data on the enzootic range of Ft in Pakistan. Since Ft has the potential to survive and persist in the environment for a longer period of time (SjÖstedt, 2007), we undertook a study to determine the prevalence of Ft in soil from eight districts of Punjab province of Pakistan followed by an evaluation of seroconversion in small and large ruminants. Besides several soil characteristics that included macro-and micronutrients, the study also examined likely risk factors that could be associated with its occurrence in soil, and therefore can contribute toward human and animal exposures. It is anticipated that the findings will be valuable to local, as well as global public health agencies for evaluating potential disease burden, identifying reservoirs, and developing strategies to prevent and control tularemia in animal and human populations.

Study and Sampling Design
A cross-sectional study was conducted in Punjab province (31.1704•N and 72.7097•E) from 2011 to 2015. The province has nine administrative divisions, 36 districts and approximate 4,883 villages. It dominates agriculture, and has the largest human and livestock populations in the country. Punjab province contains five rivers ("punj" means five and "ab" means water) which together provide one of the country's largest irrigation systems for agricultural cultivation. Besides rivers and canals, groundwater (tube-well) as well as natural rain (barani) are being used to irrigate some of the areas in the province. Though mechanical (automatic) plowing is widespread across many districts in the province, animal-based plowing (manual) is also employed at some places in the province. We used three-stage sampling design. Since an incidence rate of Ft in Pakistan is not known, we selected districts representing the main livestock production areas of the province where there exists an increased annual incidence of human and animal disease (Directorates of Human and Animal Health, Punjab Province).
Assuming 50% prevalence, 95% CI and 5% margin of error, the required number of villages was 357, however, we included 456 villages representing 10% of each of the study district to increase the validity of the results using WinEpi software (http:// www.winepi.net/uk/sample/indice.htm). From each village, we conveniently selected five sites for soil sampling; four were from livestock barns where human and animals were living in close proximity, while one represented an agricultural land only. The geographical coordinates were noted using Garmin (Dakota, U.S.A). After removing 3-5 inches of top-surface soil at each site, a total of 2,280 sample (∼250-300 gram each) were collected from 456 villages representing districts Sheikhupura (n = 295), Gujranwala (n = 360), Faisalabad (n = 370), Sargodha (n = 370), Sahiwal (n = 255), D.G. Khan (n = 215), Chakwal (n = 190), and Attock (n = 225). A brief history of each study site along with information about different categorical variables or risk factors such as distance from animal market, main road and water source, animal's density in a village, number of households in a village, number of domestic animals in a village, cover ground (vegetation) was recorded ( Table 1).

Genome Extraction and PCR Based Identification
Genomic DNA was extracted from 0.25 gram of each soil sample (PowerSoil R DNA Isolation Kit, MoBio, USA) as per manufacturer's instruction, and was subjected to real time PCR (CFX 96, BioRad, USA) using a highly sensitive and specific assay (Christensen et al., 2006) with minor modifications. The quantification of DNA was performed using the NanoDrop 1000 spectrophotometer (Thermo scientific, USA). The used real time PCR (RT-PCR) primers (Forward:  . Thermal cycling conditions included one cycle of 95 • C for 5 min followed by 45 cycles of each of denaturation at 94 • C for 5 s and annealing at 60 • C for 20 s, and then one cycle of cooling at 40 • C for 1 min. The assay was optimized and validated using the control (tul4 gene PCR products) and the proficiency testing samples that were kindly provided by the Pennsylvania State University, USA. The necessary assay controls such as dsDNA PCR product (positive control) and dDiethyl-pyrocarbonate water (negative control) were used each time. To rule out any potential contamination in processing and/or false-positivity, each sample that exhibited a positive result was gel-electrophoresed (Supplementary Material), and the process described above was repeated thrice beginning from the genome extraction.

Serum Analysis
Blood samples (∼5 mL) were collected conveniently from Ftpositive site representing goat (n = 200), sheep (n = 175), cattle (n = 179), and buffalo (n = 153) representing district Chakwal, Gujranwala, Faisalabad, Attock, Sahiwal, Sargodha, and Dera Ghazi Khan in Punjab Province, Pakistan. The separated sera were stored at −80 • C until further use. Sera (1 µL) were analyzed for anti-Ft-enzyme linked immunosorbent assay (ELISA) antibodies using an SERION ELISA classic Ft kit (Virion/Serion, Germany) according to the manufacturer's instructions. A specific secondary alkaline phosphatase antibody (ThermoFisher, USA) was used for goat (F (ab') 2-rabbit antigoat IgG H+L), sheep (F(ab')2-donkey anti-sheep IgG H+L), cattle (goat anti-bovine IgG (H+L), and buffalo (goat anti-bovine IgG H+L). Optical density of sera were read within 60 min at 405 nm against substrate blank and 655 for reference and field samples as per manufacturer's recommendations. Only valid samples with OD value of substrate blank <0.25 and variation between OD values of standard serums not higher than 20% were considered for further analysis. The positive and negative serum samples were analyzed using Microsoft R Excel-based software tool SERION activity with reference to lower cut off value 0.42 and upper cut off value 1.43 provided by manufacturer. Samples below the lower cut off value were considered as negative while samples above the upper cut off value were measured as positive.

Statistical Analysis
The RT-PCR results along with numerical (soil characteristics) and categorical variables (potential risk factors) were compiled into a Microsoft Excel spreadsheet. In an earlier study (Muhammad et al., 2017), we described association between physio-chemical characteristics of soil and the presence/absence of Ft DNA in soil samples. In this study, we further explored these data and other potential variables to quantify factors associated with detection of bacterial DNA in soil. The data on physiochemical characteristics of soil were not normally distributed (Shapiro-Wilk test, p ≤ 0.05); therefore, Mann-Whitney U-test was applied to assess the effect of those numeric variables on soil positivity. In Figure 2 the variables with p ≥ 0.2 were included in further analyses. The correlation matrix plot revealed collinearity among the selected variables, and correlated variables (r ≥ 0.3) were subjected to principal component analysis (PCA) (Abdi and Williams, 2010;Pourhoseingholi et al., 2012). PCA reduces dimensionality in the data and transforms the variables into a new set of uncorrelated variables called principal components. The Bartlett test of sphericity was significant (p < 0.05) and The Kaiser-Meyer-Olkin (KMO) measure of sampling adequacy was 0.582 (p < 0.05), indicating appropriateness of dataset for PCA. The first three principal components had an eigenvalue >1 and were used as covariates in binary logistic regression. A Chi Square test was used to evaluate association between the occurrence of Ft DNA in soil and categorical variables. A Fisher Exact test was used as alternative to Chi square test where any of its assumption was violated. The regression model included the presence/absence of Ft DNA in soil samples as a dependent variable. The independent variables in the model were moisture, P, Ni, Mn, Na, Zn, PC1, PC2, and PC3. The model also included categorical variables with p ≥ 0.2. A backward likelihood ratio method was used to select covariates in the regression. A p ≤ 0.05 was considered significant in the outcome of regression analysis (Nandi et al., 2016). The data were analyzed in R using "FactoMineR, " "factoextra, " and "corrplot" packages (Table 3).

RESULTS
The Prevalence (Figure 1). Interestingly, a total of 20 samples (4.3%, 95% CI: 2.86-6.68) representing agriculture land with no apparent human interaction were also found to be positive where there was no apparent animal and human interaction.

Relationship Between Soil Chemistry and Categorical Variables to Occurrence of Ft in Soil
In the earlier study, we described association between physiochemical characteristics of soil and presence/absence of Ft DNA in soil samples using a simple t-test and chi-square (Muhammad et al., 2017). In this study, we further explored these data and other potential variables to quantify factors associated with detection of bacterial DNA in soil. Table 1 shows eigenvalues, percentage variance and cumulative percentage of variance of principal components. The first three principals had eigenvalues > 1 and cumulatively explained 77.26% of total variance (inertia) in the dataset. Varying between + 1 and −1, the values (≥ 0.4) of the loadings represent the correlation between each variable and a principal component (  Figure 3 displays quality of representation (squared cosine, cos2) of variables and individuals (sampling sites) along the first two principal components. The angles between variables show degree of correlation, and their lengths represent importance for the respective components. Clay, Cd, and Pb were evidenced to be more important for PC1. Organic matter and N were contributed greatly in PC2. Clay and slit were negatively correlated with each other. The sampling points in the periphery had better representation (cos2 values closer to 1).

Seroconversion in Domestic Animals to Ft in the Studied Districts
The seroconversion was found in 6.22% of small and large ruminants (n = 44, 95% CI: 4.67-8.25) Spatial distribution of seropositivity in animals has been illustrated in Figure 4 where a village has been red highlighted corresponding to seroconversion in any of the study animals. A significant difference (p < 0.05) was found in the prevalence of serum anti-Ft-ELISA antibodies among cattle (

DISCUSSION
The molecular diagnostic assay (RT-PCR) had a high sensitivity and specificity with a detection limit as low as <100 genome copies, and therefore allows simultaneous examination of numerous samples with rapid turnaround time (Christensen    , 2006). Such assays are typically used for surveillance of extremely dangerous pathogens (EDPs) in the natural environment, particularly in settings such as Pakistan that lack a highly contained laboratory, trained personnel, culturing facilities, and repositories for isolating archives. With regards to detection of antibodies, the biological marker (lipopolysaccharide or LPS) employed for the detection of seroconversion does not produce any cross-reaction with any component of currently known bacteria (Schmitt et al., 2005;Jenzora et al., 2008), and hence is a suitable macromolecule for the detection of antibodies while monitoring a large number of samples originating from domestic animals and humans (Sharma et al., 2013). Despite evidence of Ft in soil, there have been no reported tularemia outbreaks or cases in Punjab province to-date. There are a number of potential reasons for why this may be the case.
Firstly, clinical cases may remain undiagnosed or misdiagnosed due to lack of laboratory based diagnostic capabilities throughout the country, as well as due to the fact that the clinical manifestations of tularemia can be fairly general. Moreover, strains isolated from the northern hemisphere are much more virulent than those isolated from Europe and Asia (SjÖstedt, 2007;Oyston, 2008). Lastly, the lack of reported cases could be due to climate and micro-and macro-nutrient characteristics of the soil supporting the persistence of Ft in the environment without any concomitant infection or outbreak (Dennis et al., 2001;Oyston, 2008;Alkhuder et al., 2010). Though it requires further evaluation of potential disease reservoirs (primarily rodents), and isolation and subsequent genomic characterization of prevailing strains, the more frequent occurrence of Ft in soil at sites close to water sources may be indicative of the prevalence of a less virulent strain (Ulu-Kilic and Doganay,  (Eliasson et al., 2006;SjÖstedt, 2007). An area in close proximity to water may have a higher probability for the presence of Ft, however the presence of mammalian host may be required for manifestation of disease (Desvars et al., 2015).
Although an increased persistence of Ft has been observed in environments lacking humidity (Wilkinson, 1966;Dennis et al., 2001), the organism cannot remain viable under desiccation conditions (Faith et al., 2012). A significant association between the occurrence of Ft in soil and the presence of animals at the same site may correspond to the presence of disease reservoirs at locations in Punjab province where Ft was identified, such as the rodents that are found frequently at livestock barns and surrounding areas. It has been observed that an increase in the human population can decrease the number of rodents in the fields and corresponding areas (Sutherst, 2004); in contrast, an increase in the animal population can lead to an increase in the rodent population and consequently an increase in the risk of Ft (Gürcan, 2014). Most of the positive samples were collected in irrigated and cropland areas (Chakwal, Gujranwala, and Faisalabad), that are interconnected with each other. In fact, compared to areas with forests and water, a previous study showed that Ft was primarily observed in croplands where rodents were found to be a common vector for Ft infections (Hightower et al., 2014). Our findings show that it is difficult to define the relationship between the characteristics of the studied analytes and the persistence of Ft in soil. However, this study provides insight between Ft and its relationship with soil texture and metals requiring that further investigations are conducted in order to ascertain more definitive outcomes. Metals including iron, zinc, copper, and manganese play a key role in biochemical reactions and are present in most of the organism as constituents of different biological molecules (Hood and Skaar, 2012). Copper, magnesium, manganese, zinc, nickel, cadmium and sodium were observed to have a role in persistence of Ft in soil, although further investigation of these elements is needed. Indeed, some studies emphasized that survival and persistence of different organisms in soil is influenced by the soil texture (Nicholson et al., 2005). For instance, Burkholderia pseudomallei and Ft has the potential to survive for 30 months and 14 weeks, respectively, in clay soil (Thomas and Forbes-Faulkner, 1981;Cieslik et al., 2018). Generally, soil with an increase in clay support enhanced growth of organisms as compared to soil with increased concentration of sand (Locatelli et al., 2013). A potential association between Ft DNA and Magnesium could be correlated with its role while maintaining negative charge in LPS of outer membrane and therefore stability of organism in unfavorable conditions (Leadbetter and Poindexter, 2013;Wu et al., 2016). Similar to Ft, intracellular bacteria including Salmonella enterica, Listeria monocytogenes, Brucella abortus, and Yersinia pestis require zinc for intracellular survival and replication (Desrosiers et al., 2010;Corbett et al., 2012;Ma et al., 2015) Cadmium is normally toxic to bacteria however it has several mechanisms to develop resistance including (1) enzymes which make cell wall impermeable, (2) mechanism of efflux, (3) binding of metal ions, and (4) enzymes which convert toxic metals to non-toxic (Ron et al., 1992;Abbas et al., 2017). The persistence of Ft in cadmium high concentration soil presents a potential mechanism of resistance that need to be explored in future studies.
A large proportion of ruminants (goat, sheep, cattle, and buffalo) were found to be seropositive for anti-Ft antibodies. All the seropositive animals were clinically healthy and were devoid of symptoms suggestive of Ft infection. This is not surprising, because both active and passive forms of Ft have previously been evidenced in humans, rodents, and wild animals, (Wobeser et al., 2009) and therefore a varying rate (4-19%) of asymptomatic tularemia could be expected (Kiliç, 2010). Indeed, the clinical form of the infection depends on the route of entry of the pathogen into the body and the immune status of the host (Gürcan, 2014). The percent prevalence of anti-Ft antibodies was significantly higher in large ruminants (cattle and buffalo) than small ruminants (sheep and goat). Interestingly, previous studies have shown that cattle and sheep are comparatively resistant to Ft and remain asymptomatic (Mörner et al., 1988;Pfahler-Jung, 1989;Winn and Koneman, 2006). Some animals are very sensitive to tularemia, and succumb to infection soon after development of specific antibodies, whereas cattle, pigs, and sheep are more resistant, despite even having a significantly high antibody titer (Hopla, 1974;Valdes and Valdes, 2018). Animals may not show clinical signs even after years of infection, possibly owing to the presence of post-exposure cellular and humoral immunity (Bevanger et al., 1994;Ericsson et al., 1994;Magnarelli et al., 2007). A recent study revealed that both small and large ruminants, independent of their genders, are frequently exposed to soil contaminated with Ft supporting the findings by previous studies showing that cattle maintain some degree of Ft resistance (Dienst Jr, 1963;Scheel et al., 1992;Jacobs, 2002). Similarly, another study showed that both male and female mice are susceptible to Ft where male mice develop sever clinical signs and high mortality as compared to female (Sunagar et al., 2016). These observations may be due to differences in immune response development in male and female mice. An example of this was shown during Streptococcus pneumonia infections, where immune system inefficiency was noted to be higher in female mice as compared to male mice (Wiemken et al., 2014).
In conclusion, this study provides a preliminary confirmation of the presence of Ft in the Punjab Province of Pakistan and its potential association with several soil characteristics (macro-and micro-nutrients) at a higher statistical resolution. Future studies involving unexplored geographical areas of the country, disease reservoirs including rodents, and evaluation of the at-risk human population (the farmers and veterinarians), are needed.

ETHICS STATEMENT
Blood samples were collected from cattle, buffalo, sheep, and goat following the guidelines of International Animal Care and Use Committee (IACUC) with prior consent of the farm's owner. All samples were analyzed after approval from the Ethical Research Board at the University of Veterinary and Animal Sciences, Lahore, Pakistan vide approval No: DR/236 dated 16th May, 2013.

AUTHOR CONTRIBUTIONS
BMJ, MZS, MR and KM: conceived and designed the study; JM, MZS, MH, ZUH, HRC, MTG, TJ, and MAA: sample collection and laboratory procedure across the study districts; JM, MZS, MN, TA and MHC: data analysis; JM, MZS, GSK and BMJ: manuscript write-up and necessary editing. All authors read and approved the manuscript.

FUNDING
The study was financially supported by the Defense Threat Reduction Agency, Basic Research Award number HDTRA1-10-1-0080 to The Pennsylvania State University, USA (Prime Awardees) and University of Veterinary and Animal Sciences, Lahore Pakistan (Subcontract Awardees).