A Comparison of Cascade Impactor Data Reduction Methods

Aerosol Science and Technology 37: 187–200 (2003) c 2003 American Association for Aerosol Research ° Published by Taylor and Francis 0278-6826/03/$12.00 + .00 DOI: 10.1080/02786820390112506

A Comparison of Cascade Impactor Data Reduction Methods Patrick T. O’Shaughnessy1 and Otto G. Raabe2 1 2

Department of Occupational and Environmental Health, The University of Iowa, Iowa City, Iowa Center for Health and the Environment, University of California, Davis, California

The mass median aerodynamic diameter, dg , and geometric standard deviation, σg , of an aerosol is typically determined from a reduction of data obtained with a cascade impactor under the simplifying assumption that each stage of the impactor has an ideal collection efficiency. Two such reduction techniques are described and compared to an “inversion” method that incorporates the actual collection efficiencies of an impactor. Theoretical comparisons were made to demonstrate the difference in the estimate of dg and σg between the 2 ideal reduction methods and the inversion technique. Results indicated that, in general, both dg and σg are overestimated when the collection efficiencies are assumed to be ideal. A spreadsheet application of the inversion method is described.

are typically calibrated in terms of the “aerodynamic” diameters of the aerosol particles impacting on each stage. Therefore in this discussion, dg will be taken to be the mass median aerodynamic diameter (MMAD) of the aerosol. Several different impactor data reduction methods exist for determining dg and σg . The purpose of this paper is to describe and compare 3 reduction methods and their application to the use of computer spreadsheet programs. The methods will be described in an order associated with an increase in their computational complexity dictated by attempts to minimize the simplifying assumptions inherent to the initial methods described. IMPACTOR DATA REDUCTION METHODS

INTRODUCTION A cascade impactor is a multistage aerosol sampling device that separates particles by size according to their inertial properties in a moving air stream (Lodge and Chan 1986). The size distribution of many aerosols resulting from an analysis, or “reduction,” of impactor data approximates the lognormal distribution characterized by a frequency distribution that is skewed toward the larger particles (Hinds 1982; Raabe 1971). Given that an aerosol’s size distribution is approximately lognormal and unimodal, the 2 statistical parameters required to completely describe such a distribution are the geometric mean of the diameters, dg , and geometric standard deviation, σg (Hinds 1982). Furthermore, for a lognormally distributed variable, the geometric mean is equivalent to the median diameter, which is that diameter associated with 50% of the cumulative distribution of the mass collected by the impactor. Cascade impactors Received 18 September 2001; accepted 15 July 2002. The authors would like to thank Dr. Stephen Hillis of the University of Iowa Statistical Consulting Center for his help in formulating the method used to evaluate the standard errors of the parameter estimates. Address correspondence to Patrick T. O’Shaughnessy, Department of Occupational and Environmental Health, The University of Iowa, 100 Oakdale Campus, 180 IREH, Iowa City, IA 52242-5000. E-mail: [email protected]

Probit Method A common technique for determining dg and σg from impactor data is to first plot the cumulative mass percent versus the related cut diameter of each stage on log probability paper, where the stage cut diameter is defined as the aerodynamic diameter of a particle collected on the stage with 50% efficiency (Baron and Heitbrink 1993; Johnson and Swift 1997; Hinds 1982). If the mass fractions derived from the impactor data reduction represent those derived from a dust with a lognormal distribution of particle diameters, the data sets will plot as a straight line on graph paper of that type. Historically, the best-fit line through impactor data plotted on log probability paper was drawn by hand. An estimate of dg can be obtained directly from the plot by visually determining the diameter related to the point where the best-fit line intersects the line associated with 50% probability. If the particles were normally distributed, the standard deviation would simply be the difference between the diameter at 84.1% subtracted from the diameter at 50%, as this range (34.1%) represents the area under the normal distribution associated with one standard deviation from the mean (Hinds 1982). However, the particle distribution is lognormal. Hence the natural logarithm of σg is ln σg = ln d84.1% − ln d50% .

[1] 187

188

P. T. O’SHAUGHNESSY AND O. G. RAABE

By the properties of logarithms this difference can be expressed as a ratio and, after taking the antilog of both sides of the equation, σg can be expressed as σg =

d84.1% . d50%

[2]

More recently, computer code has been written to perform a least-squares linear regression of the data to increase the accuracy of the estimates of dg and σg (Knutson and Lioy 1995; Hinds 1986). A least-squares linear regression involves determining the parameters, α and β, of the linear equation y(x) = α + βx

[3]

that minimizes the sum of the squared difference between each observed y value, yi , and the corresponding expected value, y(xi ), for 2 < i ≤ n total x and y values with n–2 degrees of freedom: X [yi − y(xi )]2 = min. [4] Because the relationship between stage cut diameter and mass cumulative percent is nonlinear, a numerical transformation of these parameters must be performed prior to the linear regression analysis. The linear transformation of the cut diameters is relatively simplistic in that the natural logarithm is calculated for each cut diameter because the underlying distribution is assumed to be lognormal. The more involved process of computing a linear transformation of the cumulative percent values is accomplished by first realizing that the scaling of the probability axis is derived from the integral of the standardized normal probability density function (PDF) that defines the cumulative fraction 8(z) under the Gaussian distribution curve: µ 2¶ Z z −z 1 dz, [5] 8(z) = √ exp 2 2π −∞ where z = (x − µ)/σ in which µ and σ are the mean and standard deviation of the data values, respectively, and the standardized normal variable, z, has µ = 0 and σ = 1. Establishing z as the independent variable, the resulting linear equation is ln d = α + βz.

[6]

The corresponding graph is related to the use of the linear probit scale rather than the nonlinear probability scale. Under this premise, dg is equivalent to the intercept, α, as this represents the diameter at z = 0, which is equivalent to a probability of 50%. Likewise, the slope, β, can be expressed as the change in the natural log of the diameters relative to a unit change in z: β=

ln dz=1 − ln dz=0 . 1−0

[7]

Given that z values of 1 and 0 are equivalent to probabilities of 84.1 and 50%, respectively, Equation (7) is identical to Equation (1) and therefore β is equivalent to the natural logarithm of σg , hence σg = exp(β). Spreadsheet Development. As described above, the cumulative percent calculated for each impactor stage can be transformed into its equivalent (linear) z value and plotted relative to the logarithm of the cut diameters. Current versions of most spreadsheets contain a predefined function that incorporates the numerical solution relating 8(z) to z. In spreadsheets that contain such a function, the acronym “NORMSINV” is used to name the function. Therefore the stage cumulative percents (transformed into cumulative fractions), calculated as part of the impactor data reduction, can be transformed into their corresponding z values and plotted on a linear scale versus the stage cut diameters plotted on a log scale (Figure 1). Furthermore, the built in linear regression capabilities of these spreadsheets can be utilized to determine dg and σg by taking the antilogarithm of the slope and intercept resulting from a linear regression of log dc versus z.

Inverted Probit Method As the proof given in the previous section implies, the use of linear regression analysis to determine dg and σg relies numerically on establishing z as the independent variable. This condition is counter to the physical aspect of the impactor sampling process whereby z is a product of the mass measured on each stage, which is dependent on the cut diameter of each stage. Therefore the linear relationship should be inverted to correctly account for the actual relationship between the cumulative fractions as the dependent variable and cut diameters as the independent variable. A determination of dg and σg can then be accomplished from the α and β values resulting from the linear regression analysis of the inverted variables where α ln dg = − , β ln σg =

1 . β

[8] [9]

Sample Variance Inequalities. Proper application of a leastsquares linear regression is performed under the assumption that the variance of the dependent variable is constant regardless of the magnitude of the independent variable. However, the probability scale of log-probability paper (or log-probit paper) disproportionately magnifies the effect of an error in cumulative percent near the extremes of the scale compared to values near 50% (Hinds 1982, 1986). (The same is also true when applying a logarithmic scale to the ordinate but is often ignored, as when determining log-log or semilog relationships.) In general, if a regression is performed on dependent values transformed by the function, f (yi ), rather than directly on the data values, yi , then the estimated standard deviation (or “uncertainty”) associated with each data value, syi , must be

IMPACTOR DATA REDUCTION METHODS

189

Figure 1. Spreadsheet-generated log-probit paper with regression line on points determined from impactor data reduction. modified by s f (yi) =

d f (yi ) syi dy

[10]

to compensate for the distortion of each syi caused by the scaling transformation (Bevington and Robinson 1992). For example, if the dependent variable is transformed by taking the natural logarithm of each value, the standard deviations would be modified by sln( yi) =

d(ln yi ) 1 syi = syi . dy yi

[11]

For the case described here, the cumulative mass fraction, 8(z), is transformed into a corresponding z value by the inverse of

Equation (3). With reference to Equation (11), the standard deviation of the variable z, sz , is given in terms of the observed standard deviation of each 8(z) value, s8 , by √ dz 2π s8i = s8i , szi = d8(z) exp(−z 2 /2)

[12]

where for an impactor the subscript “i” is used to index the total number of stages of the impactor. As shown in Figure 2, the value of dz/d8(z) remains relatively constant between cumulative fractions of 0.2 and 0.8, but increases dramatically below and above those levels. The curve shown in Figure 2 also demonstrates the amount of distortion to be expected from plotting on a probability axis given

Figure 2. Distortion in probability scale relative to no distortion at a cumulative percent level of 50%.

190

P. T. O’SHAUGHNESSY AND O. G. RAABE

Figure 3. The effect of scaling distortion on sample error assuming a 5% error about each point. equal s8 values. This distortion results in an enlargement of any error associated with z values plotted below or above the 20 (z < −0.85) and 80% (z > 0.85) lines, respectively (Figure 3). To compensate for differences between each sz resulting from an application of Equation (12), the method of weighted least squares (WLS) can be employed (Bevington and Robinson 1992). WLS is a generalized form of least-squares regression that allows for the incorporation of the sz that necessarily vary. In general, WLS weights the squared errors between the measurements, yi , and the regression values, y(xi ), given in Equation (4), by the reciprocal of each variance term s 2f (yi) (Bevington and Robinson 1992; Raabe 1978). The WLS method then determines the intercept, α, and slope, β, that minimizes the summation given in Equation (13). X [yi − (α + βxi )]2 X [yi − y(xi )]2 = = min. s 2f (yi) s 2f (yi)

e

z(ln di ) = α + β(ln di ).

[15]

To then determine the values of α and β that minimize Equation (13), the partial derivatives of Equation (13) with respect to both parameters are evaluated when set equal to zero (Bevington and Robinson 1992). The computational equations that satisfy that criteria are given below. µ X xi X xi yi ¶ 1 X xi2 X yi − , [16] α= 1 szi2 szi2 szi2 szi2 µ ¶ 1 X 1 X xi yi X xi X yi β= − , [17] 1 szi2 szi2 szi2 szi2 where:

[13]

Since the variance term is used to weight the squared difference between measured and predicted levels of y, it is not important to accurately determine the absolute value of each variance, but only to ensure that each has a magnitude that is correctly relative to all other variance terms. In summary, determining a minimum value of Equation (13) is accomplished after (1) transforming the cumulative fractions, 8(z), associated with each stage into an equivalent z value; (2) making the logarithmic transformation of the diameters; and (3) transforming the standard deviation of the cumulative fractions for each stage, s8i , with the use of Equation (12) to obtain X [z i − z(ln di )]2 ³√ ´2 = min, 2π s 8i 2 −z /2

where

[14]

X 1 X x 2 µX xi ¶2 i 1= − . szi2 szi2 szi2

[18]

Commercially available spreadsheets do not have an option to perform weighted least-squares regression. Therefore to determine α and β, other software must be used or extra columns must be developed in a spreadsheet to compute the summations given in Equations (16)–(18). Equations (8) and (9) can then be used to determine dg and σg , respectively, from the values of α and β determined from Equations (16)–(18). Variance Assumptions. Application of the method described above relies on knowledge of the standard deviation of the cumulative fractions for each stage, s8i . However, an evaluation of the size distribution of an aerosol with the use of a cascade impactor is often obtained with a single sample. In this case, a value for each s8i cannot be determined directly. As given in the appendix, the variance of the stage cumulative fractions, 2 2 , can be determined if the variance of the mass, smi , collected s8i

191

IMPACTOR DATA REDUCTION METHODS

on each stage is known. This likewise leads to the problem of 2 determining each smi when only one sample is taken. Rather than rely on the na¨ıve assumption that sm2 is equivalent for all stages, we propose that a more realistic assumption be that sm2 is directly proportional to the mass on each stage. Although not meant to be conclusive, to test this assumption a series of 12 trials were performed during which a pulverized grain dust was aerosolized into a 1 m3 environmental chamber. During each trial the chamber flow rate, aerosol generation rate, impactor flow rate, and placement of the impactor within the chamber were held constant. After sampling, impactor substrates were weighed in a climate-controlled room with a 6 place balance and the variance of the mass collected on each stage was determined. The mass fraction on each stage for each trial was computed.

Inversion Method The inverted probit reduction method described above reduces the error induced by making linear transformations of nonlinear variables. However, that method still incorporates the simplifying assumption that each stage of the impactor has an ideal, perfectly sharp collection efficiency. Several data “inversion” techniques have been developed to incorporate the actual stage collection efficiencies (Cooper 1993; Dzubay and Hasan 1990; Puttock 1981; Raabe 1978; Ramachandran and Vincent 1997). A least-squares fitting method for analyzing particle size distribution data was described by Kottler (1950). Raabe (1978) was the first to develop a stage-efficiencies weighted least-squares method for impactor data to calculate the maximum likelihood fitted parameters of the distribution based on the variation in the variance of the stage mass fractions while also incorporating the actual impactor stage efficiency curves to account for their deviation from the ideal. A spreadsheet application of that method is described below and applied to the reduction of data derived from the commonly-used Marple personal cascade impactor (Series 290, Anderson Inst. Inc., Smyrna, GA, referred to in subsequent text as “the impactor”). Incorporation of Stage Efficiency Curves. The regression techniques described above for determining dg and σg are permissible if ideal stage collection efficiencies are assumed because a single diameter can be associated with the mass collected on a stage. However, if the actual, nonideal stage efficiency curves are considered in the analysis of dg and σg , then each stage mass fraction, f i , is associated with all diameters spanning the stage collection efficiency curve of that stage. If one assumes that a lognormally-distributed aerosol with a given dg and σg is measured by the impactor, then the fraction of dust collected over a small size range can be calculated for a given stage if the collection efficiency for that size range is known. The total mass fraction for a particular stage can then be determined by summing the individual fractions determined over the entire range of the collection efficiency curve. Given the resulting set of computed mass fractions for each stage, f i (dg , σg ), the WLS method can be used to determine the value of dg and σg that minimizes the sum of the weighted squares of the differences

between the measured fractions, f i , and the computed fractions as given in Equation (19) (Raabe 1978). X [ f i − f i (dg , σg )]2 = min. sfi2

[19]

Unlike that given in Equation (13), the weighting factor, 1/sfi2 , of Equation (19) does not compensate for differences in scaling but is rather the reciprocal of an estimate of the variance for each measured fraction (Bevington and Robinson 1992; Raabe 1978). Again, this variance is unknown if only one sample is taken. However, the distribution of a fraction, f , is best described by the binomial distribution, which has a variance defined as f (1– f ) given one sample. As part of the experimental analysis used to verify the assumption used to estimate the variance of the mass on each stage, sm2 , the variance of the mass fractions, s 2f , on each stage was likewise computed and compared to values based on the assumption that s 2f = f (1– f ). Furthermore, as mentioned above, only the relative, not absolute, value of each variance term is needed for application in Equation (19). Spreadsheet Development. Properties of various impactors have been extensively evaluated by many investigators since its invention by May (1945). For example, Raabe et al. (1988) provide results of careful stage calibration of a circular-jet impactor. Rader et al. (1991) developed equations to fit measured efficiency curves for the stages of the Marple personal cascade impactor previously calibrated by Rubow et al. (1987). The general form of the equation used by Rader et al. to describe each stage efficiency, E i , is ·µ ¶bi ¸ dar , [20] E i (dar ) = tanh ai where dar is the aerodynamic resistance diameter first described by Raabe (1976) and ai and bi are the fitted parameters that best describe the curves for each stage. The aerodynamic diameter, dae , is related to the dar by the slip-correction factor, C, p dar = dae C(dae ) and is therefore referred to as the “slip-corrected equivalent of the aerodynamic diameter” by Rader et al. (1991). The use of explicit equations to define the collection curves of the impactor evaluated here simplified the application of this method for use in a spreadsheet. However, as emphasized by Raabe (1978), this method will work equally well given a series of interpolated values for each curve. Furthermore, impactor calibration methods other than those described by Rader et al. (1991) may be applied to obtain collection curves for impactors (e.g., Raabe et al. (1988)). As applied to the inversion spreadsheet, the efficiency values for 101 diameters between 0.32 µm and 100 µm (−0.5–2.0 log10 diameter by increments of 0.025 log diameter) and for each stage were calculated with the use of Equation (20). The resulting curves are shown in Figure 4.

192

P. T. O’SHAUGHNESSY AND O. G. RAABE

Figure 4. Impactor stage collection efficiency curves as determined by Rader et al. (1991).

Under the assumption that the aerosol sampled by the impactor is unimodal and lognormally distributed, the theoretical cumulative fraction, 8(z), less than each of the 101 diameters was determined with the use of a spreadsheet (Excel, Microsoft Corp., Seattle, WA) and the predefined function NORMDIST. That function computes the 8(z) value associated with a particular diameter given a mean diameter and standard deviation. In this case, the log10 of the diameters is normally distributed. Therefore the function was supplied with the log of each diameter as well as the log of chosen values for dg and σg . The mass fraction, f i (1d), between each of the 101 diameters was then computed by subtracting adjacent 8(z) values. These f i (1d) values were then modified to account for internal losses and sampling inlet efficiency with the use of equations supplied by Rader et al. (1991). The amount of aerosol within a certain size range that deposits onto a given stage is equivalent to the product of the amount in that size range that penetrates through the previous stage and the stage collection efficiency for that size range (Raabe 1978; Rader et al. 1991; Ramachandran and Vincent 1997). Therefore starting with the modified f i (1d) values, the percent penetrating through and depositing on each stage was determined for each of the 100 size classes analyzed. The estimated fraction deposited on each stage, f i (dg , σg ), was then calculated by determining the sum of the 100 fractions for each size range computed for each stage. These estimated fractions were then used to develop the summation defined by Equation (19). Nonlinear root-solving techniques are available to find the values of dg and σg that minimize Equation (19) (Bevington and Robinson 1992; Press et al. 1989). Software programming can be used to implement these numerical methods, however, most spreadsheets contain a precoded function to perform nonlinear root-solving automatically. The particular one utilized as part

of this analysis, “Solver” in the spreadsheet Excel, utilizes a generalized reduced gradient (GRG2) nonlinear optimization code (Frontline Systems Inc., Incline Village, NV). In this case, Solver is used to minimize the contents of the cell containing the summation defined by Equation (19) while changing the contents of the cells containing the values of dg and σg used to determine the predicted fractions, f i (dg , σg ) with the additional constraints of dg ≥ 0.01 and σg ≥ 1.01. A common difficulty that arises when attempting to find the minimum (or maximum) of a complex numerical system, such as the one described here, is the likelihood of determining a “local” minimum rather than the desired minimum. To avoid this problem here, the Probit method described above can be used to find initial values for dg and σg that will likely be near those found when using the Inversion method. Employing the Probit method is also useful for performing an initial evaluation as to whether the cumulative mass fractions fall approximately on a straight line when plotted and therefore indicate that the aerosol size distribution is lognormal and unimodal.

Confidence in the Parameter Estimates A measure of the confidence in the estimates of the parameters dg and σg can be obtained by computing their standard error (SE). This is a relatively straightforward process if dg and σg were computed using linear least squares regression as when using the Probit method. In that case, the SE values are determined from an estimate of the experimental error variance, σ 2 , by computing the error sum of squares (SSE) equivalent to the summation given in Equation (19) for the weighted least squares case. An unbiased estimate of σ 2 can then be obtained by dividing the SSE by the appropriate number of degrees of freedom, which is the sample size, n, minus the number of parameters, p. The variance

IMPACTOR DATA REDUCTION METHODS

of a single parameter estimate, say p(1), is then computed from SSE/(n– p) multiplied by a function of the independent variable (or matrix of multiple independent variables), f (x), by s 2p(1) =

SSE f (x) n– p

[21]

and the standard error of p(1), s p(1) , is the square root of s 2p(1) . (The nature of the term, f (x), depends on the parameter evaluated. A statistical text may be consulted for a complete description of the formulation of SE values for estimated parameters.) Given a SE value for a parameter estimate, the (1–α) 100% confidence limits for the true parameter value can be obtained from p(1) ± tν,α/2 S p(1) ,

[22]

where ν represents the n– p degrees of freedom. The output displayed when using a spreadsheet to perform a univariate linear regression analysis will typically contain both the SE values and corresponding confidence intervals for the calculated slope and intercept. The method described above results in an unbiased estimate of the SE of a parameter estimate obtained from a linear least squares regression analysis. This estimation is also “exact” in that the function, f (x), can be solved explicitly. However, an exact solution of f (x) cannot be obtained for a nonlinear function. Therefore the SEs of the parameter estimates cannot be solved explicitly. One method for estimating SE values for parameters associated with a nonlinear function is to utilize a Taylor series expansion to approximate the nonlinear model with linear terms (Neter et al. 1996). Application of this “Gauss–Newton” method can result in both a least-squares estimate of the parameters, dg and σg , as well as their standard errors. However, in this case, the method was used only to determine the standard errors once best estimates were determined using the spreadsheet’s precoded iterative solving routine. An explanation of this method is given in the appendix.

Method Comparison Theoretical comparisons between the 2 linear regression methods described above and the Inversion method were made under the assumption that use of the Inversion method resulted in the most accurate estimate of the true dg and σg values. Therefore rather than use the Inversion method to determine unknown dg and σg values from a set of measured stage fractions, known dg and σg values were applied to a spreadsheet containing the equations describing the stage collection efficiency curves to determine the predicted stage fractions. These stage fractions therefore represent those that would result from a perfectly lognormally distributed aerosol captured by the impactor with nonideal stage collection efficiencies. A primary difficulty associated with the use of the linear regression methods occurs when the data points do not fall on a

193

straight line. This may occur when the distribution of the aerosol deviates from that of a perfect lognormal distribution. Deviation from a straight line may also occur if the aerosol distribution is not unimodal. However, for these comparisons only the mass fractions of a unimodal, lognormal distribution was determined. Any deviations from the straight line therefore occurred as a consequence of predicting an insignificant amount of mass for stages near the top or bottom of the impactor. For example, very low mass fractions were predicted for the lower stages given an aerosol with a large dg and low σg . This caused the points associated with the lower stages to drop relative to the slope created by the upper stage points. In practice, one may perform a linear regression on only those points that form a straight line. For comparison purposes, however, only those combinations of dg and σg that resulted in linear regression through all points with an r 2 value greater than 0.995 were compared to the Inversion method. This limitation excluded the comparison of distributions with σg < 2.0 and dg > 8. Given the limitations described above, the comparisons were conducted for the 15 combinations of dg values of 1, 2, 4, 6, and 8 µm and σg values of 2.0, 2.5, and 3.0. The stage mass fractions predicted from stage collection efficiency equations were applied to the separate spreadsheets used to determine dg and σg by the linear regression methods. The percent differences in the resulting values of dg and σg , compared to those of the original set applied, were then computed for each combination of dg and σg . The application of corrections for internal losses and inlet sampling efficiency were not applied when making these comparisons as they were assumed to hold true regardless of method. Furthermore, the set of stage cut diameters applied to the linear regression methods were those reported in the Rader et al. (1991) paper rather than the commonly-used cut diameters reported in the Rubow et al.(1987) paper. RESULTS AND DISCUSSION

Variance Assumptions For the 12 dust trials performed, the variance of the mass collected on each stage and back up filter was computed. A statistical analysis (F test) determined a significant difference between the smallest (0.0002, stage 8) and largest mass variance (0.137, stage 3) ( p < 0.001). Therefore an assumption that the mass variance, sm2 , is constant when applied to the calculation of the variance of the cumulative mass fraction, s82 , is not valid. Because the variances obtained by estimation were of a different magnitude than those calculated, each set of estimated and calculated variances were normalized relative to the sum of all variances obtained for each set. As shown in Figure 5, a 2 = mi calculation of s82 for each stage using the assumption smi gives a closer approximation to the actual variance of the mass measured for each stage compared to values obtained under 2 = constant. The SSE between the measured the assumption smi

194

P. T. O’SHAUGHNESSY AND O. G. RAABE

Figure 5. Normalized stage cumulative fraction variances compared to predicted variances. variances and the estimated variances was 236 when applying the first assumption compared to 796 when applying the second. Likewise, the variance of the 12 measured mass fractions for each stage, s 2f i , was computed for each stage and back up filter. A statistical analysis (F test) determined a significant difference between the smallest (0.152, stage 8) and largest variance (2.769, stage 2) ( p < 0.001). Therefore the need to utilize weighted least squares to compensate for differences in s 2f is justified when using the Inversion method. A comparison was made between these calculated variances and estimations of the variance based on the assumption that s 2f i = f i (1– f i ) (Figure 6). After normalizing each set of values as a percentage of the sum of all values, the

SSE (145.3) was less than that obtained while assuming a linear relationship between the two (188.6). Although these analyses were not performed in such a way as to demonstrate conclusive results to support use of the variance assumptions, they do lend reasonable guidance for an estimation of these variances when only one sample has been taken.

Method Comparison As shown in Figures 7 and 8, use of the 2 linear regression reduction methods resulted in values of dg and σg greater than those supplied to the Inversion method. However, there was closer agreement between the Inversion method and the

Figure 6. Normalized stage mass fraction variances compared to predicted variances.

IMPACTOR DATA REDUCTION METHODS

195

Figure 7. Percent difference in dg (a) and σg (b) predicted by the Probit method relative to those applied to a spreadsheet that incorporated the actual stage efficiency curves. Inverted Probit method compared to its agreement with the Probit method. The increased estimate of dg by the 2 linear regression methods is primarily related to the nature of the actual stage collection efficiency curves relative to the ideal collection assumed when using these 2 methods. The sigmoidal nature of the collection curves indicates that each stage will fail to collect some particles greater than the cut diameter and will likewise collect some particles less than the cut diameter. Conversely, under ideal collection each stage will collect only those particles with sizes falling between its cut diameter and the larger cut diameter of the previous stage. Actual collection on each stage will there-

fore result in collecting particles with sizes above and below this range. If for all stages the mass of particles smaller than the stage cut diameter is greater than the mass of particles larger than the cut diameter of the previous stage, the linear-regression methods will overestimate dg . As shown in Figure 9, this is the case for the Marple personal cascade impactor regardless of the underlying size distribution. The long tail of the collection efficiency curve for stage 1 (Figure 4) has a pronounced effect on this difference in the estimate of dg . The Inversion method correctly associates the mass collected on stage 1 with the contribution of particles less than the cut diameter for that stage, whereas the linear-regression

196

P. T. O’SHAUGHNESSY AND O. G. RAABE

Figure 8. Percent difference in dg (a) and σg (b) predicted by the Inverted Probit method relative to those applied to a spreadsheet that incorporated the actual stage efficiency curves.

methods assume the mass to be contributed by all particles greater than the cut diameter, thereby adding to the inflation of the estimate of dg by that method. A sensitivity analysis was also performed in which stage fractions predicted from the equations describing the stage collection efficiency curves were applied as “measured” mass fractions to the spreadsheet containing the Inversion method. If starting values for the dg and σg were purposely chosen to be much different than those used to calculate the stage fractions, then the Solver routine would often converge on a minimum located at a combination of dg and σg that were not near the original set. However, when the Probit method was used to generate starting

values, as suggested above, the largest percent error between supplied and determined dg and σg was 7.5 × 10−5 and 5.6 × 10−5 , respectively. The 95% confidence intervals about the estimates of dg and σg shown in Figure 1 were determined using the unbiased method associated with linear regression analysis. This method resulted in confidence intervals of 1.04 µm about the estimated dg value of 4.34 µm and 1.03 about the estimated σg value of 2.64. The estimation technique described in the appendix for determining confidence intervals about parameters estimated by nonlinear least squares was also applied to the same impactor data set. This method resulted in confidence intervals of 0.25 µm about

IMPACTOR DATA REDUCTION METHODS

197

Figure 9. The ratio of all mass collected below, relative to all mass collected above, successive stage cut diameters given various combinations of dg and sg . an estimated dg value of 3.85 µm and 0.14 about the estimated σg value of 2.40. The fitted distributions using these 2 methods are plotted in Figure 10 together with a histogram of the relative stage mass percent relative to aerodynamic particle diameter, dae . As shown in Figure 10, the distribution determined with the Probit method is wider and shifted toward the larger particle diameters, a relationship consistent with the general results found when comparing these 2 methods.

CONCLUSION The spreadsheet application of 2 linear regression methods and an inversion method for the reduction of data collected by a cascade impactor were described. Application of least-squares linear regression to impactor data plotted on log-probability paper requires that the stage cut diameter be the dependent variable in order to derive the dg and σg values from the resulting intercept and slope, respectively. To correct for this condition, a linear

Figure 10. Histogram of the relative mass percent on each stage compared to lognormal distributions resulting from application of the Probit and Inversion methods.

198

P. T. O’SHAUGHNESSY AND O. G. RAABE

regression method was described that inverts the relationship between cut diameter and the associated stage cumulative percent. These linear regression methods were compared to an inversion method that incorporates the stage collection efficiency curves in the analysis of dg and σg . The method comparison demonstrated that both linear regression methods overestimated both the dg and σg within the ranges analyzed. The linear regression methods will therefore underestimate the fraction of the smaller, respirable particles. Furthermore, the increased complexity of the Inverted Probit method over that of the Probit method does not appear to be justified as this method also produces dg and σg values consistently higher than the Inversion method. Future work will involve the addition of features to the Inversion method spreadsheet to allow for the analysis of bimodal and trimodal distributions. In that regard, the extent to which the impactor can distinguish the separate median diameters of a bimodally distributed aerosol will be evaluated. An electronic copy of any, or all, spreadsheets described in this paper is available upon request to the corresponding author. NOMENCLATURE d aerosol diameter log d log base 10 of aerosol diameter slip-corrected aerosol diameter dar mass median aerodynamic diameter of the aerosol dg lognormal size distribution measured mass fraction on the ith stage of an imfi pactor f i (dg , σg ) expected mass fraction on the ith stage based on impactor stage collection efficiency curves and unimodal, lognormal size distribution defined by dg and σg f i (1d) expected mass fraction for the size interval 1d on the ith stage of an impactor estimate of the standard deviation of the mass fracsfi tion of the ith stage of an impactor estimate of the standard deviation of mass collected smi on the ith stage of an impactor estimate of the standard deviation of a dependent syi variable measured at a dependent level, xi s8i estimate of the standard deviation of the cumulative fraction of mass collected on the ith stage of an impactor s8i normalized to compensate for scaling distorszi tions when transformed to its equivalent normal standard variate, z an independent variable indexed by i from 1 to n xi a dependent variable associated with the indepenyi dent variable, xi z normally distributed variable having a mean of 0 and standard deviation of 1 z value associated with the ith stage of a cascade zi impactor

α β 8(z) σg σi2

intercept of a linear regression slope of a linear regression the cumulative distribution function of the standard normal distribution, N (0, 1) geometric standard deviation of the aerosol lognormal size distribution population variance of mass collected on the ith stage of an impactor

REFERENCES Baron, P. A., and Heitbrink, W. A. (1993). Factors Affecting Aerosol Measurement Quality. In Aerosol Measurement: Principles, Techniques, and Applications, edited by K. Willeke and P. A. Baron. Van Nostrand Reinhold, New York, pp. 146–176. Bevington, P. R., and Robinson, D. K. (1992). Data Reduction and Error Analysis for the Physical Sciences, 2nd ed., McGraw-Hill, New York. Cooper, D. W. (1993). Methods of Size Distribution Data Analysis and Presentation. In Aerosol Measurement: Principles, Techniques, and Applications, edited by K. Willeke and P. A. Baron. Van Nostrand Reinhold, New York, pp. 146–176. Dzubay, T. G., and Hasan, H. (1990). Fitting Multimodal Lognormal Size Distributions to Cascade Impactor Data, Aerosol Sci. Technol. 13:144– 150. Hinds, W. C. (1982). Aerosol Technology: Properties, Behavior, and Measurement of Airborne Particles, John Wiley & Sons, New York. Hinds, W. C. (1986). Data Analysis. In Cascade Impactor: Sampling & Data Analysis, edited by J. P. Lodge and T. L. Chan. ACGIH, Inc., Cincinnati, OH. Johnson, D., and Swift, D. (1997). Sampling and Sizing Particles. In The Occupational Environment—Its Evaluation and Control, edited by S. R. DiNardi. American Industrial Hygiene Association, Fairfax, VA. Knutson, E. O., and Lioy, P. J. (1995). Measurement and Presentation of Aerosol Size Distributions. In Air Sampling Instruments for Evaluation of Atmospheric Contaminants, 8th ed., edited by B. S. Cohen and S. V. Hering. ACGIH, Cincinnati, OH, pp. 121–137. Kottler, F. (1950). The Distribution of Particle Sizes, J Franklin Inst. 250:339– 352, 419–441. Lodge, J. P., and Chan, T. L. eds. (1986). Cascade Impactor: Sampling & Data Analysis, ACGIH, Inc., Cincinnati, OH. May, K. R. (1945). The Cascade Impactor: An Instrument for Sampling Coarse Aerosols, J. Sci. Instrum. (London) 22:187–195. Neter, J., Kutner, M. H., Nachtsheim, C. J., and Wasserman, W. (1996). Applied Linear Statistical Models, 4th ed., Irwin, IL, pp. 531–552. Press, W. H., Flannery, B. P., Teukolsky, S. A., and Vetterling, W. T. (1989). Numerical Recipes in Pascal: The Art of Scientific Computing, Cambridge University Press, New York. Puttock, J. S. (1981). Data Inversion for Cascade Impactors: Fitting Sums of Log-Normal Distributions, Atmos. Environ. 15(9):1709– 1716. Raabe, O. G. (1971). Particle Size Analysis Utilizing Grouped Data and the Log-Normal Distribution, Aerosol Sci. Technol. 2:289–303. Raabe, O. G. (1976). Aerosol Aerodynamic Size Conventions for Inertial Sampler Calibration, APCA Journal 26(9):856–860. Raabe, O. G. (1978). A General Method for Fitting Size Distributions to Multicomponent Aerosol Data Using Weighted Least-Squares, Environ. Sci. Technol. 12(10):1162–1167. Raabe, O. G., Braaten, D. A., Axelbaum, R. L., Teague, S. V., and Cahill, T. A. (1988). Calibration Studies of the DRUM Impactor, J. Aerosol Sci. 19:183– 195. Rader, D. J., Mondy, L. A., Brockmann, J. E., Lucero, D. A., and Rubow, K. L. (1991). Stage Response Calibration of the Mark III and Marple Personal Cascade Impactors, Aerosol Sci. Technol. 14:365–379.

199

IMPACTOR DATA REDUCTION METHODS Ramanchandran, G., and Vincent, J. H. (1997). Evaluation of Two Inversion Techniques for Retrieving Health-Related Aerosol Fractions from Personal Cascade Impactor Measurements, Am. Ind. Hyg. Assoc. J. 58:15–22. Rubow, K. L., Marple, V. A., Olin, J., and McCawley, M. A. (1987). A Personal Cascade Impactor: Design, Evaluation and Calibration, Am. Ind. Hyg. Assoc. J. 48(6):532–538.

Confidence in the Parameter Estimates The least squares criterion for the nonlinear response function f i (dg , σg ) was given in Equation (19) as

APPENDIX

If the 2 parameters, dg and σg , are expressed as a vector, γ , the response function can be rewritten as f i (γ ), which is the mean response for the ith stage plus filter. To minimize the least squares criterion, the Gauss–Newton, or linearization, method uses a Taylor series expansion to approximate the nonlinear regression model with linear terms and then employs ordinary least squares to estimate the parameters (Neter et al. 1996). Considering the case where the weights, 1/s 2f i , are not applied, then given starting values, gk(0) , for the 2 parameters, the mean responses for the n stages and back up filter are approximated by the linear terms in the Taylor series expansion, which, for the ith stage, is

Variance of Cumulative Fractions For a cascade impactor with n stages, with the top (largest cut diameter) stage numbered 1 and n + 1 the filter, each stage collects a mass of material m i (with standard deviation smi ) associated with particles larger then the cutoff diameter for that stage. The cumulative mass fraction, F, less than the cutoff diameter for a stage, k, is therefore the ratio of all mass collected by higher numbered stages to the total mass collected on all stages as defined below. Pn+1

v i=K +1 m i = = Fk = P n+1 u+v i=1 m i

u v

k X

¸ p−1 · ¡ (0) ¢ X ¡ ¢ ∂ f i (γ) f g + γk − gk(0) , f i (γ) ∼ = i ∂γk γ=g(0) k=0

1 1 = , +1 w

where u and v are independent of each other and given by u=

X [ f i − f i (dg , σg )]2 = min. s 2f i

where g(0) is the vector of starting values for γ . If the notation is simplified as ¡ ¢ αi(0) = f i g(0) ,

mi ,

i=1

v=

n+1 X

·

mi .

(0) = Dik

i=k+1

·

s2 s2 = u2 + v2 v u

¸· ¸2 u , v

¸ , γ=g(0)

βk(0) = γk − gk(0) ,

Since u and v are unrelated, the variance of (u/v) is given by

2 s(u/v)

∂ f i (γ) ∂γk

then (0) f i (γ) ∼ = αi +

p−1 X

(0) (0) Dik βk

k=0

where su2 =

k X

2 smi

sv2 =

i=1

n+1 X

2 smi .

i=k+1

and therefore an approximation of the nonlinear regression model is (0) fi ∼ = αi +

p−1 X

(0) (0) Dik βk + εi ,

k=0

Also, 2 sw2 = su/v

so that the variance of the cumulative fraction for stage k is · s F2 k =

· 2 ¸· ¸2 ¸ sw2 sw s2 1 2 = = w4 . [ ] F k 2 2 w w w w

where εi is the random error term for the ith case. Defining Fi(0) = f i − αi(0) , then the linear approximation can be expressed in matrix notation as F(0) ∼ = D(0) β (0) + ε,

200

P. T. O’SHAUGHNESSY AND O. G. RAABE

which has the same form as the general matrix expression for a linear regression model: Y = Xβ + ε. Given this linear approximation of the nonlinear model, an approximate variance-covariance matrix of the regression coefficients can be estimated by s2 {g} = MSE(D0 D)−1 , where SSE = MSE = n– p

P

[ f i − f i (γ)]2 . n– p

With the application of weighted least squares, the least squares criterion given in Equation (19) can be expressed as ¸ X · fi f i (γ) 2 − = min. sfi sfi An evaluation of the variance-covariance matrix can therefore be evaluated as given above, except that each of the partial derivatives, Dik , defined above must be multiplied by 1/si , the square root of the weights, wi = 1/s 2f i . Since each derivative is multiplied by every other derivative when evaluating D0 D, this is equivalent to multiplying each pair of derivatives by the weights, wi . In matrix form this involves creating an n × n matrix, W = diag[w i ], consisting of the wi along the diagonal and all zeros elsewhere, and computing the variance-covariance matrix by s2 {g} = MSEw (D0 WD)−1 ,

where 1 MSEw = = n– p

P

[ f i − f i (γ)]2 . s 2f i

Application to Impactor Inversion Routine. Once the Solver function has been used to determine the value of dg and σg that minimize the criterion function, the optimal values for dg and σg can be used to determine the matrix of derivatives, D. The derivatives are first computed while first holding σg to its optimal value and varying dg by a very small amount, 1dg (say, 10−9 ). The derivative associated with each of the stage plus filter mass fractions is then Dik =

f (dg + 1dg , σg ) − f (dg , σg ) . 1dg

The same is performed while holding dg constant and varying σg by a small amount. Assuming all 8 stages of the personal cascade impactor are used, this will then develop a 9 × 2 matrix of derivatives, D. Likewise, a 9 × 9 matrix, W, can be shown in the spreadsheet as a 9 row by 9 column table with the weights along the diagonal and all 0s elsewhere. Likewise, the transpose of D, D0 , can easily be formulated in the spreadsheet. Matrix “manipulation” functions in the spreadsheet can then be used to multiply and invert the matrices in order to evaluate the 2 × 2 matrix (D0 WD)−1 . The diagonal values of this matrix are then multiplied by MSEw to determine the standard errors in the estimates of dg and σg . Finally, a 95% confidence interval about the estimates can be obtained by multiplying the square root of the variance values obtained by the t value at α/2 and n– p degrees of freedom, which in this case is t0.975,7 = 2.365.