Aerosol Size Distributions from
Genetic Inversion
of Polar Nephelometer Data
B.
R. Lienert, J. N. Porter and S. K. Sharma
Hawaii Institute of Geophysics and Planetology,
Abstract
We show how genetic
inversions can be used to recover lognormal aerosol size distributions from
multi-angle optical scattering cross-section data measured by a polar nephelometer at a wavelength of 0.532 μm. The inversions can also
be used to recover the absolute calibration factor of the polar nephelometer. We demonstrate the method by applying it to
polar nephelometer data measured during the Shoreline Environmental Aerosol Study (SEAS) at
Introduction
Inversion of optical scattering
cross-section measurements for aerosol size distributions allows critical
optical parameters such as aerosol extinction, phase function and lidar ratio
to be predicted as a function of wavelength using Mie theory. However, due to
the oscillatory variation of the Mie scattering cross-sections with wavelength
and aerosol size, such inversions can be highly non-unique (Post, 1975).
Uniqueness can be represented by the width of resolution kernels (Backus
and Gilbert, 1970). These are localized weighted averages of a continuously
varying quantity (in this case aerosol size) that can be estimated by inverting
a finite set of measured data (in this case scattering cross-sections). As the
width of a resolution kernel is decreased, the error in its weighted average
increases, resulting in a “trade-off” between resolution and error (
In a previous paper (Lienert et al., 2001), we showed how
genetic inversion (Holland, 1975; Goldberg,
1989) could be used to recover multimodal lognormal aerosol size distributions
from multi-wavelength optical extinction data. The advantage of genetic
inversion is that it explores a wide range of size distribution parameters,
similar to the
In this paper, we use genetic inversion to investigate the inversion of multi-angle scattering, similar to that measured using sky radiances (Dubovik and King, 2000) or polar nephelometers (Barkey and Liou, 2001). Since multi-angle phase functions are also predicted by Mie theory, we suspected that similar problems of non-uniqueness could be present in the recovered distributions.
Description of the Polar Nephelometer
The
instrument design follows the system of
Numerical
Methods
In order to explicitly define the measured quantities and size distribution parameters, which are defined in a number of different ways in the literature, we have included a summary of the theory we used to calculate angular differential scattering cross section, dσi/dΩ, in the appendix. dσi/dΩi was calculated at each radius and wavelength using a program incorporating the routine of Bohren and Huffman (1983). The integral in eq. A10 was performed using quadrature integration (2000 points equally spaced in ln r from r=0.002 to 200μm). This number of points (400/decade) was the minimum found necessary to prevent the occurrence of small oscillations in calculated values of P1 and P2 (defined in the Appendix). The accuracy of the calculations was confirmed by calculating and comparing results for two models (Haze M and Cloud C-3) tabulated by Deirmendjian (1969).
To demonstrate how the Mie-scattering
phase functions vary with the type of size distribution, P1 is plotted in Fig. 2 for the three unimodal area distributions in Fig. 1 with mean radii
ranging from 0.0825 to 1.35μm.
We have plotted the size distributions as area, rather than number, as area is
used to calculate extinction. We used a value of 1.36+0i for the refractive
index of the aerosol particles, which is the value expected for a mass-weighted
sea salt-water solution (Tang et al., 1997) at 85% relative
humidity, the value measured at the time of the polar nephelometer
measurements (Masonis et al., this issue). It is apparent from Fig. 2
that the shapes of the resulting phase functions, as well as the values of the
lidar ratios, are critically dependent on particle size, particularly at
scattering angles of less than 30 degrees and greater than 130°. Note in
particular, the characteristic undulations between angles of 90° and 130°, as
well as between 140° and 180°, for the largest size distribution. The variation
between 140° and 180°, frequently termed a rainbow, is characteristic of large
spherically-shaped aerosols (Bohren and Huffman,
1983).
The phase functions calculated in Fig. 2
are for Mie scattering only. The measured scattering will be a combination of
Mie scatter and molecular scatter. For a laser polarized at right angles to the
plane of scattering, the molecular scatter, given by eq.
A14, is independent of angle. Its contribution to the total scatter will
therefore depend on the relative size of Mie scatter.
This is illustrated in Fig. 3, where we have plotted total scatter normalized
by integrated Mie scatter (i.e., the extinction for
zero absorption) for one of the models in Fig. 1 with three different lognormal peak amplitudes
giving values of scattering ranging from 1.7x10-5 to 4.4x10-4
m-1. The depth of the “valley” between 40° and 140° gives an
approximate indication of size of the aerosol extinction relative to the
molecular scatter. Since the molecular scatter is known, it should be possible
in principle to invert the combined scattering data for the polar nephelometer
calibration factor, as well as for the aerosol size distribution causing Mie scattering.
The genetic inversion
algorithm that we used is that described by Carroll (1996) and used in the previous paper by Lienert et al. (2001). Briefly, a
random number generator is used to select the parameters for a lognormal size
distribution having the form of eq. A5. The
scattering predicted by Mie theory is calculated at each measurement angle and
the "fitness" of this distribution is calculated from its mean
absolute error of fit to the measured data. Two sets of lognormal parameters
are then concatenated into binary gene strings (parents), which are used to
breed child solutions using randomly selected sections of each parent gene.
Successive generations of children are then used to maximize the fitness. In
this paper, for parameters such as the mutation fraction, population size,
number of children, etc, we used the values given in Table 1. Note that we
enabled the elitism feature, which ensures that the individual having the best
fitness is always retained in subsequent iterations. In our previous paper, we
defined the fitness as a specified standard deviation divided by the error of
fit. In this paper, we use a different definition of fitness, which is similar
to the inverse of a misfit function used in simulated annealing (e.g.,
(1)
where dPk,observed(θ)
/dΩ are n measured differential
scattering values and dPk(m,x,θ) /dΩ is given
by eq. A8. The error of fit is then
(2)
We used the logarithm of the phase
function in eq. 1 to increase the relatively low weight
of the scattering values observed at intermediate angles. We also found that by using
the mean absolute error (the M1 norm) rather than the mean sum of
squared data differences (the M2 norm used by Lienert et al., 2001), we were able to
significantly reduce the effect of large outliers on the solutions. Iterative
inversions (e.g.,
Since genetic inversion uses a random number generator to choose
individuals, it is not nearly as susceptible to becoming trapped in a local
fitness maximum as are iterative processes. However, because all possible
parameters are not explored, it is still possible that the global fitness
maximum will not be found, particularly in highly non-unique problems such as
this one. The only way to ensure that a truly global fitness maximum is
recovered is to increase the population size to include all possible parameters
which would make the calculation time prohibitive. Our approach is to repeat
the inversions multiple times using different starting seeds for the random
number generator (Lienert et al., 2001). In this way we obtain a family of
solutions representing local fitness maxima in the solution space. We also no
longer chose a convergence criterion based on fitness level (Lienert et al.,
2001) as this decreased the exploration of the parameter space. In this study,
we always ran 50 generations and used the five maximum fitness results for a
total of ten repeat inversions with different random number seeds.
As a test of the genetic inversion method,
we applied it to artificial angular scattering data generated at a wavelength
of 0.532μm and 1.25 degree
intervals. We added random noise of ±3% (peak) to the scattering calculated for
a unimodal distribution having r1=0.5 μm and μ1=2.0.
We then attempted to recover the original distribution using repeated genetic
inversions. The results are plotted in Figs. 4 and 5 as number and area
distributions, respectively. All of the recovered distributions had errors of
fit (eq. 2) that agreed to within ±0.1%. Although
there is considerable scatter in the number distributions, the area
distributions are reasonably well grouped about the original distribution
(heavy curve). The reason for the decrease in apparent scatter of the area
distributions is that the data being inverted (scattering cross sections)
depend on area rather than number. When a lognormal number distribution is
converted to area, the lognormal median radius, rm,
shifts to a larger value, rm,p,
according to equation
(3)
(e.g., Tsay et al., 1991) with p=2,
that depends on the
width, μi, of the number
distribution. Although the distribution's width is unchanged by converting to
area, different widths appear to be compensated for by changes in the peak
amplitude, Ni. with little
effect on the angular scattering. The net effect is a larger apparent non-uniqueness
in the inverted number distributions due to shifts in their central radii, Ri
that disappear when number is converted to area.
As
realistic size distributions are more adequately represented by a bimodal
distribution (Porter and Clarke, 1997), we repeated the numerical experiment
with a bimodal distribution as the starting model. We chose the initial bimodal
distribution to have similar mode radii and amplitudes to the distribution
experimentally measured by Clarke et al. (this issue) during SEAS. The resulting inverted area distributions are
shown in Fig. 6. It is clear that although the coarse mode peak (r1=1.5μm) is well
recovered, the amplitudes of the accumulation modes (r1=0.15μm) are too low, indicating that the polar nephelometer data is relatively insensitive to the latter
modes.
Inversion of
The polar nephelometer
was operated from 0400-0600 hrs GMT on
Although the theoretical scattering
values (Fig. 7) agree well with the experimental data at angles of up to 135°,
there is a systematic deviation at angles larger than this. To investigate
whether this is due to an incorrect value of the aerosol refractive index, in
Fig. 9 we have plotted the theoretical scattering for one of the inverted models
in Fig. 8 at three different values of refractive index. Fig. 9 demonstrates that
the model misfit cannot be accounted for by varying the refractive index. The
other possibility is the presence of non-spherical aerosols. Although the hydrated
coarse mode salt droplets are spherical, the accumulation mode, which Clarke et
al. (this issue) concluded was due to volcanic ash, could be causing the
misfit. Although we have demonstrated that the polar nephelometer
data is insensitive to the accumulation mode (Fig. 6), this conclusion may not
be valid if the aerosols are non-spherical. We experimented with adding an
additional spherical mode with a similar radius to that measured by Clarke et
al (this issue) having a refractive index of 1.6+0.01i, but were still unable
to fit the data at angles of greater than 135° without seriously degrading the
fit at angles of less than 30°.
Discussion
The inverted area distributions in Fig. 8 compare reasonably well in both amplitude and radius with the sizer data (heavy curve) of Clarke et al (this issue), although they do not include the accumulation mode aerosols. However, Clarke et al (this issue), in a plot of aerosol area times Mie cross-section efficiency, show that the accumulation mode contribution is less than 5% of the total scattering. Our mean calculated extinctions (Table 2) are similar to the values inferred from scanning lidar over breaking waves (Porter et al., this issue) and of inlet-corrected nephelometer measurements made during the SEAS experiment (Clarke et al., this issue). The fairly large scatter in all three of the unimodal parameters (Table 2) is due to non-uniqueness in the inversion, as the errors of fit are close to identical. This non-uniqueness can only be reduced by improving the scatter in the measured data or increasing their number. We have since reduced some of this error by using a logarithmic amplifier (Lienert et al, 2002). The lidar ratios we calculated from the different aerosol models in Table 1 ranged from 13.9 to 14.5. These are significantly lower than the lidar ratios of 20-30 measured directly at an altitude of 11 m by Masonis et al. (this issue). However, our instrument was directly measuring large salt spray droplets carried up from breaking waves on the beach which may not have been detected by the other instruments due to inlet losses (Porter and Clarke, 1997).
Appendix: Mie
Scattering Calculations
Mie’s solution (Mie, 1908; van de Hulst, 1981) for a spherical
particle of radius r microns (μm) gives the complex scattering
components 
and
, perpendicular and
parallel, respectively, to the plane of scattering where m is the complex refractive index of the particle, θ is the
scattering angle (in the plane of scatter), k=2π/λ,
l is
the wavelength (μm) and x=kr. For unit
incident radiation flux polarized either perpendicular (j=1) or parallel (j=2) to
the scattering plane, respectively, the differential angular scattering cross
sections are then
(A1)
μm 2sr-1., where Ω is the solid
angle. Note that we have used a derivative with respect to Ω to make the
distinction between total and differential scattering cross-section explicit.
Integrating eq. A1 over Ω gives the total
scattering cross-section as
(A2)
μm2. Normalizing eq. A2 by the cross-sectional area of the spherical particle, pr2 gives the (dimensionless) integrated scattering efficiency
A3)
For comparison with measurements by other instruments such as lidar, nephelometers, etc., it is also useful to calculate the extinction cross section (van de Hulst, 1981)
(A4)
μm2,
where
. For the case of zero absorption (Im(m)=0), 
.
In order to extend eqs. A1-A4 to an
assemblage of spherical particles having different radii, we assumed an M-modal lognormal distribution in ln r (Lienert et al., 2001)
(A5)
m-3, where Ni is the peak number
concentration in m-3, 
is the central radius
and 
is the spread in each
mode (we use 
rather than σi for the spread to distinguish it from
scattering cross section). We have also defined Ni as the peak, rather than the integrated number
concentration to keep it independent of σi. The advantage of the lognormal
distribution is that it has the same form in number, area and volume
distributions.
The total differential
scattering cross sections are now obtained by integrating eq.
A1 over ln r
to give
(A6)
m-1 sr-1, where the factor of 1012
results from converting k from (μm) -1
to m-1. Similarly, the
total integrated scattering becomes
(A7)
m-1. The differential scattering efficiencies
(phase functions), given by
(A8)
sr-1 are also useful to calculate. The ratio of total scattering to the scattering at 180O, termed the lidar ratio, is then
(A9)
sr, where the total extinction is given by
(A10)
m-1. When the absorption is zero, 
and
(A11)
In the case of non-polarized light, it is easily shown that the equations for the total integrated scattering efficiency, Q(m), and the differential phase function, dP(θ,m)/dΩ, are given by
(A12)
and
(A13)
The corresponding molecular
(Rayleigh) scattering equations for linear
polarization are (e.g., Coulson, 1988)
(A14)
at right angles to the
scattering plane and
(A15)
m-1 sr-1, perpendicular to the
scattering plane where ma is the atmospheric refractive index (close
to 1 for air) and N is the molecular
concentration in m-3 which is
calculated as a function of altitude using standard corrections (Coulson,
1988). Here, eqs. A14 and A15 have been separately
normalized to give unit integrals over solid angle, rather than normalizing
their sum, as is done in the unpolarized case (e.g., Coulson,
1988).
References
G. Backus and F. Gilbert, 1970. Uniqueness in the
inversion of inaccurate gross earth data, Phil.
Trans. R. Soc. London, 266,
123-130.
B. Barkey and K.N. Liou, 2001. Polar nephelometer for light-scattering
measurements of ice crystals, Optics Lett., 26,
232-234
C.F. Bohren and D.R.
Huffman, 1983. Absorption and Scattering
of Light by Small Particles, Wiley Interscience,
S. D. Billings, 1994. Simulated annealing for
earthquake location, Geophys. J. Int. 118,
680-692.
D.J. Carroll, 1996. Chemical laser modeling with
genetic algorithms, AIAA J., 34, 338-346.
A. Clarke, V. Kapustin, S. Howell, K. Moore and B.
Lienert, 2003. The Shoreline
Environment Aerosol Study (SEAS-2000): Marine aerosol measurements influenced
by a coastal environment and long range transport,
K. L. Coulson, 1988. Polarization and Intensity of Light in the
Atmosphere, DEEPAK Publishing,
D.
Deirmendjian, 1969. Electromagnetic
Scattering on Spherical Polydispersions,
O. Dubovik and M.D. King,
2000. A flexible algorithm for retrieval of aerosol optical properties from sun
and sky radiance measurements, J. Geophys. Res. 105,
20, 673-20,696.
J.E. Gentile, 1998. Random
Number Generation and Monte Carlo Methods, Springer Verlag,
D.E. Goldberg, 1989. Genetic
Algorithms in Search, Optimization and Machine Learning, Addison-Wesley,
J.H. Holland, 1975. Adaptation
in Artificial and Natural Systems, Univ.
H. C. van de Hulst, 1981. Light
Scattering by Small Particles,
D.D. Jackson, 1972. Interpretation
of inaccurate, insufficient and inconsistent data, Geophys. J. Roy. Astr. Soc., 28, 97-109.
B.R. Lienert, J.N. Porter, S.K. Sharma, N. Ahlquist and D. Harris, 2002. A 50 MHz logarithmic
amplifier for use in lidar measurements,
B. R. Lienert, J.N. Porter and S.K. Sharma, 2001,
Repetitive genetic inversion of optical extinction data, Applied Optics, 40, 3417-3427
S. J. Masonis, T.L.
Anderson, D.S. Covert , V. Kapustin, A.D. Clarke, S.
Howell, and K. Moore, 2003. A study of the
extinction-to-backscatter ratio and it's relation to other aerosol optical
properties during the Shoreline Environment Aerosol Study, with a
comparison to a polluted site and to Mie theory,
G. Mie, 1908. Beiträge zur optic trüber medien spieziell kolloidaler metallösungen, Ann. Phys., 25, 377-445.
S. Oshchepkov,
H. Isaka, J.-F. Gayet, A. Sinyuk, F. Auriol and S. Havemann, 2000.
Microphysical properties of mixed-phase & ice clouds retrieved from in in situ airborne "Polar
Nephelometer" measurements, Geophys. Res. Lett., 27 , 209-212,
Porter, J.N., S.K. Sharma, B. R.
Lienert, E. Lau and K. Horton, 2003. Aerosol scattering fields over Bellows Beach Oahu during the SEAS experiment,
J. Atmospheric and Oceanic Tech.
(this issue).
Porter, J.N., T.F. Cooney, C. Motell, 1998, Coastal aerosol phase function measurements with a custom polar nephelometer, ONR
Ocean Optics XIV Conference, Kona, Hawaii, 5.
J. Porter and A. D. Clarke, 1997. Aerosol size
distribution models based on in situ
measurements, J. Geophys.
Res., 102, 6035-6045
M.J. Post, 1975. Limitations of cloud droplet size
distribution by Backus-Gilbert inversion of optical scattering data, J. Opt. Soc. Am., 5, 483-486.
I. N. Tang, K. Tridico and
K. H. Pung, 1997. Thermodynamic and optical
properties of sea-salt aerosols, J. Geophys. Res., 102, 23,269-23,275.
S-C. Tsay, G. L. Stephens
and T. J. Greenwald, 1991. An investigation of aerosol microstructure on visual
air quality, Atmos. Environ., 25A, 1039-1053.
Table 1
|
Micro genetic algorithm enable |
YES |
|
Population size |
50 |
|
Number of parameters |
4 |
|
Mutation fraction |
0.04 |
|
Crossover probability |
0.4 |
|
Tournament selection enable |
YES |
|
Elitism enable |
YES |
|
Parameter creep enable |
NO |
|
Uniform selection enable |
YES |
|
Niche selection enable |
NO |
|
Number of children |
1 |
Table 1. Parameters used in
the genetic inversions.
Table 2
|
N1, cm-3 |
R1, μ |
Μ1 |
σex, m-1 |
L, sr |
C |
ε |
|
1.02 |
1.34 |
2.47 |
1.36x10-4 |
14.0 |
3.53x104 |
9.4% |
|
2.00 |
0.87 |
2.55 |
1.39x10-4 |
14.2 |
3.36x104 |
9.4% |
|
1.92 |
1.63 |
2.63 |
1.36x10-4 |
14.3 |
3.43x104 |
9.4% |
|
0.96 |
3.03 |
2.36 |
1.39x10-4 |
13.9 |
3.45x104 |
9.4% |
|
1.54 |
1.46 |
2.87 |
1.35x10-4 |
14.5 |
3.49x104 |
9.4% |
|
1.5±0.5 |
1.7±0.8 |
2.6±0.2 |
(1.37±0.02)x10-4 |
14.2±0.2 |
(3.45±0.06)x104 |
9.4% |
Table 2. Parameters and
their means for the five inverted unimodal
distributions shown in Fig. 8.
Figure 1
Figure 2
Figure 3
Figure 4
Figure 5
Figure 6
Figure 7
Figure 8
Figure 9.
Figure Captions
Figure 1. Three lognormal
area distributions approximately representing the accumulation (r1=0.083μm), coarse (r1=0.305μm) and large (r1=1.36μm) modes of aerosol size
distributions.
Figure 2. Combined Mie and
Molecular scatter phase functions for three different amplitudes of the coarse
mode peak in Fig. 1. The lidar ratios in sr are also shown for each curve.
Figure 3. Combined Mie and
Molecular scatter phase functions for three different amplitudes of the coarse
mode peak in Fig. 1. The theoretical extinctions (ext) in m-1
are also shown for each curve.
Figure 4. Number
distributions obtained by repeated genetic inversions of synthetic data
generated by adding ±3% random noise to the theoretical phase functions
calculated for a unimodal distribution (heavy solid
curve).
Figure 5. The size distributions in
Fig. 4 converted to cross-sectional area, dA /dln r
Figure 6. Cross-sectional area distributions obtained by
repeated genetic inversions of synthetic data generated by adding ±3% random
noise to the theoretical phase functions calculated for a bimodal distribution
(heavy solid curve). The inverted models (light curves) do not show an obvious
second mode due to the low amplitudes of the accumulation modes and their
closeness of their central radii to the coarse mode peak radii.
Figure 7. Measured data
(squares) collected at Bellows beach at 1200 hrs,
Figure 8. Cross-sectional area distribution models obtained by five repeated genetic inversions of the data in Fig. 7 (solid lines). Also shown (squares) are the aerosol sizer data of Clarke et al. (this issue) converted to cross-sectional area.
Figure 9. Theoretical scattering curves (solid curves) calculated for one of the unimodal models in Fig. 8 at a number of different refractive indices. The measured data are again shown as squares.