Selection and Characterization of Fugitive River Dust Episodes over Zhuoshui River in Taiwan

Received: July 16, 2018
Revised: November 10, 2018
Accepted: December 5, 2018

Download Citation: RIS|BibTeX|https://doi.org/10.4209/aaqr.2018.07.0259  

Cite this article:

Chang, I.C., Chiang, Y.C. and Yu, T.Y. (2019). Selection and Characterization of Fugitive River Dust Episodes over Zhuoshui River in Taiwan. Aerosol Air Qual. Res. 19: 995-1006. https://doi.org/10.4209/aaqr.2018.07.0259

HIGHLIGHTS

• Principal component analysis was applied to screen river fugitive dust episodes.
• Component scores of the PC1 for PM₁₀ levels are utilized as the suitable indicators.
• Statistical indicators related to PM₁₀ and meteorological parameters are evaluated.

Keywords: Principal component analysis; Synoptic weather; Fugitive river dust.

INTRODUCTION

In Taiwan, the base flow of certain rivers has dropped sharply because of topographic characteristics of drainage basins, climate change, water resource allocation, watershed management, and riverbank reclamation and development. After the 1999 Jiji earthquake, riverbeds were raised and exposure of downstream riverbeds increased. Furthermore, after typhoons, a large amount of debris is flushed down rivers from upstream; this occurs during the northeast monsoon season, resulting in substantially elevated levels of dust emissions. According to air quality monitoring conducted by the Taiwan Environmental Protection Administration (TEPA), dust episodes primarily occur from October to April of the following year. In addition to reducing the quality of life of and causing inconvenience for residents of downwind areas, dust episodes may affect residents’ physical and mental health. As well as long-term changes in the natural environment, possible causes of river dust include part of the riverbed gravel being exposed when a dry spell occurs during a strong northeast monsoon in winter; moreover, when farmers do not suppress dust carefully, river dust may be generated when the farmers use riverbanks during periods of fallowing, planting, and turning. Recently, numerous river dust-related grievances have been filed with the Taiwanese government in the surrounding basins of the Zhuoshui, Beinan, Daan, and Dajia Rivers, indicating that river dust emissions are a serious problem. Take Syuguang Station as an example, a river dust monitoring station along the Zhuoshui River; the maximum daily average concentration of PM₁₀ (particulate matter of 10 µm in diameter) appeared on October 22, 2013 (417 µg m^–3), during which the maximum 1-hour concentration was 1102 µg m^–3.

Screening for episodes of air pollution is essential preparation work for formulating regulations on air pollution, air quality simulation, and air quality management. In particular, simulating secondary air pollutants and related topics require great resources to process the complex physical and chemical reactions involved. Simulation calculations can be simplified and time saved by using representative air pollution episodes. Additionally, verifying the simulation results of air quality models through air quality monitoring data can facilitate analysis and understanding of the causality between numerous factors such as air quality, weather, pollution sources, and topography. Conventional screening of air pollution scenarios begins with classifying meteorological patterns and preparing relevant statistics. The classification of meteorological patterns must rely on meteorological data and the professional capability of meteorological experts. However, because the determination of classification methods and methods of expert assessment are not unified, the classification of meteorological patterns cannot objectively demonstrate the statistically extreme concentration of pollutants. The statistical methods used for meteorological patterns often adopt extreme values and the exceedances of ambient air quality as screening standards. When using extreme values, the statistical methods often cannot derive consistent results for numerous monitoring stations due to their differing characteristics.

Regarding the classification of meteorological patterns, the United States Environmental Protection Agency (U.S. EPA, 1991) selected meteorological conditions that were prone to high ozone and then selected the 3 days with the most severe ozone pollution levels in each meteorological pattern to define an ozone episode. Cassmassi (1999) used synoptic weather patterns and atmospheric stability to establish the Meteorological Potential of Atmospheric Pollution (PMCA) index, and adopted multivariate regression tools to predict PM₁₀ episodes. Perez and Reyes (2002) used a neural network to predict PM₁₀ episodes and found that temperature difference between day and night was a crucial factor in specific meteorological patterns. Regarding statistical research on monitoring values, Ames et al. (1985) selected episodes based on a frequency of not more than 10 times in 3 years (3 times a year); permissible concentration was also used to determine an episode. Meyer et al. (1997) attempted to simulate the maximum 8-hour ozone values from June 1 to August 31, 1987, using the Regional Oxidant Model, and selected meteorological conditions with high concentrations of pollutants, identifying these as ozone episodes. However, the three-dimensional grid air quality model entails high resource costs, which is unsuitable as a tool for determining an immediate pollution episode. Using the aforementioned meteorological classifications, extreme values statistical methods, and air quality models is not feasible for instantly selecting appropriate and statistically representative pollution scenarios. In terms of using statistical methods to screen pollution episodes, Yu and Chang (2000) and Yu (2013) adopted principal component analysis (PCA) to screen ozone pollution and PM₁₀ episodes, respectively. Kuebler et al. (2002) used classification and regression tree analysis to screen ozone episodes on the Swiss Plateau. Beaver and Palazoğlu (2006) adopted K-means cluster analysis to screen ozone episodes in California. Zhang et al. (2014) used PCA and a nonparametric T² control chart to predict episodes of over-standard ozone concentrations. Sun et al. (2015) predicted over-standard ozone episodes by adopting generalized linear mixed effects models (GLMMs). In the same study, the researchers compared the difference between linear regression models, generalized linear models, multilayer perceptron, and support vector machines, and found that GLMMs provided predictions with lower prediction errors and more correct over-standard ozone episodes for various weather patterns.

In the field of air pollution research, PCA is an effective and objective data analysis tool that is often used to simplify variables (Yu and Chang, 2006; Dai et al., 2015; Iodice et al., 2016; Yao et al., 2016; Chen et al., 2017), identify sources of pollution (Viana et al., 2006; Shi et al., 2009; Deka et al., 2014; Huang et al., 2015; Luo et al., 2015; Mari et al., 2016; Arhami et al., 2017; Widiana et al., 2017), classify meteorological patterns (Maheras, 1984; Maryon and Storey, 1985; Eder et al., 1994; Cheng and Lam, 2000), and evaluate model diagnostics (Eder et al., 2014; Li and Wen, 2014). Researchers can use large amounts of existing monitoring data to instantly select statistically representative air pollution scenarios and, in particular, provide consistent high concentrations of air pollutants at most stations (Du et al., 2017; Lang et al., 2017; Tian et al., 2018). These data can serve as a decision-making reference for formulating future strategies on air quality management to screen for air pollution episodes. The PCA has two advantages; first, it simplifies the relevant variables of air pollution to achieve economic effectiveness; and second, for screening pollution episodes, it exhibits more objective statistical representativeness than the meteorological classification. For screening river dust episodes, this study adopted the research method of Yu and Chang (2000) and Yu (2013). Based on air quality monitoring data, this study adopted PCA to screen pollution episodes, and then analyzed the air quality and meteorological characteristics of dust episodes. Finally, the appropriateness of the screening process was evaluated.

 
RESEARCH METHODS

 
Air Quality Data

The Taiwan Air Quality Monitoring Network was officially launched by the TEPA in September 1993. To date, 76 general stations, 63 regular stations, 4 industrial stations, 2 national park stations, 4 background stations, and 6 traffic stations have been installed. The following five stations are located along the Zhuoshui River basin: Jushan, Douliou, Lunbei, Mailiau, and Taishi. The monitoring parameters are SO₂, CO, O₃, PM₁₀, NO, NO₂, NMHC, THC, PM_2.5, wind speed, wind direction, temperature, and relative humidity. The Zhuoshui River is the longest river in Taiwan with the largest amount of sediment, and its upstream tributaries carry a large amount of sediment from the mountains, depositing into the estuary (Dadson et al., 2003). Its estuary is a wide and braided channel with numerous exposed sandbars. During the winter dry spell, these sandbars easily form a large bare area. Coupled with the prevalence of the northeast monsoon, PM₁₀ levels are elevated by wind, leading to serious river dust disasters in the downwind areas of the south bank of the Zhuoshui River. In Taiwan, river dust episodes usually occur during the dry spell of the northeast monsoon season (from October to April of the following year). In addition to general air quality stations, the TEPA has established automatic PM₁₀ monitoring equipment near rivers that commonly have fugitive river dust to understand river dust conditions and forecast air quality. High PM₁₀ concentrations are likely to be detected by the monitoring stations near riverbeds because dust accumulated on the riverbed generates dust emissions when blown by strong wind. River dust data measured by a monitoring station can only represent the influence of dust in the local area and cannot represent regional air quality. The TEPA established two stations at Yisian and Syuguang Elementary Schools. Along the Zhuoshui River basin, five parameters are monitored: PM₁₀, wind speed, wind direction, temperature, and relative humidity. To simultaneously evaluate the effect of river dust on air quality, the air quality monitoring data of this study were taken from the EPA’s general and river dust stations from 2011 to 2017. The target area was the Zhuoshui River basin, and the following seven stations were included: Jushan, Douliou, Yisian, Syuguang, Lunbei, Mailiau, and Taishi (Fig. 1).

Fig. 1. Locations of ambient air quality monitoring stations over Zhuoshui River.

This study’s primary tasks were as follows:

1. To screen river dust episodes of the Zhuoshui River: This study employed PCA of multivariate statistics to analyze data from the TEPA’s air quality and river dust monitoring stations for the Zhuoshui River basin from 2011 to 2017 (data from the river dust monitoring stations began from 2011). The river dust monitoring stations only recorded PM₁₀ data; therefore, PM₁₀ was the research target and PM₅ was excluded.
2. To mitigate the effects of dust storms and foreign pollution episodes on the screened river dust episodes: Affected by dust storms in Mongolia and foreign air pollution, the PM₁₀ concentration can increase sharply in a short period. Therefore, this study compared dates when dust storms and foreign pollution affected Taiwan and removed PM episodes that may have been caused by dust storms and foreign pollution.
3. To evaluate the appropriateness of episodes for statistical representativeness. Existing data of synoptic weather patterns, meteorology, and PM₁₀concentrations were used to evaluate the appropriateness of pollution episodes.

 
Principal Component Analysis

Monitoring data from air quality stations were analyzed
from the EPA’s air quality stations. The data period was from 2011 to 2017, and the daily average PM₁₀ value was taken as the sample value. The data dimension of PM₁₀ could be treated as a collection of time series of 7 vectors (7 stations × 2557 days). Subsequently, the monitoring data were normalized using the following equation: 

where Z_ik is the score of the k^th Z of station i; C_ik is the k^th pollutant value of station i; µ_i is the average value of station i; and S_i is the standard deviation of station i. The relationship between the standardized Z-score and the nonrotating principal component value is as follows: 

where L_ij is the factor loading of the j^th principal component of station i; and P_jk is the component score of the k^th variable in the j^th principal component. The score of the principal component can be derived using the inverse matrix of the aforementioned formula:

where λ_j is the eigenvalue of the j^th principal component and also represents the variance of the j^th principal component. All of the principal components are arranged based on the order of their explained variance from greatest to least; therefore, using the first few principal components enables the grasping of most of the variance in the entire sample space and achieves the simplification of variables.

The principal components of these linear combinations exhibit not only the smallest variance in their components, but also the maximum individual differences of the individual components. Therefore, the first principal component can explain the maximum variation in the concentration field. This study screened appropriate river dust episodes using the characteristics of the first principal component.

RESULTS AND DISCUSSIONS

 
Geographical Distribution of PM₁₀ Concentrations during Dry Season

The highest monthly average of PM₁₀ daily average values (Fig. 2) from the monitoring stations in the Zhuoshui River basin occurred in November, followed by December and March. The monthly mean PM₁₀ during the dry season at all monitoring stations was substantially higher than that during the wet season. The difference in the monthly average concentrations during the wet and dry seasons was approximately 30–45 µg m^–3. The monitoring results from the river dust monitoring stations at Yisian and Syuguang Elementary Schools were compared, as shown in Fig. 3(a). At Syuguang, the cumulative probability density distribution indicated that differences in PM₁₀ concentration during wet and dry seasons at 25%, 50%, and 75% cumulative probability were 21, 30, and 40 µg m^–3, respectively; at Yisian, they were 21, 31, and 40 µg m^–3, respectively. Notably, the difference in extremely high PM₁₀ concentrations at Syuguang was much higher than those at Yisian. According to the locations of the stations and the Zhuoshui River, both stations were located on the south bank of the river. From the difference in PM₁₀ concentrations between wet and dry seasons, the difference was consistent at less than 90% cumulative probability; however, for extremely high PM₁₀ concentrations, the river dust exhibited greater influence at Syuguang than at Yisian. For the differences in PM₁₀ concentrations between wet and dry seasons at the TEPA’s general stations, the two sites exhibiting the largest difference—Taishi and Mailiau—were selected for comparison (Fig. 3(b)). The cumulative probability density of difference in PM₁₀ concentration during the wet and dry seasons at Taishi were 19, 26, and 30 µg m^–3 at the cumulative probabilities of 25%, 50%, and 75%, respectively; those at Mailiau were 25, 33, and 37 µg m^–3, respectively. However, the cumulative probability density of 90% and 95% for PM₁₀ difference (the difference in PM₁₀ concentrations for wet and dry seasons) were 37 and 56 µg m^–3 in Taishi; and 42 and 80 µg m^–3 in Mailiau. Compared with the aforementioned seven stations, river dust was found to exhibit the greatest influence on the closest stations at Syuguang and Yisian on the south bank of the Zhuoshui riverbed, followed by the Mailiau and Taishi Stations.

Fig. 2. The monthly mean PM₁₀ concentrations of monitoring stations. (a) Syuguang and Yisian (b) Taishi and Mailiau

Fig. 3. Distribution of daily mean PM₁₀ concentrations in dry and wet seasons.

PCA

PCA was performed based on the daily average PM₁₀ concentration from the seven stations. The results (Table 1) showed that the eigenvalues of the first three principal components were 4.57, 0.89, and 0.55, respectively, and the explanation for the concentration variance was 65.3%, 12.7%,and 7.9%, respectively. Therefore, the first unrotated principal component accounted for 65.3% of the variance in PM₁₀ concentration. The factor loading value represented the correlation between the principal component and each station. From the factor loading of each station and the first principal component (Table 2), the top three stations exhibiting the highest factor loadings of the first principal component were Douliou (0.888), Lunbei (0.848), and Taishi (0.834). The factor loading of each station and the first principal component was higher than 0.72. For the absolute value of the factor loading of the second principal component, only Jushan exhibited a value higher than 0.5.

If the first unrotated principal component exhibited a high value, the concentration of all stations would collectively increase. Therefore, to understand harmful river dust episodes, the number of exceedances (in 1 day, if one station detected a PM₁₀ concentration higher than daily mean PM₁₀ standard, the number of exceedances was 1; if it was two stations, the number was 2) for PM₁₀ was selected and the component scores of the first principal component were compared. As for the standard, this study selected two criteria: 150 and 125 µg m^–3 (150 was STN1 and 125 was STN2). According to the results in Fig. 4, the higher the component score of the first principal component, the higher the numbers of exceedances there were. The results of the analysis were as follows:

1. The component score of the first principal component that was higher than 1 accounted for 13.5%, and could screen for 78.9% of exceedanceswith PM₁₀ concentrations exceeding 125µg m^–3. The average number of exceedances was 1.89. The component score could screen for 85.0% of exceedances with PM₁₀ concentrations exceeding 150 µg m^–3, and the average number of exceedances was 0.88.
2. The component score of the first principal component that was higher than 2 accounted for 3.30% and could screen for 42.3% of exceedanceswith PM₁₀ concentrations exceeding 125 µg m^–3. The average number of exceedances was 4.14. The component score could screen for 60.4% of exceedances with PM₁₀ concentrations exceeding 150 µg m^–3, and the average number of exceedances was 2.55.

Fig. 4. Relationship between PC1 score and percentage of exceedance for daily PM₁₀ standard and average numbers of exceeding daily PM₁₀standard.

In addition to assess the relationships between component scores of the first principal component and numbers of PM₁₀ episodes detected by monitoring stations, this study evaluated synoptic weather patterns (Soong et al., 2005), wind speed, wind direction, height of mixing layer (Air Quality Modeling Center, TEPA, 2018), and ventilation index (wind speed recorded in Wuqi Weather Station, height of mixing layer in CWB stations; Air Quality Modeling Center, TEPA, 2018). Therefore, according to existing monitoring data, parameters such as the PM₁₀ concentration, ventilationindex, and height of mixing layer exhibited low degrees of correlation. No moderate positive correlation (0.4–0.7) was found among theseparameters. The correlation between wind speed, wind direction, and PM₁₀ concentration was investigated. The number of exceedances for PM₁₀standard in each year was compared, and the highest number of exceedances occurred in 2013 (100 times) and 2014 (121 times). The station with the most severe PM₁₀ concentration in 2013 was Syuguang (23 times). The highest PM₁₀ concentration recorded by Syuguang was in October. Therefore, the hourly concentration, wind speed, and wind direction in October 2013 were analyzed to realize correlations between these parameters. The most severe PM₁₀ episode detected by Syuguang in October 2013 was October 21–27. During this period, the air quality was also affected by Typhoon Francisco and Typhoon Lekima. The relationship between the hourly PM₁₀ concentration and wind speed at Syuguang was analyzed during these two selected periods. The relationship between hourly PM₁₀ concentration and wind speed at Syuguang in October 2013 (Fig. 5(a)) and October 21–27, 2013 (Fig. 5(b)), demonstrated that the hourly PM₁₀ concentration and wind speed did not exhibit favorable linear or quadratic correlations. The relationship between hourly wind direction and hourly PM₁₀ concentration at Syuguang in October 2013 and the average PM₁₀ concentration at every time period are shown in Figs. 5(c) and 5(d). The hourly wind direction in Fig. 5(c) shows that in October 2013, winds with PM₁₀ concentrations higher than 200 µgm^–3 were from the north, northeast, and northwest, whereas PM₁₀ concentrations higher than 400 µg m^–3 were only from the north and northwest. Syuguang is located on the south bank of the Zhuoshui River; thus, the air quality at Syuguang was obviously affected by fugitive river dust. As for the periods of high PM concentration and high wind speed (Fig. 5(d)), they all occurred at approximately 15:00 and 16:00. Therefore, high wind speeds resulted in high PM₁₀ concentrations during river dust episodes.

Fig. 5. Relationships between hourly PM₁₀ concentrations and meteorological parameters at Synguang Station (July, 2013).

To evaluate the effectiveness of screening river dust episodes, daily mean PM₁₀ values exceeding 125 and 150 µg m^–3 were compared. The numbers of exceedances for PM₁₀ standard over 125 µg m^–3, which was more statistically representative (the number of exceedances was 497 in 2011–2017) than 150 µgm^-3. The primary stages of principles for screening river dust episodes are as follows: (1) During the monitoring period, 1 day occurs where the component score of the first principal component is higher than 3. (2) Dust storms and transboundary pollution are eliminated. (3) The component score of the first principal component is higher than 1 for 5 consecutive days.

According to the first stage, the total number of days that the score of the first principal component was higher than 3 was 18 days. The results of screened PM₁₀ concentrations at stations are shown in Table 3. The highest score of the first principal component was found on January 21, 2014, and the component score of the first principal component was 6.05. The daily average PM₁₀ concentration value was 345 µg m^–3 at Syuguang and 368 µg m^–3 at Taishi. For the second stage, the dates of dust storms and transboundary pollution were compared. During the 18 days, except for October 22, October 25, and October 13, 2017, which were affected by the three typhoons Francisco, Lekima, and Khanun, respectively, the high PM₁₀ concentrations on the remaining dates were caused by transboundary pollution. According to the third stage, river dust episodes that exhibited a component score of the first principal component higher than 1 for 5 consecutive days were selected. Therefore, only the dust episodes during October 22–26, 2013, could be selected. To provide multiple choices, the score of the first principal component was modified to be higher than 2, enabling three alternative episodes to be selected (Table 4). Synoptic weather patterns were an essential factor that affected air quality; the cross-correlation among PM₁₀concentrations, principal components, and synoptic weather patterns was explored. The monitoring results from Syuguang (Table 5) indicated the following: (1) The first three synoptic weather patterns exhibiting the highest daily PM₁₀ concentrations were northeast monsoon and peripheral circulation of the typhoon, peripheral circulation of the typhoon, and high-pressure circulation mixing with a warm sector. The average wind speeds of the three weather patterns recorded at Lunbei were 3.7, 3.0, and 1.9 m s^–1, respectively. (2) The average wind speed and PM₁₀ levels of strong northeast monsoons, standard northeast monsoons, and weak northeast monsoons were compared. The strong northeast monsoon exhibited the highest average wind speed (4.0 m s^–1) and the highest PM₁₀ concentration. The six stations all exhibited these characteristics except for Lunbei. (3) Regarding the weather pattern of the northeast monsoon and peripheral circulation of the typhoon, Syuguang, Mailiau, Taishi, and Douliou Stations exhibited the highest average PM₁₀ concentrations.

(4) The top three synoptic weather patterns exhibiting the highest component scores of the first principal component were the northeast monsoon and peripheral circulation of the typhoon, high-pressure circulation mixing with a warm sector, and high-pressure recirculation. Whether high-pressure recirculation and wind speed severely influence fugitive river dust will be examined later. (5) For the synoptic weather pattern exhibiting the highest PM₁₀ concentration, a strong northeast monsoon was determined to carry the highest PM₁₀ concentration at Mailiau, Pacific high pressure stretching westerly at Yisian, and high-pressure recirculation at Lunbei.

In addition to considering synoptic weather patterns, a cross-analysis was conducted to carefully evaluate the correlation between the scores of the first principal component, weather patterns, and wind speeds at Lunbei (wind speed was divided into high and low wind speeds at 3.5 m s^–1) (Table 6). The findings were as follows: (1) High wind speed elevated the component score of the first principal component of weather patterns including the northeast monsoon and peripheral circulation of the typhoon, peripheral circulation of the typhoon, strong northeast monsoon, and the Pacific high pressure. That is, these four synoptic weather patterns increased the PM₁₀ concentration at each station under the condition of high wind speed. (2) Comparing a standard northeast monsoon with a weak northeast monsoon, the component score of the first principal component of high wind speed was not necessarily higher than that of low wind speed. (3) In the Zhuoshui River basin, river dust episodes were easily formed in the aforementioned four synoptic weather patterns at high wind speeds. Northeast monsoon and peripheral circulation of the typhoon as well as peripheral circulation of the typhoon were directly related to typhoons. Based on previous synoptic weather patternsin Table 6, the Pacific high pressure stretching westerly exhibited higher scores in the first principal component. This Pacific high pressure stretching westerly often occurred at 1 day before or after two weather patterns of northeast monsoon and peripheral circulation of the typhoon, and peripheral circulation of the typhoon. (4) A typhoon that affected Zhuoshui River areas was usually accompanied by rainfall, and the rainfall reduced the PM concentration.

Relevant studies on fugitive river dust pollution in Taiwan have begun by separating the seasons into wet and dry. However, the present study investigated monitoring data from the Zhuoshui River basin and found that typhoons in July, August, and September still caused local river dust episodes; for example, September 15, 28, and 29, 2012, and October 16–21, 2015. In addition to transboundary pollution (Oh et al., 2015; Lai et al., 2016) and dust storm episodes, two characteristics—typhoons that did not pass over the land in Taiwan and northeast monsoons during dry season—are also likely to cause river dust episodes in the Zhuoshui River basin. Based on the past experiences, fugitive river dust episodes easily occurred during dry season at Taiwan, and the construction period of performing yearly control measures to combat fugitive river dust was then set in September. However, analytical results of this study indicated the typhoon season (often occurring July–September) (Fang et al., 2009; Cheng and You, 2010) could be another crucial factor to cause fugitive river dust pollution. Therefore, control method and construction period of mitigating fugitive river dust at Taiwan must take into account the typhoon season.

 
CONCLUSION

This study applied unrotated principal component analysis to screening fugitive river dust episodes and analyzed the PM₁₀ concentrations, synoptic patterns, and meteorological parameters for Zhuoshui River. The first component, as confirmed by its component score, is an adequate indicator of air pollution episodes. The number of exceedances at a monitoring station with regard to the PM₁₀ standard increased with the component score of the unrotated first principal component. Aside from PM₁₀ episodes resulting from transboundary pollution and dust storms, the three dominant synoptic patterns that caused fugitive river dust episodes on Zhuoshui River were the peripheral circulation of a typhoon, the northeast monsoon and the peripheral circulation of a typhoon, and a strong northeast monsoon. This study identified two phenomena favorable to fugitive river dust pollution: firstly, a typhoon that did not pass over the land in Taiwan; and secondly, a strong northeast monsoon during the dry season. Since the typhoon season is highly conducive to fugitive river dust episodes in Taiwan, as demonstrated by the results of multivariate analysis, it must be taken into account when proposing and implementing control measures for this type of pollution.

 
ACKNOWLEDGMENTS

We express our gratitude to the Ministry of Science and Technology of Taiwan (MOST-106-EPA-F-013-001) and Environmental Protection Administration of Taiwan for funding this study.