Mapping dust storm PM2.5 pollution risk using indicator kriging in northern Taiwan

For this study, dust storm data was collected on 24 25 March, 2012. This dust storm impacted the northern area of Taiwan for about two days, so our study was confined to the northern air quality region of Taiwan using indicator kriging for assessing the probability of exceeding the threhold concentration of PM2.5. Results showed that the highest levels of dust storm pollution in the study area were in the Zhongzheng District of Taipei City, with a concentration of 39.86 μg m-3. The lowest levels were 21.88 μg m-3, in Xinzhuang District, New Taipei City. Among the four administrative regions, variation was greatest in Taipei City, with concentrations between 22.74 and 39.86 μg m-3, and lowest in Keelung City, with concentrations between 23.97 and 27.56 μg m-3. We categorized pollution concentrations with probabilities greater than 0.7 into three categories to represent health risk hazards. Eighty-seven percent of areas were low risk, with PM2.5 concentrations above 25 μg m-3. Areas with PM2.5 concentrations greater than 30 μg m-3 were considered moderate risk areas, and comprised 1.2% of all area; this category appeared in the Zhongshan, Songshan, Datong, and Daan Districts of Taipei City and the Zhongli District of Taoyuan City. Concentrations greater than 35 μg m-3, labeled high risk areas, accounted for only 0.4% and were concentrated in the Zhongshan, Songshan, and Zhongli Districts of Taipei City. The methods used in this study and its results can be key references for future early warning and prevention of dust pollution by relevant authorities. Article history: Received 7 April 2019 Revised 31 July 2019 Accepted 7 November 2019


INTRODUCTION
Meteorologically, a "dust storm" is characterized by sandy and dusty weather that reduce surface visibility to less than 1 km (Goudie and Middleton 1992). Dust storms are a very active spring weather phenomenon in the arid regions of East Asia. Although the probability of dust storms directly affecting Taiwan is low, they can cause large-scale air quality deterioration in Taiwan. When dust storms affect Taiwan, what are the effects of this pollution on our living environment? In which areas does the pollution threaten health? In which areas is human health at high risk? Where should people avoid going outdoors? All these problems can be addressed using "dust storm pollution risk mapping". Such mapping is an important and urgent task because it can not only assess the probabilities of hazards, but also help to avoid and reduce resultant damage. Air pollution refers to the aerial presence of one or more substances of which the amount, nature, and time of exposure can harm human, animal, or plant life; damage property; or interfere with a comfortable living environment. Taiwan is an island located in the monsoon climate region of the western Pacific. Northeasterly winds prevail in the winter due to cold, high pressure areas in Siberia. Taiwan is often affected by dust storms during winter and spring. These dust storms are characterized by increased PM 2.5-10 concentrations, lowered relative humidity, increased wind speed, and stable wind direction. The PM 2.5 comes mainly from pollutants in aerosol phase such as ammonium sulfate and ammonium nitrate, as well as from motor vehicle exhaust. PM 2.5 poses a greater hazard to human health than other aerosols.
On 12 March 1995 in Keelung, Yilan, Xizhi, Ruifang, and other areas in northern Taiwan, vehicles were covered with a layer of reddish-brown mud after a rainfall. Reddishbrown mud also precipitated from roadside stagnant water and accumulated on the roadsides. The local people suspected that this was caused by discharge of pollutants from a nearby power plant. Subsequently, with the help of several meteorological and environmental experts, it was found that the reddish-brown muddy rain originated from the long-distance transmission of dust from the Asian mainland. Since then, Taiwan has begun to pay attention to the effects of dust from the Asian mainland on air quality in Taiwan. Although there are few cases in which dust storms in the mainland have seriously affected air quality in Taiwan, the concentration of suspended particulates in large areas of air in Taiwan may increase significantly in a short period of time. Therefore, the Taiwan EPA (Taiwan Environmental Protection Administration) attaches great importance to this issue. In past years, Taiwan EPA focused on several acts such as policy planning for protection of air quality, stationary air pollution source control, mobile air pollution source control, air quality purification zone, and indoor air quality management. In recent years, as desertification in northwestern China has become increasingly serious (Feng et al. 2015), the frequency and scale of sandstorms has increased (Chang 2002;Liu and Young 2002). This means that Taiwanese' health and environment will continue to be affected. Mainland dust storms occurred five times in the 1950s, eight times in the 1960s, 13 times in the 1970s, 14 times in the 1980s, 23 times in the 1990s, and 12 times in the 2000s. Yuan et al. (2004) investigated the Asian dust storm periods during 2001 and 2002 and showed that the Asian dust storm invaded Taiwan from either the northwest or the northeast. Increasing both PM 10 concentration and coarse particle mode in the size distribution of atmospheric aerosols validated the invasion of Asian dust storms. Chan and Ng (2011) analyzed 14 years data of air quality and mortality in Taipei between 1994 and 2007. Compared with reference days, particulate matter with aerodynamic diameter less than 10 and 2.5 μm (PM 10 and PM 2.5 ) increased statistically significant by 24.2 and 7.9 μg m -3 per dust day, respectively. There were also statistically significant increases in sulfur dioxide (SO 2 ) and ozone (O 3 ) but decreases in temperature during Asian dust storms. Zhao et al. (2011) reported that heavy dust storm recorded in the past 50 years during 21 -23 March 2010, the compositions of ions in dust aerosols showed that the multi-sources of aerosol were ubiquitous during the dust episode and highest concentrations of for PM 2.5 and PM 10 reaching 454.51 and 990.24 μg m -3 , respectively. Fu et al. (2010) also reported the dust storm event in Shanghai from 20 March to 19 April 2007, the daily 24 h average PM 10 concentration of 648 mg m -3 was observed on 2 April 2007 and PM 2.5 /PM 10 was 15.5% on that day. The result indicated that the dust passed through the East China Seas before reaching Shanghai, which is one of the typical dust pathways that lead to heavily polluted days in Shanghai due to dust transport. Miller-Schulze et al.
(2015) analyzed two remote sites in Central Asia from July 2008 to July 2009 to show mineral dust and organic matter (OM) was the dominant contributors to PM 10 and PM 2.5 . Mineral dust was a more significant contributor to the coarse PM (PM 10-2.5 ) during high event samples and explained the majority of the variance in PM concentrations, and that the major apportioned factors of PM 10 and PM 2.5 were chemically similar between sites.
In addition, in the summer Taiwan is influenced by high pressure areas in the Pacific with southwesterly monsoon winds prevailing (Chen et al. 1999). Studies of three major cities in the northern, central, and southern regions of Taiwan (Taipei, Taichung, and Kaohsiung, respectively), have shown that the ratio of PM 2.5 to PM 10 concentrations differs in all three cities. The variation between the ratios shows that the original pollutant sources in the three cities are different and change over time. Chuang et al. (2008) found that PM 2.5 and organic matter tend to accumulate in the Taipei Basin and that local vehicle exhaust emissions are the main cause of PM 2.5 and organic matter formation. Their study also showed that PM 2.5 has accumulated in the Taipei Basin due to geomorphic barriers. It showed that air quality in the Taipei basin is affected by the northeast monsoon passing through the Chinese mainland, resulting in an increase in PM 2.5 concentrations; in summer, when the monsoon does not pass through the Asian continent, air quality is relatively unaffected. In addition, Taiwan is often affected by dust storms from winter to spring each year, with the impacts felt as far as Siberia, China, South Korea, Japan, and even North America. The dust storms affecting Taiwan are characterized by an increase in PM 10-2.5 concentrations, decrease in relative humidity, increase in wind velocity, and stable wind direction. After a dust storm, PM 2.5 concentrations increase, and research indicates that this is strongly correlated with an anticyclone caused by atmospheric retention; larger particles from dust storms on the Asian mainland mainly affect areas near the source and then settle to the ground. Smaller particles, however, can be transported to higher air layers and then transported by westerly winds toward the east. In addition to affecting the Asian mainland, this smaller particulate matter affects regions further east, including Japan, South Korea, the North Pacific, Hawaii, and as far as Canada and the western coastal United States; further south, Taiwan, Hong Kong, and even the Philippines are affected. The scope of influence is therefore quite large. When air passes over the Pacific Ocean, its composition includes non-sea salt sulfates, although sea salt comprises the most abundant coarse particles. The high altitude areas in Taiwan may be affected by long-range transport of products from biomass burning in southeastern Asia (Chen et al. 2013). Wang et al. (2012) applied the April 2001 dust storm episode over the trans-Pacific domain using CMAQ model developed by US Environmental Protection Agency (EPA) to demonstrate the long-range transport of Asian pollutants can enhance the surface concentrations of gases by up to 3% and aerosol species by up to 20% in the Western US. According to a study conducted in Hanoi, Vietnam, the main local PM 2.5 sources were local cooking of foods, traffic activities, and industrial emissions. In addition to sea salt, the bulk of PM 2.5 was composed of soil, road dust, and matter from construction, with construction being the main source (Hai and Kim Oanh 2013). Moreover, analysis of air pollution distribution in metropolitan areas due to residential and transportation activities has shown that these activities have increased the impact of PM 2.5 on climate change (Diaz and Dominguez 2009).
So-called dust storms are caused by strong winds that pick up many sand and dust, resulting in dusty weather that limits visibility. They are basically the result of arid desert environments. It is therefore distinct from "strong wind", "scattered dust", or "floating dust" that drifts in from other places during windy weather. Dust storms in Asia come from the deserts of northwestern mainland China. Despite the fact that dust storms on the mainland have seriously affected air quality in Taiwan over a long timeframe, once a dust storm arrives in Taiwan, it may cause a rapid concentration of aerosols in a large volume of air within a short time period. Consequently, the Taiwan EPA pays considerable attention to these phenomena. In general, the arrival of dust storms is characterized by a rapid increase in PM 10 concentrations, an increase in wind speed, and a decrease in relative humidity. Derivative particles are the main substance in urban air pollutants and have a significant effect on the environment and human health. In December 2011, Taiwan incorporated PM 2.5 into its air quality regulations, adopting the same standard values as the United States (adopted 2006) and Japan (adopted 2009). These were revised and released on 14 May 2012, with a 24-hr value (i.e., 24-hr continuous sampling) of 35 μg m -3 and a mean annual value (arithmetic mean for the daily mean for each day of the year) of 15 μg m -3 , which is comparatively strict national standard that in the world (Environmental Protection Agency, Executive Yuan of Taiwan 2020, https://airtw.epa.gov.tw). The results of this study showed that, when dust storms occurred, the concentration of PM 2.5 in Taiwan's northern air quality region exceeded the annual average standard of 15 μg m -3 , and even the 24-hr standard value exceeded 35 μg m -3 in the center of Taipei City. Liu et al. (2016) used MATLAB spatial interpolation tools to predict the average concentration of PM in Beijing during 2013 and 2014. The results showed the main factors of PM pollution were dust storms and strong winds in spring and autumn, rainfall and the warm wet climate in summer, and cold fronts and snowfall in winter. Pollution characteristics in the Beijing area were higher in the south and lower in the north, and the pollution sources might be regional transport as well as local anthropogenic sources. Cao et al. (2015) developed a new approach to identify sand and dust storm (SDS) source areas in Iran using a combination of nine relat-ed datasets, namely drought events, temperature, precipitation, location of sandy soils, SDS frequency, human-induced soil degradation (HISD), human influence index (HII), rain use efficiency (RUE), and net primary productivity (NPP) loss to identify the source areas in Iran using kriging interpolation technique.
In order to estimate and map air pollution risk during development of a dust storm, geographical statistics which offer optimal linear unbiased estimation, in combination with indicator kriging approach was first attempted in this study. In addition to estimating the air pollution concentration in areas without measuring stations, it is also successful to delineate to construct an air pollution risk map and identify the probability of harm to health. The main purpose of this study was to use the indicator kriging method to conduct risk mapping for dust storm pollution. For this purpose, data on dust storm events affecting the Taiwan area were collected from the website of the Bureau of Environmental Protection of the Executive Yuan's dust database. The database showed that there were no dust impact events from 2013 to 2015. There were only two dust impact events in 2012, which were a dust storm event from 24 to 25 March 2012, mainly affecting the northern region, and another occurring from 31 March to 1 April 2012, which had only an extremely slight island-wide impact (URL: https://taqm.epa. gov.tw/dust/tw/Database.aspx). This study focuses on the 26 air quality monitoring stations in the air quality region of northern Taiwan. It uses geostatistics, which are well-suited for optimal linear unbiased estimation, combined with the indicator kriging method to estimate air pollution risks during dust storms and construct a map of these risks. In this study, the spatial variability of PM 2.5 concentrations during dust storms was first explored, then pollution risk maps for the entire region were constructed. Finally, the probability of hazardous concentrations of PM 2.5 in the northern air quality region was estimated.

Research Stations
In August 2005, 76 air quality monitoring networks across Taiwan were completed, including placement of PM 2.5 automatic monitoring equipment in major metropolitan areas and eastern and outer islands. The 76 air quality monitoring stations in Taiwan have been rationally distributed across all counties and cities, and monitoring data is also regionally representative. In 1997, the EPA divided Taiwan into seven air quality regions according to regional characteristics and since then has conducted characterization analyses, implemented a grading management system for control zones in 2002, and incorporated Taiwan's offshore islands zone. The northern air quality region is the northern yellow block and its largest boundary area is New Taipei City. This block includes Keelung City on the east side and Taipei City on the west side. The block on the southwest side is Taoyuan City. Except for Keelung, all these cities are among the six largest municipalities of Taiwan. The population of this air quality region is estimated 9.2 million which accounts for 39.16% of Taiwan's total population of 23571227. The northern air quality region with four administrative districts depicted in Fig. 1 Fig. 3 which clearly identify the first 12 hr exhibits the higher PM 2.5 concentration with highest value of 57.8 μg m -3 .

Ordinary Kriging
Many natural phenomena are spatially or temporally variable; kriging usually aims to analyze a given point along a spatial or temporal scale. These quantities of these phenomena were often described using random variables. They may vary both spatially and temporally and form a random field. The random variables in different spatial or temporal fields are not independent and tend to be correlated. These physical quantities, represented by random variables, are called regionalized variables. Geostatistics, based on the theory of regionalized variables, analyzes spatially distributed, random, and structural natural phenomenon. In other words, geostatistics contains two analyses of spatial structure and interpolation estimation: the first aims to build spatial correlations by utilizing a semi-variogram for analysis of the random field, while the second aims to interpolate the data using a correlated spatial structure (Cheng et al. 2003;Chen et al. 2008). The detailed theory of ordinary kriging is described in the appendix. In particular, the semi-variogram obtained from the observations of Z(x) is an experimental semi-variogram, formed by several discrete points. It can be used for the kriging method only after fitting a continuous model. The continuous model is defined as a theoretical semi-variogram and must satisfy the conditional semipositive definite. Four theoretical semi-variograms are often applied, as follows Eqs. (1) -(4): a. power model In Eqs.
(1) -(4), ~ and a are model parameters, and h is the distance between any two points in space. The semi-variogram is the essential core of the kriging method, because it represents spatial variation and continuity. ( ) h c will increase as the distance h increases and approach a sill. At this distance, h will be called the influence range. The value of the influence range for the spherical, exponential, and Gaussian models is a, 3a, and a 3 , respectively. In addition, the fit of a theoretical semi-variogram is selected to optimize model parameters using numerical methods. In order to assess the suitability of the selected model and the reasonableness of the assumption for kriging estimation, the cross-validation scheme is used to check the criteria of unbiasedness and coherence. If the average error between estimation and observation (Kriging Average Error, KAE) approaches zero, this means the estimator is unbiased and can be expressed as in Eq. (5). The ratio of expected value of square for kriging error and kriging variance (Kriging Reduced Mean Square Error, KRMSE) approaching 1 means the model is even better; this coherence is expressed in Eq. (6).

Indicator Kriging
Indicator kriging is a method of nonlinear geostatistics method, which is also known as nonparametric geostatistics. Journel (1983) proposed an indicator function transforming data to a binary distribution; for example, data lower (higher) than or equal to a certain cutoff value would be transformed to 1 (0), as expressed by Eq. (7). The transformed data will be distribution-free and consequently, can be used to solve the extreme values problem that commonly emerges in natural phenomenon. The conditional probability of each  class can be estimated by simple kriging, and the estimate for the value of a random variable can be expressed by medium or maximum probability and expected value (Bierkens and Burrough 1993;Chiang et al. 2014).
z c is a cutoff value in the above equation. For example, if z c is a critical value for pollution and non-pollution, then the indicator value can be set as 1 for concentrations lower than or equal to z c and set as 0 for concentrations larger than z c . The probability of finding a value lower than z c in space P Z Z c # , using kriging, can be expressed as Eq. (8).

Dust Storm Semi-Variogram Analysis
The geostatistical model is a method for investigating the characteristics of a spatial structure. It uses the amount of variability between data points, calculates the weights between known observation stations and unknown points, and uses observations from known stations. Values and weights can be used to derive estimates for unknown points. The semi-variogram function is the core of geostatistics. Because the critical variation value of the semi-variogram function represents the variation of the random variable [Z(x)], the dimensionless variable must be used if the data has a proportional effect. The semi-variogram function divides the fine aerosol concentration at each station by the standard deviation of the PM 2.5 concentration at all stations [σz], i.e., [Z(x)/σz]. In order to calculate the semi-variogram function, this study isolated the impact and spatial variability of PM 2.5 concentrations in different regions. Average of 48 hr data of 26 stations in northern air quality region on 24 -25 March 2012 was used to compute the semi-variogram. Table 1 shows that the theoretical semi-variogram obtained by the exponential model had the smallest sum of squared error (1.084) and an R 2 coefficient of 0.316. For the Gaussian model, the sum of squared error was 1.088, and the coefficient of determination was 0.312. The theoretical semivariogram obtained using the power model had the largest sum of squared error (1.356) and the lowest R 2 coefficient (0.145). In addition, the cross-validation method was used to test the suitability of each semi-variogram model. Estimated average error between estimation and observation (Kriging Average Error, KAE) are used for unbiased verification of the model. KAE approaching zero indicates that the nonbiased effect is better. Optimization of the model is judged by the expected value (KRMSE) of the ratio between the estimated error and the kriging variation. The closer KRMSE is to one, the better the optimization effect of the model is.
The results of the cross-validation of the four different models are shown in Table 2.

Estimation and Analysis of Dust Storm Spatial Distribution
In order to obtain spatial distribution data for PM 2.5 concentrations in Taiwan's northern air quality region, this study performed an overall regional estimation using the geostatistics method of ordinary kriging for the region including Taipei City, New Taipei City, Keelung City, and Taoyuan City. To do so, the northern air quality region was first divided into 250 m grids: 4326 grids in Taipei City, 32742 grids in New Taipei City, 2059 grids in Keelung, and 19283 in Taoyuan, for a total of 58410 uniformly distributed grid points. The index model was then used to estimate PM 2.5 concentrations during dust storms. Since PM 2.5 has an extremely small particle size, it easily enters the body via respiration and affects human health. Therefore, in order to improve environmental quality and maintain human health, the Taiwan EPA had amended air quality standards, adopting 24-hr average for PM 2.5 air quality standards. The 24-hr average and annual averages are 35 and 15 μg m -3 , respectively. Therefore, this study will classify estimated PM 2.5 concentration values in the northern air quality region into 6 levels: level 1, PM 2.5 < 15 μg m -3 ; level 2, 15 -20 μg m -3 ; level 3, 20 -25 μg m -3 ; level 4, 25 -30 μg m -3 ; level 5, 30 -35 μg m -3 ; and level 6, > 35 μg m -3 .
Estimated results are shown in Fig. 4. The box-and-whisker plot of the estimated PM 2.5 concentrations in the administrative districts (Fig. 5) shows that the mean values for the four administrative districts are largely similar, with the highest concentrations occurring in Taipei City (27.95 μg m -3 ) and the lowest in Keelung City (26.91 μg m -3 ).

Risk Analysis of Pollution from Dust Storms
When dust storms affect Taiwan, what are the pollution conditions in our living environment? In which areas does the pollution reach health-threatening levels? When and where does one need to avoid going outdoors? In order to estimate and map air pollution risk during development of a dust storm, this study will use geographical statistics, which offer optimal linear unbiased estimation, in combination with indicators used in kriging technology. In addition to estimating the air pollution concentration in areas without measuring stations, we will also be able to estimate the probability of PM 2.5 exceeding certain threshold and construct an air pollution risk map.
In order to assess the risk of health hazards from dust storm pollution, this study uses indicator kriging methods to estimate the probability of exceeding certain pollution concentrations. During dust storms in the northern air quality region, pollution concentrations exceeded the annual average standard value of 15 μg m -3 . This means that the pollution concentration values in this study can be placed into three categories according to their level of risk: low risk when the critical value of PM 2.5 is greater than 25 μg m -3 ; moderate risk at concentrations greater than 30 μg m -3 ; and high risk at levels greater than 35 μg m -3 . In addition, in this study any area with a risk probability greater than 0.7 was considered to be a health hazards risk area. The estimation results are shown in Fig. 6.

Dust Storm Semi-Variogram Analysis
From Table 2, The KAE values were the best in the power model, then Gaussian, then exponential, and finally spherical model. However, the KRMSE values did not vary much between the models, although the spherical model was least optimal. The exponential model produced the better cross-validation results than other models. Overall, the Gaussian model was selected as the representative theoretical variation of PM 2.5 concentrations during a dust storm.

Estimation and Analysis of Dust Storm Spatial Distribution
The results showed that the minimum PM 2.5 value was 21.88 μg m -3 , the maximum value was 39.86 μg m -3 , the average was 27.24 μg m -3 , and the standard deviation was 1.01 μg m -3 . Dust storms refer to both sand and dust storms, in which strong winds lift large quantities of dust from the ground into the air, causing the air to be particularly turbid and horizontal visibility to be less than 1 km. Therefore, during the occurrence of dust storms, air quality is relatively poor due to high PM 2.5 concentrations. The PM 2.5 concentrations estimated by this study in the northern air quality region exceeded 21.88 μg m -3 , so no areas were assigned to Levels 1 and 2. Most of the areas (96.58%) were assigned to Level 4, and 1.32%, 1.84%, and 0.26% of the areas were assigned to Levels 3, 5, and 6, respectively. Among the four administrative districts, the greatest variation was in Taipei City, with 17.12 μg m -3 of total concentration value, followed by 13.50 μg m -3 in Taoyuan City, 11.33 μg m -3 in New Taipei City, and 3.59 μg m -3 in Keelung City. In addition, the maximum dust pollution value in the region occurred in the Zhongzheng District of Taipei City, where the PM 2.5 concentration was 39.86 μg m -3 ; the minimum value (21.88 μg m -3 ) was in the Xinzhuang District of New Taipei City.

Risk Analysis of Pollution from Dust Storms
The study showed that the majority of areas (86.95% of the total area) in the northern air quality region were areas with low risk of health hazards. This included 74.39% of the area of Taipei City, 85.10% in New Taipei City, 52.21% in Keelung City, and 96.62% in Taoyuan City. Areas with moderate risk accounted for only 1.18% of the total area: 11.14% of Taipei City, 0.11% of New Taipei City, and 0.89% of Taoyuan City. There were no moderate risk areas in Keelung City. Finally, only 0.42% of the regions had a high risk of health hazards. These occurred in Taipei City (3.01%) and Taoyuan City (0.60%). The study also found that the high risk areas were concentrated in Taipei City's Zhongshan and Songshan Districts and Taoyuan City's Zhongyuan District of Taoyuan.

CONCLUSION
This study collected data on PM 2.5 concentrations in Taiwan's northern air quality region associated with the dust storm events of 24 -25 March 2012. The study also aimed to analyze spatial distribution characteristics of daily average PM 2.5 concentrations using geostatistical methods. The results showed that when dust storms occurred, concentrations of PM 2.5 in the northern air quality region exceeded the annual average concentration of 15 μg m -3 . Thus, on a six-level scale, no areas fell within the first and second levels, and most areas (96.58%) were categorized as level 4. One point three-two percent of the areas were level 3, 1.84% were level 5, and 0.26% were level 6. The level 6 area appeared only in Zhongshan District of Taipei City and in Zhongmu District of Taoyuan City. In order to assess the risk of health hazards from dust storm pollution, this study used indicator kriging methods to estimate the probability of exceeding certain pollution concentrations. The results showed that most of Taiwan's northern air quality region was an area with low risk of air pollution health hazards; low-risk areas accounted for 86.95% of the total area. Most low-risk areas were located in Taipei City. Only 1.18% of the areas were medium-risk areas for health hazards: 11.4% of Taipei City, 0.11% of New Taipei City, and 0.89% Taoyuan City. High-risk areas occupied only 0.42% of the total area, including 3.01% of Taipei City and 0.60% of Taoyuan City, concentrated in Taipei City's Zhongshan and Songshan Districts and Taoyuan City's Zhongyuan District. This study examining the spatial distribution characteristics of PM 2.5 provides a risk evaluation approach using indicator kriging for the pollution of dust storm and also successful to delineate to construct an air pollution risk map and identify the probability of harm to health. The future work lies combing the seeking medical, place of residence information associated with the citizens in northern air quality region because the living and work area are not necessary associated with medical region. Due to very limited dust storm events in the recent years, the result of this study can only be used as the reference of spatial pollution information of dust storm exposure. However, by integrating the weather and atmospheric information, air quality regional forecasting still can be performed to provide public protection and enhance early warning systems for vulnerable groups, provide a reference for health risk and exposure assessment, and establish relevant standards for development of regulatory strategies.
By partially differentiating i m and n and letting the differentiation form be zero in Eq. (A6), the kriging system equation is obtained in Eq. (A7):  This can also be expressed as Eq. (A9):