Aerosol extinction coefficient and backscattering coefficient are the main optical parameters of aerosols, which can be used to further invert aerosol parameters, such as aerosol PSD, optical refractive index and so on. In this paper, the Fernald method is used to solve the extinction coefficient: α A z = − S A S M α M z + z 2 P z exp 2 S A S M − 1 ∫ z z c α M z ′ d z ′ z 2 P z α A z c + S A S M α M z + 2 ∫ z z c z ′ 2 P z ′ exp 2 S A S M − 1 ∫ z z c α M z ″ d z ″ d z ′ (1) where, z represents the detection height, α A z is the aerosol extinction coefficient at height z. α M z is the extinction coefficient of atmospheric molecules and may be determined by meteorological data or by standard atmospheric models. S M is the ratio of atmospheric molecular extinction coefficient and backscattering coefficient, and S A is the ratio of aerosol extinction coefficient and backscattering coefficient (namely, lidar ratio). In this paper the standard atmospheric model is used to calculate the extinction coefficient of atmospheric molecules.

In the Fernald method, the lidar ratio must be assumed. Recent observations suggest that lidar ratios range between 40 and 55 Sr In fact, Yinchuan area is located in the mid latitude region, in using the Fernald or Klett methods, in this paper, the lidar ratio of different wavelengths is selected as the constant value, namely, 40, 50 and 50 Sr for the wavelengths of 1,064, 532 and 355 nm, respectively (Mao et al., 2016b).

The color ratio can reflect the size characteristics of the aerosol particles and refers to the ratio of the backscattering coefficient of 1064 nm to the total backscattering coefficient of 532 nm, which is defined as: R C R z = β 1064 z β 532 T z (2) where, β _532T(z) represents the total backscattering coefficient, namely, the sum of the β 532 ⊥ z and β 532 / / z . The color ratio is positively related to the particles size, that is, the increase of the ratio with the particles size increase. Generally, for dust aerosols, the color ratio is between 0.3 and 1.5, and mostly gathers around 0.8, while for soot-type aerosols, the color ratio of the is mostly clustered around 0.35 (Twomey, 1977).

The relationship between the aerosol extinction coefficient and the PSD can be expressed by: α λ , z = ∫ r 0 r M K α r , m , λ , s n r , z d r (3) where n r , z represents the PSD with particle radius r at height z, which can be described by number, surface area and volume distribution, respectively. In Eq. 3, K α r , m , λ , s is the kernel function of the extinction efficiency factor, which depends on the particle radius r, the complex refractive index m, the wavelength λ, and the particle shape s. According to Mie scattering theory, aerosols can be approximately considered to be spherical particles, The kernel function K α r , m , λ , s is given by: K α r , m , λ , s = π r 2 Q m , λ , s (4) where, Q m , λ , s is the extinction efficiency factor.

Similarity, the relationship between the aerosol backscattering coefficient and the PSD can also be expressed by: β λ , z = ∫ r 0 r M K β r , m , λ , s n r , z d r (5) K β r , m , λ , s = π r 2 Q β m , λ , s (6) where, K β r , m , λ , s is the kernel function of the backscattering efficiency factor. For the kernel function, as r increases, the efficiency factor Q β is decayed in the oscillation and eventually tends to be a constant. Its oscillatory properties are closely related to the complex refractive index and have obvious differences under three different complex refractive indexes for cities, biomass combustion and desert aerosols, as shown in Figure 1. Figure 2 shows the relationship between the extinction efficiency factor and the particle radius at different wavelengths. Figure 3 shows the relationship between the backscattering efficiency factor and the particle radius at different wavelengths. From Figure 2 and Figure 3, the first main peak lies in between 0.1 and 1 μm, and as the wavelength λ increases, the position of the first main peak moves in the same direction with the particle radius r. Since the Yinchuan area is surrounded by four deserts, with high altitude, dry and dusty climate, it belongs typical semi-arid climate in northwest China, the complex refractive indexes m is chosen 1.55–0.01i to retrieve (Mao and Li, 2014). The extinction efficiency factors and backscattering efficiency factors.

The scattering that occurs when the size of the aerosol particles is close to or greater than the wavelength of the incident laser when the laser is transmitted in the atmosphere is called Mie scattering. The corresponding scattering and extinction efficiency factors for Mie scattering are given by: Q = 2 x 2 ∑ n = 1 ∞ 2 n + 1 R e a n + b n (7) Q β = 1 x 2 ∑ n = 1 ∞ 2 n + 1 − 1 n a n − b n 2 (8) where the n can be calculated by the following equation: n max = x + 4 x 1 3 + 2 (9)

According to the size and production mechanism, aerosols can be divided into two types: fine ( 2 r < 2 μ m ) and coarse ( 2 r ≥ 2 μ m ) particle patterns. In some functions, n r representing the number size distribution. Because the lognormal distribution can represent the size distribution of the dust grains quite well, a combination of two lognormal distributions is selected for describing the PSD: d n r d ⁡ ln r = ∑ i = f , c N t i 2 π 1 / 2 ⁡ ln ⁡ σ i exp − ln ⁡ r − ln ⁡ r i n 2 2 ln ⁡ σ i 2 (10) where N t i is the total particle number concentration of the ith mode, r i n is the median radius, and σ i is the geometric standard deviation. The subscript i = f , c refers to the fine mode and coarse mode, respectively. Therefore, two modes have a total of six unknown parameters. Among them, r i n and σ i have a small range, but the scope of N t i is particularly large, its search scope needs to be narrowed.

Another parameter, the aerosols particle effective radius, can be used to describe the physical properties of the aerosol particle, that is, judging the proportion of coarse and fine particles increases or decreases, which is defined as: r e f f = ∫ r min r max r 3 n r d r ∫ r min r max r 2 n r d r (11)

GA can find optimal solutions by simulating the “survival of the fittest” principle during natural evolution. In practice, the simple genetic algorithm (SGA) has been popular for the optimization problem application. However, SGA has many disadvantages. For example, SGA is easy to sink into local optima, and its search efficiency is very low. SGA has only one population and all operations aim at this single population.

In contrast, MPGA can use multiple subpopulations to replace a single population of SGA and perform their respective search processes in feasible solution domains. Each subpopulation can evolve in parallel according to different evolutionary strategies and genetic operations. Compared with SGA, MPGA has some obvious advantages. For example, to maintain the evolutionary stability of the best individuals, MPGA can accelerate the evolutionary speed and outcome the premature convergence of SGA to some extent. However, because of the large range of N t f , and N t c in Eq. (12), the MPGA algorithm cannot avoid sinking into local optima.

As mentioned above, the difference between the theoretical optical parameters at each wavelength and the optical parameters obtained by the experimental data can be defined as objective functions and the total sum of errors is as a objective function for minimization. The advantage of using the total error as the objective function is that the algorithm is simple, but the disadvantage is that the impact of the error value at each wavelength on the total error cannot be taken into account. This problem exists in SGA and MPGA. Therefore, if the error at each wavelength is taken as an objective function and the contribution of each error to the objective function is considered, the multi-objective genetic algorithm is undoubtedly the most appropriate. In this paper, the multi-objective genetic algorithm and successive approximation method are selected, which can avoid premature convergence and have significant superiority in optimizing multiple-objective functions.

For multi objective genetic algorithm, the existence of multiple objectives in a problem will in principle produce a group of optimal solutions, rather than a single objective optimal solution. In this paper, the NSGA-II is used to invert the aerosol PSD. The unknown parameters are randomly obtained in a large solution space, which yields different n r , and the Mie theory is used to calculate the aerosol extinction coefficient and backscattering coefficient. The difference between the theoretical value and the measured value is calculated, and the square of the difference is taken as the objective function value of the NSGA-II algorithm. After multiple generations of evolution, the optimal solution was calculated. Because the parameters range of N c and N f is too large, to improve the accuracy of the NSGA-II algorithm, the smaller ranges can be obtained by the successive approximation method. The detail steps can be given as follows:

(1) The parameter ranges of the three models are used to retrieve the PSD. The 40-bit binary string is used to replace the five real parameters, and the binary string is decoded to obtain the parameters within the value range. These unknown parameters σ f , N t f , r f , N t c , and r c and a certain parameter ln σ c are substituted into Eq. 6 to obtain n r .

Figure 4 shows the flow chart of the aerosol PSD inversion algorithm based on the NSGA-II algorithm.

The Fernald method is used to retrieve aerosol extinction coefficient and backscattering coefficient at 355 and 532 nm and backscattering coefficient at 1064 nm.

The value range of the six parameters σ c , σ f , N t f , r f n , N t c , r c n that need to be inverted are listed in Table 1, where N t i is the total particle number concentration of the ith mode, r i n is the median radius, and σ i is the geometric standard deviation, and the subscript i = f , c refers to the fine mode and coarse mode. Because σ c has been identified in the three modes, so we only need to reverse the five parameters.

Moreover, the initial parameter setting of the multi-objective genetic algorithm are listed in Table 2.

To verify the feasibility of NSGA-II algorithm for retrieval of aerosol PSD, some simulations, in which the different errors are superimposed on the input optical data, were carried out. A group of data from the parameter ranges of the three models listed in Table 1 was selected to obtain the PSD, which acted as a real value. The group of data is substituted into Eqs 3, 5 to calculate the aerosol extinction coefficient and backscattering coefficient. While a certain error is added to the calculated extinction coefficient and backscattering coefficient to invert the PSD, which is acted as a calculated value, respectively. The inversion average error of the PSD is defined as: δ = ∑ r = 0.01 5 n a r − n I r n a r 500 (12) where r is the aerosol particle radius, n a r is the real PSD value, and n I r is the calculated PSD value.

As mentioned above, Yinchuan area is located in the northwest of China, surrounded by four deserts. In the west, north and east of Yinchuan area, there are Tengger Desert, Badain Jaran Desert, Ulanbuhe Desert and Maowusu Desert. Yinchuan area belong to the typical continental arid and semi-arid climate, there are four distinct seasons, and there are many sand and dusty weather in spring, with an average annual rainfall of 200 mm. A Multi-wavelength lidar has been devolved and located in Teaching Building 17 at North Minzu University (38°29N, 106°06E) in Yinchuan, China (Mao et al., 2016b). Figure 9 shows the schematic diagram of multi-wavelength lidar system. Table 5 list the parameters of the lidar. To verify the feasibility of using the proposed algorithm, some experiments have been conducted since March 2021.

Figure 10 shows the extinction coefficients at wavelengths of 1,064, 532, and 355 nm at 21:22 Beijing time on 7 April 2021, as well as the Ångström exponent distributions, aerosol PSDs inverted from the NSGA-II algorithm, effective radius profile of the aerosol particles. It was a sunny on that day.

In Figure 10C, the Ångström index at 532/355, 1,064/355, and 1,064/532 from 2 to 2.6 km, 3.3–4 km all increase, so under the lower troposphere, the number of fine particlesincreases with height, and the number of coarse particles decreases with height.

In Figure 10D, when the particle radius is 0.03 < r < 0.1 μm, the aerosol PSD at 2.6 km described by the green curve above the red curve indicates that the number of fine particles is greater than that at 2 km, which is consistent with the conclusion obtained according to the Ångström index above. When the particle radius is 0.07 < r < 0.3 μ m, the aerosol PSD at 4 km of the orange curve above the black curve is indicated the number of fine particles is greater than that at 3.56 km. At r > 1 μ m, the PSDs at 3.56 km and 4 km are less than those at 2, 2.4, 2.6, and 2.8 km, respectively, indicating that the number of coarse particles at 2, 2.4, 2.6 and 2.8 km are all greater than those at 3.56 km and 4 km, which is also consistent with the conclusion.

In Figure 10E, it is obvious that the effective radius of the aerosol particles from 2 km to 4 km is constantly decreasing, which is consistent with the increasing trend of the Ångström index.

Figure 11 shows the extinction coefficients at wavelengths of 1,064, 532, and 355 nm at 21:43 Beijing time on 22 August 2022, as well as the Ångström exponent distributions, aerosol PSDs inverted from the NSGA-II algorithm, effective radius profile of the aerosol particles and depolarization ratio profile at different heights. It was cloudy on that day.

The type of cloud can be preliminarily judged according to the depolarization ratio. Sassen et al. found that the depolarization ratio of water clouds was less than 0.15, and the depolarization ratio of ice-water mixed clouds was between 0.15 and 0.5 (Sassen et al., 1992). From the depolarization ratio profile shown in Figure 11F, there are water clouds at 2.8–4 km. Figure 11C also shows the increasing 1,064/355 Ångström exponent at 2.8–4 km, indicating the coarse particles increase. This further confirms the presence of water clouds at 2.8–4 km. In Figure 11A, the extinction coefficient at 3.8 km is the largest, and in Figure 11D, the number concentration of aerosol at this height with radius greater than 0.2 μm is also the largest (similar to 2 km). From Figure 11E, the effective radius of aerosol particles decreases firstly and then increases significantly after the cloud layer. The depolarization ratio profile, extinction profile and Ångström exponent profile all reflect the accuracy of the inversion trend of the effective radius of aerosol particles.

Figure 12 shows the extinction coefficient, Ångström exponent, PSDs, and effective radius of aerosol particle under dusty weather conditions on 15 March 2021. Due to the dusty weather, in Figure 12C, the Ångström exponent of 1,064/355 and 1,064/532 is decreasing at height 3–3.5 km, while in Figure 12E the aerosol particle effective radius at same height is increasing, which is consistent with Ångström exponent.

Figure 13 shows the extinction coefficient, Ångström exponent, PSD, effective radius and color ratio of aerosol particles under dust weather at 02:20 on 16 March 2021. By comparison, the backscattering coefficient of 1,064 nm in Figure 13 is greater than that in Figure 11, and the number concentration of aerosol particles in Figure 13 is much more than that in Figure 11 Moreover, in Figure 13D, the color ratio above 1.6 km ranges from 0.4 to 0.8, consistent with the range of 0.35–1.5 observed for dust aerosol. In Figure 13F there is a small rise at the height of 2.8–3.2 km. In Figure 13C, there is a significant decrease for 1,064/532 Ångström exponent. In Figure 13E, the effective radius of the aerosol particles at the height of 2.8–3.2 km increases slightly.

If the five objective functions are linearly weighted and a comprehensive effect function is used to represent the objective of the overall optimization, the solution corresponding to the optimal effect function is considered the optimal solution of the problem, so that the multi-objective optimization problem is transformed into a single objective optimization problem. The MPGA and SGA inversion algorithms use this model. The feature of this model is that the RE of the overall objective function is small, but it is not guaranteed that the relative REs of the five optical parameters are small at the same time. In the actual inversion process, the RE of extinction coefficient at 1,064 nm is always too high. However, the NSGA-II algorithm can balance the RE of the five optical parameters to a certain extent.

For better comparison, the best settings of the three intelligent algorithms listed in Table 2 are used. Figure 14A shows the PSDs at 1.6 km retrieved by the three intelligent algorithms on 8 April 2021. In Figure 14A, the number concentration inverted from the three algorithms are close to each other at a radius greater than 1 μm, but there are great differences at a radius from 0.02 to 1 μm. In Figure 14B, the RE of the 1,064 nm extinction coefficient is relatively large, while the REs of the extinction coefficient of 355, 532 nm and the backscattering coefficient of 355 nm are particularly small. This inversion result is not what we need. Generally, if the RE of all objective functions are close to each other, that the inversion result will be closer to the real value.

In this paper, the dispersion degree of error is taken as the main evaluation reference, including the objective function J (itself is the absolute error with the experimental observation value), the RE, the range of the RE (ED), the standard deviation (SD) and the dispersion coefficient (COV). J = ∑ i = 355,532,1064 3 α T i − α A i + ∑ i = 355,532 2 β T i − β A i (13) R E = J / ∑ i = 355,532,1064 3 α T i + ∑ i = 355,532 2 β T i (14) R E i = J i / α T i , i = 355,532,1064 J i / β T i , i = 355,532 (15) A r r a y = α T i − α A i , α T i − α A i , α T i − α A i , α T i − α A i , α T i − α A i (16) E D = max A r r a y − min A r r a y (17) S D = ∑ i n R E i − R E i ¯ 2 n , n = 5 (18) C O V = S D / R E i ¯ (19) where, α T i and α A i are the actual extinction coefficient measured by lidar and the theoretical values calculated by Mie theory, respectively, β T i and β A i are the actual backscattering coefficient measured by lidar and the theoretical value calculated by Mie scattering theory, respectively.

Here, the ED, SD, and COV are all used to evaluate the dispersion degree of the overall distribution of the RE. Table 6 lists the REs of the five optical parameters for the three inversion algorithms. The dispersion degree of the RE of MPGA and SGA is much larger than that of NSGA-II in terms of ED, SD, and COV. The REs of the five optical parameters of NSGA-II inversion are relatively closer. Table 7 lists the RE of the total objective function linearly accumulated by the five objective functions, it can be seen that the total RE of the three inversion algorithms is not much different, compared with SGA and MPGA, NSGA-II algorithm has a smaller RE. Therefore, NSGA-II has certain advantages over SGA and MPGA in the inversion of aerosol PSD.