Airborne Lidar Observations of a Spring Phytoplankton Bloom in the Western Arctic Ocean

: One of the most notable effects of climate change is the decrease in sea ice in the Arctic Ocean. This is expected to affect the distribution of phytoplankton as the ice retreats earlier. We were interested in the vertical and horizontal distribution of phytoplankton in the Chukchi Sea in May. Measurements were made with an airborne proﬁling lidar that allowed us to cover large areas. The lidar proﬁles showed a uniform distribution of attenuation and scattering from the surface to the limit of lidar penetration at a depth of about 30 m. Both parameters were greater in open water than under the ice. Depolarization of the lidar decreased as attenuation and scattering increased. A cluster analysis of the 2019 data revealed four distinct clusters based on depolarization and lidar ratio. One cluster was associated with open water, one with pack ice, one with the waters along the land-fast ice, and one that appeared to be scattered throughout the region. The ﬁrst three were likely the result of different assemblages of phytoplankton, while the last may have been an artifact of thin fog in the atmosphere.


Introduction
One of the most obvious effects of climate change is the rapidly decreasing sea ice cover in the Arctic Ocean. Perennial sea ice cover decreased by 4.0% per decade between 1978 and 2010, but by 8.3% per decade over the latter years (1996)(1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010) of this time period [1]. Another study found that September sea ice cover decreased by 4.8% per decade between 1979 and 1996, but by 18.6% per decade between 1997 and 2014 [2]. The decrease in the winter was less, with a decrease of 2.4% per decade from 1979 to 1999 and a decrease of 3.4% per decade from 2000 to 2017 [3]. The rate of decrease in spring and autumn lie between the summer and winter values [4]. Arctic sea ice is also getting thinner, with an average reduction of 1.75 m between 1980 and 2009 [5]. This decline in thickness also means that ice is drifting faster [6].
This decline in sea ice has been recognized as a major driver of ecological changes in the Arctic [7]. Thinner sea ice lets more light through and increases the productivity of algae growing on the bottom of the ice [8]. Loss of sea ice is also increasing phytoplankton productivity [9,10]. It has been estimated that 30% of the increase is due to a larger area of ice-free water and the remaining 70% to a longer growing season. That increase in growing season has created a new phytoplankton bloom in the fall [11]. A model suggests that the increased absorption caused by increased phytoplankton is contributing to a warming of the surface waters in the Arctic Ocean [12].
One other biological effect is an increase in large phytoplankton blooms under the ice. The first detailed report was of a massive bloom in the Chukchi Sea in July of 2011 [13], apparently fed by upwelling at the shelf break [14]. The area of these summer blooms seems to be increasing at a rate of about 2% per year [15]. One important driver is increased Remote Sens. 2021, 13, 2512 3 of 13

Materials and Methods
The primary instrument was the NOAA oceanographic lidar [18,32,33]. The source was a linearly polarized, flashlamp-pumped, Q-switched, frequency-doubled Nd:YAG laser producing 100 mJ, 12 ns pulses of 532 nm light at a repetition rate of 30 Hz. The beam was expanded to 17 mrad full angle and steered to overlap the receiver telescopes at a distance of 300 m.
Two receiver telescopes (Kepler) were used to collect light at two orthogonal polarizations. Telescope diameters were 6 cm for the co-polarized receiver and 15 cm for the cross-polarized receiver. A sheet of polaroid material on the front of each telescope was used to select the polarization. An aperture at the focus was used to set the field of view of each receiver to be slightly larger than the transmitter divergence. A 3 nm wide interference filter was used to reduce the amount of background light reaching the photomultiplier-tube (pmt) photodetectors. The detector outputs were fed through a logarithmic amplifier to increase dynamic range and digitized with an eight-bit, giga-sample-per-second digitizer. A computer recorded 2000 samples of each pulse, global-positioning system (GPS) information, and the photomultiplier gains.
The lidar was mounted over the camera port in a NOAA Twin Otter aircraft pointed about 15 • forward from vertical to reduce specular reflections from the surface. The aircraft was generally flown at an altitude of 300 m and a speed of 60 m s −1 . Flights were made on seven days between 19 and 28 May 2019 over the Chukchi and Beaufort Seas of the Arctic Ocean ( Figure 1) out of Utqiaġvik, Alaska, USA. The objective of the flights was to make multiple passes over the ice edge, but fog and low clouds were constraints. The figure also shows the satellite-derived ice edge (ice fraction = 0.15) for 23 May 2019 and 23 July 2014 [18], each in the middle of this and of an earlier lidar survey. In the study area shown, the position of the ice edges is similar, and there were actually 9500 km 2 more of open water in this area in May 2019 than in July 2014.
The lidar data were processed in several steps. First, the data were visually inspected to remove those data files contaminated by fog and low clouds. Next, lidar pulses hitting ice at the surface were identified by a saturation of the return signal in the cross-polarized channel. The ice fraction was calculated for each~1 km section of flight track by taking the ratio of the number of these pulses to the total number of pulses in the segment, which was fixed at 500. For pulses without ice at the surface, the pmt gains and response functions of the log-amps and digitizer were used to convert digitizer count to photocathode current for each sample. The last 100 samples, which correspond to depths well below the laser penetration depth, of these pulses were averaged to provide an estimate of the background light signal, and this was subtracted from the rest of the samples in the pulse. The surface sample for each pulse was taken to be the one with the largest signal, and depths were measured from this sample. Data were converted from photocathode current to attenuated volume backscatter coefficient using calculated calibration coefficients at 300 m (422 m −1 sr −1 A −1 for the co-polarized channel and 79.2 m −1 sr −1 A −1 for the cross-polarized channel) and correcting for the difference in the range-squared loss between the actual altitude and 300 m using the GPS altitude. The co-and cross-polarized channels were added together to get the unpolarized return, and a linear regression of the logarithm of this signal against depth was performed over the depth range of 3-5 m. The volume scattering coefficient at the lidar scattering angle of 180 • , β(π), is the exponential of the intercept of this regression. We converted this to an estimate of backscattering coefficient, b b , by multiplying by 2π. This is equivalent to assuming uniform scattering in the backward direction [33]. The lidar attenuation coefficient, α, is −0.5 times the slope of the regression, since the signal decays as exp(−2αz). For the conditions of these surveys, α is very nearly equal to the diffuse attenuation coefficient, K d [34,35], and we will refer to K d from here on. Both K d and b b were averaged over 1 km flight segments; segments were eliminated if the variability was too great (standard deviation of b b values > 0.7 times the mean). made on seven days between 19 and 28 May, 2019 over the Chukchi and Beaufort Seas of the Arctic Ocean ( Figure 1) out of Utqiaġvik, Alaska, USA. The objective of the flights was to make multiple passes over the ice edge, but fog and low clouds were constraints. The figure also shows the satellite-derived ice edge (ice fraction = 0.15) for 23 May, 2019 and 23 July, 2014 [18], each in the middle of this and of an earlier lidar survey. In the study area shown, the position of the ice edges is similar, and there were actually 9500 km 2 more of open water in this area in May, 2019 than in July, 2014.  To convert b b into the particulate backscattering, b bp , the pure sea water value was subtracted. The pure water value of 1.08 × 10 −3 m −1 was used throughout, based on an approximation [36] to temperature and salinity dependent measurements [37]. A typical value of 32 PSU was used for salinity. A variation of ±1 PSU would produce an error of 0.9%. The measured sea-surface temperature was near 0 • for all flights, and this value was used for the calculation. A variation of ±1 • , which is larger than what we observed, would produce an error of 0.6%.
All data, including the measured depolarization ratio, δ = β × (π)/β (π), and the lidar ratio, LR = K d /β(π), were averaged into 0.3 • longitude by 0.1 • latitude bins for visualization. Here, β × (π) refers to the value that is cross-polarized with respect to the transmitted light and β (π) refers to the co-polarized value. All values, including δ and LR, were based on the 3-5 m regression. At 71 • latitude, the averaging bins are nearly square (10.9 km by 11.1 km). Reported statistical properties were calculated before binning.
To check for clusters of lidar ratio and depolarization, the latter was plotted as a function of the former. The plotted data were observed to be separated into four distinct clusters. These were defined by: cluster C1 with LR < 30 and δ < 0.08, cluster C2 with LR < 30 and δ > 0.08, cluster C3 with 30 < LR < 180, and cluster C4 with LR > 180. It was clear from the figure that the separation between clusters was not complete. For example, there was more overlap between clusters C2 and C3 than others. In an attempt to quantify this, we assume that the distribution of points within each cluster is Gaussian and calculate the probability that a point from that distribution actually lies within a different cluster.

Results
Lidar measurements of K d and b bp showed no vertical structure in the profiles down to the limits of lidar penetration at about 30 m. The subsurface plankton layers observed in July 2014 were completely absent in this data set.
The binned ice fractions from the lidar ( Figure 2) show a pattern similar to the satellite ice edge in Figure 1. There is a large region of mostly clear water west of Utqiaġvik, with ice farther west and to the north. No lidar measurements were made in the land-fast ice along the shore, because there were no holes in this ice through with the lidar could penetrate. Note that there are some differences between the satellite-derived ice edge, which was observed on one day in the middle of the survey period, and the lidar ice cover, which was spread out over several days. Lidar measurements of Kd and bbp showed no vertical structure in the profiles down to the limits of lidar penetration at about 30 m. The subsurface plankton layers observed in July, 2014 were completely absent in this data set.
The binned ice fractions from the lidar ( Figure 2) show a pattern similar to the satellite ice edge in Figure 1. There is a large region of mostly clear water west of Utqiaġvik, with ice farther west and to the north. No lidar measurements were made in the land-fast ice along the shore, because there were no holes in this ice through with the lidar could penetrate. Note that there are some differences between the satellite-derived ice edge, which was observed on one day in the middle of the survey period, and the lidar ice cover, which was spread out over several days.  The particulate backscattering coefficient ( Figure 4) has a pattern similar to Kd. In fact, the correlation between them was R = 0.66 (P < 10 −99 ). The correlation between bbp and ice cover (R = −0.06, P = 0.03) was less than the correlation between Kd and ice cover. Lidar measurements of Kd and bbp showed no vertical structure in the profiles down to the limits of lidar penetration at about 30 m. The subsurface plankton layers observed in July, 2014 were completely absent in this data set.
The binned ice fractions from the lidar ( Figure 2) show a pattern similar to the satellite ice edge in Figure 1. There is a large region of mostly clear water west of Utqiaġvik, with ice farther west and to the north. No lidar measurements were made in the land-fast ice along the shore, because there were no holes in this ice through with the lidar could penetrate. Note that there are some differences between the satellite-derived ice edge, which was observed on one day in the middle of the survey period, and the lidar ice cover, which was spread out over several days.  The particulate backscattering coefficient ( Figure 4) has a pattern similar to Kd. In fact, the correlation between them was R = 0.66 (P < 10 −99 ). The correlation between bbp and ice cover (R = −0.06, P = 0.03) was less than the correlation between Kd and ice cover. The particulate backscattering coefficient ( Figure 4) has a pattern similar to K d . In fact, the correlation between them was R = 0.66 (P < 10 −99 ). The correlation between b bp and ice cover (R = −0.06, P = 0.03) was less than the correlation between K d and ice cover. Remote Sens. 2021, 13, x FOR PEER REVIEW 6 of 13  The lidar ratio ( Figure 6) was negatively correlated with ice fraction (R = −0.08, P = 0.005). It was also negatively correlated with Kd (R = −0.42, P < 10 −53 ) and bbp (R = −0.52, P < 10 −86 ). However, it was positively correlated with the depolarization ratio (R = 0.18, P < 10 −9 ). Since the lidar ratio is proportional to attenuation, it seems strange that these two quantities should be negatively correlated. The explanation is that attenuation and scattering are highly correlated, and the range of scattering values is larger than the range of attenuation values. When the attenuation increases, the scattering tends to increase by an even larger amount, so the lidar ratio actually decreases.   The lidar ratio ( Figure 6) was negatively correlated with ice fraction (R = −0.08, P = 0.005). It was also negatively correlated with Kd (R = −0.42, P < 10 −53 ) and bbp (R = −0.52, P < 10 −86 ). However, it was positively correlated with the depolarization ratio (R = 0.18, P < 10 −9 ). Since the lidar ratio is proportional to attenuation, it seems strange that these two quantities should be negatively correlated. The explanation is that attenuation and scattering are highly correlated, and the range of scattering values is larger than the range of attenuation values. When the attenuation increases, the scattering tends to increase by an even larger amount, so the lidar ratio actually decreases. The lidar ratio ( Figure 6) was negatively correlated with ice fraction (R = −0.08, P = 0.005). It was also negatively correlated with K d (R = −0.42, P < 10 −53 ) and b bp (R = −0.52, P < 10 −86 ). However, it was positively correlated with the depolarization ratio (R = 0.18, P < 10 −9 ). Since the lidar ratio is proportional to attenuation, it seems strange that these two quantities should be negatively correlated. The explanation is that attenuation and scattering are highly correlated, and the range of scattering values is larger than the range of attenuation values. When the attenuation increases, the scattering tends to increase by an even larger amount, so the lidar ratio actually decreases.
The scatter plot of depolarization vs. lidar ratio (Figure 7) shows the four clusters that were identified visually. These were defined by: cluster C1 with LR < 30 and δ < 0.08, cluster C2 with LR < 30 and δ > 0.08, cluster C3 with 30 < LR < 180, and cluster C4 with LR > 180. Each cluster is plotted in a different color. The heavy black lines denote the one standard deviation contour. The mean and standard deviation of each cluster lidar ratio and depolarization were: cluster C1-(13.7 ± 5.2, 0.046 ± 0.012), cluster C2-(15.3 ± 6.5, 0.126 ± 0.036), cluster C3-(61.6 ± 24.6, 0.109 ±0.035), and cluster C4-(96.53 ± 53, 0.110 ± 0.029).  The scatter plot of depolarization vs. lidar ratio (Figure 7) shows the four clusters that were identified visually. These were defined by: cluster C1 with LR < 30 and δ < 0.08, cluster C2 with LR < 30 and δ > 0.08, cluster C3 with 30 < LR < 180, and cluster C4 with LR > 180. Each cluster is plotted in a different color. The heavy black lines denote the one standard deviation contour. The mean and standard deviation of each cluster lidar ratio and depolarization were: cluster C1 -(13.7 ± 5.2, 0.046 ± 0.012), cluster C2 -(15.3 ± 6.5, 0.126 ± 0.036), cluster C3 -(61.6 ± 24.6, 0.109 ±0.035), and cluster C4 -(96.53 ± 53, 0.110 ± 0.029). It is clear from the figure that the separation between clusters was not complete. For example, there was more overlap between clusters C2 and C3 than others. In an attempt to quantify this, we assumed that the distribution of points within each cluster is Gaussian and calculated the probability that a point from that distribution actually lies within a different cluster. Based on this analysis, C1 was the tightest cluster, with a probability of extending into C2 and C3 of 3.1 × 10 −5 and 4.6 × 10 −6 , respectively. C2 overlapped with C1 and C3 with probabilities of 0.035 and 6.9 × 10 −4 . C3 overlapped with C2 and C4 with probabilities of 0.033 and 2.8 × 10 −12 . Finally, C4 overlapped with C3 with a probability of  The scatter plot of depolarization vs. lidar ratio (Figure 7) shows the four clusters that were identified visually. These were defined by: cluster C1 with LR < 30 and δ < 0.08, cluster C2 with LR < 30 and δ > 0.08, cluster C3 with 30 < LR < 180, and cluster C4 with LR > 180. Each cluster is plotted in a different color. The heavy black lines denote the one standard deviation contour. The mean and standard deviation of each cluster lidar ratio and depolarization were: cluster C1 -(13.7 ± 5.2, 0.046 ± 0.012), cluster C2 -(15.3 ± 6.5, 0.126 ± 0.036), cluster C3 -(61.6 ± 24.6, 0.109 ±0.035), and cluster C4 -(96.53 ± 53, 0.110 ± 0.029). It is clear from the figure that the separation between clusters was not complete. For example, there was more overlap between clusters C2 and C3 than others. In an attempt to quantify this, we assumed that the distribution of points within each cluster is Gaussian and calculated the probability that a point from that distribution actually lies within a different cluster. Based on this analysis, C1 was the tightest cluster, with a probability of extending into C2 and C3 of 3.1 × 10 −5 and 4.6 × 10 −6 , respectively. C2 overlapped with C1 and C3 with probabilities of 0.035 and 6.9 × 10 −4 . C3 overlapped with C2 and C4 with probabilities of 0.033 and 2.8 × 10 −12 . Finally, C4 overlapped with C3 with a probability of It is clear from the figure that the separation between clusters was not complete. For example, there was more overlap between clusters C2 and C3 than others. In an attempt to quantify this, we assumed that the distribution of points within each cluster is Gaussian and calculated the probability that a point from that distribution actually lies within a different cluster. Based on this analysis, C1 was the tightest cluster, with a probability of extending into C2 and C3 of 3.1 × 10 −5 and 4.6 × 10 −6 , respectively. C2 overlapped with C1 and C3 with probabilities of 0.035 and 6.9 × 10 −4 . C3 overlapped with C2 and C4 with probabilities of 0.033 and 2.8 × 10 −12 . Finally, C4 overlapped with C3 with a probability of 2.8 × 10 −12 . The number points in each cluster was 331, 202, 478, and 220, so C1 and C3 were the most populated clusters.
The spatial distribution of clusters ( Figure 8) shows that cluster C1 was dominant in the broken ice to the west and north, cluster C2 was dominant along the edge of the land-fast ice, cluster C3 was dominant in open water, and cluster C4 was scattered geographically.
The spatial distribution of clusters ( Figure 8) shows that cluster C1 was dominant in the broken ice to the west and north, cluster C2 was dominant along the edge of the landfast ice, cluster C3 was dominant in open water, and cluster C4 was scattered geographically.

Discussion
One of the most surprising features of this investigation was that the ice extent in the survey area was actually less in May, 2019 than it was in July, 2014 (see Figure 1). This is consistent with an overall decrease in Arctic sea ice between 2014 and 2019. For example, the minimum ice extents for those two years were 4.58 × 10 6 km 2 and 3.91 × 10 6 km 2 , which is a decrease of 14.6% in two years [38]. This decrease is much greater than the long-term trends discussed in the Introduction, and it reflects inter-annual variability more than a long-term trend.
With similar ice extent in the study area, one might expect that the stratification of the water column might be similar in July, 2014 and May, 2019. In fact, in July 2014, the water column was stratified and the depths of the subsurface plankton layers observed in the lidar corresponded to the depth of stratification [18]. The situation in May, 2019 was different. An Alamo float in the Chukchi Sea (float 9234) reported no significant salinity stratification until about 15 June, 2019 [39]. Recall that, at temperatures near freezing, density stratification closely follows salinity stratification.
Kd and bbp are important optical parameters, but they are also important because of their relationships to phytoplankton. For example, a common model relating Kd at 532 nm to the chlorophyll concentration C in mg m −3 is [30,40] 0.67 0.0452 0.0474 .
The distribution of chlorophyll concentration (Figure 9) is similar to the distribution of Kd in Figure 3. Due to the nonlinear nature of the relationship between Kd and C, the correlation coefficient between C and ice fraction (R = −0.16, P < 10 −8 ) is slightly less than that between Kd and ice fraction.

Discussion
One of the most surprising features of this investigation was that the ice extent in the survey area was actually less in May 2019 than it was in July 2014 (see Figure 1). This is consistent with an overall decrease in Arctic sea ice between 2014 and 2019. For example, the minimum ice extents for those two years were 4.58 × 10 6 km 2 and 3.91 × 10 6 km 2 , which is a decrease of 14.6% in two years [38]. This decrease is much greater than the long-term trends discussed in the Introduction, and it reflects inter-annual variability more than a long-term trend.
With similar ice extent in the study area, one might expect that the stratification of the water column might be similar in July 2014 and May 2019. In fact, in July 2014, the water column was stratified and the depths of the subsurface plankton layers observed in the lidar corresponded to the depth of stratification [18]. The situation in May 2019 was different. An Alamo float in the Chukchi Sea (float 9234) reported no significant salinity stratification until about 15 June 2019 [39]. Recall that, at temperatures near freezing, density stratification closely follows salinity stratification.
K d and b bp are important optical parameters, but they are also important because of their relationships to phytoplankton. For example, a common model relating K d at 532 nm to the chlorophyll concentration C in mg m −3 is [30,40] K d = 0.0452 + 0.0474C 0.67 . (1) The distribution of chlorophyll concentration (Figure 9) is similar to the distribution of K d in Figure 3. Due to the nonlinear nature of the relationship between K d and C, the correlation coefficient between C and ice fraction (R = −0.16, P < 10 −8 ) is slightly less than that between K d and ice fraction.
Backscattering is more commonly related to Particulate Organic Carbon (POC) [41], and a relationship has been developed specific to the Alaskan Arctic [42], This relationship was developed for a wavelength of 550 nm, but the difference between that and 532 nm is small. Phytoplankton, of course, are the drivers of primary productivity, and models of productivity have been developed [43][44][45][46][47][48]. Calculation of these parameters is beyond the scope of this paper.
Our data show relatively weak, negative correlations between K d and ice fraction and between b bp and ice fraction. This suggests the strong spring phytoplankton bloom observed in open water is accompanied by an under-ice bloom that is only slightly weaker. This is consistent with a recent modeling effort that estimated that 63% of the total primary production in Arctic waters occurs in water with 50% or more ice fraction [49]. That work also estimated that 41% occurs in water with 85% or more ice fraction. Backscattering is more commonly related to Particulate Organic Carbon (POC) [41], and a relationship has been developed specific to the Alaskan Arctic [42], This relationship was developed for a wavelength of 550 nm, but the difference between that and 532 nm is small. Phytoplankton, of course, are the drivers of primary productivity, and models of productivity have been developed [43][44][45][46][47][48]. Calculation of these parameters is beyond the scope of this paper.
Our data show relatively weak, negative correlations between Kd and ice fraction and between bbp and ice fraction. This suggests the strong spring phytoplankton bloom observed in open water is accompanied by an under-ice bloom that is only slightly weaker. This is consistent with a recent modeling effort that estimated that 63% of the total primary production in Arctic waters occurs in water with 50% or more ice fraction [49]. That work also estimated that 41% occurs in water with 85% or more ice fraction.
Our data also show that depolarization is positively correlated with ice fraction. This is somewhat surprising, given that the negative correlation of scattering and ice fraction. It suggests that the composition of scattering particles is different in open water than under the ice.
Since depolarization and lidar ratio depend more on the characteristics of the particulate scattering and less on the number density of particles, we expect that the different clusters represent different types of particles. Neeley et al. [50] showed that there were very distinct phytoplankton groups associated with sea ice and open water. This would seem to explain the difference in scattering properties between clusters C1 and C3 that are strongly associated with sea ice and open water, respectively. Cluster C2 is strongly associated with the edge of the land-fast ice along the northwest coast of Alaska. It is characterized by low lidar ratio and high depolarization, which suggests that this region has yet another group of phytoplankton species. The characteristics of this cluster could also be explained by a higher concentration of mineral particles than other areas.
The final cluster, C4, is spread throughout the area. It is characterized by very high lidar ratios, which might reflect another group of phytoplankton species, but might be the result of atmospheric attenuation. For example, a thin fog would reduce the energy reaching the water, thus reducing the amount of backscattered light. It would not affect the attenuation, so it would result in high attenuation-to backscatter ratios. Table 1 presents the correlations using the full data set and those without including cluster C4. The differences are small, and do not affect the overall conclusions. The relationships between ice fraction and the optical parameters Kd and bbp are slightly strengthened by removing the data in cluster C4. Our data also show that depolarization is positively correlated with ice fraction. This is somewhat surprising, given that the negative correlation of scattering and ice fraction. It suggests that the composition of scattering particles is different in open water than under the ice.
Since depolarization and lidar ratio depend more on the characteristics of the particulate scattering and less on the number density of particles, we expect that the different clusters represent different types of particles. Neeley et al. [50] showed that there were very distinct phytoplankton groups associated with sea ice and open water. This would seem to explain the difference in scattering properties between clusters C1 and C3 that are strongly associated with sea ice and open water, respectively. Cluster C2 is strongly associated with the edge of the land-fast ice along the northwest coast of Alaska. It is characterized by low lidar ratio and high depolarization, which suggests that this region has yet another group of phytoplankton species. The characteristics of this cluster could also be explained by a higher concentration of mineral particles than other areas.
The final cluster, C4, is spread throughout the area. It is characterized by very high lidar ratios, which might reflect another group of phytoplankton species, but might be the result of atmospheric attenuation. For example, a thin fog would reduce the energy reaching the water, thus reducing the amount of backscattered light. It would not affect the attenuation, so it would result in high attenuation-to backscatter ratios. Table 1 presents the correlations using the full data set and those without including cluster C4. The differences are small, and do not affect the overall conclusions. The relationships between ice fraction and the optical parameters K d and b bp are slightly strengthened by removing the data in cluster C4.  18 1 There has been previous work using depolarization to distinguish between different types of scattering particles. A ship-based lidar was used to show a strong relationship between depolarization and the particulate backscattering ratio, b bp /b p , where b p is the total particulate scattering coefficient [51]. The backscattering ratio has been independently used to infer properties of scattering particles in the ocean [52,53]. Depolarization from that same lidar has been used to identify a coccolithophore bloom [54]. Another attempt used the ratio δ/b bp instead of δ to discriminate different types of scattering particles [31]. The argument for this parameter was that multiple scattering contributes to the depolarization, and dividing by b bp compensates for this effect, while retaining the differences in depolarization due to differences in particle characteristics. While this worked for the case presented, it is probably not generally valid, as the multiple-scattering component depends on the depth of the measurement and on the single-scattering albedo [29].
We did not expect the negative correlation between depolarization and b bp that was observed. If changes in b bp are a result of changes in the number density of the same particle type, we would expect this correlation to be zero. The correlation was calculated separately for each of the four clusters. In three of these, the correlation was not statistically significant. Only in cluster C1, along the ice edge, was the correlation significant (R = −0.62, P = 0.0012). This suggests that the composition of scatterers along the ice edge is changing, even within the cluster C1. On the other hand, the correlation between depolarization and K d was significant for all of the clusters except for C2, along the edge of the landfast ice. Where the correlation with K d is significant, but the correlation with b bp is not (e.g., C3 in open water), it is tempting to infer that the difference is due to variations in absorption by dissolved organic compounds (DOC). As the phytoplankton assemblage becomes more depolarizing, the concentration of DOC decreases. This seems to be more important at larger values of K d . For K d < 0.25 m −1 , for example, the correlation between depolarization and K d is actually positive (R = 0.35, P < 10 −4 ).
The results of this study suggest that airborne lidar can be very useful for measurements of the spring bloom in the Arctic Ocean. It can measure ice fraction, diffuse attenuation coefficient and particulate backscatter, and the relationships between these parameters. With the lidar-specific parameters of depolarization and lidar ratio, it can identify regions of distinct assemblages of scattering particles. It can cover large areas of the ocean quickly, penetrate into the water in broken ice, and is unaffected by weather except for low clouds and fog. Previous work has validated lidar measurements of K d [55] and b bp , [33] so future work in the Arctic can concentrate on converting these optical parameters to C and POC. These, in turn, can provide important information on Arctic Ocean ecosystems.

Conclusions
The ice edge in the Chukchi Sea in May 2019 was about the same place it was in July 2014. The vertical distribution of phytoplankton was very different, however, with strong subsurface layers in July 2014, but a uniform bloom right to the surface in May 2019. In both cases, scattering and attenuation of the lidar were greater in open water than under the ice. A cluster analysis of the 2019 data revealed four distinct clusters based on depolarization and lidar ratio. One cluster was associated with open water, one with pack ice, one with the waters along the land-fast ice, and one that appeared to be scattered throughout the region. The first three were likely the result of different assemblages of phytoplankton, but we cannot determine from the lidar return what the makeup of these assemblages might have been. The last cluster may have been an artifact of thin fog in the atmosphere, but we showed that this would not affect the main conclusions of the study.

Data Availability Statement:
The lidar data used in this paper can be found on the NOAA Chemical Sciences Laboratory at https://csl.noaa.gov/groups/csl3/measurements/2019ArcticStratification/, accessed on 25 May 2021.