Dust Particle Size Distribution Inversion Based on the Multi Population Genetic Algorithm

The aerosol number size distribution is the main parameter for characterizing aerosol optical properties and physical properties, it has a major influence on radiation forcing. With regard to some disadvantages in the traditional methods, a method based on the multi population genetic algorithm (MPGA) is proposed and employed to retrieve the aerosol size distribution of dust particles. The MPGA principles and design are presented in detail. The MPGA has better performance compared with conventional methods. In order to verify the feasibility of the inversion method, the measured aerosol optical thickness (AOT) data of dust particles taken by a sun photometer are used and a series of comparisons between the simple genetic algorithm (SGA) and MPGA are carried out. The results show that the MPGA presents better properties when compared with the SGA with smaller inversion errors, smaller population size and fewer generation numbers to retrieve the aerosol size distribution. The MPGA inversion method is analyzed using the background day, dust storm event and seasonal size distribution. The method proposed in this study has important applications and reference value for aerosol particle size distribution inversion.


INTRODUCTION
Aerosols have significant influences on the radiation budget of the Earth's atmospheric system.On the one hand, aerosols can directly perturb the Earth's radiation balance through scattering and absorbing solar radiation.On the other hand, by acting as cloud condensation nuclei, aerosols can change the optical properties, cloud cover and cloud lifetime, influence the cloud formation process and have indirect effects on climate change.By absorbing solar radiation aerosols can affect and change the structure and stability of the atmospheric temperature, which is a semi-direct effect on climate change.The aerosol concentration and size distribution are the main parameters for the characterization of aerosol optical properties, physical properties and important physical radiation quantity.The particle concentration can intuitively reflect the temporal and spatial distribution characteristics of aerosols.The size distribution can affect the scattering phase function, optical thickness, single scattering albedo and other optical parameters, and therefore further affect aerosol radiative forcing (Zhou 2012).
In the sun photometer technique the aerosol optical thickness (AOT) is taken from solar direct irradiance, and using AOT inversion the columnar size distribution of the aerosol can be inferred.The relationship between the AOT and size distribution can be described using the first kind of Fredholm integral equation, but the equation is ill-posed due to the oscillatory nature of the Mie kernel.In fact, for the first kind of Fredholm inversion integral equation, the Phillips-Twomey inversion algorithm, which was proposed originally by Philips and Tikhonov and perfected gradually by Twomey (Twomey 1977), has been widely applied in particle size distribution estimation by many researchers (Heintzenberg et al. 1981;Pahlow et al. 2006;Arias and Frontini 2007).For example, King et al. (1978) employed the Phillips-Twomey method to retrieve the particle size distribution through measuring the AOT.The AOT retrieval from direct sun measurements is relatively simple.However this inversion method has some disadvantages, for example, the aerosol size distribution is not a smooth function and the change in the distribution shape is very large.If the smooth constraint is unreasonable, the Phillips-Twomey method will lead to failure.The Phillips-Twomey method requires a priori assumptions about the nature of the aerosol size distribution.The uniqueness of the solution based on the Phillips-Twomey method is also questionable due to the non-linearity of the problem and the local nature of linearized solutions (Lienert et al. 2001).
Also the sun photometer uses diffuse sky radiation measurement data to retrieve the aerosol size distribution.Nakajima et al. (1996) used the sky diffuse almucantar or principal plane measurement to retrieve single scattering albedo, phase function, refractive index and volume distributions.Li and Mao (1989) proposed a new approach for retrieving the size distribution by measuring the direct solar radiation and sky light on the solar almucantar.Qiu et al. (1983) combined solar spectral extinction and forward angular scattering observations to obtain the size distribution.
Other methods were also developed to retrieve the particle size distribution, including least squares optimization algorithm (Qiao and Zhang 2004;Su et al. 2005), Powell optimization algorithm (Liu and Zheng 1998), neural network algorithm (Guardani et al. 2002), dynamic system identification (Zheng et al. 2002), and simulated annealing algorithm.These methods have been widely used in the measurement of air particles and liquid suspensions.
In 1967, the genetic algorithm (GA) concept was first proposed by Bagley (1967).In 1975, systematic research into the GA theory and method was implemented by J. H. Holland of Michigan University (1975).The GA was previously used to illustrate the self-adaptation process of natural and artificial systems.At present the GA has been widely applied to many research fields, such as machine learning, adaptive control, pattern recognition, image processing, and social sciences.
The GA was also was used in aerosol size distribution inversion (Lienert et al. 2001(Lienert et al. , 2003)).Lienert et al. (2001) described a genetic method for deriving aerosol size distributions from the multi-wavelength extinction measurements and employed the genetic inversion method to obtain particle size distributions from multi-angle optical scattering cross-section data measured by a polar nephelometer at a wavelength of 0.532 mm in 2003.Lienert et al. (2003) used the simple genetic algorithm (SGA) to complete their work.The SGA has some disadvantages, for example, SGA easily produces convergence premature problem, and the search efficiency is very low in the later stage of the evolution process.In view of this, the novel multi population genetic algorithm (MPGA) is proposed to replace the conventional SGA.
In this paper the MPGA is employed to retrieve the dust particles size distribution using AOT data taken by a sun photometer CE-318.The aerosol size distribution principle based on the MPGA is discussed.Some experiments, including comparisons with SGA method, were carried out during a background day, dust storm event and different seasons, and analyzed to verify the feasibility of the MPGA inversion method.

Aerosol Size Distribution
According to the relationship between aerosol size distribution and AOT, the AOT is defined as (Li et al. 1995) ( , ) ( , , ) ( , , ) ( ) where ( , , ) K r m m is the extinction efficiency factor, r is particle radius, r 0 and r M are the lower and upper limit of particle radius, respectively, n(r, z) is the aerosol size distribution at height z, ( ) ( , ) n r n r z dr is the aerosol size distribution of entire atmosphere column, and z M is the top height of atmosphere column.
In Eq. ( 1) the kernel function is πr 2 n(r, z)K (r, λ, m).According to the Mie scattering theory, with the increase in r, the extinction factor K decays oscillation and finally tends to constant 2. The oscillation characteristic is closely related to the refractive index.
In fact the size range of atmospheric aerosol particles is very large and it covers five orders of magnitude from 0.001 -100 μm.Generally, according to their sizes and formation process, aerosol particles are categorized into three types, coarse particles with diameters larger than 2.5 μm and smaller than 10 μm, the fine particles with diameters between 0.1 -2.5 μm, and ultrafine particles with diameters smaller than 0.1 μm.In this study the coarse particles include mainly dust produced from deserts surrounding the Yinchuan area, dust stirred up by vehicles traveling on roads and dust from the soil surface brought into the air by the wind.The fine particle sources include vehicular emissions, fossil fuel combustion from residential cooking and heating.
The size distribution n(r) can usually be described by several functions, such as the Junger distribution, Diermedjean distribution, lognormal distribution and log quadratic curve distribution and so on.In fact, the real aerosol distributions can be well represented by a combination of three lognormal distributions, which can be written by (Lienert et al. 2001 where N i , r mi , and i v are integrated number concentration, median diameter and geometric stand deviation, respectively, subscript i represents the three aerosol modes, respectively.According to the number concentration, the area distribution can be also given by ( )

Principle of MPGA
As an artificial intelligent algorithm the GA simulates the selection, crossover and mutation process of survival of the fittest discipline in the natural world.Its basic principle is as follows: starting from an initial population including some random individuals; after a series of processes, such as random selection, crossover and mutation operations, new individuals that are more suitable for the environment are produced.This permits the population to enter a better region of the search space.Through continuous evolution the solutions will converge to a final population that is the most suitable for the environment, where, the optimal solution will be obtained.
Compared with SGA the MPGA uses multiple sub populations to replace a single population.Each sub population uses its own different evolution strategies and genetic operands to perform parallel evolution.For example, each sub population has its own crossover probability and mutation probability.These different sub populations evolve independently.When all sub populations perform some evolution generations, the current most optimal individual will be assigned into each sub population in order to promote the evolution of each sub population.Under the condition of maintaining the stability of the best individual, this process can select and retain the best individual of each sub population and speed up the evolution process.The MPGA can also avoid convergence premature of single population evolutionary process (Deng et al. 2008).The MPGA has several characteristics as follows (Shi et al. 2011): (1) The MPGA breaks through the SGA framework in which only one population performs the genetic evolution and introduces multiple populations to perform optimal search in parallel.Different control parameters are used in each population to obtain different search objectives.
(2) By using the immigration operator these populations can link to each other to achieve the co-evolution of multiple populations.Therefore the most optimal solution is the combined result of multiple co-evolving populations.
(3) The artificial selection operator preserves the best individual in each population and uses this as the convergence criteria.
Figure 1 shows a schematic diagram of the MPGA algorithm structure (Shi et al. 2011).

MPGA Design
In the MPGA design chromosome encoding is the first problem to be considered.In Eq. ( 2), N i , r mi , and i v (i = 1, 2, 3) are all unknown parameters.Therefore there are nine unknown parameters that need to be determined by the MPGA.The binary encoding is a conventional method.In this paper a binary string of 10-bits in length is selected to represent one parameter, therefore a binary string 90-bits in length can denote the nine unknown parameters mentioned above.For example, for a binary string of length of 90-bits (0 1 0 0 1 0 1 0 1 1 0 1 1 0 1 1 1 1 1 1......... 0 0 1 0 0 1 1 0 1 0), the first 10-bits represents the N 1 , the second 10-bits represents r m1 , the 3 rd 10-bits represents 1 v , and so on, the 9 th 10-bits represents 3 v .In general, in order to convert a binary number to decimal number, decode formula can be given by The AOT taken by the sun photometer CE-318 is the measured value Ti x .By substituting Eq. ( 2) into (1), the computed AOT value Ai x can be obtained.The design goal of the MPGA is to decrease the error between the computed AOT value and measured AOT value.Therefore the error sum of square is selected as the objective function, which is written by ( ) , ,..., where the subscript i represents the measurement wavelengths of the sun photometer CE-318, which are 1640, 1020, 870, 670, 500, 440, 380, and 340 nm, respectively.Therefore, the fitness function can be defined as Figure 2 shows a flow chart of the aerosol particles size distribution inversion algorithm based on the MPGA.In this algorithm the random values of unknown parameters N i , r mi , and i v are given first, respectively, and according to these random values the aerosol particles size distribution and calculated AOT values are retrieved.After calculating the error between the calculated AOT and measured AOT, the fitness value is obtained from the calculation error.The MPGA then makes use of the fitness value to optimize the unknown parameters.The unknown parameters are iteratively updated continuously, until the evolution generation reaches the set value.The optimization calculation process is terminated and the final optimized values for the unknown parameters are derived.The final optimized values are used to obtain the aerosol particles size distribution.In the inversion algorithm the measured AOT are used as the goal value.Through repeated iteration operation the calculated AOT will be consistently close to the goal value by optimizing unknown parameters.

THE EXPERIMENTS AND RESULTS ANALYSIS
Yinchuan City is a medium sized city that lies in eastern region of Northwest China.It has the typical temperate continental climate.There are four large deserts, the Badanjilin desert, Wulanbuhe desert, Tenggeli desert and Maowuso desert, located in the northwest, west and east of the Yinchuan area, respectively.Yinchuan is an important transportation channel and the source of desert dust in China.Beifang University of Nationalities (106°06'E, 38°29'N) is located west of Yinchuan City.West and North of the campus are the famous Three-North Shelterbelt and Gobi deserts at the eastern foot of Mt.Helan.The geographic details of Yinchuan and the surrounding regions were introduced and can be seen in Mao et al. (2014).A sun photometer CE-318 was installed on teaching building No. 17 of the campus and operated from September 2012.Every day the sun photometer automatically measured the direct solar radiation from 8:00 -17:00 Beijing time (except in situations when the sun was obscured by clouds, rain, snow or dust).The AOT temporal and spatial evolution data and aerosol size distribution for two years in the Yinchuan area were obtained.These data have important scientific value for quantitatively investigating the aerosol influence on atmospheric physics and chemistry and monitoring the dust and air pollution in the Yinchuan area.
A combination of lognormal distributions for three aerosol modes can reasonably describe the real aerosol size distribution.There are nine unknown lognormal distribution parameters that need to be determined by MPGA.Table 1 shows the variation range for the nine unknown input parameters.
Some comparisons between the SGA and MPGA were performed.In order to facilitate comparison the two algorithms were setup with the most optimal status, namely the setting options and parameters values for the MPGA and SGA, such as population size, crossover probability, mutation probability, fitness function and unknown parameters range were determined experimentally.Table 2 lists the setting options and parameters for the MPGA and SGA.It should be pointed out that the complex refractive index used is 1.5 -0.01 i, which is a reasonable value for dust particles.

No
Figure 3 shows the fitness function variation value for the best chromosome under different evolution generation times based on the SGA and MPGA, respectively.It is clear that the maximum MPGA fitness value is far larger than that of SGA with fewer generation times, which shows that the MPGA can solve the convergence premature problem in the single population SGA evolution to a large extent.
Figure 4a shows the size distributions retrieved by the SGA and MPGA, respectively.Figure 4b shows the inversion error between the computed AOT value and measured value at different wavelengths using the SGA and MPGA, respectively.In Fig. 4, the measured AOT value is (0.1650, 0.194, 0.200, 0.218, 0.270, 0.309, 0.359, 0.380), the computed AOT values retrieved by SGA and MPGA are (0.1763, 0.1998, 0.1977, 0.2169, 0.2753, 0.3061, 0.3462, 0.3819) and (0.1644, 0.1944, 0.1994, 0.2207, 0.2742, 0.3067, 0.3503, 0.3883), respectively.It is evident that the MPGA inversion error is less than that for the SGA.Moreover, Fig. 4 shows that there is one obvious peak at radius 1.0 μm. Figure 5 shows the area distributions using the data in Fig. 4 based on the MPGA.It can be seen that two obvious peaks are observed where the particle size is at 0.1 and 1 μm, respectively.Figure 6 shows a set of size distributions which employ the same measured AOT in Fig. 4 as retrieved by the MPGA.It is clear that the MPGA also has a disadvantage, namely, the inversion solutions are not unique, which is the general problem of many size distribution inversion methods.However the differences in inversion solutions are not large.In the radius range less than 0.5 μm the differences are very small.
In order to further verify the retrieval property of the two methods, five groups of measured AOT data were randomly selected to compute the inversion errors for different wavelengths and standard deviations using the SGA and MPGA, respectively.The results are listed in Table 3.It is obvious that both the MPGA errors and standard deviations are less than those for the SGA.Also from the previous comparison in Fig. 4 and Table 2, the MPGA has other better attributes than the SGA, such as smaller population size and fewer generation times.Figure 7 shows the correla-tion between the measured AOT value and computed AOT value obtained by the MPGA.There are hundred groups of data used in the correlation calculation.The correlation coefficient is up to 0.9998.
Figure 8 shows the AOT and particle size distribution based on the MPGA on 31 October 2012.On that day the weather was sunny and cloudless.In the morning the AOT value of 1640 nm is only 0.0665, which can be viewed as a background value for this area.With the increase of human activities, the aerosol content in the atmosphere gradually accumulated.At 16:20 the AOT value reached the maximum value.After 16:30 the AOT value started to gradually decrease.The size distributions were obtained from the AOT value at five different times.In the morning, the number concentration was lower.With the accumulation of aerosol particles, the number concentration gradually increased.After arriving at the maximum the number concentration started to decrease.However in the radius range from 0.1 -1 μm the number concentration variation trend was consistent with the AOT value.This suggests that the aerosol particles with radius ranging from 0.1 -1 μm are the major factors contributing to the variation in AOT value.In this radius range both dust particles and anthropogenic aerosols can cause the increase in AOT value.As the weather was clear and calm and the dust particle contribution was very small, the afternoon AOT value increase was caused mainly by anthropogenic aerosol particles, produced mainly from vehicular emissions during rush hour, fossil fuel burning from residential cooking and heating and dust stirred up by vehicles traveling on roads.
Figure 9 shows the AOT and particle size distributions based on the MPGA inversion method during a dust storm event (on 8 February 2013).From morning, the AOT value gradually increased until 11:25 and the AOT value at 500 nm increased to 1.07 when the dust storm event occurred.Before the dust storm, both the AOT value and aerosol particle concentration remained at a relatively low level.With the increase in AOT value the aerosol particle concentration started to increase.After the dust storm event both the particle concentrations at 11:25 and 12:10 were significantly higher than before the dust storm.In this case, the study result shows that the concentration of coarse particles increased due to the dust storm event and the ultrafine particle concentration decreased due to coagulation with the coarse particles.This result is in agreement with the results of Wu et al. (2009) and Hussein et al. (2011).Past field measurement studies indicated that anthropogenic aerosols can be transported along with dust particles during a dust storm, resulting in a bi-modal distribution (Wu et al. 2009).As can be seen from Fig. 10 the particle size distribution was also bimodal during the dust storm.For example, there were two peaks at radius 0.1 and 0.7 μm, respectively, for the curve at 12:10.This was consistent with the result obtained by Wu et al. (2009).
It should be pointed out that for the very big AOT     value obtained in the dust storm event there are some inversion differences and limitations in using the MPGA.For example in the boundary point calculation sometimes some distortions exist, resulting in the concentration being near zero at the 5 μm radius.However, in another dust storm case no distortions existed, such as on 8 February 2013.That suggests that for the number size distribution inversion in dust storm the MPGA algorithm has some disadvantages and further research needs to be made.
Figure 10 shows the season average AOT and particle size distribution and area distribution at different seasons from September 2012 to August 2013.It is obvious that the AOT is maximum in spring and minimum in autumn.As the weather is very dry in spring in the Yinchuan area, the soil is very loose, eolian dust is easily brought into the air under the influence of strong northwesterly winds.Dust storms, which are originated from the deserts located northwest of the Yinchuan area, also occurs occasionally.When a dust storm passes through the Yinchuan area the dust particle concentration increases significantly.
In Fig. 10b, the particle number concentration in spring is far larger than that in other seasons with the dust particle radius ranging from the 0.2 -2 μm.There is an evident peak located at 1 μm.This shows that the eolian dust particles with radii ranging from 0.2 -2 μm are the main aerosol component in spring, where the particle number concentration is significantly higher than in other seasons.It also shows that the dust particles in this radius range are the main factors that cause the increase in AOT value in the spring.
In addition, for area distribution the maximum peaks for all seasons are located in the range greater than 1 μm.In spring the area distribution is unimodal and its peak is located at 1 μm.In other seasons the area distributions are biomodal.These peaks are located at the range from 1 -3 μm and from 0.07 -0.2 μm, respectively.

CONCLUSION
In aerosol size distribution inversion the traditional Phillips-Twomey inversion algorithms have some disadvantages such as requiring very reasonable smooth constraints and a priori assumptions about the nature of the aerosol size distribution.As an artificial intelligent optimization algorithm the genetic method can overcome the Phillips-Twomey algorithm disadvantages to a large extent.This study used the measured AOT data taken by a sun photometer and the MPGA algorithm to retrieve dust particles size distribution.Compared with the SGA approach MPGA has better properties, such as smaller inversion errors, smaller population size and fewer generations.The MPGA inversion method successfully retrieved and analyzed the particle size distribution during background days and dust storm events.The method proposed in this paper has important reference value in particle size distribution inversion.

Fig. 1 .
Fig. 1.Schematic diagram of MPGA algorithm structure. b Fig. 2. Flow chart of aerosol particles size distribution inversion algorithm based on MPGA.

Fig. 3 .
Fig. 3. Variation of fitness function value of the best chromosome with evolution generation times based on SGA and MPGA, respectively.

Fig. 4 .Fig. 7 .
Fig. 4. (a) Particle size distributions retrieved by SGA and MPGA.(b) Inversion error between the computed AOT value and measured value at different wavelength using SGA and MPGA.

Fig. 8 .
Fig. 8. (a) The AOT, (b) number concentration based on MPGA during the background day on 31 October of 2012.

Fig. 9 .
Fig. 9. (a) The AOT, (b) number concentration and area distributions based on MPGA during the dust storm event on 8 February 2013.

Table 1 .
Variation range of the unknown parameters.

Table 2 .
Setting options and parameters of MPGA and SGA.

Table 3 .
Inversion errors of different wavelength and standard deviations by using of SGA and MPGA.