Abstract
The ability to estimate mean particle size using simple, low-power optical instruments promises to greatly expand coverage of particle size measurements in the ocean and advance understanding of myriad processes from sediment transport to biological carbon sequestration. Here we present a method for estimating the mean diameter of particles in suspension from high-resolution time series of simple optical measurements, such as beam attenuation or optical backscattering. Validation results from a laboratory clay aggregation experiment show a good fit with independent mean particle diameter estimates in the 10–80 μm diameter range, with relative biases of 17%–38% and relative root mean square errors of 10%–24%. In the 80–200 μm range, quantitative validation data were not available, but our mean diameter estimates correlated strongly with particle settling rates.
© 2013 Optical Society of America
1. Introduction
Size is a fundamental characteristic of particles in the ocean, playing a critical role in how particles interact with each other and their environment. For example, organism size is a major determinant of ecosystem structure through its effects on metabolic rates and trophic interactions [1,2]. Particle size is also fundamental to biogeochemical cycling through its effect on sinking rate; large organic particles such as phytoplankton aggregates and fecal pellets sink rapidly, transporting carbon to the deep ocean and facilitating atmospheric drawdown [3–7].
Particle size is also critical to the understanding of sediment transport and deposition, primarily as a driver of sinking and aggregation rates within the benthic boundary layer. These rates significantly affect dispersal patterns of sediments and associated nutrients, carbon, and contaminants in coastal systems. Improvement in understanding of these sediment processes has been a major aim of scientists and engineers hoping to improve prediction of the effects of dredging, contaminant release, and seasonal riverine and other physical processes on water quality in coastal regions [8–11].
Widespread measurements of particle size could therefore greatly advance our understanding of multiple oceanic processes, from phytoplankton ecology to drawdown to sediment contaminant transport. Many instruments, both bench-top and in situ, are used for measuring size of marine particles, including microscopes, electrozone (Coulter) counters, flow cytometers, cameras, and laser diffraction systems. However, the size, cost, ease of use, power, and data requirements of these technologies limit the systematic measurement of particle size in the ocean.
Shifrin [12,13] developed an inversion method for estimating mean particle size in a water sample using fluctuations in beam transmissometer measurements caused by random distribution of particles. The low cost, moderate size, and widespread deployment of beam transmissometers make this an appealing method for expanding the coverage of particle size measurements in the ocean. However, Shifrin’s inversion has received little attention in the oceanographic literature and has not, to our knowledge, been tested using commercially available equipment. Here we revisit Shifrin’s inversion and present a similar but simpler inversion based on the variance-to-mean ratio (VMR; see Table 1 for a list of all abbreviations used in this paper) for estimating mean particle size using not only beam transmissometers, but also, unlike the Shifrin inversion, using smaller, low-power in situ backscattering and fluorescence sensors. Such sensors are already widely deployed on platforms such as moorings and autonomous floats and gliders, so the VMR inversion has the potential to greatly increase the coverage of in situ estimates of particle size at little or no added cost.
In Section 2 we present the VMR inversion, and in Section 3 we present validation data from a laboratory experiment, in which the particulate beam attenuation and particulate backscattering (via scattering at 124°) of aggregating and settling clay particles were measured over time [14]. We use both the VMR and Shifrin inversions to convert the variability in to mean particle size and apply the VMR inversion to as well. The estimated mean particle sizes vary by more than an order of magnitude during the experiment and agree closely with both independent estimates of size and measures of particle settling throughout most of this size range, with the VMR inversion performing slightly more consistently than the Shifrin inversion. In Section 4, we discuss the accuracy of the VMR inversion and its potential for field application.
2. Methods
A. VMR Inversion Equation
The VMR inversion is defined in this section for a time series of the particulate beam attenuation coefficient but can equally be applied to other additive optical properties, such as . The VMR inversion is used to estimate a weighted mean cross-sectional area of particles in suspension. is formally defined in Eq. (1), where is time, is the mean number of particles of class within a sensor’s sample volume and is the cross-sectional area of a particle of class . The VMR inversion is given in Eq. (2), where is the VMR inversion estimate of , and are the variance and mean of , is the sensor sample volume, is the attenuation efficiency of a particle (also known as the extinction efficiency ), is a correction for water movement during the integration time of a single measurement (), defined in Eq. (3), and is the ratio of particle residence time within the sample volume to [Eq. (4)]. The input parameters of Eq. (2) are described in detail in Subsections 2.A.1–2.A.4. Equations (1)–(4) are derived in Appendices A–C, given the assumptions that ambient particle concentration is constant over the period of a single VMR inversion calculation, particles are randomly distributed in space, all particles share the same (see Subsection 4.E.1 and Appendix A), and random fluctuations of the number and size of particles in are the only source of variance in . For ease of interpretation, we convert estimates of into an equivalent diameter by assuming a circular cross section [Eq. (5)]:
1. Time Series Statistics, and
The variance and mean are calculated over a finite period of time, during which changes in particle size and concentration should be minimal. If present, variance from sources other than the random distribution of particles in space should be removed to avoid a positive bias in the VMR inversion. Variance from low-frequency sources may be removed using a high-pass filter, as discussed in Subsection 2.C.2.
2. Sample Volume,
If the response of the sensor to a particle within its sample volume is the same regardless of particle location (homogeneous), then is simply the sensor’s geometric sample volume. We define for the nonhomogeneous case in Appendix B [Eq. (B7)].
3. Attenuation Efficiency,
The attenuation efficiency is defined as the fraction of light incident on the cross-sectional area of a particle that is scattered away from the detector or absorbed. For the evaluation dataset presented here, we use two different values as estimates of : first, literature values of to simulate field applications of the VMR inversion in which independent size measurements would not be available; and second, we calculate empirically, via Eq. (6) in order to isolate other potential sources of bias. Note that the VMR inversion assumes that all particles share the same , so that the mean in Eq. (6) applies to all particles. Implications of variable are discussed in Subsection 4.E.1:
4. Correction for Water Movement,
Movement of the sample relative to the sample volume during will increase the volume sampled, requiring a further correction factor , given in Eq. (3) (derived in Appendix C). For the validation dataset presented here, movement of the sample was slow relative to for the ac-9 and ECO BB, so was close to one (Table 2). However, the longer averaging time settings of the LISST, especially during Exp. 2 cause to drop as low as 0.35, becoming an important term in the VMR inversion.
B. Evaluation Dataset
The VMR inversion was tested using previously collected data from two laboratory clay aggregation experiments (Exp. 1 and Exp. 2) [14]. Setup and results of the two experiments were similar; however, during Exp. 1, a computer malfunction created a 10 min data gap, so the data presented in this section come from Exp. 2, unless noted otherwise. Initially, 4.2 g of bentonite clay were added to 120 L of reverse osmosis water in a 46-cm deep, wide black laboratory sink. Subsequently, 47.6 g of salt dissolved in a small volume of water was added to the sink to induce aggregation; the water in the sink was stirred for and not stirred thereafter. However, acoustic Doppler velocimeter measurements recorded particle speeds for the remainder of the experiment (data not shown), most likely due to convection. During the course of the 4 h experiments, the particle size distribution (PSD) was estimated using a Sequoia Scientific LISST-100 (type B; standard inversion; size range 1.4–230 μm diameter). Equation (1) was applied to the LISST PSD to calculate , which was then converted to an equivalent diameter using Eq. (5). The volume scattering function was measured at 124° [15] and 660 nm by a WET Labs ECO BBRT in Exp. 1 and by a WET Labs ECO BB3 in Exp. 2. The beam attenuation coefficient was measured with an acceptance angle of 0.93° and a path length of 10 cm by a WET Labs ac-9 (nine wavelengths between 412 and 715 nm; data reported here are at 650 nm for comparison with the ECO BB data unless otherwise noted). The LISST also measured , but with a smaller acceptance angle of 0.0269° and a shorter path length of 5 cm (670 nm only). The LISST was configured to sample at 0.25 Hz (Exp. 1) and 0.1 Hz (Exp. 2), the ac-9 at 6 Hz, and the ECO BB sensors at 1.1 Hz (Exp. 1) and 1.0 Hz (Exp. 2). Blank measurements were collected in the sink with both the ac-9 and the ECO BB sensors before addition of the clay. Mean blanks were subtracted from all subsequent measurements to yield and the particulate scattering coefficient . The equation was used to derive from at 124°, using [15].
The first of the experiment was dominated by aggregation (referred to hereafter as the “aggregation period”); increased by a factor of , accompanied by increased variability in and on time scales , while mean and on scales of minutes changed by less than 25% (Fig. 1). At 1.6 h (Fig. 1, gray dashed line), aggregates were recorded in the LISST’s largest size bin in significant numbers (10% of total aggregate area concentration), so after this point some aggregates likely exceed the LISST’s upper size limit and are not included in . The remainder of the experiment after 1.6 h was a period of rapid settling (referred to hereafter as the “settling period”), with decreases in mean and exceeding 70% (Fig. 1).
C. Data Processing for VMR Inversion
We calculated at 3-min intervals over the course of the experiment using Eq. (2). We similarly calculated by substituting for and for in Eq. (2). In this section we describe in detail the steps taken to calculate and from the time series of and .
1. Bin Width of VMR Calculation
The VMR inversion was calculated in 3-min bins to minimize the change in over a single bin while still obtaining an acceptable VMR precision ( relative standard deviation). If rare particles () do not dominate the or signal, then the VMR precision can be approximated by the precision of the sample variance of a normal distribution [Eq. (7)], where is the number of independent samples and is the standard deviation of the sample VMR [16]. Here can be approximated by the ratio of the bin width to or (whichever is larger) yielding of 32–45 for all sensors except the LISST in Exp. 2, whose 10-s integration time reduced to 18. Theoretical minimum relative VMR precisions calculated via Eq. (7) and converted to relative or precisions are given in Table 3. In addition to providing sufficient precision, the bin width should also be greater than the filter width (see Subsection 2.C.2) so as not to reduce variance more than intended:
2. Detrending Filter
To minimize variance from sources other than particle size, (e.g., small scale patchiness or nonrandom distribution of particles), and time series were detrended by subtracting a filter [e.g., Fig. 2(a), red line] consisting of a running median followed by a running mean, both with windows of 50 s. This filter was chosen after experimentation with various filters on simulated time series. The 50-s window is approximately 10 times (see Subsection 2.C.3). In simulated time series, this detrending filter with window 10 times either or the sample interval (whichever is larger) preserved of the variance due to particle size, while minimizing bias due to any variability in concentration on scales of minutes and longer (data not shown). Note that in this study is calculated using autocorrelations in the detrended data (after applying the filter), so a guess or alternative estimate of must be used when deciding the initial filter width. If the subsequent autocorrelation-based estimates (Subsection 2.C.3) are very different from the initial guess, then a new detrending filter can be applied and a new calculated until the filter width is close to 10 times .
3. Estimation of Residence Time
To estimate , temporal autocorrelations were calculated on the detrended time series at 15 min intervals using Eq. (8), where is lag time. Strong positive autocorrelations () were found at for both ac-9 and measurements, and residence time was estimated as twice the lag time corresponding to . The median estimated residence times across all 15 min bins for each of the two experiments were 4.3 and 4.6 s for ac-9 and 5.6 and 4.5 s for (Table 2). We estimate the mean path lengths of particles passing through the sampling volumes of the ECO BB and ac-9 as 1 cm and 8 mm, respectively, so these residence time estimates are consistent with acoustic Doppler velocimeter measurements (measured in Exp. 2 only) indicating particle movement throughout the course of the experiment (data not shown). The LISST time series were not strongly autocorrelated due to lower sampling frequency, preventing estimation of via Eq. (8). Instead, ac-9 was scaled to the LISST’s smaller weighted beam width (7 mm versus 8 mm: ) to yield estimates of 3.8 and 4.0 for Exps. 1 and 2, respectively.
4. Outlier Removal
To decrease uncertainty in the VMR, large outliers in the detrended time series [Fig. 2(c), red circles] with absolute values greater than a threshold of 4 times the interquartile range of each 3 min bin [Fig. 2(c), green lines] were eliminated from both the raw and detrended time series ( of data points removed from each time series). Note that very large particles can cause positive outliers, which, if present in sufficient numbers, may be analyzed separately [17–19].
5. Calculation of Mean and Variance
After outlier removal, mean and were calculated for each bin [Fig. 2(b)] from the unbinned time series [Fig. 2(a), blue line]. Variances [Fig. 2(d)] were calculated from the detrended time series [Fig. 2(c), blue line] for the same bins.
6. Estimation of and
Empirical estimates of and over the course of the experiment (Fig. 3) were obtained by dividing mean values of and by the total cross-sectional area concentration derived from the LISST PSD [Eq. (6)]. Patterns in and were similar during both experiments, but was consistently higher during Exp. 2, most likely due to uncertainty in the calibrations of the two ECO BB sensors [20]. Divergence between LISST and ac-9 after the first hour is likely due to the ac-9’s larger acceptance angle, which causes it to underestimate somewhat the attenuation due to larger particles [21]. During the settling period the largest aggregates likely surpassed the maximum LISST size bin, leading to an under-estimate of total area concentration and a positive bias in and estimates. During the aggregation period (Fig. 3, left of vertical dashed lines), changes in estimated and were small () despite a large increase in (Fig. 1), so the median and from the aggregation period (Table 2) were applied to the entire time series. The use of a single also ensures that any trends in and are entirely independent from trends in the LISST PSD. Note, however, that the VMR inversion is most valuable in applications where PSD measurements are not available, in which case literature values of and must be used in lieu of these empirical estimates, adding uncertainty to the VMR inversion.
In order to simulate an application in which could not be determined empirically, we chose the theoretical asymptotic value for [22]. is less well constrained in the ocean, but we can set a likely range of based on literature values of marine particulate backscattering ratios (0.005–0.035) [23,24] and single-scattering albedos (0.8–1.0) [25]. In this case, a choice of (; ) constrains additional uncertainty in within a factor of 3, equivalent to within a factor of 1.7 (Table 2). In this study we evaluate overall VMR inversion performance using these literature values of and (Figs. 4 and 5), but we use our empirical estimates of and for the remaining analyses (Figs. 6–8) in order to explore other possible sources of bias.
7. Instrument Sampling Volume and Integration Time, and
We used an ac-9 sample volume of 5.0 ml based on manufacturer specifications (beam diameter: 0.8 cm; path length: 10 cm). The ECO BB and LISST sample volumes (0.62 and 1.9 ml, respectively) were calculated using Eq. (B7) to account for inhomogeneous sensing efficiency within the sample volume (see Appendix B). The ac-9 at each wavelength is 3 ms (WET Labs ac-9 User’s Guide revision T), and we approximate the ECO BBRT and BB3 as the inverse of the sampling frequency minus a 50 ms overhead time: 0.87 and 0.95 s in Exp. 1 and Exp. 2, respectively. LISST is simply the averaging period: 4 and 10 s in Exps. 1 and 2, respectively. These values of , combined with the values of reported in Subsection 2.C.3, yield correction factors of 0.9998 for the ac-9, and 0.95 and 0.93 for the ECO BBRT and ECO BB3, and 0.65 and 0.35 for the LISST in Exps. 1 and 2, respectively (Table 2).
3. Results
A. VMR Inversion and Comparison with LISST
Patterns of (Fig. 4, black and red lines) and (Fig. 4, blue line) were similar. At the start of both experiments, there was a factor of discrepancy between and both and . This discrepancy disappeared gradually over the first 0.75 h. This “period of initial bias” is discussed in Subsection 4.A.2 and is excluded from the validation regressions in Fig. 5. During the “validation period” from hour 0.75 (Fig. 4, solid gray line) to the start of the settling period (Fig. 4, dashed gray line), agreed closely with both and . Validation regressions of both and against from this period (Fig. 5, bold colors) were highly linear () and slopes were close to . Mean relative biases of and ranged from to (Table 3). When literature and were replaced by empirical estimates (Table 2), relative biases decreased to (Table 3). During the settling period beginning at 1.6 h (dashed gray lines in Fig. 4), and continued to increase to a peak of approximately 200 μm at 2.2 h and then declined to by the end of the experiment (Fig. 4), presumably as the largest aggregates settled out of suspension. remained lower and more constant during the settling period, never exceeding 100 μm, likely due the LISST’s inability to detect particle sizes . During the settling period, the LISST area PSD was sharply peaked at the largest size bin, strongly suggesting the presence of larger aggregates, but a “stranded” population of particles maintained below 100 μm [14].
B. Spectral Dependence of
A power law of the form (where is wavelength) was fit to the ac-9 time series, and the resulting spectral slope varied between 1.35 and 1.6 during the experiments [14]. Both and the empirical (calculated separately at each wavelength) have the same spectral slope: . If the VMR inversion assumptions hold, the spectral slope of should be twice . These spectral dependencies should cancel out in Eq. (2), so , the spectral slope of , should be zero (i.e., is accurate at all wavelengths), regardless of . However, during the first 0.75 h of the experiment, coincident with the disagreement between and , we observed a strong spectral dependence of (Fig. 6). The spectral slope was steepest () at the start of the experiments, increased to zero (no spectral dependence) in the first hour, and remained within for the remainder of the experiments (Fig. 6). The low during the first 0.75 indicates a violation of at least one VMR inversion assumptions and provides further justification for rejecting the early VMR inversion results. We discuss possible sources of both the low and the early disagreement between and in Subsection 4.A.2.
C. Relationship Between and Clearance Rate
The rate at which water at the depth of the sensors is cleared of particles is a function of sinking rate (as well as vertical concentration gradient and convective mixing rate). The rate of relative decrease in per minute, [Eq. (9)] is a measure of particle clearance rate. Because particle size is a primary driver of sinking rate [26], we expect that will be positively related to ; Fig. 7(a) shows that it is. The time series was smoothed using a second-order polynomial spline with a 70 min window before the calculation of . The first derivative of the smoothed time series was approximated at time (in min) using a 6 min window and normalized to to yield [Eq. (9)]:
was further averaged into 9-min bins to reduce variability. There was a strong positive relationship between and all estimates of during the aggregation phase [Fig. 7(a), plus sign]. During the settling phase [Fig. 7(a), open circles], the relationships of with and change little while the relationship with changes sharply. Figure 7(a) supports the idea that the LISST PSD misses large particles during the settling period and provides qualitative support for the validity of when the LISST could not be used for quantitative validation. For comparison, the spectral slopes of and , two other proposed proxies for particle size [27,28], did not maintain strong, consistent relationships between either [14] or [Fig. 7(b)] during this experiment. The relationship between and is discussed further in Subsection 4.A.3.
D. Comparison with Shifrin’s Inversion
The Shifrin inversion is given in Eq. (10), adapted from Shifrin [13] (his Eq. 6.28), where is the transmittance, defined as the ratio of transmitted light in the presence of particles to transmitted light in the absence of particles, is the attenuance, defined as , is the cross section of the transmissometer beam, and is a function of , given in Eq. (11) (adapted from Shifrin [13]; his Eqs. 6.18, 6.21) and is an integration variable in Eq. (11):
Note that Shifrin’s [13] Eq. 6.18 contains a misprint, but the correct equation is in an earlier paper [12] (his Eq. 3).
As with the VMR inversion, we detrended and removed outliers before calculating the variance and calculated the Shifrin inversion for 3 min bins. The integral in Eq. (11) was numerically approximated using a quadrature algorithm (MATLAB, “integral.m”). As before, was converted to using Eq. (5).
The magnitude and patterns in and were very similar, including the initial bias and the disagreement with during the settling period, but the two methods diverged at higher values of . Both and were calculated from all ac-9 wavelengths during the validation period using wavelength-specific empirical estimates and divided by corresponding values of to estimate bias. Mean relative biases were 1.09 for and 0.91 for . Type I linear regressions against show a significant negative correlation between and bias in but not between and bias, causing the two inversions to converge as approached zero (Fig. 8). This result is discussed in Subsection 4.C.
4. Discussion
A. VMR Accuracy
The strong agreement between , , and throughout the experiment (Fig. 4) indicates that the VMR inversion is reliable and robust to differing principles of instrument operation, differing sampling volumes, differing bin averaging, and even differing transmissometer acceptance angles (even though the larger acceptance angle of the ac-9 causes it to underestimate due to the larger particles [14]). However, the VMR inversion does not always agree with . We asses the accuracy of the VMR inversion (and the LISST) during three distinct periods.
1. Validation Period
During the validation period, and were strongly correlated with over nearly an order of magnitude (10–80 μm; —log scale) with slopes close to and mean biases of 17%–38% (Fig. 5; Table 3). The largest single cause of bias is uncertainty in and ; use of direct estimates (Table 2) reduced biases to . The causes of these remaining small but significant mean biases are unknown, but may include uncertainty in , slightly nonrandom particle distribution, a small size-dependence of and , instrument noise proportional to the signal, and/or a slight low bias in due to particles . Given the uncertainties involved in estimating particle size by any method, these results are encouraging. For comparison, mean discrepancies of 11% were found between median particle diameters measured by a LISST and a Coulter counter across a wide variety of particle types [29]. Biases of up to 38% may often be acceptable in natural systems where varies by a factor of 2 or more. Furthermore, the strong relative agreement between , , and (regression slopes on log scale: 1.02–1.16) supports the use of the VMR inversion as a relative size indicator, regardless of the choice of or .
2. Initial Bias
The assumption that a single value of or applies to all particles may introduce some bias in the VMR inversion, but it is unlikely that variability in efficiencies can explain the large (factor of 6) and ephemeral bias, observed at the start of the experiments, given the stability of and estimates throughout the aggregation phase (Fig. 3), indicating minimal relationship between or and size.
The discrepancy between and both and during the first 0.75 h is most likely caused by an additional source of variance other than particle size, leading to a positive bias in the VMR inversion. This extra variance apparently decreased over most of this early period [e.g., Fig. 2(d)], causing corresponding decreases in and (despite increasing particle size). The strong negative associated with this extra variance indicates that its spectral slope is likely flat or even negative (higher at higher wavelengths and opposite to ). The cause of this extra variance is not clear. Nonrandom (patchy) distribution of particles at the start of the experiment would increase variance. However, such clustering would cause the same relative bias at all wavelengths and cannot explain the negative (Fig. 6). The initial mixing of salt into the sink would have caused inhomogeneity in index of refraction (“schlieren”), which could lead to fluctuations in , but this effect would also increase mean and interfere with the LISST PSD inversion [30]. The source of the variance could be electronic noise, but it is unclear why it would affect all sensors in both experiments in such similar ways.
Regardless of its source, the relative effect of added variance is greatest at low , when variance due to particle size is small; after 0.75 h, the rapid increase in variance due to increasing likely overwhelmed any remaining bias.
3. Settling Period
During the settling period, the presence of aggregates too large for the LISST inversion to constrain prevented quantitative assessment of VMR inversion accuracy. However, the strong positive relationships between , , and during this period [Fig. 7(a)] provide qualitative validation that and continue to track particle size. Perhaps more importantly, this relationship demonstrates the potential of the VMR inversion as a tool for studying the dynamics of particle settling in the ocean.
B. VMR Precision
The time series was much smoother than the and time series (Fig. 4), despite much shorter averaging times (10 s versus 3 min), suggesting much higher precision. The root mean square errors (RMSEs) of the validation regressions should therefore be close to the theoretical precision of the VMR, plus the effect of any nonlinearity in the regressions. The RMSEs of the validation regressions correspond to relative precisions of 10% and 12%, close to the theoretical relative precision of 11% (Table 3). and relative RMSEs were higher than theoretical values by 2–8 percentage points, perhaps due to added uncertainty from inhomogeneous sample volumes. These results support the use of Eq. (7) to provide an initial estimate of VMR precision for planning sampling strategy, although sensor-specific details may decrease precision somewhat, in which case wider bins may be needed to achieve target precision.
C. Comparison with Shifrin
The performance of the Shifrin inversion was similar to that of the VMR inversion; agreement between the two was usually within 20% (Fig. 8) and improved at longer wavelengths and lower . However, the disagreement at high values of is a result of differing inversion assumptions, which are discussed in this section.
While the VMR inversion attributes entirely to random fluctuations in the number of particles in , the Shifrin inversion includes an additional source of variance: fluctuations in the degree of particle-particle shading in the beam (see Appendix D). Consequently, Shifrin predicts a higher for a given and , or a lower for a given and . The difference between and is primarily a function of attenuance , because at higher , particle-particle shading is greater and fluctuations in shading have a larger potential effect on . If variable shading indeed contributes to measured , we expect to be more robust than to changes in . Instead, we found that was significantly correlated with bias in but not (Fig. 8).
This result suggests that the effect on of variations in shading within the sample beam was significantly smaller than that predicted by Shifrin. The irregular, elongated shapes [14] and possible semi-transparency of the aggregates in this study may account for some of the reduction in variability versus Shifrin’s assumption of solid, circular particle cross sections. Additionally, the fluctuations that were produced by shading, especially of the smallest particles, may have occurred at frequencies too high for the ac-9 to resolve. Therefore, while the VMR inversion appears to be a more consistent method for estimating from fluctuations in this experiment, the Shifrin inversion may perform better in conditions of more solid particles, slower relative particle speed, or shorter sample integration times. High variance on time scales much shorter than is expected in such cases, however, and the VMR inversion may still be used if this higher-frequency variance is identified and filtered out.
D. Comparison with Spectral Slope
The spectral slopes of and have been proposed as in situ optical proxies for particle size [27,28]. Spectral slopes provide certain advantages over the VMR method: no single must be assumed, can be calculated over shorter temporal and spatial windows than the VMR, is not biased by noise or other variance, and ocean color measurements may allow estimation of from space [28,31]. On the other hand, requires more complex, multiwavelength sensors and conversion of to particle size requires assumptions of Jungian PSD shape and homogeneous spherical particles, often violated in natural particle assemblages. During this study, the VMR inversion performed qualitatively better during the settling period and as a whole over the entire experiment, whereas (although not ) performed qualitatively better during the initial 0.75 h. Ultimately, our results show that each method can be sensitive to factors other than mean particle size, and more work needs to be done to understand and quantify these other effects. Then, the combination of both methods may provide enhanced information about different parts of the PSD or perhaps about particle composition.
E. Considerations for Field Deployment
1. Choice of
The literature values of and presented in this study (Table 2) are intended as first guesses to apply to particles of unknown composition sampled using visible wavelengths. If local particle optical properties are better constrained, site-specific values of may increase VMR inversion accuracy. Alternatively, if the ratio of or to another physical particle measurement (e.g., suspended particulate mass or particulate organic carbon concentration) is approximately constant, this ratio may be substituted for in Eq. (2), yielding an estimate of weighted mean particle mass, for example, in place of .
It is important to note that even if bulk is well constrained, a strong size dependence of within a single sample violates an assumption of the VMR inversion. In this case, it is more accurate to use Eq. (A11), yielding , the mean weighted by optical cross section and not physical cross section (see Appendix A). still provides useful bulk particle size information, but it is more difficult to validate using PSD measurements and may not correlate as strongly with particle clearance rate. In this study, neither nor varied strongly with size, so both and , estimates were accurate within 40% or better. Some work to date has shown that appears well constrained in situ, in accordance with theory [32], but is likely more variable, especially for mixed populations of mineral and organic particles. Ultimately, in situ field studies are needed to test the accuracy of the VMR inversion in natural particle populations.
2. Uncertainty in
In sampling regimes in which is smaller than the interval between consecutive samples (as with the LISST in Exp. 2), temporal autocorrelation will be minimal, preventing the estimation of via autocorrelation (Subsection 2.C.3). If possible, therefore, sampling frequency should be increased or platform speed decreased until a temporal autocorrelation is observed. If this is not possible, but the velocity of the sample relative to the sensor is known, may be calculated based on this velocity and sensor geometry. Otherwise, uncertainty in will propagate to uncertainty in or . However, as long as is constant, VMR inversion relative size estimates will still be useful.
3. Uncertainty Due to
In order to apply the VMR inversion to vertical profile data containing strong gradients in and , profiling speed should be slow enough and sampling frequency high enough that within a layer of interest to achieve a theoretical relative VMR precision of via Eq. (7). Alternatively, the results of multiple profiles through the same water mass may be combined to increase precision.
4. Unconstrained Variance and the Limits of Detectable Size
The dimensions of set a hard upper limit on the size of particles detectable by the VMR inversion. Particles whose diameters approach (or exceed) the smallest dimension of will often (or always) occur only partially within , “appearing” as smaller particles to the VMR inversion. In practice, however, the largest size of particle included in the VMR inversion may be set by the likelihood of encountering a particle of a given size in during the course of a single VMR calculation. In this case, large particles too rare to be sampled reliably should be explicitly excluded using an outlier filter, as explained in Subsection 2.C.4.
The lower limit of the size of particles that contribute to scattering (and therefore to VMR inversion estimates of ) is set by the relationship between and (or . In general, decreases with decreasing for (i.e., for visible wavelengths), although the precise lower size limit will also depend on particle shape and index of refraction [33]. In addition, for any given mean signal, there will be a lower limit of detectable , set by the variance due to sources other than particle size (e.g., patchy distribution of particles, schlieren [30], or instrument noise). Such “extra” variance, along with uncertainty in , has the potential to cause significant bias in field application of the VMR inversion. During the first 0.75 h of this experiment, unconstrained noise set the lower limit of detectable at . Further studies will be necessary to better understand and constrain potential sources of variance and determine this lower limit in field applications. If variance from a particular source is independent of size and well characterized, it can be subtracted from or before calculating the VMR inversion, lowering the limit of detectable .
5. Optimal Sample Volume
In general, larger sample volumes will be more effective at capturing larger particles and smaller sample volumes will be better at distinguishing smaller sizes. An ideal would therefore be just large enough to reliably capture the largest particles of interest, both physically and statistically (see previous section). In some cases, a pair of sensors with different sample volumes may yield extra information about the PSD, as the smaller volume will exclude the larger particles and be more sensitive to smaller particles. However, other factors, such as sensor sensitivity and signal-to-noise ratio, sample duration, internal averaging and platform speed may be just as important as for determining sensitivity to small particles and should all be taken into account when choosing an appropriate sampling strategy for use with the VMR inversion.
5. Conclusions
We find that the VMR inversion consistently produced accurate estimates of in the range of 10–80 μm during two laboratory clay aggregation experiments. At higher sizes (up to 200 μm), the VMR inversion yielded plausible results that retained the expected correlation with clearance rates, but these higher estimates could not be quantitatively validated because the larger aggregates likely exceeded the size range of the LISST inversion. For , VMR estimates were biased high, likely due to an additional unconstrained source of variance at the start of the experiment.
Testing the VMR inversion in situ will be important to determine the effect of the added complexity of marine environments on inversion accuracy. However, in light of the large range of particle sizes that contribute to optical backscattering and beam attenuation in the ocean, from submicrometer bacteria to aggregates of several millimeters or greater, we believe that the VMR inversion is likely to provide useful mean particle size information in a variety of marine environments, even where sources of bias exist. If so, this inversion, combined with the widespread current and future deployment of small optical sensors, promises to greatly expand the spatial and temporal coverage of marine particle size estimates at little extra cost, leading to increased understanding of particle size in the ocean and the impact of size on the many processes that particles mediate.
Appendix A: Derivation for an Ideal Sensor
The VMR inversion can be applied to any additive optical property. Here, as in Section 2, is used as an example. The VMR inversion rests on the assumption that particles are randomly distributed in space. For the derivation in this section, let us also assume an ideal sensor, which can measure instantaneously the particulate attenuation at time within the sample volume . Let be the portion of due to a single class of identical particles. Equation (A1) defines [22] as a function of , the number of particles of class within , and the attenuation cross-section of a single particle:
In the case of multiple particle classes, is simply the sum of from all particle classes:
Similarly, the expected value (mean) of , , can be expressed in terms of :
The equality (where is constant), and the equality (where and are independent, random variables) allow us to similarly express the variance in terms of :
The assumption that particles are randomly distributed in space allows us to approximate as a Poisson random variable, in which case
After substituting Eq. (A5) into Eq. (A4), the VMR becomes
where is the mean , weighted by the total attenuation cross-section of each class:Particle cross-sectional area can be related to through the attenuation efficiency factor , defined as the fraction of light incident on a particle that is attenuated:
In the event that is constant across all particle classes in suspension, we can rewrite Eq. (A7) as
where is the mean cross-sectional area, weighted by the total cross-sectional area of each size class [Eq. (1) in text]. We obtain , the VMR inversion estimate of by substituting Eq. (A9) into Eq. (A6) and solving for :Equation (2) in the text is derived from Eq. (A10), combined with the correction of for water movement derived in Appendix C. In the event that varies strongly with particle class we may replace Eq. (A10) with
where is the mean , weighted by and is a weighted average of that depends on a nonlinear interaction between the PSD and but must fall between the minimum and maximum in suspension.Appendix B: Correction of for Inhomogeneities in Sensing Efficiency
In practice, the contribution of a particle within to a measurement of may depend not only on and , but also on the position within . This may result from inhomogeneous illumination of , e.g., by a noncollimated or Gaussian beam source, or by variation in the efficiency of the detection of light scattered from differing positions, e.g., as a function of distance from the detector. Such deviations, combined with the random distribution of particles, will add random variability to measurements of , introducing a positive bias in . This bias can be corrected by modifying as follows.
Let be the full geometric sample volume of a sensor, and let a homogeneous scattering medium fill (see Table 4 for a list of additional terms introduced in the Appendices). Let be the signal per unit volume coming from point , , , normalized to the total signal. It follows that the spatial integral of within the limits of is equal to one:
As a particle moves around the sample volume of a nonhomogeneous sensor, its will appear to change. If we ignore any variations in over the spatial extent of a single particle, then we can define , the “apparent” at point , , , as follows:
We can now treat each particle position as a separate particle class and replace in Eq. (A6) with , yielding
where is defined asSubstituting Eqs. (B2), (B1), and (A7) into Eq. (B4) we obtain
Substituting Eq. (B5) into Eq. (B3), we obtain
Therefore, in the case of a nonhomogeneous sample volume we can simply define as a “weighted sample volume”, given by
and continue to use Eq. (A6) and all of its derivatives, including the VMR inversion [Eq. (A10)].For the ac-9 we assume homogeneous sensing efficiency and use (beam diameter: 0.8 cm; path length: 10 cm). We calculate LISST of 1.9 ml via Eq. (B7) using manufacturer specifications: a Gaussian beam with a waist radius of 3.5 mm, cut off at a radius of 5 mm, and a path length of 5 cm. Assuming minimal change in beam width along the path length, we calculate as a function of radius from the beam center:
and substitute into Eq. (B7).For the ECO BB sensors, we obtained an estimate of using Eq. (B7), with defined as
where is the power received by the detector of a scattering meter from elemental volume , centered at point .We calculated following Sullivan et al. [15], using the Sullivan and Twardowski [34] phase function and the following sensor geometry: angles between center of source or detector cone and sensor face: 49° and 73°; source and detector cone half angles: 17.5° and 42.5°; distance between source and detector: 1.7 cm (Mike Twardowski, pers. comm.). The calculation was modified from Sullivan et al. [15] to correct an error in the original equation, replacing the solid angle of viewed from the detector with the solid angle of the detector viewed from . We also replaced the beam attenuation term with the absorption coefficient of water ( at 660 nm), because sensitivity analysis showed that attenuation by particles along the path of light in the ECO BB sensors will simultaneously decrease both and apparent , with minimal net effect on the VMR inversion.
Appendix C: Correction for Water Movement
Each measurement taken by an optical sensor must be integrated over a finite duration . Movement of the water relative to the sensor causes , the total geometric volume of water sampled during , to increase beyond . We may derive , the weighted volume sampled over , from using Eq. (B7). In this section we find the ratio of to in terms of and the residence time of a particle in for the simplified case in which the sample moves at a constant velocity in direction relative to and the cross section of in the plane is constant along the -dimension. Let us also assume that the sensor response within is homogeneous () so that the relative contribution of a point along to the total signal is determined only by the length of time that point spends in during . In this case we can rewrite Eq. (B7) as
where is the length of in the -direction. We can express in terms of as follows: and we can divide into the length of in the -direction and the length that the sensor moves relative to the sample during . Let us now define a correction factor in terms of , , and as follows:There are two solutions for Eq. (C3), corresponding to and . At the two solutions converge.
For , we can divide into three sections, defined in relation to and :
The first section represents the water that is within at the start of a sample. The second section represents water that is at or beyond the leading edge of at the start of the sample and passes through the entire length of during . The third section represents water that is within at the end of a sample. We can obtain the following solutions for the two integrals in Eq. (C3) for after substituting Eq. (C4):
andSubstituting Eq. (C5) and Eq. (C6) into Eq. (C3) we obtain the solution
We can find in a similar fashion for as well. In this case, becomes
Equation (C8) is identical to Eq. (C3), except that all instances of and have been reversed, and replaces . Substitution into Eq. (C3) and integration in this case yields
Let us now define as the ratio of to and note that must also equal the ratio of to [Eq. (4) in text]:
We can now combine Eqs. (C7) and (C9) and express them in terms of :
Note that is continuous for .
In order to obtain the full VMR inversion equation, corrected for water movement, we must replace with in Eq. (A10). The expression may then be substituted for [Eq. (C3)] to obtain Eq. (2) from the text.
Appendix D: Effect of Particle-Particle Shading
Implicit in Eq. (A1) and all derivations that follow is the assumption that the contribution of each particle to is not changed by the presence of other particles. This is a fundamental assumption of any measurement of an inherent optical property such as or . However, violation of this assumption will affect the VMR inversion differently than it affects or , so a brief discussion is warranted.
In a beam transmissometer, the light attenuated by a particle is detected by the reduction in source light that reaches the sensor. If this particle is in the shadow of another particle, this reduction in detected light is decreased. If particles are randomly distributed in space, the calculation of takes into account the mean effect of shading on the detected signal. However, random variations around this mean—not accounted for by the VMR inversion—will cause a positive bias when present by adding to the variance. The magnitude of variations in caused by random particle-particle shading will depend on the concentration, size, shape, and transparency of particles in suspension. The time scale of these variations will depend on the particles’ size and relative speed. Shifrin’s (1988) particle size inversion relies on the same fundamental assumptions as the VMR inversion, but additionally incorporates the effect of variation in particle-particle shading for completely opaque spherical particles (i.e., the maximum possible effect). Therefore, if all other assumptions of both size inversions are met, but if particles are not completely opaque and spherical or if some of the fluctuations caused by particle-particle shading occur too rapidly to be detected by the transmissometer, then the true mean particle size will fall somewhere in between the estimates of the Shifrin inversion and the VMR inversion. The effect of particle-particle shading on the validation dataset presented in this paper is discussed in Subsection 4.C.
A scattering meter detects light that travels from the light source to a particle, where it is scattered and then to the sensor. Absorption by other particles along this path will decrease the scattered light that is detected. At very high particle concentrations, or for scattering meters with long path lengths, an absorption correction factor may be necessary to correct the mean scattering measurement. If uncorrected scattering measurements are used, absorption will induce a low bias in the VMR inversion. If an absorption correction is applied to scattering measurements before the VMR inversion, accounting for mean shading, there may be a small high bias in the VMR inversion due to variability in this shading. While absorption was not measured during this study, any particle-particle shading effect was likely minimal due to the short path length of the backscattering instruments.
Inversion development and analysis were funded by a NASA fellowship (grant NNX-10AP29H to Briggs), the NASA Ocean Biology and Biogeochemistry (OBB) program (grant NNX-13AC42G), and the National Science Foundation (grants OCE-0628107 and OCE-0628379); laboratory experiments were funded by the Office of Naval Research (contract N00014-10-1-0508). The authors thank the two reviewers, Michael Twardowski and Giorgio Dall’Olmo, for their thorough and insightful comments, which greatly improved this manuscript. Thank you also to George Jackson and Ivona Cetinić for helpful discussions and comments and to Clementina Russo and Jim Loftin for help with the laboratory experiments.
References
1. L. Dickie, S. Kerr, and P. Boudreau, “Size-dependent processes underlying regularities in ecosystem structure,” Ecol. Monogr. 57, 233–250 (1987). [CrossRef]
2. S. Jennings and K. J. Warr, “Smaller predator-prey body size ratios in longer food chains,” Proc. R. Soc. B 270, 1413–1417 (2003). [CrossRef]
3. V. S. Smetacek, “Role of sinking in diatom life-history cycles: ecological, evolutionary and geological significance,” Mar. Biol. 84, 239–251 (1985). [CrossRef]
4. A. L. Alldredge and M. W. Silver, “Characteristics, dynamics and significance of marine snow,” Prog. Oceanogr. 20, 41–82 (1988). [CrossRef]
5. V. L. Asper and W. O. Smith, “Abundance, distribution and sinking rates of aggregates in the Ross Sea, Antarctica,” Deep Sea Res., Part I 50, 131–150 (2003).
6. G. A. Jackson, “Effect of mixed layer depth on phytoplankton removal by coagulation and on the critical depth concept,” Deep Sea Res., Part I 55, 766–776 (2008).
7. G. Fischer and G. Karakas, “Sinking rates and ballast composition of particles in the Atlantic Ocean: implications for the organic carbon fluxes to the deep ocean,” Biogeosciences 6, 85–102 (2009).
8. J. C. Winterwerp and W. G. M. van Kesteren, Introduction to the Physics of Cohesive Sediment Dynamics in the Marine Environment (Elsevier, 2004).
9. G. C. Chang, T. D. Dickey, and A. J. Williams, “Sediment resuspension over a continental shelf during Hurricanes Edouard and Hortense,” J. Geophys. Res. 106, 9517–9531 (2001). [CrossRef]
10. J. R. Gray and J. W. Gartner, “Technological advances in suspended-sediment surrogate monitoring,” Water Resour. Res.45 (2009). [CrossRef]
11. K. J. Curran, P. S. Hill, and T. G. Milligan, “Fine-grained suspended sediment dynamics in the Eel River flood plume,” Cont. Shelf Res. 22, 2537–2550 (2002). [CrossRef]
12. K. S. Shifrin, B. Z. Moroz, and A. N. Sakharov, “Determination of characteristics of a dispersion medium based on its transparency data,” Dokl. Akad. Nauk SSSR 199, 589–591 (1971).
13. K. S. Shifrin, “The method of fluctuations,” in Physical Optics of Ocean Water (American Institute of Physics, 1988), pp. 224–232.
14. W. Slade, E. Boss, and C. Russo, “Effects of particle aggregation and disaggregation on their inherent optical properties,” Opt. Express 19, 7945–7959 (2011). [CrossRef]
15. J. M. Sullivan, M. Twardowski, J. R. Zaneveld, and C. Moore, “Measuring optical backscattering in water,” in Light Scattering Reviews 7, A. A. Kokhanovsky, ed. (Springer, 2013), pp. 189–224.
16. J. E. Freund, Mathematical Statistics, 2nd ed. (Prentice-Hall, 1971), pp. 212–216.
17. W. D. Gardner, M. J. Richardson, and W. O. Smith, “Seasonal patterns of water column particulate organic carbon and fluxes in the Ross Sea, Antarctica,” Deep Sea Res., Part II 47, 3423–3449 (2000).
18. J. K. B. Bishop and T. J. Wood, “Particulate matter chemistry and dynamics in the twilight zone at VERTIGO ALOHA and K2 sites,” Deep Sea Res., Part I 55, 1684–1706 (2008).
19. N. Briggs, M. J. Perry, I. Cetinić, C. Lee, E. D’Asaro, A. M. Gray, and E. Rehm, “High-resolution observations of aggregate flux during a sub-polar North Atlantic spring bloom,” Deep Sea Res., Part I 58, 1031–1039 (2011).
20. G. Dall’Olmo, E. Boss, M. J. Behrenfeld, and T. K. Westberry, “Particulate optical scattering coefficients along an Atlantic Meridional Transect,” Opt. Express 20, 21532–21551 (2012). [CrossRef]
21. E. Boss, W. H. Slade, M. Behrenfeld, and G. Dall’Olmo, “Acceptance angle effects on the beam attenuation in the ocean,” Opt. Express 17, 1535–1550 (2009). [CrossRef]
22. C. F. Bohren and D. R. Huffman, Absorption and Scattering of Light by Small Particles (Wiley, 1983).
23. E. Boss, W. S. Pegau, M. Lee, M. Twardowski, E. Shybanov, G. Korotaev, and F. Baratange, “Particulate backscattering ratio at LEO 15 and its use to study particle composition and distribution,” J. Geophys. Res. 109, C01014 (2004). [CrossRef]
24. M. S. Twardowski, E. Boss, J. B. Macdonald, W. S. Pegau, A. H. Barnard, and J. R. V. Zaneveld, “A model for estimating bulk refractive index from the optical backscattering ratio and the implications for understanding particle composition in case I and case II waters,” J. Geophys. Res. 106, 14129–14142 (2001). [CrossRef]
25. M. Babin, A. Morel, V. Fournier-Sicre, F. Fell, and D. Stramski, “Light scattering properties of marine particles in coastal and open ocean waters as related to the particle mass concentration,” Limnol. Oceanogr. 48, 843–859 (2003). [CrossRef]
26. W. Dietrich, “Settling velocity of natural particles,” Water Resour. Res. 18, 1615–1626 (1982). [CrossRef]
27. E. Boss, M. S. Twardowski, and S. Herring, “Shape of the particulate beam attenuation spectrum and its inversion to obtain the shape of the particulate size distribution,” Appl. Opt. 40, 4885–4893 (2001). [CrossRef]
28. T. S. Kostadinov, D. A. Siegel, and S. Maritorena, “Retrieval of the particle size distribution from satellite ocean color observations,” J. Geophys. Res. 114, C09015 (2009). [CrossRef]
29. R. A. Reynolds, D. Stramski, V. M. Wright, and S. B. Woźniak, “Measurements and characterization of particle size distributions in coastal waters,” J. Geophys. Res. 115, 1–19 (2010). [CrossRef]
30. O. Mikkelsen, T. Milligan, and P. Hill, “The influence of schlieren on in situ optical measurements used for particle characterization,” Limnol. Oceanogr. 6, 133–143 (2008).
31. H. Loisel, J. M. Nicolas, A. Sciandra, D. Stramski, and A. Poteau, “Spectral dependency of optical backscattering by marine particles from satellite remote sensing of the global ocean,” J. Geophys. Res. 111, C09024 (2006). [CrossRef]
32. M. J. Behrenfeld and E. Boss, “Beam attenuation and chlorophyll concentration as alternative optical indices of phytoplankton biomass,” J. Mar. Res. 64, 431–451 (2006). [CrossRef]
33. W. Clavano, E. Boss, and L. Karp-Boss, “Inherent optical properties of non-spherical marine-like particles—from theory to observation,” Oceanogr. Mar. Biol. 45, 1–38 (2007). [CrossRef]
34. J. M. Sullivan and M. S. Twardowski, “Angular shape of the oceanic particulate volume scattering function in the backward direction,” Appl. Opt. 48, 6811–6819 (2009). [CrossRef]