In this study, we used satellite data from MODIS Aqua (Moderate Resolution Imaging Spectroradiometer on the Aqua satellite), chosen for its consistent long-term data record. We acquired daily Level 2 MODIS Aqua remote sensing reflectance [R_rs(λ)] data, with the 1-km resolution at nadir, from the National Aeronautics and Space Administration (NASA) Ocean Biology Distributed Active Archive Center (OB.DAAC; https://oceancolor.gsfc.nasa.gov/). We mapped each daily Level 2 scene employing the l2bin and l3mapgen tools available in the SeaDAS-OCSSW software version 8.3 (the NASA Ocean Biology Processing Group Science Software).

The GoM is an optically complex environment, and consequently, the atmospheric correction process can over-correct blue wavelengths, resulting in non-physical negative R_rs(λ) values or unrealistic spectral shapes. Thus, to ensure the quality of the R_rs(λ) dataset, we applied two quality control (QC) approaches. First, we removed all R_rs(λ) spectra that were negative at any wavelength. Second, we applied the Quality Water Index Polynomial (QWIP) method of Dierssen et al. (2022). QWIP is a quantitative metric approach and was developed to assess the quality of R_rs(λ) in the visible wavelength spectral range of 400–700 nm. The QWIP score is used to identify spectra that fall outside the general trends observed in aquatic optics for optically deep waters. For this study, R_rs(λ) spectra that fall outside the QWIP score of |0.2| were flagged to eliminate poorer quality spectra, as suggested by Dierssen et al. (2022).

The quality-controlled satellite R_rs(λ) were converted to enhanced RGB (eRGB) images using a standardized linear contrast image stretching method and gamma correction (McCarry et al., 2023). Briefly, the eRGB image processing procedure standardizes the RGB bands by applying a standardized stretch technique and a gamma correction specifically for the blue band. The standardized stretch range had a minimum of zero and a maximum taken as the 90th percentile of each band in a global in situ dataset of R_rs(λ) (Valente et al., 2022). This approach was implemented to resolve visual inconsistencies between individual satellite images. When creating an RGB image, typically each color band is scaled based on the maximum and minimum values present in the given image, resulting in different scalings being used for each image (McCarry et al., 2023). By standardizing the range over which each band is scaled, different images are directly comparable. This is crucial for accurate analysis and comparison of satellite data, especially when dealing with a large number of satellite images. Following McCarry et al. (2023), we generated eRGB images using three R_rs(λ) wavebands from the quality-controlled satellite R_rs(λ) data: 443 nm for the blue channel, 488 nm for the green channel, and 555 nm for the red channel. This approach provides standardized eRGB data that are consistent across the whole dataset.

To verify the presence of surface patches in eRGB satellite imagery as C. finmarchicus, we utilized in situ abundance data of C. finmarchicus stage V–VI from the Continuous Plankton Recorder (CPR) collected in the GoM, sourced from a CPR database maintained by the National Oceanic and Atmospheric Administration (NOAA). The CPR is towed horizontally behind ships of opportunity, enabling continuous sampling of surface zooplankton (Richardson et al., 2006), making it an ideal dataset to compare with surface satellite observations. The CPR program has generated a long-term (1961–2013) and standardized database of marine zooplankton counts, which includes C. finmarchicus (Parent et al., 2012). The GoM CPR transect spans 450 km across the GoM, extending from Massachusetts to Nova Scotia in a west–east direction.

The initial identification of C. finmarchicus patches in eRGB satellite imagery involved comparing in situ CPR C. finmarchicus stage V–VI and eRGB satellite data from 2003 to 2013. We selected this timeframe due to the overlap of data availability from both CPR and MODIS Aqua satellite data. We focused initially on identifying high-abundance C. finmarchicus patches by selecting the highest 1% of the CPR values. From a dataset of 1,143 values, this resulted in 85 data points above a C. finmarchicus abundance of 10,000 individuals m⁻³ (ind m⁻³). We chose this threshold to enhance computational efficiency while ensuring the initial identification of C. finmarchicus presence. Based on the time and location of the top 1% C. finmarchicus abundance observations, we identified potential patches in the satellite eRGB imagery. We checked imagery from the same day of the CPR observation in the region surrounding the location of the observation.

In this study, we aim to refine the method presented by McCarry et al. (2023) to estimate the surface abundance of C. finmarchicus. This approach relies on satellite ocean color observations; thus, due to the diel vertical migration undergone by C. finmarchicus, the method presented provides a snapshot of the surface community in the middle of the day, rather than an estimate of the full community. Here, we give a brief overview of the McCarry et al. (2023) approach (which we denote as MC23 for brevity) and highlight the areas for refinement.

The core of the MC23 remote sensing method is based on two look-up tables (LUTs). These LUTs were created using a Case 2 (Morel and Prieur, 1977) bio-optical model based on Lo Prejato et al. (2020) and using specific inherent optical properties from Bengil et al. (2016) and Basedow et al. (2019). In the first LUT, referred to as the “Case 2” LUT, the bio-optical model contained the OSCs of phytoplankton chlorophyll-a and algal biogenic detritus (denoted as CHL), CDOM, and inorganic mineral suspended solids (MSS). The second LUT is the same as the first, except that it includes C. finmarchicus and is referred to as the “Case 2 + C. finmarchicus” LUT. For each LUT, typical ranges of concentrations for each of the OSCs were used, and R_rs(λ) was calculated using the radiative transfer software EcoLight v5.3 (Mobley and Sundman, 2016). Finally, R_rs(λ) spectra were converted to eRGB values in the same manner as described above for the satellite data. Thus, the final LUTs contained OSC concentrations with corresponding eRGB values.

In MC23, for one eRGB satellite image, every pixel was compared with the eRGB values in the LUTs (Case 2 and Case 2 + C. finmarchicus). The comparison between the satellite and LUTs eRGB was quantified using the Delta E 2000 ( Δ E 2000 ) algorithm, which is based on the CIE2000 model developed by the International Commission on Illumination (CIE). Δ E 2000   is a standard way of quantifying the difference between two colors in the L * a * b (where L = lightness, a = red–green axis, and b = yellow–blue axis) color space. Δ E 2000 values range from 0 to 100, where low Δ E 2000 value indicates that the two colors are similar, while a higher value indicates a high difference. Therefore, in this study and MC23, a “best match” between a satellite pixel and the LUT was taken as the lowest Δ E 2000 value ( ΔE 2000 m i n ). For a given satellite image, a map of ΔE 2000 m i n acts as an anomaly map: low ΔE 2000 m i n values indicate the satellite pixels are well described by the bio-optical model (i.e., there is a low anomaly between the bio-optical model and the satellite observation), and conversely, higher ΔE 2000 m i n   values indicate a less well described (more anomalous) system. MC23 showed that using the Case 2 + C. finmarchicus LUT reduced the anomaly (i.e., the ΔE 2000 m i n values were lower) compared to the Case 2 LUT in a dense patch of C. finmarchicus. Further, the C. finmarchicus concentration was then estimated from the LUT value that corresponded to the ΔE 2000 m i n for every pixel.

In this study, we used an updated version of the Case 2 + C. finmarchicus LUT from MC23, in which MSS values are varied, instead of the constant MSS value of 0.03 g m⁻³. The updated version of the Case 2 + C. finmarchicus LUT has 4,410 entries, an increase from 1,176 entries (in MC23). In this paper, all future mentions of the Case 2 + C. finmarchicus LUT refer to this updated LUT (rather than the original in MC23). With a larger LUT, the process of comparing every pixel with every eRGB value in the LUT becomes computationally more expensive. In addition, there is the potential for even further increases in computational time and resources as we start to explore larger numbers of satellite images or create a finer-resolution LUT (i.e., more concentrations for the OSC). Thus, we explored the effectiveness of using a threshold-based approach. Here, we compared a given satellite image to the Case 2 LUT (which has a total of 1,008 entries, as in MC23). The ΔE 2000 m i n   anomaly map then highlights the regions where the Case 2 model does not sufficiently describe the system. However, where the anomaly is low, the Case 2 model is considered sufficient. Therefore, only pixels with a Case 2 ΔE 2000 m i n   value exceeding the implemented threshold had the updated Case 2 + C. finmarchicus LUT applied. Thus, with this approach, we only applied the Case 2 + C. finmarchicus LUT to a subset of the pixels, not the whole image, and for those pixels, we were able to estimate the C. finmarchicus concentration based on the value from the LUT that corresponds to the ΔE 2000 m i n , as in MC23.

We detected potential patches of C. finmarchicus in the eRGB imagery based on the highest 1% of the C. finmarchicus abundances in the CPR data. Figure 1 illustrates an example: the detection of a potential C. finmarchicus patch in eRGB imagery taken on June 17, 2009, identified through comparison with in situ CPR C. finmarchicus data. The eRGB image shows a region of red/brown pixels, potentially indicative of a C. finmarchicus patch (McCarry et al., 2023). The recorded CPR C. finmarchicus abundance at 19:33 local time (EST) was 21,402 ind m⁻³ at the location depicted by a white star in Figure 1 . Thus, the high CPR observation coincides with the potential patch of C. finmarchicus visible in the satellite imagery. The match between the eRGB satellite imagery and the CPR C. finmarchicus data provides evidence supporting the results of McCarry et al. (2023): dense patches of C. finmarchicus can be detected in eRGB imagery. Further, this demonstrates the applicability of the MC23 approach in a new region.

To investigate the Δ E 2000 m i n threshold limit to use for the threshold-based approach, we binned the Case 2 Δ E 2000 m i n values (see Figure 2 ). Different threshold values are used for defining a color difference (i.e., the Δ E 2000 ) depending on the application. In the printing industry, colors are considered visually consistent within a Δ E 2000   threshold of 4 (Liu et al., 2012). Therefore, to test this threshold for this particular application, we binned Δ E 2000 m i n values of 0–4 within the first bin, and then we binned higher Δ E 2000 m i n values in steps of 2.

In Figure 2 , regions in the GoM and off the coast of Norway identified as potential C. finmarchicus patches show Δ E 2000 m i n values larger than 4 (see Figure 1 and figure 10 of McCarry et al., 2023). These higher Δ E 2000 m i n values indicate a greater deviation between the “best-matched” satellite and Case 2 LUT eRGB values. This deviation suggests that the pixels remain anomalous and are not optically resolved.

For the GoM image, we applied the Case 2 + C. finmarchicus LUT to the pixels where the Case 2 Δ E 2000 m i n   ≥4.

Figure 3 shows a refined version of the ΔE 2000 m i n anomaly map, which combines ΔE 2000 m i n values from the Case 2 LUT and the Case 2 + C. finmarchicus LUT. This refined binned ΔE 2000 m i n   anomaly map shows a lower ΔE 2000 m i n across the entire image. The reduced ΔE 2000 m i n values suggest that incorporating the C. finmarchicus component into the LUT resolved the anomalies in the regions where the Case 2 Δ E 2000 m i n values were ≥4. Thus, the binned Δ E 2000 m i n   anomaly maps indicate that a Case 2 Δ E 2000 m i n threshold of ≥4 typically identifies pixels that are anomalous due to C. finmarchicus.

Surface C. finmarchicus concentrations were quantified using the threshold-based approach, where each satellite pixel was first compared to the Case 2 LUT, and for any pixels where ΔE 2000 m i n was above the threshold (≥4), the Case 2 + C. finmarchicus LUT was used. The surface concentration of C. finmarchicus was then estimated only for those pixels where the Case 2 + C. finmarchicus LUT was used, and the C. finmarchicus concentration was taken from the combination of OSCs that produced the identified Δ E 2000 m i n   value. Using this approach on the MODIS Aqua eRGB imagery shown in Figure 1 , Gulf-wide surface C. finmarchicus concentrations ( Figure 4 ) were estimated, with peak concentrations of 150,000 ind m⁻³ within the patch.

The GoM is known to have a wide variety of zooplankton species, many of which have the pigment astaxanthin, for example, Centropages and krill (Bandaranayake and Gentien, 1982; Koomyart et al., 2017). As the MC23 method is based on the detection of astaxanthin, other astaxanthin or red-colored species may mistakenly be detected as C. finmarchicus. To explore this potential misidentification, we examined satellite images at times of the year when C. finmarchicus is not typically present in the surface waters of the GoM but other astaxanthin or red-colored species are. Here, we present two case studies, highlighting the effect on C. finmarchicus abundance estimates from other astaxanthin or red-colored species.

In the first case, we observed an apparent high abundance of surface C. finmarchicus in the late fall, as shown in Figure 5 . Typically, C. finmarchicus is expected to be in diapause in the late fall; thus, it is less likely they will be present in surface waters (Häfker et al., 2018). In this case, the Case 2 ΔE 2000 m i n anomaly map (see Figure 5B ) shows high ΔE 2000 m i n values across much of the region. The refined ΔE 2000 m i n map shows low ΔE 2000 m i n values in these anomalous regions, indicating that these pixels appear to be well resolved by the Case 2 + C. finmarchicus LUT ( Figure 5C ). The estimated surface C. finmarchicus concentration (see Figure 5D ) corresponding to the refined   ΔE 2000 m i n anomaly map is high in these regions that were anomalous when using the Case 2 LUT only, but not when using the Case 2 + C. finmarchicus LUT. In this case, local knowledge suggests that it is less likely to be C. finmarchicus present at the surface during late fall and more likely to be Centropages, a species known to be in the surface waters of the GoM in November (Stegert et al., 2012). If so, this set of results could indicate that the ΔE 2000 m i n anomaly method identifies astaxanthin pigment, irrespective of species, and the high surface concentrations of C. finmarchicus estimated during the late fall may actually be attributed to other astaxanthin containing zooplankton species.

The second case study illustrates how dark red/brown phytoplankton blooms may be misidentified as C. finmarchicus. There was a large bloom of Tripos muelleri, a brown dinoflagellate, in the GoM in the spring and summer of 2023 (Ray, 2023). Thus, to explore the impact of a different red/brown organism on this ocean color remote sensing approach, a single scene from spring 2023 was acquired and processed in the same way as previously described ( Figure 6 ). In this case, the anomalous pixels in the Case 2 ΔE 2000 m i n anomaly map are high and extensive across the GoM ( Figure 6B ). While the ΔE 2000 m i n anomaly is reduced using the refined, threshold-based approach (see Figure 6C ), the anomalies are still high. The presence of anomalous pixels in the refined ΔE 2000 m i n anomaly map indicates that the pixels are not well resolved by the Case 2 + C. finmarchicus LUT; i.e., the addition of C. finmarchicus did not fully remove the anomaly. The estimated surface C. finmarchicus concentration (see Figure 6D ) is high across the GoM. When C. finmarchicus swarm at the surface, they tend to be in small patches, and the Gulf-wide feature in Figure 6D is unusual. Therefore, the unresolved refined ΔE 2000 m i n anomaly ( Figure 6C ), the spatial extent of the C. finmarchicus patch, and the known presence of T. muelleri suggest a misidentification C. finmarchicus.

We applied the eRGB image processing technique from MC23 to identify potential patches of surface C. finmarchicus in the GoM. While the eRGB standardization technique developed by McCarry et al. (2023) was initially intended for the Norwegian Sea, the range was derived from a comprehensive global dataset (Valente et al., 2022). Thus, the direct application of the MC23 eRGB approach is equally effective in standardizing satellite RGBs for the GoM. The implementation of the eRGB image standardizing technique is particularly crucial when working with large time-scale datasets, as it ensures consistency across the eRGB images, thereby allowing direct comparison between them.

The anomaly detection process is a test of optical closure. The refined anomaly detection used a threshold-based approach, where a Case 2 ΔE 2000 m i n ≥ 4 was identified as anomalous (see Figure 2 ). In the case of the June 2009 image ( Figures 2A , 3 ), we observed high ΔE 2000 m i n values (≥4) in the Case 2 ΔE 2000 m i n anomaly map ( Figure 2A ), particularly in regions of potential C. finmarchicus. This suggests that the Case 2 bio-optical model in these areas is incomplete, and the potential C. finmarchicus patch is anomalous without the inclusion of the C. finmarchicus constituent. In the refined, threshold-based approach, the Case 2 + C. finmarchicus LUT was applied on the anomalous pixels, resulting in significantly lower ΔE 2000 m i n values, indicating that the anomaly in the optical model was resolved with the addition of C. finmarchicus ( Figure 3 ). Therefore, in this case, the addition of C. finmarchicus component into the Case 2 standard bio-optical model in the anomalous region resolves the anomalies between the bio-optical model and the satellite imagery, as evidenced by the lower ΔE 2000 m i n values in the refined ΔE 2000 m i n anomaly.

Both case studies, illustrated in Figures 5 , 6 , highlight factors that impact the accuracy of the ocean color anomaly detection method for estimating C. finmarchicus concentrations. The method presented here, and in McCarry et al. (2023), relied on the detection of astaxanthin to quantify C. finmarchicus. Astaxanthin is a pigment common in various zooplankton species, and there are red/brown phytoplankton species, with many of these plankton and phytoplankton species present in the GoM. We found an unusually high C. finmarchicus concentration in late fall, which is likely to be another astaxanthin-rich species such as Centropages. Furthermore, a bloom of T. muelleri was misidentified as a large swarm of C. finmarchicus. These instances suggest that the method can be used to detect the presence of red/brown/astaxanthin-rich species, but without a priori knowledge of the system, a specific species cannot be identified. In the case of the T. muelleri bloom, the anomalous pixels in the refined anomaly map (high ΔE 2000 m i n values in Figure 6C ) indicate that the Case 2 + C. finmarchicus LUT did not fully resolve the anomaly. This, along with the estimation of high C. finmarchicus concentrations (concentrations at the maximum value included in the LUT, see Figure 6D ), suggests that even higher concentrations of C. finmarchicus are needed in the LUT to resolve the anomaly. Alternatively, the failure to eliminate the anomaly using the existing C. finmarchicus range in the LUT could be used to identify a further mismatch, and therefore, the resulting C. finmarchicus estimates could be eliminated as being unrealistic. Here, we estimated C. finmarchicus concentrations (in individuals m⁻³) based on the absorption coefficient for individual C. finmarchicus. However, other red/brown/astaxanthin-rich plankton package their pigments differently, and therefore, each species will absorb different amounts per individual. Thus, the C. finmarchicus concentrations shown in Figures 5D , 6D cannot simply be taken as concentrations of Centropages or T. muelleri, and appropriate species and pigment-specific absorption measurements are required. This package effect is likely to happen within a given species, related to the size of the individuals, stage of life cycle, or physiological effects; thus, the variability in species-specific absorption also needs to be quantified. Consequently, we conclude that the refined, threshold-based anomaly detection method has the potential for detecting surface astaxanthin-rich or red-colored species. These could then be converted to species concentration with local knowledge of the species present and an appropriate, species-specific absorption measurement. To obtain species-specific absorption, it is important to measure the in situ absorption of astaxanthin-rich zooplankton under varying conditions, including different locations, life stages, and concentrations.

Optical closure may fail because either there is a problem in our bio-optical model or because there is an issue with the satellite data. Given the optically complex nature of the GoM, there can be issues with the atmospheric correction, particularly near the coast. While we have implemented a QC approach to remove spectra that have atmospheric correction issues, it is likely this approach is not perfect. The bio-optical model we used was built from Ligurian Sea data, and its application in the GoM may have caused misclassifications, although we are encouraged by the degree of optical closure that has been achieved. Implementing a GoM bio-optical model (including additional OSCs) and improving satellite atmospheric correction may further improve optical closure and the accuracy in identifying and quantifying C. finmarchicus for this region.

The in situ CPR C. finmarchicus population in the GoM recorded at approximately 19:33 local time (EST) on June 17, 2009, was 21,402 ind m⁻³ at the location depicted by a white star in Figures 1 , 4 . The refined anomaly detection method estimated the surface concentration of C. finmarchicus at the same location to be 80,000 ind m⁻³, obtained from satellite imagery captured earlier that day at 13:12 local time (EST). Although these in situ and satellite measurements were taken at different times on the same day, other factors contributing to this discrepancy should also be considered. However, the discrepancy between in situ and satellite measurements is consistent with the findings of McCarry et al. (2023). A potential factor for the discrepancy between satellite and in situ data could be partially attributed to the self-propelling movement of C. finmarchicus. Some studies have shown the ability of C. finmarchicus to avoid capture in net sampling, highlighting factors that contribute to the underestimation of population densities, including copepods’ ability to detect nets in daylight and their tendency to avoid them (Fleminger and Clutter, 1965; Giering et al., 2022). This behavior, often linked to their sensitivity to the noise and disturbances created by ship movements and the CPR sampling process, could lead to a significant underestimation of their true population numbers. Furthermore, a comparative study by Baumgartner (2003) on C. finmarchicus abundance estimates, utilizing net samples and optical plankton counters in the GoM, revealed that in patchy distribution environments, substantial concentrations of C. finmarchicus may not be adequately captured by nets, CPR, or instruments. It has also been shown that the CPR in particular tends to underestimate zooplankton populations due to its small aperture, which makes evasion relatively easy (Dippner and Krause, 2013). Further, the in situ CPR data only include stages V and VI C. finmarchicus, whereas the satellite measurements will include all life stages that contain astaxanthin. Consequently, the in situ values would be lower than the satellite estimates, as seen in our results. Collectively, these insights underline the complexities in accurately gauging the population densities of such copepods. The limitations inherent in both traditional zooplankton nets and optical plankton counters pose challenges for ecological modeling like NARW modeling, which depends on these data sources. However, for modeling NARW populations, thresholds of C. finmarchicus concentration become important. Ross et al. (2023) found that C. finmarchicus concentrations exceeding 10,000 ind m⁻³ are indicative of NARW presence. Further, McCarry et al. (2024) found that the minimum detectable level of C. finmarchicus via the MC23 approach was 1,300 ind m⁻³ in Case 1 waters and 8,700 ind m⁻³ in Case 2 waters. Given the potential for underestimation by in situ sampling techniques, it is likely that the Ross et al. threshold would translate into an even larger number for remote sensing observation. Thus, in the optically complex, Case 2 waters of the GoM, the minimum detectable C. finmarchicus concentration from the satellite is worse than the concentration indicative of NARW presence and may in fact be significantly lower.

One potential way to improve the C. finmarchicus quantification would be adopting a spectral matching approach that utilizes the full spectral signatures, moving beyond the three RGB bands used here. This solution, however, has challenges set by the limited availability of spectral bands on heritage ocean color satellites (SeaWiFS, MODIS, and VIIRS) that are capable of effectively differentiating spectral signatures between these species. The new NASA PACE mission will provide global-scale hyperspectral data and could consequently enhance methods for quantifying C. finmarchicus and other optically significant zooplankton. The spectral differences in zooplankton due to different pigments could be detectable via PACE due to its higher spectral resolution; however, McCarry et al. (2024) found that reflectance spectra from astaxanthin-rich and non-astaxanthin-rich zooplankton were similar. Nevertheless, there only exist a limited number of measurements of zooplankton optical properties. Thus, more in situ measurements of zooplankton optical properties, coupled with in situ measurements of zooplankton concentrations (ideally from multiple sampling approaches, e.g., nets, optical plankton counters, and acoustics), are important for developing and refining approaches for ocean color-based zooplankton detection.

The bio-optical model used in this study did not include a component for the backscattering coefficient of C. finmarchicus. Although the backscattering signal of C. finmarchicus is likely minimal and deemed negligible due to their larger particle size (Davies et al., 2021), the development of an instrument to measure large particle backscattering would help to measure the backscattering component, thus potentially enhancing the accuracy of the C. finmarchicus bio-optical model. Additionally, considering the varying concentrations of astaxanthin pigment in different astaxanthin-rich species, which affect their absorption and backscattering properties, and understanding the species-specific astaxanthin pigment concentration could help to differentiate these species and improve the accuracy of species concentration estimates. Thus, it is important to build up a database of in situ zooplankton absorption measurements and quantify the variability in these observations.