Intrinsic and extrinsic drivers of transmission dynamics of hemorrhagic fever with renal syndrome caused by Seoul hantavirus.

Seoul hantavirus (SEOV) has recently raised concern by causing geographic range expansion of hemorrhagic fever with renal syndrome (HFRS). SEOV infections in humans are significantly underestimated worldwide and epidemic dynamics of SEOV-related HFRS are poorly understood because of a lack of field data and empirically validated models. Here, we use mathematical models to examine both intrinsic and extrinsic drivers of disease transmission from animal (the Norway rat) to humans in a SEOV-endemic area in China. We found that rat eradication schemes and vaccination campaigns, but below the local elimination threshold, could diminish the amplitude of the HFRS epidemic but did not modify its seasonality. Models demonstrate population dynamics of the rodent host were insensitive to climate variations in urban settings, while relative humidity had a negative effect on the seasonality in transmission. Our study contributes to a better understanding of the epidemiology of SEOV-related HFRS, demonstrates asynchronies between rodent population dynamics and transmission rate, and identifies potential drivers of the SEOV seasonality.

Introduction

There are 60,000–100,000 cases of hemorrhagic fever with renal syndrome (HFRS), a rodent-borne zoonosis, reported annually, globally [1, 2]. HFRS is caused by subtypes of hantavirus, each of which is associated with a distinct rodent species [3]. Humans are usually infected by inhaling air contaminated with saliva, feces, or urine from infected rodents [4]. China is the main epidemic region for HFRS, and accounts for 90% of cases, globally [1]. The two predominant hantavirus strains circulating in endemic areas of China are Hantaan virus (HTNV), which is carried by striped field mice (Apodemus agrarius), and Seoul virus (SEOV), which is carried by Norway rat (Rattus norvegicus) [5].

HFRS epidemics vary significantly across seasons and are influenced by extrinsic factors across regions, like climate [6-10]. Following the Hantavirus Pulmonary Syndrome outbreaks that occurred in the Four Corners region of the southwestern United States in 1993, the correlation between hantavirus infections and climatic conditions was described using a cascade hypothesis [11]. The hypothesis posited that favorable climate conditions would lead to more available food, and to greater rodent population sizes, thereby enhancing the risk of hantavirus infections. Similarly, the 2007 HFRS outbreak in temperate southern Europe may have been caused by increased population density of the bank vole, the vector of Puumala virus (PUUV), which, in turn, may have resulted from abundant food due to preceding warmer than usual autumn and winter weather [12, 13]. In central China, higher temperature and precipitation in the previous summer led to favorable food conditions for the striped field mouse, the rodent host of HTNV, which led to an autumn peak in incidence of HFRS [14].

Hantaviruses carried by species of wild rodents, such as HTNV and PUUV, have been extensively studied, but SEOV is vectored by the Norway rat, a domestic rodent that is found in urban environments and is associated with humans [15-17]. SEOV infections, have milder symptoms and lower fatality rates (1–2% vs. 5–10%) than HTNV infections, so less attention is focused on SEOV [16, 17]. Besides, SEOV infections have a higher asymptomatic infection rate of 8–20%, compared with HTNV at 1–4% [18]. In the past two decades, SEOV-related HFRS cases have been reported in the United Kingdom (2012) [19], France (2014) [20], the United States, and Canada (2017) [21] and most provinces of China [22-26]. It has been estimated that more than 25% of hantavirus infections in South Korea and China are caused by SEOV [27, 28]. To better control the potential threat, it is important to understand how critical factors impact SEOV-related HFRS transmission dynamics.

Here, we explore the roles of related intrinsic and extrinsic factors in transmission dynamics of HFRS caused by SEOV in the pre-vaccination era. Intrinsic factors refer mainly to herd immunity or host population dynamics, and extrinsic factors include external forces such as climatic or environmental factors. We collected detailed epidemiological and rodent trapping data from a typical SEOV-endemic area of China, Huludao City, where the Norway rat is the dominant rodent species in residential areas [29]. Additionally, we investigated potential causal relationships among the environmental variables, the population dynamics of the rodent reservoir, and the incidence of HFRS. We subsequently constructed a mathematical model to quantify intrinsic transmission dynamics and extrinsic effects.

Materials and methods

Data

Huludao City (40°56′N, 120°28′E) in the Liaoning Province of China (S1 Fig) has a temperate monsoon climate with a hot and rainy summer (mean temperature, 23 °C; mean monthly precipitation, 123 mm) and a cold and rainless winter (mean temperature, -5.7 °C; mean monthly precipitation, 2 mm). The city covers 10,434 square km and has a population of 2.6 million people. HFRS has always posed a severe threat to human health in Huludao City since 1984 [30]. To control HFRS transmission, mass vaccination campaigns against hantavirus were conducted beginning in 2005, and combined with rat extermination programs during 2005 and 2012 [31].

Local demographic data were collected from the Liaoning Statistical Yearbooks. Land cover changes from 1998 to 2015 in Huludao City were derived from the annual European Space Agency (ESA) Climate Change Initiative (CCI) land cover maps, with a 300 m spatial resolution. Meteorological data for Huludao City from 1998 to 2015 were obtained from the Chinese Bureau of Meteorology, and then processed into monthly climate data, including monthly mean maximum temperature, mean temperature, mean minimum temperature, cumulative precipitation, relative humidity, and absolute humidity. Absolute humidity is calculated from temperature and relative humidity [32].

Data of human HFRS cases in Huludao City from 1998 to 2015, were obtained from the Huludao City Center for Disease Control and Prevention, China. All HFRS cases were confirmed according to the standard diagnosis set out by the Ministry of Health of the People’s Republic of China, and by detecting antibodies against hantavirus in serum samples. Serum samples were tested for immunoglobulin (Ig) G and IgM antibodies against HTNV and SEOV by indirect immunofluorescent assay (IFA) [33]. Serum-epidemiological surveys on human hantavirus infection were conducted annually between 1998 and 2015. Anonymous (non-personal) information was used.

Rodent population density was investigated indoors and outdoors by the powder-trace method on a monthly basis, obtained from the historical literature [34]. According to the criteria, 400 powdered panels are placed in the special sites for spot checks across the city every month, of which 240 panels are set in restaurant, hotel, station and other public places where rodents are commonly found, and 160 panels are set in general household. Two powdered panels were placed in each room (about 15 m2) in selected spots. Outdoors, each panel was placed along a wall at 5–10 m intervals. All panels were set at night and recovered before sunrise. Rodent population density was calculated as the number of positive powdered panels (with footprints or tail tracks from rodents) divided by the number of effective powdered panels used.

Ethics statement

The study protocol was reviewed by the institutional review board of the Huludao Municipal CDC and ethics approval was not required. We have received consent from home/land owners to collect rodent data on private land and in private homes. The Animal Ethics Committee of the Huludao Municipal CDC also waived approval for this study.

Wavelet analysis

Wavelet analysis is widely used in ecology and epidemiology studies to explore the variety in the periodicity of a time series through its decomposition properties [35]. Here, the Morlet wavelet was used to detect the non-stationary characteristics of the incidence of HFRS fluctuations over time, and the bootstrapping method [36] was performed to test the statistical significance of the results. In the significance test, 1000 surrogate datasets were simulated by bootstrapping to test the null hypothesis, where a P value of < 0.05 was considered to be statistically significant.

Convergent cross mapping

Convergent cross mapping (CCM) is used to detect nonlinear causal relationships between time series of HFRS incidence (Y) and environmental factors (X), which is designed to measure causality in nonlinear dynamical systems. CCM can help to identify bidirectional causality (i.e. X and Y are mutually coupled) or unidirectional causality (e.g., X time series variable in a system has a causal influence on Y, but not vice versa) [37, 38]. In this study, CCM was used to identify time-delayed effects of a causal interaction between time series of environmental variability, the population dynamics of rodent hosts, and HFRS incidence based on nonlinear state space reconstruction [39]. In a system where x causes y, Sugihara et al. purposed that the state of x(t) can be estimated from the reconstruction of y(t) using the nearest-neighbor forecasting method, also called cross mapping [40]. Pearson’s correlation between the estimated x(t) and observed x(t) reflects the causal effect of x on y, called “cross map skill” [39]. In analysis, the embedding dimension (E) were set according to the simplex projection results [40]. Time-delayed effect (tp) were set as 0–6 months for detecting the time lags [41].

Given the seasonality of HFRS incidence, climate variables, and rodent population density, the seasonal component and the response of HFRS risk and rodent population density to the anomalies of other variables were examined to avoid spurious correlations. For each variable, 1500 surrogate time series having the same degree of shared seasonality, but with randomized anomalies, were generated for seasonal surrogate test. For a specific variable, the month of year average (seasonal cycle) was calculated first, and the seasonal anomaly was represented by the difference between the observed value and the seasonal cycle. Then the random shuffling of the time series data of seasonal anomalies was added back to the season average [42]. The above analysis was implemented in the “rEDM” package in R.

The SEOV-related HFRS transmission model

We estimated the epidemiological parameters for HFRS epidemics in Huludao City by fitting the time series from the observed monthly incidence and rodent population density to a discrete-time susceptible-infection model in the Bayesian framework (Table 1) [43].

Table 1

The goodness of fit for the candidate models for Norway rat population and HFRS dynamics.

Climate variable	R2	DIC
Norway rat population dynamic model
AHt-3	0.50	9.32
TMINt-4	0.50	9.17
AHt-3, TMINt-4	0.47	11.41
none	0.50	9.17
HFRS dynamic model
RHt-1	0.75	7.57
none	0.72	7.89

Norway rat population dynamic model

First, we modelled the Norway rat population dynamic model based on Verhulst-Pearl logistic growth model [44], without consideration of SEOV infection, and identified the effect of climate on population fluctuation. Population fluctuations are regulated by an identical seasonal birth rate (b), natural mortality (m), environmental carrying capacity (K), and potential extrinsic effect (θ) in Eq (1). There is no significant sex difference in the Norway rat population in our study area [45]. To reduce the dimensionality of the model, age classes were omitted.

Here, 1 month is used as the time interval in the time-discrete model; b represents an array consisting of 12 values, representing the reproduction rates of the rodent population over 12 months. As Norway rat breeds throughout the year in the regions studied [45], all of the reproduction rates were set at >0. m denotes the natural mortality of the Norway rat population, whose mean life span is generally not more than 1 year [46]. K represents the environmental carrying capacity and is used to model the environmental pressure from intraspecific competition for food availability or other resources. For identifying the extrinsic effect on Norway rat population fluctuation, we constructed θ to reflect the potential climate effect. Based on the results of the correlation analyses and CCM, minimum temperature and absolute humidity were potential drivers for rodent population dynamics. We constructed equation θ = δTMIN+δAH+ε, where TMIN and AH refer to monthly mean minimum temperature and absolute humidity, respectively, and the parameters δ and δ quantify the degree of climate effect and ε is an intercept.

Integrated SEOV-related HFRS transmission model

The mathematical model for host rodent population outlined in Eq (1) can be extended to incorporate rodent-human SEOV transmission. SEOV infections were modeled as occurring among rodents and from rodents to humans in Eqs (2)–(4). The rodent host population (R) was divided into two classes: susceptible (R) and infected (R). A susceptible rodent could be infected when bitten by an infectious rodent or by inhaling air contaminated with SEOV in the excrement from infectious rodents [47].

In the above equations, all newborn rats were regarded as being susceptible to infection because there is no vertical transmission of SEOV [47]. Given that there is evidence that SEOV infection has no evident impact on death and fertility in Norway rat [48], we did not differentiate mortality (m) and birthrate (b) between R and R β is the transmission rate between rodents and α allows for non-linear contacts from rodents to humans. H represents the whole human population and H represents susceptible human population. Observed cases Y are assumed to be correlated to the unobserved number of cases according to a binomial observation process with constant observation rate ρ, due to non-reporting from mild symptoms overlooked by the health system and asymptomatic infections. The prior information of ρ (0.006, 0.002–0.009) was assessed from our epidemiological survey. β represents SEOV transmission, which contains seasonal contact rates between rodents and humans β and climate-driven transmission potential β, β = ββ. Based on the results of the correlation analyses and CCM, a causal relationship between relative humidity and HFRS incidence was detected. Therefore, we constructed that β = δRH γ t-1, where RH was relative humidity with 1-month lag, and δ was used to quantify the effect of relative humidity on HFRS transmission dynamics and γ allowed for the non-linearity of the effect. The results of the estimation for γ and δ were the mainly focused where the direction of the climate effect depends on whether γ (δ>0) is positive or negative. All parameter settings are shown in Table 2.

Table 2

The prior setting and posterior probability for parameters in the optimal model.

Description	Par.	Prior/Range	Posterior mean	Std
Norway rat population dynamic model
Environment carrying ability	K	12 [5, 25]	8.24	1.11
Natural mortality	m	0.07 [0.04, 0.5]	0.16	0.08
Initial rodent population	RN,0	1.4 [0.5, 3]	1.31	0.24
Reproduction rate in Jan.	b(1)	1 [0.01, 5]	0.83	0.16
Feb.	b(2)	1 [0.01, 5]	0.90	0.16
Mar.	b(3)	1 [0.01, 5]	1.06	0.17
Apr.	b(4)	1.5 [0.01, 5]	0.93	0.16
May.	b(5)	1.5 [0.01, 5]	0.77	0.15
Jun.	b(6)	1 [0.01, 5]	1.24	0.19
Jul.	b(7)	1 [0.01, 5]	1.14	0.17
Aug.	b(8)	1 [0.01, 5]	0.85	0.17
Sep.	b(9)	1.2 [0.01, 5]	0.93	0.17
Oct.	b(10)	1.2 [0.01, 5]	1.22	0.18
Nov.	b(11)	1 [0.01, 5]	1.15	0.18
Dec.	b(12)	1 [0.01, 5]	0.98	0.18
HFRS dynamic model
Scale factor	δ3	0 [0, +∞]	17.68	4.58
Forcing factor	γ	0 [-∞, +∞]	-1.23	0.20
Initial virus-carrying rodent	RI,0	0.1 [0, 1.5]	0.02	0.01
Infectious rate (RI, RS)	βR	1 [0.01, 1]	0.14	0.02
Error term of rodent population	τ	0 [0, 3]	2.24	0.21
Observation rate	ρ	0.006 [0.002, 0.009]	0.006	0.002
Nonlinear effect of rodent population	α	1[0, 100]	9.07	1.15
Seasonal contact rate (RI, human) in Jan.	βsea (1)	1 [0.01, 100]	47.12	3.45
Feb.	βsea (2)	1 [0.01, 100]	89.26	3.86
Mar.	βsea (3)	5 [0.01, 100]	84.24	3.77
Apr.	βsea (4)	5 [0.01, 100]	87.09	3.72
May.	βsea (5)	5 [0.01, 100]	54.16	3.16
Jun.	βsea (6)	1 [0.01, 100]	65.05	3.42
Jul.	βsea (7)	1 [0.01, 100]	41.96	5.71
Aug.	βsea (8)	1 [0.01, 100]	32.56	5.13
Sep.	βsea (9)	3 [0.01, 100]	40.22	5.25
Oct.	βsea (10)	3 [0.01, 100]	40.75	4.41
Nov.	βsea (11)	3 [0.01, 100]	42.73	4.33
Dec.	βsea (12)	5 [0.01, 100]	45.98	4.46

Std, standard deviation of the sample mean.

The Markov chain Monte Carlo (MCMC) sampling method was used for parameter estimation with the MATLAB toolbox delayed rejection adaptive metropolis algorithm [49]. The values of the parameters were estimated by sampling from the convergent posterior distributions of the Markov chains. Gelman-Rubin-Brooks MCMC convergence diagnostic was used to test convergence by identifying whether the series came from a stable distribution or not [50]. Prior distributions of the parameters were mainly set according to scientific literature and values of biological relevance, but those parameters that lacked prior information were set to a large variance (100) with no boundaries (Table 2). The chains were set to perform 1 million iterations with burn-in of 20,000. Deviance information criterion (DIC) was used to measure the goodness of fit for the models with various framework, and a difference value >10 was regarded to be significant. R2 was used to measure the interpretability of the model for the observed data.

Results

HFRS epidemics in Huludao City

During 1998 to 2015, 6796 HFRS cases were reported in Huludao City, with an annual incidence from 2.50 (1/100,000, minimum in 2008) to 26.99 (1/100,000, maximum in 2005). The pattern of HFRS incidence showed a significant seasonality with an annual peak in spring from March to May (Fig 1, S2 Fig). The epidemic variability could be divided into three periods: period I (1998–2004), characterized by a high-level of HFRS incidence of over 17.6 (1/100,000) and a stable periodicity without limited intervention; period II (2005–2012), vaccination campaign and rat extermination program were implemented simultaneously and the local HFRS incidence decreased dramatically to a relatively low level at 10.0 (1/100,000); period III (2013–2015), the incidence of HFRS gradually rebounded. By using wavelet analysis, we found changes in endemic periodicity after 2006 (Fig 1B), possibly due to the rat extermination program and mass vaccination campaign that were conducted beginning in 2005. Here, to reduce the effect of human interventions on disease dynamics, we focused our analysis on the study period between 1998 and 2004 (period I). In period I, a total of 3914 HFRS cases were reported. The local rodent community mainly comprised Norway rat and the house mouse, with the former accounting for 94.52% of the total. The population density of the rodents ranged from 4.89% to 7.28%, with two relatively low values in January of 1998 and 1999.

Fig 1

Seasonality and periodicity of hemorrhagic fever with renal syndrome (HFRS) incidence in Huludao City, 1998 to 2015.

(A) HFRS incidence per month. Red bars represent disease incidence in spring (March to May), while blue bars represent the other months. (B) Periodic variety of HFRS incidence. Colors indicate the power of the wavelet, where red to blue represent strong to weak power and the black line indicates the maximum power. The white line represents statistical significance (P < 0.05).

Causal relationships between environmental factors, Norway rat population density and HFRS transmission risk

We first examined the impact of land use and land cover changes on HFRS transmission. Land use changed very slightly in period I (S3 Fig), so the effect of habitat and land use change on rodent hosts and disease transmission were minimal during this time period. We then conducted a cross-correlation analysis of the climate variables, rodent population density, and HFRS incidence. The Norway rat population was found to be positively correlated with mean temperature (Pearson’s r = 0.23, P < 0.05), maximum temperature (r = 0.25, P < 0.05), minimum temperature (r = 0.23, P < 0.05), cumulative precipitation (r = 0.28, P < 0.05), and absolute humidity (r = 0.24, P < 0.05), with the maximum correlation coefficients occurring at 3, 3, 4, 4, and 4-month time lags, respectively (S4A Fig). HFRS incidence was strongly correlated with the climate variables. All the tested climate variables showed significant negative correlations with HFRS incidence, and maximum cross-correlations between temperature and HFRS incidence occurred at a lag of 3 months and other climate variables at a lag of 2 months (S4B Fig). Given the survival time of hantavirus outside the host and the incubation period for HFRS [51, 52], the longest time lag was set as 6 months.

Convergent cross mapping (CCM) was used to detect causality between these time series. However, rodent data had weak cross mapping skills for the climate variables. Only a marginally significant, weak causal effect was found between minimum temperature (with a 4-month lag, tp = -4) and absolute humidity (with a 3-month lag, tp = -3) on rodent population density, with cross map skills of 0.18 and 0.12, respectively (S5 Fig). To test the non-stationary (transient) correlation between Norway rat and climate variation, we conducted a further test of wavelet coherence analysis. The results showed scattered and small-sized distribution of significant areas with inconsistent phase differences (S6 Fig). In all, the combined results indicate that there might be a weak causal relationship between Norway rat and climate variables. Only relative humidity with a 1-month lag was identified as a significant causal factor for HFRS transmission (cross map skill, 0.86) by the CCM method (S7 Fig). The significance test distinguished anomalies of climate effects from shared seasonal cycle effects (S7B Fig).

Model validation

To quantify the effect of climatic factors on rodent population density, we tested the models (see Materials and methods for details) with θ containing (i) minimum temperature, (ii) absolute humidity, (iii) both variables, and (iv) none of these climate variables. The goodness of fit for all the candidate models is provided in Table 1. The model without the climate influence was regarded as optimal to infer the population dynamics of Norway rat in a local urban setting (DIC = 9.17, R2 = 0.50, Fig 2A). The reproduction rates of the local Norway rat population, as estimated by our model, peaked in June, July, and October, and reached the lowest level in May and August (Table 2). The average life span of the Norway rat is estimated to be about 6 months, and the local environmental carrying capacity is estimated to be 8.24. The trace plot of Markov chains and the posterior distribution for the parameters in the optimal rodent model is shown in S8 Fig. All the parameter traces have passed the convergence test.

Fig 2

Fitting results of the dynamic model.

(A) Time series of Norway rat population density and (B) HFRS cases. Black points indicate the real observations (obs) and lines indicate the simulated (sim) time series. Shaded areas indicate the 95% credible interval. (C) The season pattern of Norway rat population density and (D) HFRS cases. The 5% and 95% of simulated data are shown in shaded area.

Based on the causal relationship tested by CCM, we constructed a HFRS dynamic model containing the infected rodent population and seasonal contact rate to quantify the specific effect of relative humidity on risk of HFRS (see Materials and methods for details). We constructed two models: one containing relative humidity and the other without. The model with relative humidity had a slightly better fit when the simulated HFRS cases correlated well with the observations and 75% of the variance was explained (R2 = 0.75, DIC = 5.57) (Fig 2B), relative to the model without relative humidity (R2 = 0.72, DIC = 7.89). The optimal models also captured the main seasonal patterns of both the rodent population density and reported cases (Fig 2C and 2D). The parameter γ in β = δRH γ t-1, was estimated to be −1.23 (95% CI, -1.62– -0.83) (Fig 3A), indicating that relative humidity had a negative impact on HFRS transmission (Fig 3B). Additionally, during winter–spring, the estimated seasonal contact rates peaked from February to April (Fig 4), which may result from the behaviors of human and rodent. Most of the patients are farmers (77%) and it’s the slack time of winter-spring when people stay at home. Besides, due to the severe weather outside during that time, Norway rats would live and feed in closer proximity to human residence for favorable living conditions and available food, which could be a reason for the high contact rate. The trace plot of Markov chains and the posterior distribution for the parameters used in the optimal HFRS transmission model are shown in S9 Fig.

Fig 3

Force of infection and transmission rate of the seasonality of HFRS risk.

(A) Posterior distribution of climate-driven transmission potential (β), represented by forcing factor (γ). The power function with negative exponent (γ) means a negative relationship between transmission and relative humidity. (B) The effect of seasonal changes in relative humidity on climate-driven transmission potential estimated from Eq (4).

Fig 4

Season epidemics, rodent population dynamics and estimated contact rate.

Error bars show the 95% credible intervals.

Discussion

We used a mathematical model to quantify the impact of intrinsic and extrinsic factors on risk of HFRS, and to test the mechanistic understanding of HFRS transmission dynamics in a SEOV endemic area of China. By analyzing time-series data of Norway rat populations and SEOV infections, we found that the population dynamics of Norway rat in urban settings were relatively well-predicted by a simple model which excluded climatic conditions. Unusually, we demonstrated that a potential biotic driver (relative humidity) could enhance predictive ability of the HFRS transmission analysis.

Our results highlight the crucial role of mass immunization campaigns and rat extermination program in HFRS transmission control. It should be noted that SEOV infections began to rise again after 2013 in the post-vaccination era. In HFRS endemic areas where hosts are mainly wild rodents, fluctuations in the population density of wild rodents, such as striped field mouse or bank vole, are influenced by climatic factors and can mediate the effect of climate on the risk of HFRS transmission, mainly because of the close association between climate and host food supply [53]. However, in the SEOV-endemic area, the CCM results showed that the change in population density of the Norway rat was insensitive to climate variability. Our results differed from those of a previous study which found a relationship among the virus-carrying index, climate variables, and HFRS incidence by using structural equation modeling (SEM). This study reported that the Norway rat served as an important mediator of disease transmission when climate variability was found to influence the risk of HFRS [54]. This inconsistency may be due to a difference in sampling frequency as we surveyed on a monthly basis, while Guan et al. (2009) sampled quarterly. Additionally, the SEM analysis used a latent variable to account for measurement error but did not evaluate the contribution of specific variables to the transmission dynamics.

Norway rat is a domestic rodent whose lifestyle differs greatly from those of wild rodents. Norway rat lives in residential areas and feeds on various food items from humans instead of field crops [55]. Food resources and habitats are relatively stable and change only slightly with climate variability. Additionally, Norway rats breed throughout the year in Huludao City, where the winter temperatures rarely drop below −10 °C and residential areas have warmer conditions because of the winter heating policy in China. Winter breeding of wild rodents is limited because of low temperatures [56].

Our findings indicated that an abiotic factor, relative humidity, was a critical indicator of HFRS risk in the SEOV-endemic area of Huludao City, which is consistent with previous studies in which SEOV was the main virus type [57, 58]. We assumed that relative humidity may be associated with the survival or infectivity of SEOV in the environment, the activity of the rodent, and human hosts or the transmission process, although the underlying mechanism is not clear. We do not know of any relevant experiments on how environmental conditions affect survival of SEOV outside the host, while PUUV and HTNV have been confirmed to have longer infective periods at low temperatures and high humidities in experimental environments [59]. In our model, interannual variation in SEOV infections was partly explained by changes in relative humidity. A less intense outbreak in 2000 might be associated with a lower relative humidity-driven transmission potential (Fig 2B).

Some limitations in our study should be mentioned. First, rodent population density was extremely low in the springs of 1998 and 1999, which reduced the overall explanatory ability of our model. Second, prevalence of SEOV infection in Norway rats was assessed by mathematical model due to lack of available data. Third, to ameliorate the influence of vaccination campaign and other interventions, only the selected time range of data was used to detect the relationship between intrinsic/extrinsic drivers and SEOV transmission. Further efforts should be made to improve model reliability through efforts such as addition of data to the available dataset. More attention should be paid to the association between economic development, such as infrastructure improvement, the spatial distribution of Norway rat, and induced SEOV infections across districts in subsequent studies.

The SEOV-related HFRS is posing a health threat to an increasing number of people. Although the Norway rat is found everywhere, SEOV cases are not. We suspect this may mainly be due to increased surveillance and attention focused on this pathogen, as SEOV previously caused many asymptomatic human infections previously. In conclusion, based on the longitudinal and complete dataset, our study yielded a proper framework for understanding the intrinsic transmission dynamics and extrinsic effects on HFRS risk caused by SEOV, especially at the human-animal-environment interface. The framework can be flexibly adjusted when more information is available. Furthermore, we provided clues to potential environmental drivers on SEOV transmission dynamics, which would be useful for further research related to public health issues.

The geographical location of Huludao City.

(A) Huludao City (orange) is located in northwest of Liaoning Province, adjacent to the Bohai sea. (B) Liaoning Province is highlighted in red in China map.

(TIF)

Click here for additional data file.

Seasonal distribution of HFRS cases over three periods.

(TIF)

Click here for additional data file.

Land cover composition in Huludao City from 1998 to 2015.

The land cover classification is according to IPCC land categories. Other includes shrubland, sparse vegetation, bare area and water.

(TIF)

Click here for additional data file.

The relationship between hemorrhagic fever with renal syndrome (HFRS) incidence, Norway rat population density and environmental variables with time lags.

(A) Cross-correlations between the population density of Norway rats and the climate variables. (B) Cross-correlations between HFRS and climate variables. TM, monthly mean temperature; TMAX, monthly mean maximum temperature; TMIN, monthly mean minimum temperature; PREC, monthly cumulative precipitation; RH, monthly mean relative humidity; AH, monthly mean absolutely humidity; RN.(t), Norway rat population density at time t; IN.(t), HFRS incidence at time t.

(TIF)

Click here for additional data file.

The causal relationship between climate variables and rodent host population, and between climate variables and HFRS incidence.

(A) The convergent cross mapping results for Norway rat population density and climate variables with time lags (tp) and for HFRS incidence (B). The causal relationship was regarded as significant (P < 0.05) when the solid line (cross map skill of the observed data) exceeded the dashed line (surrogate data). TM, monthly mean temperature; TMAX, monthly mean maximum temperature; TMIN, monthly mean minimum temperature; PREC, monthly cumulative precipitation; RH, monthly mean relative humidity; AH, monthly mean absolutely humidity.

(TIF)

Click here for additional data file.

Wavelet coherence between Norway rat population density and climate variables from 1998 to 2004.

Arrows showed phase differences between Norway rat population density and climate variables. The straight down represents climate variables lead Norway rat population density 1/4 period. Legend is wavelet coherence level. The parts inside white contour shows a significant relationship. The parts outside the transparent cone have affected by edge effect. TM, monthly mean temperature; TMAX, monthly mean maximum temperature; TMIN, monthly mean minimum temperature; PREC, monthly cumulative precipitation; RH, monthly mean relative humidity; AH, monthly mean absolutely humidity; RN, Norway rat population density.

(TIF)

Click here for additional data file.

Detecting cross-map causality beyond shared seasonality of environmental drivers on hemorrhagic fever with renal syndrome.

(A) Time series of hemorrhagic fever with renal syndrome (HFRS) incidence (orange) and relative humidity (green). (B) The convergent cross mapping results at different time series length for HFRS incidence and relative humidity (RH) and surrogate RH (surr_RH) with 1-month lag (tp = -1). The causal relationship was regarded as significant (P < 0.05) when the solid line (cross map skill of the observed data) exceeded the dashed line (surrogate data).

(TIF)

Click here for additional data file.

Estimated parameters in Norway rat population dynamic model.

(A) The trace plot of markov chains of the parameters after 1 million burn-in. (B) The posterior distribution of the parameters after 1 million burn-in. All the estimations have passed the Gelman-Rubin-Brooks MCMC convergence diagnostic.

(TIF)

Click here for additional data file.

Estimated parameters in HFRS transmission dynamic model.

(A) The trace plot of markov chains of the parameters after 1 million burn-in. (B) The posterior distribution of the parameters after 1 million burn-in. All the estimations have passed the Gelman-Rubin-Brooks MCMC convergence diagnostic.

(TIF)

Click here for additional data file.