The Effect of Climate Variables on Dengue Burden in Indonesia: A Case Study from Medan City ()
1. Introduction
Dengue is the most important vector-bone diseases which are transmitted to human by mosquitoes and it has grown severely in recent decades worldwide in tropical and sub-tropical region. World Health Organization (WHO) estimated that 2.5 billion people are lived at risk from dengue (WHO, 2011) and most of the severe dengue cases were occurred predominantly in infants and children (Hammond et al., 2005) . Urbanization, demographic change, sanitation and warming temperature are associated with the spread of Aedes aegypti and Aedes albopictus as the main vectors of dengue (Murray et al., 2013) . Furthermore, more than 125 countries have Ae. aegypti and Ae. albopictus vectors and it is estimated to be ubiquitous year around in the tropics with the highest risk region in America and Asia (Ebi & Nealon, 2016) . In addition, Kraemer et al. predicted that Ae. aegypti concentration mainly exists in the Northern Brazil and Southeast Asia countries (Kremer et al., 2015) .
Dengue is endemic with the transmission throughout the year in Indonesia, (Karyanti & Hadinegoro, 2009) . Since 1968-2009, dengue incidence rate increased with the peak cases every 10 years in 33 provinces (Karyanti & Hadinegoro, 2009) . In spite of great effort in Aedes control in Indonesia, the outbreak had occurred with higher frequency and intensity from 2006-2009 where 97% of the region were infected (Ministry of Health Republic of Indonesia, 2010) . WHO reported that Indonesia is the second largest with dengue cases in the world after Brazil (WHO, 2012) . The most prevalent of dengue cases in Indonesia outside Java Island was North Sumatra Province (Ministry of Health Republic of Indonesia, 2018) where dengue has widely spread with high morbidity and mortality rate. According to Health Ministry of North Sumatra, Medan city is the region with the highest cases of dengue (Medan City Public Health Agency, 2011) . There were more than 1300 annually dengue cases in the past 5 years in this city (Medan City Public Health Agency, 2011) . Medan is largest city outside Java in Indonesia and in the past 30 years, the urban population in Medan city has increased by 718 thousand people (53.3%) and in the next 30 years it is expected to increase by a further 1.4 million people (132.1%) (BPS Medan, 2015) . This condition creates Medan more vulnerable to dengue prevalence because high population density provides human virus reservoirs which allow rapid dengue transmission (Fullerton et al., 2014) . Despite the socio-economic condition, climate factor also gives significant impact to the rate of infectious vector-borne disease including dengue (Husain & Chaudhary, 2008) .
According to the previous studies, suitable local temperature, precipitation, relative humidity and water vapor are strongly associated with dengue risk (Bhatt et al., 2013) . Temperature, for example, is an important element of biting rate, egg and immature mosquito development, virus development time in the mosquito and survival of the mosquito life cycle (Christoper, 1960) . Furthermore, precipitation provides habitat for mosquito life cycle and highly influences the vector distribution (Morin et al., 2013) . It also can compose the size, population and behavior of Ae. aegypti (Morin et al., 2013) . Other weather variables, such as humidity and evaporation rate influence the longevity of Ae. aegypti which can lead the infected female Aedes to complete more than one cycle of the virus replication (Anderson & Ricco-Hesse, 2006) . In Indonesia, climate variability study and dengue cases were conducted by using bivariate analysis (Sintorini, 2007) and descriptive statistic (Dini et al., 2010; Wirayoga, 2013) in different site location. They stated that there was no significant correlation between the climate factors and the incidence rate of dengue. Moreover, those researchers employed simple linear regression known as the most widely used to analyze the relation between response variables (i.e. dengue cases) and explanatory variables (i.e. climate variables). The common situations we deal with human health, environmental and climate data is that the relationship between response variable and explanatory variables is non-linear. One of the approaches to cope with nonlinearity in regression issue is Generalized Additive Models (GAM) which fits a smoothing curve through the data and keeps the conditions of normal errors, constant variance and independence. Thus, we performed GAM in this analysis to represent the data driven rather than the model driven.
As we show later in this paper, there was an observed trend of increase in temperature and dengue cases in the period between 1980-2015 and 2002-2015 (BMKG, 2016; BPS Medan, 2002-2015) . In this paper, we estimated the association between the number of dengue patients with several weather factors. This was performed through GAM, considering the number of dengue patients as response variable and climate factor such precipitation, minimum temperature, maximum temperature, average temperature and relative humidity as predictor variables. One of the aims of the study was to estimate the effects of climate variables on the evidence of dengue. Another objective was to use the model to predict vulnerable area for dengue cases.
2. Method
The overall research framework on impact of climate variables on dengue was summarized in Figure 1. There was three main process which conducted in this paper, namely compile data base, import data to R and assess GAM, design vulnerability index and generate in geographic information system (GIS). To assess the impact of climate variables on dengue, GAM was applied. More explanation about GAM was described in Section 2.3. Furthermore, the blue color indicates predictor variables to generate vulnerability index map.
2.1. Study Area
Medan is located in North Sumatra, Indonesia at approximately between 3˚27'N and 3˚47'N latitude and 98˚35'E and 98˚44'E longitude as shown in Figure 2. Generally, climate in the study area is characterized by high temperature, high relative humidity and huge quantity of rainfall. The climate condition of study area is strongly influenced by Asia-Australia monsoon system namely Northeast Monsoon (December, January, February) and Southwest Monsoon (June, July, August) (Wrytky, 1961) . Moreover, climate system is also influenced by year to year El-Nino/Southern Oscillation (ENSO), Indian Ocean Dipole Mode (IOD) and local topography (As-syakur et al., 2013) . The ENSO and IOD event caused extreme weather condition, resulting in incidence of flood and forest fire. In
![]()
Figure 2. The study area in Medan City, North Sumatra.
addition, this study area is located in low altitude and most of the land use is dominated by residential area (BIG, 2016) .
2.2. Climate and Epidemiological Data
The study period consisted of 35-year period of daily weather data and 13-years period of annual epidemiological data (i.e. number of Dengue Hemorrhagic Fever (DHF) patients in each sub-district). The weather data was downloaded online in BMKG (i.e. Indonesia Meteorological Agency) portal website (BMKG, 2016) and it consists of geographical information, date and the magnitude of climate data such as daily minimum temperature, daily maximum temperature, daily average temperature, daily relative humidity and daily precipitation. The epidemiological data were collected from Medan in Figures publication from 2002 to 2016 and it was downloaded online in Medan statistical agency portal website (BPS Medan, 2002-2015) .
2.3. Generalized Additive Model
The GAM is a semi-parametric extension of a generalized linear model, which has the smooth components of the predictor variables (Wood, 2006) . The additional value of GAM is their ability to convey highly nonlinear and non-monotonic relationships between the response and predictor variables (Wood, 2006) . GAM models were used in this study to assess the influence of climate variables on the number of dengue patients. This statistical method has been commonly used for ecological studies (Guisan et al., 2002; Setiawati et al., 2015) , but has rarely applied for health risk assessment. GAM model was built in R version 3.3.0 software, using the gam function of the mgcv package (Wood, 2006) , with the number of DHF patient as a response variable and precipitation, minimum temperature (Tmin), maximum temperature (Tmax), average temperature (Tavg) and relative humidity (RH) as predictor variables. GAM models were explained in Equation (1). In addition, GAM model also was employed to assess the trend of climate variability in annual and monthly scale, where time was used as predictor variables and climate data as response variables.
(1)
where
is response variables, fi is smoothing function, n is number of predictor variables and Xi is predictor variables. The model was evaluated based on the significance level of predictors (P-value), deviance explained (DE), and the Akaike information criterion (AIC) values. In addition, DE was the same role as determination coefficient (R2) in linear regression analysis. The smallest value of AIC and the highest value of DE were selected as the best model.
2.4. Vulnerability Index
Vulnerability map was constructed by combining climate variables, population density and population under 10 years old. Weight calculation in each variable was estimated by using ranking method (Stillwell et al., 1981) as explained in Equation (2).
(2)
where n is the number of predictor variables, rj is rank position, wj is weight normalized. Climate was known as the main factor which triggers dengue cases (Ebi & Nealon, 2016) . Hence, we put climate variables as the highest rank, followed by population density and population under 10 years old (Fullerton et al., 2014) . We used raster calculator function in the spatial analysis tools in ArcGIS 10.3 to processed vulnerability index. Vulnerability index was conducted by combining the predictor variables and weight factor (wj) as shown in Equation (3)
(3)
where Xj is the predictor variables for vulnerability index (i.e. climate variables, population density and population under 10 years old). Xj was divided into two value; 1 as vulnerable area and 0 as non-vulnerable area. The threshold of them was calculated by GAM (i.e. climate data) and using literature review for population density and population under 10 years old (Reid et al., 2009) . The vulnerability index was ranged from 0 to 1 where 0 indicates less vulnerable and 1 indicates most vulnerable.
3. Results
3.1. Time Series of Variables
Variation of climate influences the existence of Aedes in several mechanisms (Morin et al., 2013) . Figure 3 shows the annual trend of climate data from 1980 to 2015 and the annual trend of DHF patient from 2002 to 2015.
![]()
Figure 3. Annual trend of daily maximum temperature, daily minimum temperature, daily average temperature, daily relative humidity, daily precipitation and DHF patient. Y axis represent as “residual” which means the difference between GAM predictive value and the average value (i.e. residual means GAM predictive value minus the average value). The thick line represents the mean value while the dotted line represents the variance.
Figure 3 was built by using GAM analysis, where year as predictor variable, climate data and DHF patient as response variables. Y axis was defined as residual where the positive value means higher than the average value and vice versa. According to Figure 3, minimum temperature and precipitation tend to increase annually from 1980 to 2015. Furthermore, maximum temperature, average temperature and relative humidity were fluctuated. However, since 1992 maximum and average temperature were increased consistently. In addition, RH was consistently above the average from 1990-2003. This result was basically consistent to Alexander et.al, 2006 which projected that Sumatra Island temperature and precipitation trend was tend to increase (Alexander et al., 2006) .
The higher cases of dengue were occurred in the years of 2005 to 2010, with the highest dengue cases in 2009 (2090 cases). In the period between 2005 to 2009, the maximum rainfall volume recorded was 466.7 mm (in May 2009); the highest maximum temperature and minimum temperature was of 34.14˚C (May 2010) and 25.03˚C (May 2010). The highest relative humidity in the period was 94% (in September 2006, July and December 2007 and September and October 2009).
Figure 4 shows monthly trend of climate variables from 1980-2015. The highest temperature mostly found in transition month from Northeast monsoon to southwest monsoon (i.e. April, May, June) and lowest temperature was found in Northeast Monsoon (i.e. December, January, February). In addition, the highest rate of precipitation and RH occurred between September and December. According to Karyanti and Hadinegoro, dengue incident rate increased between January and March (Karyanti & Hadinegoro, 2009) . It was probably because from September to December, large amount of precipitation provides
![]()
Figure 4. Monthly pattern of maximum temperature, minimum temperature, average temperature, relative humidity and precipitation. The variables were transformed from daily to monthly data. Y axis represent as “residual” which means the difference between GAM predictive value and the average value (i.e. residual means GAM predictive value minus the average value). The thick line represents the mean value while the dotted line represents the variance.
mosquito habitat which is required for survival, egg laying, increased reproduction, and increased the viral replication (Morin et al., 2013) . In addition, ideal temperature according laboratory studies for all life phases of the Aedes is between 20˚C and 30˚C (Tun-Lin et al., 2000) which is fitted with Medan City temperature range (i.e. 20˚C - 35˚C) (BMKG, 2016) . As a result, the transmission time occurred after peak season of precipitation.
3.2. Climate Variables on Dengue by GAM
Table 1 shows the model, variable, P-value, DE and AIC for some GAM models. Data are presented only on peak season of dengue incident (i.e. January, February, March) (Karyanti & Hadinegoro, 2009) . Higher time resolution of response variable is better. However, only annual DHF cases can be obtained from local government. Hence, we assumed that number of DHF patient lied on the peak season of dengue incident as response variable, followed by predictor variables.
The variables which were most statistically significant (P < 0.05) were precipitation and maximum temperature and so the assessment was conducted with only these variables. Moreover, the deviance explained (DE) from the best model were 66.1%. DE has the same definition as the determination value in the linear regression. So that, maximum temperature and precipitation influence the number of DHF patient as 66.1%. The effect of climate conditions, deduced from GAMs, indicated precipitation was more important than maximum temperature. This was indicated by precipitation has the lowest Pvalue among other parameter. Figure 5 shows GAM plots developed to interpret the individual effect of each predictor variables on the number of DHF patients. The positive effect of maximum temperature on DHF patient was ranged from 31.8˚C to 33˚C. It means that when the maximum temperature is higher than 31.8˚C, it increased the DHF cases. Furthermore, the positive effect of precipitation on DHF patient was fluctuated. Even though, the threshold was difficult to find, higher precipitation pattern increased the probability of DHF risk.
![]()
Table 1. GAM models and obtained value.
Note: *Statistically significant (P value is less than 0.05). AIC defines the residual error. The highest value of DE and the smallest value of AIC were selected as the best model if all the predictor variables are significant.
![]()
Figure 5. Effect on maximum temperature and precipitation on the number of DHF patient. Y axis represent as “residual” which means the difference between GAM predictive value and the average value (i.e. residual means GAM predictive value minus the average value). The thick line represents the mean value while the dotted line represents the variance. Positive residual value defined the tendency of an increased case of dengue fever.
3.3. Vulnerability Map
The climate pattern in Medan district was generally similar, particularly on precipitation variables. Thus, we added socio demographic variables to build the vulnerability map. The vulnerability map was constructed by using four main variables; precipitation (i.e. rank 1), maximum temperature (i.e. rank 2), population under 10 years old (i.e. rank 3) and population density (i.e. rank 4). Climate is the main driver for dengue risk (Bhatt et al., 2013) . Thus, we put climate factors at the top rank of vulnerability analysis. Population under 10 years old gave more priority than population density because population under 10 years old was the main patient of DHF (i.e. 35% of total DHF patient) in Medan City (Fitriani, 2011) .
The description of vulnerability map was shown in Figure 6. The vulnerability score was range from less than 0.5 to 1 with the red colour was indicated as the highest risk of dengue. In addition, number of DHF patient was expressed with different sized of dark circle symbol (i.e. bigger circle means high number of DHF patient). Figure 6 stated that higher risk of dengue was located in the main city of Medan (i.e. 10 sub-district) which indicated with the red colour. It was due to higher population density, higher number of sensitive people to get infected by dengue (i.e. kids) and more frequent heavy rain and higher temperature (i.e. urban heat island effect). One hypothesis for urban climate is that urban heat island effect creates warmer temperature and generates the unstable air which leads more rain. Thus, in urban area, the problems such as, poor water and waste management, rapid urbanization, high population density, and changing climate conditions can increase water associated diseases and vector borne diseases including DHF.
Furthermore, we also checked the correlation between vulnerability index and the average DHF patient/ha as shown in Figure 7. The figure stated that there
![]()
Figure 7. The correlation between vulnerability index and DHF patient/ha.
was strong positive correlation between vulnerability index and density of DHF patient. This information is necessary to develop diseases early warning system in regional scale.
4. Discussion
Our results indicate the climatic characteristics of Medan; more rainy in September to December and less rainy in January to August, general condition of tropical climate. However, we have to keep in mind that it rains over the year and temperatures are also high over the year. The analysis of dengue cases and climate variables showed statistically significant between the disease and the maximum temperatures and precipitation. The model was also being able to estimate the magnitude of climate influence to dengue cases which showed by the DE value. However, we only have annual scale of DHF cases data so that the time lag between climate data and DHF incidence cannot be evaluated.
Temperature was identified as the main trigger for affecting dengue fever epidemic (Bhatt et al., 2013; Ebi & Nealon, 2016) . Though previous study, some investigators found that maximum temperature is more closely with dengue fever vector growth (Eastin et al., 2014) . According to this paper, it was found that maximum temperature is the most proper temperature to predict the DHF outbreak in the study area. We found that increasing daily maximum temperature can rise the number of DHF patient in the study area. Our result is also consistent with the Kristie and Nealon, 2016 which stated that at the higher temperature feeding behaviour of Ae. aegypti is more frequent which can increase the dengue risk (Ebi & Nealon, 2016) . Moreover, Ae. aegypti can accept a wider range of temperature, theoretically by utilizing habitats in the urban area with suitable temperature (Ebi & Nealon, 2016) . Some previous empirical relationship also reported that there were significant delays between temperature and dengue incidence in south Taiwan (Chien & Yu, 2014) . The estimated temporal lags can range from one month to three months (Eastin et al., 2014) .
Precipitation is well recognized as an important factor for dengue fever (Choi et al., 2016; Corwin et al., 2001) , though it also depends on the human activities about water storage and activities. This study identifies the relationship among three month accumulated rainfall during the DHF outbreak (i.e. January, February, March) and DHF patient. Precipitation explained a large portion in the variance in the study area with the DE more than 50%. It might happen because precipitation can rise the water level which also increases the mosquito survival rate and the development stage from larva to adult will be faster. Cheng et al., 2016 also stated that precipitation played a major role in forming the outbreak in Guangzou China in 2014 (Cheng et al., 2016) . However, precipitation has a complex relation to DHF cases because heavy rainfall also can destroy the breeding sites (Cheng et al., 2016) . Thus in Figure 5, the relation between precipitation and DHF cases are fluctuated and non linier. This result is also consistent with the finding from Cheng et al., 2016 which stated that the relationship between number of dengue cases and precipitation was not linear (Cheng et al., 2016) .
It is critical to consider that, in recent publication, Haryanto, 2018 did not explain clearly any significant correlation between climate parameters (temperature and precipitation) and dengue cases in Indonesia and his findings also did not include climate parameters for vulnerability assessment (Haryanto, 2018) . Kasertyaningsih et al., 2018 also stated that precipitation and temperature were not significant to the DHF cases in Yogjakarta, Indonesia (Kesetyaningsih et al., 2018) . It is because those authors applied a simple linear statistical model where in the real world linear relationships are very hard to find. Thus, it is important to show some non-linear statistical approach for analyzing the relationship between climate factors and DHF cases. One of the approaches is by using GAM Models. Finally, the dengue epidemic and its increase can be monitored by GAM Model, because it showed that the precipitation and maximum temperature has the DE value about 66.1% and the rest was influence by other factors.
Vulnerability map is important to develop in order to understand which area is more vulnerable for dengue, so it can be used by local authorities to set-up the priority area for dengue prevention and mitigation. The vulnerability assessment was built by integrating climate factor (i.e. parameter was chosen by using GAM) and non-climate factors. The result stated that the highest vulnerability index was located in the central of the city. Since vulnerability index was developed using weather and socio demographic data therefore our results can serve valuable reference for developing early warning system for local governmental agencies.
5. Conclusion
The relationship between climate variables and dengue cases was examined by using Generalized Additive Models. The result stated that climate variation influenced the number of DHF patients as 66.1%, particularly precipitation and maximum temperature. The highest risk of dengue was located in the main city of Medan, where the kid population and population density were high and the climate condition was favorable for all stages of mosquito life cycle. However, several limitations were recognized in this study; for example, a low temporal resolution of DHF in-situ data and the analysis was at the sub-district level rather than village level. For further analysis, higher temporal resolution data of DHF (i.e. weekly or monthly) are necessary to understand the time lag between peak season of precipitation and dengue cases. For further application, this model can be adopted as an alternative to develop “an early outbreak warning system” for dengue. Moreover, this initial study provides more information about local vulnerabilities, so it can be used to help decision maker to identify which area should be prioritized for dengue eradication programs.
Acknowledgements
We would like to thank to Indonesia Meteorological Agency (BMKG) and Medan Statistical Agency (BPS Medan) for providing the climate data and DHF data in Medan.