The Association Between Climate Factors and Scrub Typhus Incidence in South Korea in 2001–2019

DOI: https://doi.org/10.21203/rs.3.rs-1218122/v1

Abstract

Background: Global warming has compromised human health by increasing the incidence of infectious diseases. Scientific evidence is required to expand the knowledge of the association between meteorological factors and the incidence of infectious diseases. Our study focused on meteorological factors such as ambient temperature (AT), humidity, and scrub typhus incidence.

Objective: We aimed to investigate the long-term effects of AT and relative humidity (RH) elevation on scrub typhus incidence in South Korea.

Method: Meteorological variables were combined with scrub typhus cases reported from 2001 to 2019 in South Korea. A generalized additive model was used to explore the nonlinearity of the lagged association between meteorological variables and weekly scrub typhus incidence. To explore the long-term association between meteorological factors and scrub typhus incidence, the difference between annual mean AT or the annual number of heatwave days during 2001 to 2019 and those of the period 1971–2000 were linearly regressed on annual cumulative scrub typhus incidence.

Results: Association between weekly scrub typhus incidence and mean temperature or relative humidity with 15 weeks lags showed V- or U- shaped relationship. Above the threshold temperature (14.9°C to 17.0°C), scrub typhus incidence increased by 1.85% (95% CI: 1.5, 2.2) per 0.1°C elevation in mean temperature with 15 weeks lags. Scrub typhus incidence increased by 3.7% (95% CI: 2.7%, 4.8%) and 2.5% (95% CI: 1.6%, 3.3%) per 0.1°C increase in annual mean temperature and per one day increase in heat wave days.

Conclusion: Mean AT elevation and RH increase in summer were associated with an increased incidence of scrub typhus in the fall. Increases in annual mean AT and annual number of heatwave days were associated with an increase in the incidence rate during 2001–2019. 

Introduction

Global mean ambient temperature (AT) has increased by approximately 1°C as compared with the pre-industrial era. The Intergovernmental Panel on Climate Change (IPCC) has warned of more frequent flooding, reductions in crop production, and disruptions of ecosystem functions 1. Climate change alters environmental conditions and, consequently, increases the incidence of infectious diseases, including vector-borne diseases 2,3. Increased rainfall may induce an increase in vector density, and high temperatures may lead to an extended incubation period for the vectors. AT and relative humidity (RH) 4 may be related to determining vector survival, such as affecting the distribution and life cycle of vectors 5–7. Ticks carry Rickettsia as an endosymbiont and transmit these pathogens to humans through human–vector interactions.

Mites are the most common vectors of infectious diseases worldwide 8,9. Among these, scrub typhus, denoted as “tsutsugamushi disease” in Japan, has become a critical concern in endemic regions as well as in other areas as traffic to and from the endemic regions increases in numbers. 10,11.

AT and RH are known to be key factors in the development of Leptotrombidium, the vector mites of scrub typhus 12. The response of chiggers (“mites”) to climate change is associated with scrub typhus incidence. In a previous study, climate factors such as AT, RH, mean wind speed, precipitation, and sunshine duration ahead of the outbreak season were associated with scrub typhus incidence during the peak periods 13,14. The geographical distribution of scrub typhus vectors across the Republic of Korea (ROK) 15, which has been studied since 1980, has been dominated by two species, Leptotrombidium pallidum in the northeastern provinces and L. scutellare in the southwestern provinces 16. Japanese studies found that L. scutellare grew from mites to adults between April and August, whereas L. pallidum emerged from June and peaked in late July, August, and early September 17.

Scrub typhus usually causes symptoms such as fever, malaise, skin rash (e.g., eschar), and conjunctival injection. However, severe illnesses such as secondary acute meningoencephalitis with neurological sequelae or even fatalities can occur in rare cases 18–21. To date, the geographical spread of the disease and its vectors has expanded 10, and climate change and global warming are suspected to be the cause of the spread of scrub typhus 2,3,22.

Considering the geographical distribution of vectors in the country, we assumed that the incidence of scrub typhus in the ROK might be influenced by climatic conditions. Hence, we aimed to explore the association between climatologically preferable conditions and scrub typhus incidence from 2001 to 2019 in the ROK, with a focus on AT and RH.

Results

1. Seasonal variation of climate factors and scrub typhus

Weekly scrub typhus incidence and meteorological variables showed seasonal variations from January 1, 2001 to December 31, 2019 (Figure 1). For the first few years from the beginning of the online system activation, low numbers of scrub typhus cases have been reported in literature. Nonetheless, the annual peak of scrub typhus incidence remained constant, and the incidence tended to increase up to 2019, with a few exceptions in some years (2016–2017). Between epidemiological weeks 43 (mid-October) and 47 (mid-November), weekly cases reached the 50th percentile or above in most regions, independent of latitude. Meteorological variables such as the minimum, mean, and maximum temperatures, mRH, and weekly cumulative precipitation showed seasonal trends, with maximum observation values in the summer (June, July, and August) and minimum values observed in the winter (December, January, and February). The magnitudes of the seasonal variation for weekly cumulative sunshine duration, DTR, and MWS were comparatively small (Table S1, Figure 1, Figure S1).

Descriptive statistics

Weekly cumulative scrub typhus cases peaked in the epidemiologic weeks between 39 and 52/53 week (mean: 27 cases, IQR: [2, 32]). In addition, the weekly case reporting was less than one, if averaged. The medians of the minimum, mean, and maximum temperatures (°C) across the country were 20.8°C (IQR: [18.5, 23.0]), 24.0°C (IQR: [22.2, 26.1]), and 28.5°C (IQR: [26.7, 30.5]) in the summer (from June to October, epidemiologic weeks 22 to 35). In winter, these values were -2.9°C (IQR: [-5.6, -0.1]), 1.5°C (IQR: [-1.2, 4.2]), and 6.4°C (IQR: [3.7, 9.3]), respectively (from December to February, epidemiologic weeks 1–8 and 48–52). The DTR during winter was slightly wider than that during summer (median: 9.3°C vs. 7.8°C). The mean RH during the summer reached 98.6% (4th quartile), whereas it decreased to 22.8% (1st quartile) during winter (Table 1). Likewise, the weekly cumulative precipitation, especially during the monsoon season in the summer, rose to 695.5 mm (4th quartile). The weekly cumulative precipitation in winter was approximately one-seventh of the value in summer (mean: 7.0 mm vs. 47.7 mm). Weekly cumulative sunshine duration and MWS had the least variance compared with other meteorological variables (mean ± SD in summer vs. winter: 35.5 ± 19.1 vs. 35.4 ± 14.6 for weekly cumulative sunshine duration; 2.1 ± 0.7 vs. 2.3 ± 1.0 for MWS).

Table 1

Distribution of mean ambient temperature (°C) and reported incident numbers of scrub typhus by location (i.e., Seoul, Busan, or Jeju-do).

Province

Mean temperature (°C)

Mean relative humidity (%)

Scrub typhus cases

(cases per week)

Min

Mean ± SD

Max

Min

Mean ± SD

Max

Min

Mean ± SD

Max

Seoul

-10.1

13.0 ± 10.3

32.3

25.6

60.2 ± 11.0

95.6

0

4.0 ± 10.4

79

Busan

-1.8

15.1 ± 7.9

31.7

22.8

62.3 ± 15.1

98.6

0

9.7 ± 27.5

200

Daegu

-4.8

14.6 ± 9.3

31.7

26.1

58.3 ± 13.1

97.1

0

3.8 ± 12.4

161

Incheon

-9.7

12.5 ± 9.6

29.5

36.5

69.5 ± 11.7

98.5

0

1.3 ± 3.5

35

Gwangju

-4.1

14.4 ± 9.2

30.8

33.1

67.8 ± 10.5

98.0

0

4.2 ± 12.0

117

Daejeon

-7.2

13.3 ± 9.8

31.8

34.1

67.4 ± 11.0

96.4

0

4.7 ± 14.2

118

Ulsan

-3.4

14.6 ± 8.4

31.1

23.6

63.5 ± 14.9

97.9

0

5.5 ± 17.2

164

Kyeonggi-do (KK)

-10.4

11.8 ± 10.4

30.7

30.9

67.5 ± 10.6

96.6

0

10.6 ± 29.0

223

Kangwon-do (KW)

-9.4

11.3 ± 9.7

28.9

30.0

65.9 ± 11.8

94.4

0

3.2 ± 12.3

168

Chungcheongbuk-do (CB)

-8.8

11.9 ± 10.1

29.9

29.7

65.9 ± 10.4

90.8

0

4.5 ± 13.6

104

Chungcheongnam-do (CN)

-7.4

12.4 ± 9.7

30.0

38.9

71.3 ± 8.6

93.8

0

14.2 ± 39.6

281

Jeollabuk-do (JB)

-6.6

12.8 ± 9.6

29.6

49.1

71.6 ± 8.1

93.3

0

14.9 ± 40.1

287

Jeollanam-do (JN)

-3.5

13.7 ± 8.6

29.9

44.3

72.6 ± 9.3

93.4

0

16.4 ± 42.0

396

Gyeongsangbuk-do (GB)

-6.0

12.6 ± 9.2

29.6

30.5

65.6 ± 12.1

90.6

0

8.6 ± 24.2

177

Gyeongsangnam-do (GN)

-3.3

14.0 ± 8.8

31.0

35.8

65.7 ± 11.8

96.6

0

19.3 ± 53.6

562

Jeju-do (JJ)

1.8

16.2 ± 7.3

30.0

44.0

71.6 ± 10.2

96.1

0

1.1 ± 3.0

25

 

2. Association between climate factors and scrub typhus incidence (threshold analysis)

Weekly scrub typhus incidence was correlated with weekly mean temperature, mRH, and weekly cumulative precipitation with lags of 10–15 weeks, and with Pearson’s correlation coefficients ranging from 0.22 to 0.44, 0.22, 0.28, and 0.08 to 0.28, respectively. DTR, weekly cumulative sunshine duration, and MWS were correlated with weekly scrub typhus incidence, with a lag of approximately 20 weeks. The linear association between scrub typhus incidence for each province and meteorological variables above or below the threshold point for all provinces from GAM are shown in Figure S2.

The GAM for the association between mean AT or mean RH and weekly scrub typhus incidence per 100,000 persons in Busan is shown in Figure 2. Association plots between climate factors (i.e., mean AT and mean RH) and weekly scrub typhus incidence (per 100,000 persons) with a simple lag of 15 weeks showed a nonlinear relationship. In fact, mean AT and mean RH had a delayed impact (15-week lag) on the increase in scrub typhus incidence above a threshold point (Figure S6, S7). From Seoul to Jeju-do, the thresholds of mean ATs ranged from 14.9°C (Jeju-do) to 17.4°C (Incheon). The threshold points for mean RH for each province ranged from 57.3% (GB) to 67.7% (JN); however, in 9 out of 16 provinces, we observed a linear association between the mean RH and scrub typhus incidence (Table 2).

Table 2

Thresholds of mean ambient temperature (°C) and mean relative humidity of each location (i.e., Seoul, Busan, or Jeju-do) from which linear association with scrub typhus incidence start

 

Mean temperature (°C)

Mean relative humidity (%)

Threshold

Beta coefficient

(95% CI)

p-value

Threshold

Beta coefficient

(95% CI)

p-value

Seoul

17.0

0.117 (0.085, 0.150)

< 0.001

NA

0.056 (0.043, 0.069)

< 0.001

Busan

16.7

0.204 (0.154, 0.255)

< 0.001

NA

0.058 (0.045, 0.070)

< 0.001

Daegu

16.1

0.309 (0.139, 0.480)

< 0.001

NA

0.064 (0.046, 0.082)

< 0.001

Incheon

17.4

0.103 (0.061, 0.144)

< 0.001

69.6

0.066 (0.048, 0.084)

< 0.001

Gwangju

16.2

0.225 (0.124, 0.326)

< 0.001

NA

0.054 (0.032, 0.077)

0.002

Daejeon

15.2

0.254 (0.136, 0.372)

< 0.001

NA

0.050 (0.028, 0.072)

< 0.001

Ulsan

14.8

0.283 (0.213, 0.353

< 0.001

NA

0.041 (0.024, 0.058)

0.539

GG

15.7

0.151 (0.117, 0.185)

< 0.001

55.9

0.042 (0.027, 0.057)

< 0.001

GW

16.1

0.080 (0.054, 0.107)

< 0.001

56.8

0.015 (-0.001, 0.030)

0.070

CB

15.0

0.279 (0.086, 0.472)

< 0.001

NA

0.035 (-0.007, 0.077)

0.102

CN

15.5

0.242 (0.171, 0.313)

< 0.001

NA

0.032 (0.007, 0.058)

0.013

JB

16.1

0.240 (0.172, 0.308)

< 0.001

NA

0.016 (-0.010, 0.043)

< 0.225

JN

16.3

0.105 (0.081, 0.129)

< 0.001

67.7

0.045 (0.029, 0.061)

< 0.001

GB

15.0

0.229 (0.167, 0.292)

< 0.001

57.3

0.043 (0.024, 0.062)

0.002

GN

16.2

0.172 (0.136, 0.208)

< 0.001

57.5

0.042 (0.026, 0.058)

< 0.001

JJ

14.9

0.196 (0.136, 0.256)

< 0.001

57.8

0.001 (-0.021, 0.023)

0.918

 

Above a threshold point, weekly scrub typhus incidence (per 100,000 persons) constantly increased as the mean AT or mean RH increased. This tendency was consistently observed in all the provinces (Figure S6 and S7). GLM adopted for the linear section above a threshold mean AT resulted in a statistically significant association (beta, the slope of a linear section) for all provinces (Table 2). The beta coefficient for a 1°C increase in mean AT ranged from 0.103 (Incheon) to 0.309 (Daegu). Similarly, the beta coefficient of a 1% increase in mean RH ranged from 0.001 (Jeju, 95% CI: -0.021, 0.023) to 0.066 (Incheon, 95% CI: 0.048, 0.084) (Table 2). When pooled by meta-analysis, percent change in weekly scrub typhus incidence per 100,000 persons with a simple lag of 15 weeks per 0.1°C increase in mean AT was 1.85% (95% CI: 1.50, 2.20) across the country. Correspondingly, the percent change in weekly scrub typhus incidence per 100,000 persons with a simple lag of 15 weeks per 1% increase in mean RH in the country was 4.33% (95% CI: 3.41, 5.26) by pooling provincial data in the same manner (Figure 3A, 3C).

For sensitivity analysis, to smoothen the association plots between mean AT or mean RH and scrub typhus incidence, we adjusted for seasons by applying harmonics. To view the robustness of the applied GAM formula, we drew different GAM plots with or without different harmonics (2, 4, 6, 8, and 10). Adopting harmonics showed a smoother association. When the harmonics were given a value of six, the smoothed GAM plots remained unchanged, which supported the robustness of the model (Figure S8). After the simple lag was confined to 15 weeks, GAM and GLM were repeatedly performed by replacing the df values for a natural cubic spline with different values (Figure S9). Different df values in the natural cubic spline estimated a nearly unchanged beta coefficient. To estimate an increase in scrub typhus incidence following a unit increase in mean AT or mean RH, the GLM adopted covariates from several meteorological variables. Model I included two meteorological variables (mean AT or mean RH reciprocally and weekly cumulative precipitation). Model II was additionally adjusted for weekly sunshine duration and mean wind speed. The two models resulted in similar outcomes in the estimation of scrub typhus incidence (Figure 3B, 3D).

3. Long-term association between AT change and scrub typhus incidence

Compared with the previous 30 years, the annual mean AT difference tended to increase between 2001 and 2019. The number of heatwave days per year increased in a similar manner. When visualized in five-year terms, the long-term mean AT difference showed an increasing trend. In Figure 5, we plotted the long-term trends of meteorological variables (e.g., annual mean AT, annual mean RH), differences in annual mean AT or yearly total number of heatwave days, and cumulative scrub typhus incidence per year. Through LOESS smoothing of the trends, the plots showed that weekly sunshine duration and the number of heatwave days had gradually increased. Cumulative scrub typhus incidence, annual mean AT, and the different values of mean AT or the number of heatwave days slightly increased (Figure 4 and 5). The fitted linear regression model estimated an annual mean AT elevation of 1.17°C between 2001 and 2019 (95% CI: 1.02°C, 1.33°C). Likewise, the number of heat wave days per year compared with the past 30 years increased by approximately 2.82 days (95% CI: 2.04, 3.60). Simple linear regression analysis revealed that “year” had a stronger association than the other two variables (difference in annual mean AT or the number of heatwave days per year) with total outbreak cases. In a simple linear regression, a 1°C elevation of annual mean AT compared with the previous 30 years was associated with an additional 453.7 total outbreak cases (95% CI: 400.6, 506.8). In addition, a one-day increase in the number of heatwave days per year over the first 19 years (2001–2019) was associated with 429.5 more total outbreak cases (95% CI: 378.7, 480.4). In GLME, an annual mean AT elevation of 1°C explained a 36.9% increase in annual cumulative scrub typhus incidence per 100,000 persons in a country. In addition, one more heatwave day in a year was associated with a 2.5% increase in annual cumulative scrub typhus incidence per 100,000 persons at the country level.

Discussion

Weekly scrub typhus incidence was associated with mean AT and mean RH with a simple lag of 15 weeks in the ROK. This association was non-linear, namely U-shape or V-shape; thus, we approximated a threshold point above or below which scrub typhus incidence was linearly associated with mean AT or mean RH. The difference in the annual mean AT and the annual number of heatwave days in the period 1971–2000 increased over time from 2001 to 2019, and the annual cumulative scrub typhus incidence rate (per 100,000) was positively associated with the difference in annual mean AT or annual number of heatwave days in the past 30 years (1971–2000).

The Infectious Disease Control and Prevention Act in Korea has classified scrub typhus as a class 2 infectious disease that requires mandatory reporting since 1994. All confirmed or clinically suspected scrub typhus cases are reported to local health bureaus and then to the Korea Centre for Disease Control and Prevention 23. This means that the total number of reported cases to the health agency matched our target population. The ROK allows a single national health insurance service to be operated by a public agency under the Ministry of Health and Welfare (MOHW). Health claims data might provide corresponding information on scrub typhus incidence, considering that the Korea Standard Code of Disease 10 (KCD 10) was used when the claim data were initially collected. Compared with the national mandatory reporting system of infectious diseases, health claims data relies on health professionals’ inputs, mostly based on clinical diagnosis. Therefore, serologically unidentified clinical cases could be included in the health claims data. In comparison, surveillance data comprised of clinically confirmed cases, which did not allow for the overestimation of the real numbers. Hence, for clarity and completeness of the data, we used surveillance data for scrub typhus from the National Infectious Disease Portal.

The Responsible Division for Zoonotic Infectious Diseases, previously denoted as the Division of Vectors and Parasitic Diseases of the Korea National Institute of Health (KNIH), has annually conducted vector surveillance for scrub typhus in spring and autumn. This was conducted at evenly distributed locations by reporting the ratio of infected chiggers 24–26. At least 32 species of trombiculid mites transfer Orientia tsutsugamushi to hosts, of which nine representative Leptotrombidium species distributed across the Asia-Pacific region, are mainly reported in the ROK and Japan 27,28. The detection of infected chiggers and case reports in Seoul indicate that scrub typhus has already been urbanized across the country. Seven Trombiculidae species are domestically well-known vectors, with the dominant species being L. pallidum in the northeastern provinces and L. scutellare in the southwestern provinces. The reason for why the nonlinear curves obtained from GAM for mean AT or mean RH differed for each province (Figure S5 and S6) might be due to the proportion of dominant vector species and the difference in dominant vectors’ favorable climate conditions.

Consistent with previous studies on the incidence of mosquito-borne diseases 29–31, scrub typhus incidence has a nonlinear association with mean AT or mean RH in the ROK 32. For scrub typhus, ranges of mean AT or mean RH, in which weekly scrub typhus incidence increased accordingly, were above specific threshold points 30. Modelling approaches such as transmission rate calculations using temperature-dependent entomological parameters and mosquito biting rates proved that there was an upper limit of optimum temperatures 33. Furthermore, other infectious disease incidences (e.g., herpes zoster) or elderly mortality had a similar correlation pattern. The incidence of infectious disease or chronic disease incidence increased linearly 34,35.

AT affects the development of Leptotrombidium deliense (a vector species of scrub typhus) during its life cycle. Previous studies have found critical temperatures for the development and reproduction of L. delicense, as well as optimum AT ranges for larval activity and disease transmission. Under suitable AT conditions, egg hatching and chigger development speed up as AT increases 12,36−38.

Warm humid summers enhance the transmission of diseases such as Lyme borreliosis or tick-borne encephalitis via an increased number of vectors and their exposure to humans 39–41. High summer temperatures are likely to restrict larval activity in mid-summer and promote nymphal activity in late autumn or early spring of the following year 42,43. In addition, low humidity depresses nymphal and adult tick activity 42,44,45. Our study, which is associated with summer meteorological factors (i.e., mean AT and mean RH) with the later incidence of scrub typhus in autumn, provides supporting evidence for previous study findings or suggestions accordingly.

One of our study’s main outcomes was the association between meteorological factors (mean AT and mean RH) and scrub typhus incidence. Previous studies assumed a linear relationship between scrub typhus incidence and meteorological factors to estimate the increase in scrub typhus incidence per unit change in each meteorological variable 46. To the best of our knowledge, this study is the first to approximate the threshold points for U- or V-shaped nonlinear associations. Furthermore, the annual cumulative incidence of scrub typhus showed a linear relationship with the difference between the annual mean AT or the annual number of heatwave days and those of the period 1971–2000.

Our study results are in line with previous experiences regarding the impact of climate change, especially for global warming 2,47−51. The countrywide annual mean AT gradually increased over the past 19 years, from 2001 to 2019. Extreme hot weather, depicted by the number of heatwave days, increased in parallel. We intended to prove that warmer temperatures or extreme hot summer weather influenced the outbreak curve in the following autumn season, which was quantified by the cumulative outbreak cases. Milder and shorter winters are known to induce longer tick-activity seasons 52. In contrast, an extremely cold winter in Germany in 2012 was followed by a significantly low rate of host-seeking activity in Ixodes ricinus nymphs 53. Retrospective data from Europe proves that the development and activity of ticks in summer are influenced by high temperatures, and global warming is presumably related 40. Similarly, extreme weather conditions (i.e., extreme cold or hot temperature events) are related to the increase in non-communicable diseases 54 and to an increase in the incidence of chronic infectious diseases 55.

To provide an early alarm to the public as well as the public health sector, predictive models and projectile tools such as the Early Warning and Response System 56–60 should be created in the near future 56–60. In addition, we may need to thoroughly acknowledge the interactions between environments, pathogens, vectors, and hosts in order to prepare pre-emptively. A strength of our study is that it provides useful data for future research.

The goal of this study was to relate components of an ecosystem, which consists of environmental factors, pathogens, vectors, and hosts. To do so, we repeatedly ran GAMs by applying several different models, including simple and moving average lags in weeks. We might be able to better demonstrate the real situation if the cumulative effects of climate factors preceding outbreaks were combined accordingly. However, the increased infection risk of scrub typhus due to changes in climate factors can only be quantified by drawing on the association between them.

We approximated the increased incidence of scrub typhus in the ROK in recent years (2001–2019) when mean AT, mean RH, annual mean AT, and annual number of heatwave days increased. Compared with the previous 30 years, the annual mean AT increased by just 1°C from 2001 to 2019. The number of heatwave days also increased by almost three days during the same period. Global warming has influenced the population dynamics of arthropod vectors and the transmission dynamics of the diseases they carry 61. Our results support this idea by presenting an increased incidence rate. Long-term trends consistently showed a positive association between mean AT or heatwave days and scrub typhus incidence. Future climate change is likely to extend the transmission period of vector-borne diseases and expand the geographical range of vector habitats, which may cause new emergence or re-emergence of the diseases 62. Further studies to estimate the impact of climate variability on the incidence of vector-borne diseases, such as scrub typhus, are needed, as some of these diseases might threaten vulnerable populations at risk or continue to influence human lives 63.

Methods

Study population

The number of weekly scrub typhus cases was obtained from the National Infectious Disease Portal, provided by the Korea Disease Control and Prevention Agency (KDCA), and the period of interest ranged from 2001 to 2019 since data acquisition from this platform was only available from 2001. Surveillance system of scrub typhus was based on mandatory reports of all suspected and confirmed cases district-wise. We obtained weekly reported suspected and confirmed cases in seven metropolitan cities and nine provinces. Suspected cases were defined by clinical manifestations of scrub typhus, such as acute febrile illness, lymphadenopathy and/or skin eschar, and epidemiological association. Confirmed cases were defined based on clinical symptoms and laboratory confirmation. The year 2020 was excluded to avoid the influence of the coronavirus disease in 2019. Our study only used publicly available data; therefore, the requirement for informed consent was waived off (IRB No. 2003-131-1110).

Meteorological data acquisition

We obtained meteorological data from the Weather Data Open Portal provided by the Korea Meteorological Administration. Meteorological data such as AT, RH, precipitation, wind speed, and sunshine duration on a minute-to-minute basis were collected from a total of 102 monitoring stations across the country. We computed the daily average values of these meteorological variables by city and province.

We obtained information on the daily mean, maximum, and minimum AT, mean RH, mean wind speed, sunshine duration, and precipitation. The diurnal temperature range (DTR) was defined as the difference between the daily maximum and minimum temperatures. We then computed the weekly mean AT (mean, minimum, and maximum), mean RH, DTR, mean wind speed, and weekly cumulative sunshine duration and weekly cumulative precipitation.

We also computed the difference between the annual mean AT each year (2001–2019) and the mean AT for the period 1971–2000. Likewise, we computed the difference between the number of heatwave days each year (2001–2019) and the mean number of heatwave days during the years 1971–2000. The number of heatwave days per year was defined by the number of days in which the daily maximum AT was ≥33°C.

Statistical analysis

1. Descriptive analysis

We adopted locally estimated scatterplot smoothing (LOESS) curve fitting to determine whether scrub typhus incidence showed seasonality 64,65. Time series plots of climate factors and weekly incidence of scrub typhus from January 2001 to December 2019 are shown in Figure 1. We then conducted a bivariate analysis between the meteorological variables and weekly scrub typhus incidence. As there were lags of up to 7 months between the annual peaks of the mean AT and the peaks of scrub typhus incidence, we decided to use lags from 0 to 28 weeks (Figure S1).

The correlation between each meteorological variable (mean, maximum, and minimum AT, mean RH, mean wind speed, weekly cumulative sunshine duration, weekly cumulative precipitation, and DTR) and weekly scrub typhus incidence was explored by applying a generalized additive model (GAM). Minimum, mean, and maximum AT, mean RH, and weekly cumulative precipitation were correlated with weekly scrub typhus incidence. Since the minimum, mean, and maximum AT were correlated, we only selected mean RH and weekly cumulative precipitation as covariates for GAM.

We constructed GAM for the association between AT and weekly scrub typhus incidence by adjusting for mean RH and precipitation, using a quasi-Poisson function, and simple lags and moving average lags from 0 to 28 weeks. Meteorological variables and weekly scrub typhus incidence showed seasonality between 2001 and 2019 (Figure 1 and S1). According to previous literature, Fourier terms up to six harmonics per year were included to reflect seasonality in the GAM 66,67. The GAM formula used is as follows:

Ln [E(Y)] = β0 + β1Xt + Σ Si (Xj) + sin (2πtj /N) + cos (2πtj/N) + year,

where E(Y) denotes expected scrub typhus incidence; β1 is the coefficient (slope) of a meteorological variable for week t; Si(Xj) is a smoothing function for an individual covariate i; tj indicates the week of the year j (j = 1, 2, ... ,52); N is a harmonic for the year (N = 1, 2, ..., 6); and year depicts the year when the scrub typhus cases were reported on an ordinal scale.

Associations between AT and weekly scrub typhus incidence were explored for the seven cities and nine provinces with different lags, from 10 to 15 weeks (Figure S5). Association plots obtained by applying GAM showed that there was a nonlinear relationship between a meteorological factor 10–15 weeks ahead and weekly scrub typhus incidence (Figure S5). Compared with other lags, such as simple lags of 10–14 weeks or moving average lags of 10–15 weeks, GAM plots with a simple lag of 15 weeks consistently showed V or U shape in most of the cities or provinces (Figure S5).

2. The association between weekly AT and scrub typhus incidence (threshold analysis)

To estimate the association between the increase of scrub typhus incidence and AT or RH, which showed “V” or “U” shape, we attempted to find a threshold where the slope of the graph changes from negative to positive (Figure S6 and S7). The R HEAT package uses a generalized linear model to approximate a TP if the beta coefficient (slope) is converted from negative (positive) to positive (negative) 68. In particular, the threshold function in the HEAT package uses piecewise linear regression to approximate a threshold. For instance, β1 was calculated below a threshold, where weekly scrub typhus incidence decreased as mean AT increased and β2 was calculated above a threshold, where scrub typhus incidence increased as mean AT increased. After visual approximation of a threshold for mean AT and mean RH from GAM, we applied the threshold () function to quantify the percent changes in scrub typhus incidence following a unit change of each meteorological variable above or below the threshold. The GLM formula used is as follows:

Ln [E(Y)] = β0 + β1Xt + Σ Si (Xj) + sin (2πtj /N) + cos (2πtj/N) + year, if Xt <TP.

Ln [E(Y)] = β0 + β2Xt + Σ Si (Xj) + sin (2πtj /N) + cos (2πtj/N) + year, if Xt ≥TP,

where E(Y) denotes the expected scrub typhus incidence; β1 is the coefficient (slope) of a meteorological variable for week t below the threshold point (TP) of Xt and β2 is the coefficient (slope) of a meteorological variable for week t above the threshold point (TP) of Xt; Si(Xj) is a smoothing function for an individual covariate i; tj indicates the week of the year j (j = 1, 2, ... ,52); N is a harmonic for the year (N = 1, 2, ..., 6); and year depicts the year when scrub typhus cases were reported.

Since data distributions for meteorological variables and scrub typhus incidence were skewed, we adopted quasi-Poisson regression in the GAM and GLM models. The smoothing parameter was determined by natural cubic spline by selecting an optimal degree of freedom 60, where the unbiased risk estimator (UBRE) 69,70 in the GAM and Akaike information criterion (AIC) 71,72 values in the GLM for each province approached their minimum. Considering all provinces’ UBRE values following increasing degrees of freedom (df), our GAM model used 5 as its df.

For the sensitivity analysis, we used different dfs for each city and province. We also compared GAM with and without harmonics. Harmonics of 2, 4, 6, 8, and 10 were used for further comparison. Furthermore, we compared the GLM by additionally adjusting for mean wind speed, weekly precipitation, and weekly sunshine duration.

3. Long-term association between meteorological factors and scrub typhus incidence

We investigated the long-term association between meteorological factors and scrub typhus incidence using annual cumulative scrub typhus cases. The long-term mean AT difference was defined as the difference between the mean annual AT each year (2001–2019) and the mean AT of the previous 30 years (1971–2000). The difference in the number of heat wave days was defined as the difference between the number of heat wave days each year (2001–2019) and the mean number of heat wave days during the years 1971–2000.

We plotted the long-term trends of meteorological variables from 2001 to 2019 using the LOESS regression. The difference in annual mean AT or annual number of heatwave days from the mean AT or mean number of heatwave days between 1971 and 2000 was grouped by five-year terms to detect long-term trends (Figure 4A, 4B). A simple linear regression model overlaid the point plots to determine if a linear correlation existed between the variables. The fitted model slope, p-value, and coefficient of determination (R2) were calculated accordingly. Although the time variable year best explained increases in cumulative scrub typhus cases in autumn outbreaks in the model selection process, we intended to correlate cumulative scrub typhus cases in autumn outbreaks with differences in the annual mean AT or the number of heatwave days from those of the years from 1971 to 2000. Therefore, point plots, linear regression lines, 95% confidence intervals, slopes, p-values, and R2 values were obtained in the same manner (Figure 4).

To analyze the association between the difference in annual mean AT of each year and mean AT of the past 30 years (1971–2000), we used a generalized linear mixed-effects model (GLME). Covariates were selected from the directed acyclic graph (DAG) (www.dagitty.net) (Figure S10). We considered changes in population and land use or greenness over the last 19 years; however, we did not include them as covariates since population changes were not directly associated with AT, and land use was a causal mediator rather than a confounder (Figure S10). Instead, we used the incidence rate (annual incidence of scrub typhus per 100,000 population) to consider population factors in the regression model. We put the province as a random intercept in the glmer() to diminish spatial interference in the estimation. The analysis was repeated after replacing the difference in annual mean AT by the difference in the number of heatwave days.

All statistical analyses were conducted using the statistical software packages SAS (version 9.4) and RStudio software (ver. 1.2.5042).

Declarations

Ethical statement

All the methods used in this study followed relevant guidelines and regulations. Our study only used publicly available data; therefore, the requirement for informed consent was waived off.

Data availability

The dataset generated and/or analyzed during the current study are available from the corresponding author on reasonable request.

Author contributions

D.S. and Y.C. contributed equally as first authors. D.S. and Y.C. collected and analyzed data, and wrote the first draft. All authors contributed to the conceptualization of study design, writing and approval of the final version of manuscript.

Additional Information

The authors declare no competing interests.

Funding

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.

References

  1. Field, C. B., V.R. Barros, K.J. Mach, M.D. Mastrandrea, M. van Aalst, W.N. Adger, D.J. Arent, J. Barnett, R. Betts, T.E. Bilir, J. Birkmann, J. Carmin, D.D. Chadee, A.J. Challinor, M. Chatterjee, W. Cramer, D.J. Davidson, Y.O. Estrada, J.-P. Gattuso, Y. Hijioka, O. Hoegh-Guldberg, H.Q. Huang, G.E. Insarov, R.N. Jones, R.S. Kovats, P. Romero-Lankao, J.N. Larsen, I.J. Losada, J.A. Marengo, R.F. McLean, L.O. Mearns, R. Mechler, J.F. Morton, I. Niang, T. Oki, J.M. Olwoch, M. Opondo, E.S. Poloczanska, H.-O. Pörtner, M.H. Redsteer, A. Reisinger, A. Revi, D.N. Schmidt, M.R. Shaw, W. Solecki, D.A. Stone, J.M.R. Stone, K.M. Strzepek, A.G. Suarez, P. Tschakert, R. Valentini, S. Vicuña, A. Villamizar, K.E. Vincent, R. Warren, L.L. White, T.J. Wilbanks, P.P. Wong, and G.W. Yohe. Technical summary. In: Climate Change 2014: Impacts, Adaptation, and Vulnerability. Part A: Global and Sectoral Aspects. Contribution of Working Group II to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change. (Cambridge University Press, 2014).
  2. Khasnis, A. A. & Nettleman, M. D. Global warming and infectious disease. Arch Med Res 36, 689-696, doi:10.1016/j.arcmed.2005.03.041 (2005).
  3. Antonella Rossati. Global Warming and Its Health Impact. Occup Environ med 8, 7-20, doi:10.15171/ijoem.2017.963 (2017).
  4. Farhang-Azad, A. J., W. B.; Liat, L. B.; Lund, M.; Mateva, M.; Msangi, A.; Myllymaki, A.; Shuyler, H. R.; Arata, A. A.; Drummond, D. C.; Meylan, A. Ecology and control of rodents of public health importance. Report of a WHO Scientific Group. Report No. 0512-3054 (Print) 0512-3054, 1-42 (1974).
  5. Waits, A., Emelyanova, A., Oksanen, A., Abass, K. & Rautio, A. Human infectious diseases and the changing climate in the Arctic. Environ Int 121, 703-713, doi:10.1016/j.envint.2018.09.042 (2018).
  6. Parola, P. et al. Warmer weather linked to tick attack and emergence of severe rickettsioses. PLoS neglected tropical diseases 2, e338, doi:10.1371/journal.pntd.0000338 (2008).
  7. Fouque, F. & Reeder, J. C. Impact of past and on-going changes on climate and weather on vector-borne diseases transmission: a look at the evidence. Infect Dis Poverty 8, 51, doi:10.1186/s40249-019-0565-1 (2019).
  8. Romi, R. Arthropod borne diseases in Italy: from a neglected matter to an emerging health problem. Ann Ist Super Sanità 46, 436-443, doi:10.4415/ANN_10_04_12 (2010).
  9. Adem, P. V. Emerging and re-emerging rickettsial infections. Semin Diagn Pathol 36, 146-151, doi:10.1053/j.semdp.2019.04.005 (2019).
  10. Jiang, J. & Richards, A. L. Scrub Typhus: No Longer Restricted to the Tsutsugamushi Triangle. Trop Med Infect Dis 3, doi:10.3390/tropicalmed3010011 (2018).
  11. Beugnet, F. & Marie, J. L. Emerging arthropod-borne diseases of companion animals in Europe. Vet Parasitol 163, 298-305, doi:10.1016/j.vetpar.2009.03.028 (2009).
  12. Lv, Y., Guo, X. G. & Jin, D. C. Research Progress on Leptotrombidium deliense. Korean J Parasitol 56, 313-324, doi:10.3347/kjp.2018.56.4.313 (2018).
  13. Zheng, C., Jiang, D., Ding, F., Fu, J. & Hao, M. Spatiotemporal Patterns and Risk Factors for Scrub Typhus From 2007 to 2017 in Southern China. Clin Infect Dis 69, 1205-1211, doi:10.1093/cid/ciy1050 (2019).
  14. Seto, J. et al. Meteorological factors affecting scrub typhus occurrence: a retrospective study of Yamagata Prefecture, Japan, 1984-2014. Epidemiology and infection 145, 462-470, doi:10.1017/s0950268816002430 (2017).
  15. Grokhovskaia, I. M. F. A. U. N. S. K. & Nguen Suan, K. Features of ecology and distribution of Trombiculidae mites of Vietnam. Meditsinskaia parazitologiia i parazitarnye bolezni 38, 72-78 (1969).
  16. Choi, Y. J. et al. Geographical distribution of Orientia tsutsugamushi strains in chiggers from three provinces in Korea. Microbiol Immunol 62, 547-553, doi:10.1111/1348-0421.12639 (2018).
  17. Mamoru Takahashi, H. M., Hirofumi Matsuzawa, Koji Morita, Osamu Tsuju, Eitaro Hori, Akiyoshi Kawamura, Jr., Hiroshi Tanaka. Seasonal development of Leptotrombidium scutellare (Acari: Trombiculidae) observed by experimental rearing in field conditions. Jpn J Saint Zool 45, 113-120 (1994).
  18. Wang, C.-Y., Chang, W.-H. & Su, Y.-J. Acute encephalitis syndrome caused by Orientia tsutsugamushi. Health Technology 3, 5-5, doi:10.21037/ht.2019.05.01 (2019).
  19. Sharma, S. R., Masaraf, H., Lynrah, K. G. & Lyngdoh, M. Tsutsugamushi Disease (Scrub Typhus) Meningoencephalitis in North Eastern India: A Prospective Study. Ann Med Health Sci Res 5, 163-167, doi:10.4103/2141-9248.157486 (2015).
  20. Kim YO, J. H., Cho SG, Yoon SA, Son HS, Oh SH, Chae HS, Kim KH, Chae JS, Lee CD, Lee BS. The role of hypoalbuminemia as a marker of the severity of disease in patients with tsutsugamushi disease. Korean J Med 59, 516-521 (2000).
  21. Jain, P. et al. Emergence of Orientia tsutsugamushi as an important cause of Acute Encephalitis Syndrome in India. PLoS neglected tropical diseases 12, e0006346, doi:10.1371/journal.pntd.0006346 (2018).
  22. Patz, J. A., Epstein, P. R., Burke, T. A. & Balbus, J. M. Global Climate Change and Emerging Infectious Diseases. JAMA 275, 217-223, doi:10.1001/jama.1996.03530270057032 (1996).
  23. Jeong HJ, L. J., Gill BC , Lee YS, , Kim DM, Kim YS. Molecular epidemiology of Orientia tsutsugamushi from scrub typhus patients in South Korea, 2014-2015. PHWR 9, 882-888 (2016).
  24. Park WI, N. J., Cho SH. Surveillance of Chigger Mites as the vector of scrub typhus, 2016. PHWR 10, 1058-1061 (2017).
  25. Noh JY, S. E., Joo YR. Nationwide Surveillance of Chigger Mites, as the Vector of Scrub Typhus. PHWR 7, 1146-1148 (2014).
  26. Noh, J. Characteristics and geographical distribution of chigger mites as the vector of Orientia tsutsugamushi in Korea. PHWR 5, 234-238 (2012).
  27. Kim, G. et al. Diversification of Orientia tsutsugamushi genotypes by intragenic recombination and their potential expansion in endemic areas. PLoS neglected tropical diseases 11, e0005408, doi:10.1371/journal.pntd.0005408 (2017).
  28. Elliott, I. et al. Scrub typhus ecology: a systematic review of Orientia in vectors and hosts. Parasit Vectors 12, 513, doi:10.1186/s13071-019-3751-x (2019).
  29. YL Hii, H. Z., N Ng, LC Ng. Forecast of Dengue Incidence Using Temperature and Rainfall. PLoS neglected tropical diseases, doi:10.1371/journal.pntd.0001908 (2012).
  30. Wu, X. et al. Non-linear effects of mean temperature and relative humidity on dengue incidence in Guangzhou, China. Sci Total Environ 628-629, 766-771, doi:10.1016/j.scitotenv.2018.02.136 (2018).
  31. Teklehaimanot, H. D., Lipsitch, M., Teklehaimanot, A. & Schwartz, J. Weather-based prediction of Plasmodium falciparum malaria in epidemic-prone regions of Ethiopia I. Patterns of lagged weather effects reflect biological mechanisms. Malar J 3, 41, doi:10.1186/1475-2875-3-41 (2004).
  32. Bhopdhornangkul, B. et al. Non-linear effect of different humidity types on scrub typhus occurrence in endemic provinces, Thailand. Heliyon 7, e06095, doi:10.1016/j.heliyon.2021.e06095 (2021).
  33. Chen, S. C. & Hsieh, M. H. Modeling the transmission dynamics of dengue fever: implications of temperature effects. Sci Total Environ 431, 385-391, doi:10.1016/j.scitotenv.2012.05.012 (2012).
  34. Ha, J., Kim, H. & Hajat, S. Effect of previous-winter mortality on the association between summer temperature and mortality in South Korea. Environ Health Perspect 119, 542-546, doi:10.1289/ehp.1002080 (2011).
  35. Choi, Y.-J., Lim, Y.-H., Lee, K.-S. & Hong, Y.-C. Elevation of ambient temperature is associated with an increased risk of herpes zoster: a time-series analysis. Scientific Reports 9, 1-7, doi:doi:10.1038/s41598-019-48673-5 (2019).
  36. Womersley, H. The scrub-typhus and scrub-itch mites (Trombiculidae, Acarina) of the Asiatic-Pacific region. Records of the South Australian Museum 10, 1-435 (1952).
  37. Traub, R. & Wisseman, C. L., Jr. Ecological considerations in scrub typhus. 3. Methods of area control. Bulletin of the World Health Organization 39, 231-237 (1968).
  38. Hsin-Tao C, H. P., Tzu-Chun L. The relationship between the temperature and the development of the eggs of T. akamushi var. Deliensis. Acta Zoologica Sinica, 02 (1963).
  39. Lindgren, E. J., T. G. T. Lyme Borreliosis in Europe: Influences of Climate and Climate Change, Epidemiology, Ecology and Adaptation Measures. (WHO Regional Office for Europe, Copenhagen, 2006).
  40. Gray, J. S., Dautel, H., Estrada-Pena, A., Kahl, O. & Lindgren, E. Effects of climate change on ticks and tick-borne diseases in europe. Interdiscip Perspect Infect Dis 2009, 593232, doi:10.1155/2009/593232 (2009).
  41. Bennet, L., Halling, A. & Berglund, J. Increased incidence of Lyme borreliosis in southern Sweden following mild winters and during warm, humid summers. Eur J Clin Microbiol Infect Dis 25, 426-432, doi:10.1007/s10096-006-0167-2 (2006).
  42. Schulze, T. L., Bowen, G. S., Lakat, M. F., Parkin, W. E. & Shisler, J. K. Seasonal Abundance and Hosts of Ixodes Dammini(Acari: Ixodidae) and Other Ixodid Ticks from an Endemic Lyme Disease Focus in New Jersey, USA. Journal of medical entomology 23, 105-109, doi:10.1093/jmedent/23.1.105 (1986).
  43. Schulze, T. L., Bowen, G. S., Lamat, M. F., Parkin, W. E. & Shisler, J. K. The Role of Adult Ixodes Dammini (Acari: Ixodidae) in the Transmission of Lyme Disease in New Jersey, USA. Journal of medical entomology 22, 88-93, doi:10.1093/jmedent/22.1.88 (1985).
  44. J. L. Perret, E. G., O. Rais, L. Gern. Influence of saturation deficit and temperature on Ixodes ricinus tick questing activity in a Lyme borreliosis-endemic area (Switzerland). Parasitol Res 86, 554-557 (2000).
  45. Gray, J. S. Ixodes ricinus seasonal activity: Implications of global warming indicated by revisiting tick and weather data. International Journal of Medical Microbiology 298, 19-24, doi:10.1016/j.ijmm.2007.09.005 (2008).
  46. Li, T., Yang, Z., Dong, Z. & Wang, M. Meteorological factors and risk of scrub typhus in Guangzhou, southern China, 2006–2012. BMC Infectious Diseases 14, 139, doi:10.1186/1471-2334-14-139 (2014).
  47. Singh, A. & Purohit, B. M. Public Health Impacts of Global Warming and Climate Change. Peace Review 26, 112-120, doi:10.1080/10402659.2014.876326 (2014).
  48. Singer, M. Respiratory health and ecosyndemics in a time of global warming. Health Sociology Review 22, 98-111, doi:10.5172/hesr.2013.22.1.98 (2014).
  49. Dushana Yoganathan, W. N. R. Medical Aspects of Global Warming. Am. J. Ind. Med. 40, 199-210 (2001).
  50. Cook, G. C. Effect of global warming on the distribution of parasitic and other infectious diseases: a review. Journal of the Royal Society of Medicine 85 (1992).
  51. Chan, A. W. et al. The effects of global warming on allergic diseases. Hong Kong Med J 24, 277-284, doi:10.12809/hkmj177046 (2018).
  52. Jan C Semenza, B. M. Climate change and infectious diseases in Europe. The Lancet Infectious Diseases 9, 365-375, doi:10.1016/S1473-3099(09)70104-5 (2009).
  53. Dautel, H., Kämmer, D. & Kahl, O. in Ecology and prevention of Lyme borreliosis Ecology and Control of Vector-borne Diseases 335-349 (2016).
  54. Liu, X. et al. Association between extreme temperature and acute myocardial infarction hospital admissions in Beijing, China: 2013-2016. PLoS One 13, e0204706, doi:10.1371/journal.pone.0204706 (2018).
  55. Onozuka, D. & Hagihara, A. The association of extreme temperatures and the incidence of tuberculosis in Japan. Int J Biometeorol 59, 1107-1114, doi:10.1007/s00484-014-0924-3 (2015).
  56. Reis, B. Y., Kohane, I. S. & Mandl, K. D. An epidemiological network model for disease outbreak detection. PLoS Med 4, e210, doi:10.1371/journal.pmed.0040210 (2007).
  57. Kwak, J. et al. Scrub Typhus Incidence Modeling with Meteorological Factors in South Korea. Int J Environ Res Public Health 12, 7254-7273, doi:10.3390/ijerph120707254 (2015).
  58. Keeling, M. J. & Rohani, P. Modeling infectious diseases in humans and animals. (Princeton university press, 2011).
  59. Hussain-Alkhateeb L, K. A., Olliaro P, Rocklöv J, Sewe MO, Tejeda G, Benitez D, Gill B, Hakim SL, Gomes Carvalho R, Bowman L. Early warning and response system (EWARS) for dengue outbreaks: Recent advancements towards widespread applications in critical settings. PLoS ONE, doi:10.1371/journal.pone.0196811 (2018).
  60. Halide, H. & Ridd, P. A predictive model for Dengue Hemorrhagic Fever epidemics. Int J Environ Health Res 18, 253-265, doi:10.1080/09603120801966043 (2008).
  61. Wilke, A. B. B., Beier, J. C. & Benelli, G. Complexity of the relationship between global warming and urbanization – an obscure future for predicting increases in vector-borne infectious diseases. Current Opinion in Insect Science 35, 1-9, doi:https://doi.org/10.1016/j.cois.2019.06.002 (2019).
  62. Caminade, C., McIntyre, K. M. & Jones, A. E. Impact of recent and future climate change on vector-borne diseases. Ann N Y Acad Sci 1436, 157-173, doi:10.1111/nyas.13950 (2019).
  63. Rogers, D. J. & Randolph, S. E. in Advances in Parasitology Vol. 62 (eds Simon I. Hay, Alastair Graham, & David J. Rogers) 345-381 (Academic Press, 2006).
  64. Park JH, K. S., Youn SK, Park KS, Gwack J. Epidemiology of Scrub Typhus and the Eschars Patterns in South Korea from 2008 to 2012. Jpn J Infect Dis 67, 458-463 (2014).
  65. Bhargava, A. et al. Scrub typhus in Uttarakhand & adjoining Uttar Pradesh: Seasonality, clinical presentations & predictors of mortality. Indian J Med Res 144, 901-909, doi:10.4103/ijmr.IJMR_1764_15 (2016).
  66. Kim, Y. M., Park, J. W. & Cheong, H. K. Estimated effect of climatic variables on the transmission of Plasmodium vivax malaria in the Republic of Korea. Environ Health Perspect 120, 1314-1319, doi:10.1289/ehp.1104577 (2012).
  67. Hashizume, M. et al. The effect of rainfall on the incidence of cholera in Bangladesh. Epidemiology 19, 103-110, doi:10.1097/EDE.0b013e31815c09ea (2008).
  68. Lim YH, K. H. Manual for HEAT package. (2012).
  69. Russell, R. C. Vectors vs. humans in Australia--who is on top down under? An update on vector-borne disease and research on vectors in Australia. Journal of vector ecology : journal of the Society for Vector Ecology 23, 1-46 (1998).
  70. Peter Craven, G. W. Smoothing Noisy Data with Spline Functions. Numer. Math. 31, 377-403 (1979).
  71. Chuang, Y. H. et al. Generalized linear mixed models in time series studies of air pollution. Atmospheric Pollution Research 2, 428-435, doi:https://doi.org/10.5094/APR.2011.049 (2011).
  72. Akaike, H. A new look at the statistical model identification. IEEE Transactions on Automatic Control 19, 716-723, doi:10.1109/TAC.1974.1100705 (1974).