Particulate organic carbon in the deep-water region of the Gulf of Mexico

Ocean eddies play a major role in lateral and vertical mixing processes of particulate organic carbon (POC), as well as in the transport of heat, salinity, and biogeochemical tracers. In the open waters of the Gulf of Mexico (GoM), however, there are limited observations to characterize how these mesoscale structures affect the spatial distribution of POC in the upper water column, which is important for organic matter cycling and export. We present the distribution patterns of POC relative to mesoscale features throughout the water column in the deep-water region of the GoM during three oceanographic cruises held during the summer months of 2015, 2016, and 2017. Samples were collected under well-stratified upper ocean conditions, which allowed us to assess the spatial and temporal distribution of POC as a function of non-steric sea surface height, density, apparent oxygen utilization, and chlorophyll fluorescence. We further explored the variability of integrated surface layer POC concentrations at stations located within the cores and the edges of cyclonic and anticyclonic eddies, and those collected outside these structures. Although our results indicate mesoscale eddies modulate several important physical and biogeochemical variables and POC concentrations in the upper ocean, these features do not fully explain the spatial distribution of POC concentrations throughout the deep-water region of the GoM. Relatively lower POC concentrations were observed in the border of the cyclonic and the center of the anticyclonic eddies, in contrast to the relatively higher POC concentrations at the center of the cyclonic and the border of anticyclonic eddies. We observed high variability in POC concentration variability outside mesoscale structures, which may be attributed to other processes such as upwelling over the shelves, and the contribution by rivers during the summer especially in the northern and southern GoM.


Introduction
The reservoir of particulate organic carbon (POC) in the ocean, while considerably smaller than the dissolved organic carbon (DOC) pool (with an approximate POC : DOC ratio of 1:38, Denman et al., 2007), is an essential constituent of sequestered carbon, in addition to the dissolved and particulate inorganic carbon pools (Sarmiento and Gruber, 2006). POC sinks through the water column across isopycnals, scavenging other particles and transporting carbon and associated elements to deep waters, where some of them are remineralized or resuspended, and a minor fraction is eventually sequestered in the sediments (Eppley and Peterson, 1979;Stramska, 2009).
The POC in oceanic regions consists mainly of living organisms and marine snow, which is an aggregate of phytoplankton and small zooplankton fragments, fecal pellets, and other particles also known as organic aggregates, which includes non-living material and uncharacterized components from biological, mineral, and anthropogenic sources (Riley, 1963;Knauer et al., 1979;Dunbar and Berger, 1981;Hedges et al., 2000;Lutz et al., 2002;Stramski et al., 2008;Kharbush et al., 2020). Near active natural oil seeps or oil spills, POC can incorporate "dead carbon" or hydrocarbonderived carbon, in a process that has been described as Marine Oil Snow Sedimentation and Flocculent Accumulation. This process has been shown to be an important pathway for the distribution and destination of spilled oil, which in the Gulf of Mexico accounts for up to 14% of the total oil released as a result of the Deepwater Horizon oil spill in 2010 (Daly et al., 2016).
Particulate organic carbon concentration in the euphotic layer is thought to reflect primary production, a process during which photosynthesis fuels the formation of organic carbon from dissolved inorganic carbon and nutrients, thereby decreasing the partial pressure of carbon dioxide in surface waters (Stramska, 2009). Many of these processes are sensitive to fluctuations in temperature, light availability, ocean circulation, vertical mixing, and nutrient inputs into the surface ocean's euphotic layer (Dawson et al., 2018;Frenger et al., 2018;Dobashi et al., 2022).
In oligotrophic regions, the mixed layer is generally considered as a quasi-homogeneous with little to no variability in density and is usually depleted in nutrients due to uptake by primary producers (Falkowski et al., 1992;Velaśquez-Aristizabal et al., 2022). The deepening of this layer allows for subsurface mixing of nutrients toward the surface, fueling biological productivity and consequently increasing POC concentrations (Gardner et al., 1999). These changes can be further driven by vertical mixing from the isopycnal displacement by anticyclonic and cyclonic eddies, especially under strongly stratified upper ocean conditions during summer (Gaube et al., 2014;Gaube et al., 2019).
Mesoscale processes are responsible for the mixing and transport of heat, salt, and biogeochemical tracers across the global oceans and may act as a moderating factor in global climate change (Faghmous et al., 2015). Thus, understanding ocean eddy dynamics and their role in influencing various oceanographic and biological phenomena is of keen scientific interest (Falkowski, 1998;Faghmous et al., 2015). The productivity enhancing effects of eddies is particularly important in low-nutrient environments, where mesoscale processes can regulate the net upward flux of limiting nutrients as a result of the undulation of nutrient isosurfaces: where the shoaling of isopycnal surfaces tends to bring nutrients into the euphotic zone thereby enhancing productivity, whereas their deepening leads to a lack of change or even a decrease in primary productivity (Gruber et al., 2011).
Remotely sensed sea surface heights (SSH) yield important information on the spatial distribution and intensity of ocean eddies and is strongly related to circulation patterns in the ocean, lending insight into nutrient transport crucial to the understanding of the spatial variability of biological production (Gruber et al., 2011). The measurements performed with various in situ and remote sensing platforms provide an effective tool for studying biogeochemical constituents of oceanic waters like chlorophyll-a, particulate inorganic carbon, and particulate organic carbon (Stramski et al., 2008;Ocean Biology Processing Group, G. M, 2022).
The GoM is a semi-enclosed sea linked to the Atlantic Ocean through the Florida Straits and to the Caribbean Sea through the Yucatan Channel (Oey et al., 2005). The deep-water region of the GoM is characteristic for its oligotrophic and nutrient-limited surface waters, which are relatively isolated from coastal eutrophic waters (Martıńez-López and Zavala-Hidalgo, 2009;Pasqueron de Fommervault et al., 2017). The Loop Current (LC) brings warm water into the GoM, which leads to high stratification in its area of influence, especially during the summer. Higher intrusion of the LC into the GoM is observed more frequently during the summer, which along with the upper ocean warming, influences the circulation in the interior of the gulf (Sturges and Leben, 2000). In the GoM, positive SSH anomalies are indicative of anticyclonic eddies that are characterized by a deep, warm water layer (Müller-Karger et al., 2015;Dobashi et al., 2022).
The LC and detached LC eddies transport Caribbean Subtropical Underwater into the GoM, which is clearly distinguishable from the GoM common water (Vidal et al., 1992). This common water is formed by vertical convective mixing during the winter season, or by the mixing induced by the collision of detached LC eddies against the western GoM boundary (Schmitz et al., 2005). These mesoscale features drive hydrodynamic processes throughout the central GoM that fuel mixing processes and affect biological productivity patterns (Peŕez-Brunius et al., 2012;Müller-Karger et al., 2015;Cervantes-Dıáz et al., 2022).
The most energetic mesoscale events are the detachment of large anticyclonic eddies from the LC that drift to the west (Sturges and Leben, 2000), which are dynamically linked with cyclonic and anticyclonic eddies in the Gulf's interior (Jouanno et al., 2016). The periods between detachments are highly variable (ranging between 0.5 and 18.5 months), implying that there are some years with no eddy detachments from the LC, while during other years these can happen up to three times per year (Sturges and Leben, 2000;Leben, 2005;Vukovich, 2007;Vukovich, 2012;Delgado et al., 2019).
In addition to the quasi periodic shedding of LC eddies, other processes and structures influence the productivity in the southern Bay of Campeche, like the semi-permanent cyclonic eddy in the southern Gulf (Peŕez-Brunius et al., 2012), the confluence of seasonal along-shelf currents (Martıńez-López and Zavala-Hidalgo, 2009), and upwelling processes along the eastern margin of the Bay of Campeche and off the northeast shelf of the Yucatan Peninsula (Zavala-Hidalgo et al., 2006).
The mixed layer depth's (MLD) seasonal cycle is driven by the surface ocean cooling during winter months, fueled by northern winds that lead to the erosion of the stratification of its surface waters and a deepening of the mixed layer (Müller-Karger et al., 1991). This wintertime convection of the MLD due to strong winter winds or "nortes" ends during the late spring, when the warming cycle of the surface waters begins; the MLD becomes shallower throughout the summer into the early fall, when the stratification of the upper layers is strongest (Müller-Karger et al., 1991;Martıńez-Loṕez and Zavala-Hidalgo, 2009). These cycles of winter mixing and summer stratification are thought to be the main controlling factors of the surface chlorophyll-a concentration in GoM (Müller-Karger et al., 2015).
Modeling efforts suggest the higher integrated content of chlorophyll-a (Chl-a) over 350 m occurs in winter whereas lower ocean concentrations happen in summer, while April-May and October-November are transitional periods (Damien et al., 2018). In addition, float data show that mesoscale activity is most likely the main source of variability for the Chl-a concentration in the deepwater region of the GoM, with higher concentrations observed in cyclones compared to anticyclones (Pasqueron de Fommervault et al., 2017;Linacre et al., 2019).
Additionally, cyclonic and anticyclonic eddies drive vertical mixing through the modulation of nutrient transport to the euphotic zone, thereby exerting influence on Chl-a concentrations in the euphotic layer (McGillicuddy, 2016;Honda et al., 2018) as well as POC concentrations (McGillicuddy, 2016). Still, there are limited field observations available to characterize how mesoscale physical processes may affect the POC of the deep-water region of the GoM (Müller-Karger et al., 2015).
POC in the GoM has been mostly studied over the continental shelf, particularly in the north and northeast. The transport of freshwater plumes and filaments from the Mississippi and the Atchafalaya Rivers advect higher nutrient concentration waters to the northern shelf region, especially to the continental shelves of Louisiana and Texas (Biggs, 1992;Trefry et al., 1994;Goñi et al., 1997;Bianchi et al., 1997;Wang et al., 2004). However, the role of mesoscale structures on POC concentrations in the GOM's deepwater region, especially in the southern gulf and the Bay of Campeche, known for its relatively higher productivity within the gulf (Martıńez-Lopez and Zavala-Hidalgo, 2009;Linacre et al., 2015), has yet to be examined. Hence, our study aims to understand the role of mesoscale eddies on POC concentrations in the upper layer of the southern GoM´s deep-water region.

Sampling
For the sampling, 319 discrete seawater samples were collected during three cruises covering the deep-water region of the southern GoM on board the R/V Justo Sierra using Niskin bottles mounted on a rosette (Table 1). At each station, continuous profiles of temperature, pressure, conductivity, and chlorophyll fluorescence were measured with a Sea-Bird ® 9 plus CTD and sensors. In Figure 1, we present a map of the stations and the bathymetry of the deep-water region of GoM. All stations were > 1000 m depth; casts at some stations were done to 1000 m regardless of depth (hereafter referred to as shallow cast stations), and some to about 40 m off the bottom (deep cast stations).

Sampling processing and chemical analysis
Given the low concentration of POC in the oligotrophic waters of the GoM, it was necessary to pool water samples collected at different depths of the water column for analysis (Supplementary Table 1 and Supplementary Figure 1). At the shallow cast stations (<1000 m), the samples comprised depth filtrations of seawater from three discrete depths ranging from 10-50, 75-250, and 300-1000 m (i.e., samples were pooled, and data were obtained for these three broad layers). For the deep cast stations, samples were pooled from 10 m to the depth of maximum fluorescence, between 150-800 m, and from 1000 m to a few meters above the seafloor, the sampling strategy was agreed among all the different participating groups from the beginning of the CIGoM's project. The basic rationale behind the pooling of several depths was the need to filter large volumes of water to obtain the necessary amount of carbon for the analysis. Discrete sampling depths were generally consistent, except for the targeted sampling during each cast of the chlorophyll maximum and oxygen minimum (identified based on the CTD profiles). For visualization purposes, POC values are plotted at the average depth of all the discrete samples collected within any given layer. POC sampling and processing followed the Joint Global Ocean Flux Study methodology of Knap et al. (1996) with some modifications. Specifically, the seawater was filtered through a pre-combusted (500°C for 4 hours) Whatman ® glass microfiber filter (GF/F) with a diameter of 45 mm and 0.7 µm pore size. Typically, between 7.5 L and 60 L of water was filtered to get a minimum weight per determination of 100 mg for POC. Filtered samples were placed in aluminum foil cleaned with ethanol and Kimwipes ® , stored frozen during the cruise and transported to the laboratory for processing. In the laboratory, the samples were freeze-dried and then treated by vaporization with chlorohydric acid (1M) in a desiccator bell to remove the particulate inorganic carbon. Filters were then scraped with an X-Acto ® knife to obtain the POC sample, while minimizing the glass fiber from the filters, and encapsulated in tin capsules.
The concentrations were determined on a Costech ® CHN analyzer at the Stable Isotope Laboratory at CICESE and are reported as µmol/L. The international references used for POC were USGS 40 and IAEA-CH-6. POC concentrations were corrected for blanks (unused filters that had the same treatment as the rest) and calibrated using different weighted samples of sulfanilamide, the detection limit of this method is 4.16 µmolC (50 mgC); during the measurement period the standard deviation of the references was below 0.06 µgC.

Data analysis
2.3.1 Satellite data POC data were not generated for surface layer during our cruises because water samples were used for other analyses. Hence, we used the satellite-derived surface concentrations to complement the POC stations profiles(https://oceancolor.gsfc.nasa.gov/data/aqua/). We obtained daily estimates of POC concentration (mg m -3 ) and Chl-a surface concentrations (mg m -3 ) from the Moderate-Resolution Imaging Spectroradiometer (MODIS-Aqua Satellite) Level 3 mapped product with a spatial resolution of 0.04°(~4.4 Km) for 2015, 2016(NASA Goddard Space Flight Center, et al., 2018. Because of the lack of data at some stations locations due to clouds, sun glint, low light levels, etc., we filled the gaps using the average value of the nearest data points by considering radial increments of 0.04°around the location, up to a fourth iteration when needed (for a maximum radius of 0.16°around a station). If surface POC concentrations were unavailable within maximum, the station was eliminated from further analysis.
The surface POC concentration is based on an empirical relationship between in situ measurements of POC and blue-togreen band ratios of remote sensing reflectance. The algorithms have yielded good performance over vast areas of the oceans, including different hyper oligotrophic and oligotrophic provinces with a determination coefficient between 0.7 to 0.9 (Stramski et al., 2008).

Hydrographic and data and biogeochemical proxies
The in situ density profiles were calculated from conservative temperature, absolute salinity and pressure profiles using GSW Oceanographic Toolbox from The International Thermodynamic Equation of Seawater (TEOS-10) (McDougall and Barker, 2011).
Apparent oxygen utilization (AOU) profiles were calculated as the difference between oxygen gas solubility and the measured oxygen concentration in water using TEOS-10 (McDougall and Barker, 2011).
In the absence of direct turbulent dissipation measurements, the MLD is commonly derived from oceanic profile data using threshold, integral, least squares regression, or other proxies (Thomson and Fine, 2003). In our case the MLD was defined as the depth of the maximum Brunt-Väisälä Frequency, calculated from absolute salinity, conservative temperature, pressure, and latitude using TEOS-10 (McDougall and Barker, 2011).
The maximum fluorescence depth (MFD) was defined as the highest fluorescence value obtained from the vertical profile of each station.

POC data
Quantitative estimates of integrated POC in the surface layer vary depending on the choice of the depth to which the integration is made (Stramska, 2009). We used three different depths to estimate integrated POC: (1) to a constant depth of 100 m (POC 100 ), (2) within the mixed layer at each station, with a boundary defined by the Brunt Väisällä frequency since it provides a link between the physical structure of the water column and the depth of mixing, and (3) to the depth of the fluorescence maximum, as indicative of the depth with the highest chlorophyll fluorescence (Supplementary Material-POC data). This was done to evaluate which criteria better explained the spatial patterns of POC in terms of the physical and biological variables used in this study. To examine the main controls on the POC concentration in the upper layer of the GoM, POC data for each station were interpolated with a Piecewise Cubic Hermite Interpolating Polynomial method with a 1 m resolution (Supplementary Figure 2), which preserves the shape of the data and respects monotonicity (Fritsch and Carlson, 1980).

SSH ns data
To find out how the SSH ns may influence the concentrations of POC in the water column, we correlated the surface POC concentrations against the SSH ns .
We generated two data sets of SSH ns : • Daily sea surface height, which corresponds to SSH ns for the day when the samples were taken. • Monthly sea surface height, which corresponds to the average of the SSH ns during the period of the cruise. Kelly et al. (2021) examined the contributions and mechanisms of the carbon budgets in the GoM and concluded that lateral transport of organic matter is substantial in oligotrophic ocean regions of the northeastern Gulf in the vicinity of the Loop Current and may be crucial to multiple trophic levels in the GoM. To explore the possibility of a lateral transport between the production and formation of POC, we also correlated POC concentrations with SSH ns considered a radius of 0.25 degrees around the stations (daily and monthly SSH (0.25°)).
We also correlated the surface POC concentration of each cruise vs. the surface Chl-a concentrations, MFD, and MLD, to evaluate whether there was a spatial relationship.

Identification and classification of mesoscale structures
Mesoscale structures are mainly generated by ocean large-scale circulation instabilities due to wind or topographic obstacles, creating variability around the ocean´s mean state, and it helps us to understand the role of mesoscale structures on ocean dynamics (including POC production and export) and to discriminate the effect of eddies from other processes (Mason et al., 2014;Pegliasco et al., 2022).
To understand the relationship between mesoscale features and POC concentrations we used "The altimetric Mesoscale Eddy Trajectories Atlas (META3.2 DT allsat) products, processed by CNES/CLS in the DUACS system and distributed by AVISO+ (https://aviso.altimetry.fr)". The algorithm used for these products is derived from Mason et al. (2014) and further described in Pegliasco et al. (2022).
This method follows four steps: filtering, detection, estimation of eddy characteristics, and tracking. For further description about the algorithm check the manual¨Mesoscale Eddy Trajectory Atlas Product Handboook¨(Aviso+ Altimetry, 2022).
After we downloaded the altimetric mesoscale eddy trajectories dataset, we filtered it by date, longitude, and latitude. We used this to identify daily mesoscale structures within the GoM during the three cruises, obtaining a map for each day of the cruise, in total 21, 15, and 25 maps were obtained for 2015, 2016 and 2017 respectively (Table 1). Figure 2 presents an example of the approach followed in our selection of centers and borders in cyclones and anticyclones. To aid in the understanding of the results obtained from eddy identification, we also examined geostrophic flows (downloaded from https://cds.climate.copernicus.eu/#!/home).
We classified eddies based on their rotational direction in the Northern Hemisphere, as either cyclonic (counterclockwise) or anticyclonic (clockwise). The results from the Mesoscale Eddy Trajectories Atlas allow us to visually determine whether our sampling stations were located near the center of a structure or near the border (within a radius less than~0.08°from the contour). We classified stations based on whether they were located near the center of a cyclonic eddy (CCE), near the border of a cyclonic eddy (BCE), near the center of an anticyclonic eddy (CAE), or near the border of an anticyclonic eddy (BAE). Stations that fell outside these features were classified as no eddy (NE). We grouped the stations according to their categories and made composite vertical profiles of fluorescence, density, AOU and POC using pooled data for the three cruises.

Statistical analysis
A two-sided Wilcoxon rank sum test (equivalent to a Mann-Whitney U-test) was used to compare the medians of MLD and MFD between years. A Kruskal-Wallis analysis was used to test for differences in the MLD and MFD between cruises and in surface POC concentrations and integrated POC 100 in each category of mesoscale feature (a =0.05).
Principal component analysis (PCA) was applied to examine associations between the oceanographic variables and the POC concentrations among the station's classifications. When two variables were highly correlated (Pearson's correlation r>0.7), one was removed prior to PCA. The analysis included density at 10 m depth, MLD, MFD, SSH ns , surface Chl-a, the fluorescence maximum, the depth at which the slope of the AOU changes towards positive values, surface salinity, sea surface temperature, the average temperature in the first 100 m and integrated POC 100 concentration.  Figure 3). During this cruise we further observed the lowest SSH ns values in the Bay of Campeche region. In 2016, a large anticyclone called Poseidon (3) had been recently detached from the LC, and the rest of the GoM presented higher SSH ns values than during the 2015 cruise. In 2017, the LC showed a high level of intrusion into the GoM, but there were no large anticyclones within the gulf, and the SSH ns were generally higher than during the 2016 cruise but lower than in 2015.

Oceanographic conditions and POC concentrations
A Kruskal-Wallis test indicated that the median MLD differed between years (p<<0.001). The MLD maps presented in Figure 4 indicate that during 2015 the MLD was deeper in the northwestern part of the study area (~58 m), while the shallowest depths (~20 m) were in the east. MLDs were markedly shallower during the 2016 cruise (Figure 4-middle panel) with an overall mean of 27 m. During that cruise, the Bay of Campeche showed the shallowest MLD, around~15 m. During 2017, the MLD was deeper than in Absolute dynamic topography anomaly map depicting cyclonic and anticyclonic mesoscale structures observed on June 24, 2016. The blue/red line and dot represents the cyclone/anticyclone contour and center; data extracted from the altimetric Mesoscale Eddy Trajectories Atlas at AVISO+. The areas delimited by the dotted lines represent the border and center classification zone for the sampling stations that were in mesoscale structures. The median MFD differed significantly between years, as indicated by the Kruskal-Wallis test (p = 0.012). The MFD maps in Figure 4B (cruise´s average is in Table 2) show that the highest values were consistently located over the central gulf (abyssal region) during the three cruises, with a mean depth during 2015 of 73 ± 24 m, of 88 ± 15 m during 2016, and 86 ± 21 m during 2017.

Characterization of POC profiles and the integrated POC 100
In Figure 5, POC profiles show a pattern consistent with the wellknown exponential decay with depth in the water column (Suess, 1980), with higher concentrations at the surface (between 2 to 7 µmol L -1 ) and decreasing with depth to values less than 1 µmol L -1 below 700m. The sinking velocity of POC aggregates is mainly determined by their size, internal microstructure (porosity and heterogeneous composition), density and the amounts of biogenic calcium carbonate and opal structures that act as ballast minerals (Iversen and Ploug, 2010;Armstrong et al., 2002;De La Rocha and Passow, 2007). In addition, POC's sinking velocity and flux is hypothesized to be modulated by microbial remineralization and by zooplankton grazing (Guidi et al., 2009).
To understand the control of the surface ocean conditions on POC concentrations in the GoM, we considered different integrated depths sampled during the 2015, 2016, and 2017 cruises; the maps are shown in Figure 6 and Supplementary Figure 3, and the average, standard deviations and coefficients of variation are in Table 3 and  Supplementary Table 2.
In Figure

Influence of mesoscale structures on POC concentrations
The interpolated POC, density, fluorescence, and AOU profiles for the stations in each category (CCE, CAE, BCE, BAE, and NE) were pooled across years to examine the influence of mesoscale structures; means and standard deviations were calculated to provide a visual representation of the variability (Figure 7, and Supplementary Table 3).
The CCE (yellow) density profile presented the shallowest MLD, with a shallow pycnocline at~20 m ( Figure 7A -top panel). The MFD was located close to 75 m, with increasing AOU's values below~50 m. POC showed mean concentration of 3.6 µmol L -1 near the surface, decreasing to 0.7 µmol L -1 at 250 m depth. Profiles from BCE stations (purple) showed well-stratified waters, with a relatively deeper mixed layer depth~35 m, and deeper MFD than stations near the center of cyclonic eddies (~75 m depth with a secondary peak at about 100m). The AOU profiles were similar at CCE and BCE stations, showing negative or close to 0 values from the surface to~60 m. Below this depth, consumption and recycling of organic carbon was more important than production by primary productivity, as indicated by increasingly positive AOU values.
Mean POC concentrations were higher at stations near the center of cyclonic eddies (CCE) than near their borders (BCE). CCE stations showed the shallowest of the density profiles, which is consistent with the shallowest pycnocline near the core of cyclonic eddies. Regarding the fluorescence profiles (CCE and BCE), showed a similar maximum fluorescence depth (~75m), however, consistently higher fluorescence values were observed at all depths at stations near the CCE than in the BCE, suggesting higher chlorophyll concentrations throughout the water column. As for AOU, profiles from the CCE showed changes in the slope to positive values at shallower depths than stations in BCE ( Figure 7A -top panel).
The stations located near the core of anticyclones (CAE, blue) showed a weak pycnocline and a MFD between 80 m to 125 m. The AOU profiles showed changes towards positive values at~100 m. In contrast, for BAE stations (orange), the pycnocline was located at Surface and subsurface particulate organic carbon (POC) concentrations at different depths in the GoM for 2015, 2016, and 2017 cruises. POC concentrations at depths below the surface reflect the mean of the three discrete water sampling depths from which water samples were pooled. Integrated particulate organic carbon concentrations for the first 100 m (POC 100 ).  Near-surface POC concentrations at CAE stations were slightly lower (between~1.9 and~2.3 µmol L -1 ) than stations near BAE (between~2.1 and~3 µmol L -1 ) up to~40 m depth (among~0.8 and~2.1 µmol L -1 ); below that depth POC values were very similar (with values ranging between~0.1 and~1.7 µmol L -1 ). Comparison of CAE and BAE density profiles indicated a higher density at a given depth below~30 m, indicative of a deepening of the isopycnals at the core of the anticyclonic eddies. CAE had lower fluorescence values down to 100 m, and a deeper fluorescence maximum than BAE stations. The AOU profiles for CAE show the change in slope to positive values almost 25 m deeper than BAE ( Figure 7B -middle panel).
The mean profiles of density, fluorescence, AOU and POC for NE stations are in Figure 7C. The pycnocline is observed in the first upper~25 m implying a relatively shallow mixed layer, and a relatively deep fluorescence maximum at~80 m. The AOU profile indicates a close balance between production and consumption for the first 50 m, below which respiration begins to dominate as AOU

Correlations between POC and oceanographic variables
There were weak negative and non-significant Pearson correlation coefficients between integrated POC concentrations and SSH ns (r= -0.10 to -0.22), even when SSH was estimated over a broader area around each station (daily SSH and daily SSH (0.25°) as well as monthly SSH and monthly SSH (0.25°)) in Supplementary Tables 5, 6.
The correlation between integrated POC and the MLD, MFD and Chl-a are shown in Table 4 (Supplementary Table 7). Only for the 2015 cruise showed significant correlations between the POC 100 and MLD, MFD and Chl-a. The statistically significant negative correlations most likely result from higher (lower) POC integrated values associated with shallower (deeper) mixed layer depths and a shoaling (deepening) of the maximum fluorescence depth.
The correlations between SSH ns and MLD, MFD and Chl-a are shown in Table 5. Satellite-derived SSH ns and MLD show statistically significant positive correlations during 2015 and 2017. There were also significant positive correlations between SSH ns and. MFD for 2015 and 2016; when SSH ns was negative (positive), the maximum fluorescence depth was found at a shallower (deeper) depth. Negative SSH ns values are associated with divergence in cyclonic structures, which leads to a shallower MLD and MFD, while positive SSH ns values are associated to convergence and anticyclonic structures with the opposite pattern in the MLD and MFD.
Surface Chl-a vs. MFD and MLD consistently showed negative and statistically significant correlations for the three cruises. The negative correlation coefficients imply a deeper (shallower) maximum fluorescence depth associated with lower (higher) concentrations of a Chl-a ( Table 5).
The PCA indicated that components 1 and 2 explained 52.1% of the variability (Component 1 = 33.2% and Component 2 = 18.9%; Figure 8). The eigenvector results are presented in Supplementary Table 8.
The PCA results grouped stations based on their locations at the centers and edges of mesoscale structures, showing the largest separation between stations at the CCE and CAE, and an overlap between BCE and BAE. The CCE stations (yellow) were positively correlated with Chl-a, Fluorescence, salinity and POC 100 , and negatively correlated with the average temperature at 100m, MLD, MFD, SSH ns and to a lesser extent the density. CAE stations (blue) were negatively correlated with Sal, Fluor, Chl-a and POC 100 , and positively correlated with average temperature at 100m, MLD, MFD, SSH ns and density.
Chl-a, fluorescence, integrated POC 100 and salinity cluster in the positive region of component 1 while MFD, average temperature at 100m, SSH ns density and MLD are grouped in the negative region. Most of the stations sampled during 2015 (circles) were correlated with sea surface temperature, MLD, salinity, and fluorescence. In contrast, during 2016 (triangles) the stations were correlated with the average temperature at 100m, density and the depth where the AOU changes towards positive values. The stations sampled in 2017 (squares) fell all over the multivariate space.

Discussion
The mechanisms that control POC concentrations in the ocean differ depending on the region and the season (Dobashi et al., 2022). Sources of organic carbon in the GoM can include terrestrial and  -Pacheco et al. 10.3389/fmars.2023.1095212 Frontiers in Marine Science frontiersin.org marine sources (Bianchi et al., 1997). The contribution of terrestrial organic carbon from terrestrial sources, including estuaries, is particularly important in coastal and shelf waters and constitutes an important component of global organic carbon (Hedges, 1992). Specifically, river inputs play a significant role in the distribution of POC (Meybeck, 1982). The Mississippi-Atchafalaya River and the Grijalva-Usumacinta River Systems are the largest sources of freshwater inflow to the northern and southern GoM, respectively (Kemp et al., 2016). The Mississippi-Atchafalaya River System has been shown to influence the distribution of POC and plays an important role in the organic carbon transport to the shelf and northern oceanic region (Trefry et al., 1994;Rosenheim et al., 2013), while in the southern Gulf the Usumacinta-Grijalva River System is the main source of POC (Cuevas-Lara et al., 2021). In addition, POC fluxes are further influenced by mineral ballast, such as calcium carbonate, the concentrations of which are closely related to POC concentrations (Klaas and Archer, 2002). The concentration profiles in the GoM's deep-water region showed values that range between 2 and 7µM in the upper layer during the summer season. These concentrations are comparable with the values (4.48 ± 5.22 µM) previously reported for samples obtained under oligotrophic conditions in South and North Pacific Gyres in 2004 and 2012 (Stramski et al., 2022; Table 6). Also, POC concentrations in the Equatorial Pacific Ocean collected between 30°N-15°S and 180-160°W in 2012 range between 0.75 and 8.8 µM; the smallest values were from oligotrophic regions, and the highest for the high-nutrient Equatorial Pacific region (Graff et al., 2015). Stramski et al. (2022) sampled different oceanographic conditions ranging from oligotrophic environments like the South Atlantic Gyre to higher productivity regions of the Atlantic Ocean. They report concentrations of 4.77 ± 1.68 µM along north-to-south transects covered in 2005 and 2010. In another study, POC concentrations ranging between 1.67 and 6.33 µM were reported from samples collected in the North Atlantic Drift, North Atlantic Subtropical Gyre, North Atlantic Tropical Gyre, Western Tropical Atlantic, the South Atlantic Gyre, and the South Subtropical Convergence during 2014 (Strubinger-Sandoval et al., 2022).

Contreras
Reports from comparable oligotrophic regions like the Sargasso Sea showed POC concentrations between 1.18 and 2.96 µM from 1988 to 2004 (Duforet-Gaurier et al., 2010; Table 6), while the Hawaii Ocean time-series showed a range between 2 and 5 µM (Sẃirgońand Stramska, 2015). In the southern deep-water region of the GoM, we found POC concentrations that are comparable with those reported for the northern open waters of GoM by Liu and Xue (2020), who measured concentrations of 3 ± 1 µM in May between 2010 and 2013 along the coast of Louisiana in the north GoM. However, they also reported much higher concentrations in the surface waters influenced by the plumes near the Mississippi River Delta (values around 14 ± 6 µM). These high concentrations were attributed to the discharge of the Mississippi and Atchafalaya Rivers.
Hence, POC concentrations in the southern deep-water region of the GoM fall within the range observed in oligotrophic regions like the Sargasso Sea, and the subtropical gyres in the Atlantic and Pacific oceans during the summer season. An analysis of the global distribution of integrated POC 100 from 1995 to 2012 showed a range of values from~20 to~900 mmol m -2 (Balch et al., 2018). The comparison of global estimates with our results is consistent with the oligotrophic nature of the southern deep-water region of GoM during the summer season.
POC concentrations in the water column showed the consistent well-known exponential decay with depth, with relatively high variability in the euphotic layer and decreasing concentrations below 700 m (<1 µM). To understand the mechanisms that could control this upper layer variability in the southern region of GoM, POC concentrations were integrated over different depths (100 m, MLD, MFD) and correlated with different variables. Integrated POC  within the mixed layer depth and the fluorescence maximum depth did not show significant correlations with any of the oceanographic variables for any of the 3 years (Supplementary Table 7). However, when integrating over the euphotic layer (with depths of~100m; Linacre et al., 2019;Kelly et al., 2021;Stukel et al., 2021;Stukel et al., 2022) our results showed a negative correlation with the MFD and the MLD and a positive one with the surface Chl-a concentration for the 2015 cruise. The POC 100 showed similar average integrated concentrations of~276 mmol m -2 for 2015 and 2016, although the variability was much higher in 2015 than in 2016. In contrast, during the 2017 cruise, POC 100 was lower, with a mean value of 207 mmol m -2 . Differences in POC 100 between years were significant due to the high variability between stations. These data suggest that MFD and MLD may occasionally play a role in modulating POC concentrations in the upper water column although it is not a consistent control over time and that there may be other mechanisms contributing to this intermittent relationship.
It is important to mention that although the three cruises took place during the summer season (Table 1), the level of mesoscale activity differed (Figure 3). During 2015 and 2016 mesoscale activity was higher in the central GoM than in 2017. The observed significant correlations between POC 100 and the MLD, MFD and surface Chl-a concentrations for the 2015 cruise may be related to the two large anticyclones that were present in the central deepwater region; these were the remnants of Nautilus eddy and the recently released Olympus eddy (Figure 3). Thus, the higher variability in integrated POC concentrations found for 2015 (Table 3), may be attributed to an intensification of mixing processes associated with the edges of the anticyclonic eddies and consequently with the vertical transport of nutrients fueling the POC production in the euphotic layer. According to the results from the algorithm of AVISO, the average number of anticyclones during the cruise in 2015 was 6.7, in 2016 it was 5.4 and in 2017 was 4.7. Hence, the number of mesoscale structures may be another contributing control to the variability of integrated POC 100, especially during the summer season characteristic for strongly stratified upper ocean conditions.
Large and highly energetic mesoscale eddies exert a strong influence on the biogeochemical variables of the upper water column (Gaube et al., 2019). Eddies are associated with vertical mixing, which affects nutrient transport and the Chl-a concentration in the euphotic layer (McGillicuddy, 2016), and our results in Figure 7 show that the POC upper layer concentrations may be modulated by their location within cyclonic and anticyclonic eddies. However, although we observed significant differences in the medians considering the first 5 m of the water column, we did not find significant differences in the integrated POC 100 medians.
Higher POC concentrations, associated with the shoaling of density profiles and the deep fluorescence maximum are consistently observed at stations near the CCE than at stations near the BCE. These features are consistent with a shallowing (deepening) of the pycnocline near the core (border) of cyclonic eddies. Cyclonic eddies show an anti-clockwise flow of currents with negative SSH ns at its center and counterclockwise geostrophic velocities, a shallower thermocline and a conspicuous shoaling of isopycnals bringing colder, nutrient-rich subsurface waters closer to the euphotic zone at its center (Biggs et al., 1988;Seki et al., 2001), all of them enhancing phytoplankton productivity, higher chlorophyll and POC concentrations in the euphotic layer (Bidigare et al., 2003;McGillicuddy et al., 2003;Bakun, 2006).
The BCE shows consistently lower POC concentrations than stations in the CCE, with close to~0.5 µmol L -1 difference in their mean values (Figure 7). The Chl-a concentrations at the center of a cyclonic eddy are consistently higher than at the borders, most likely as a result of the shoaling of the pycnocline allowing nutrients to reach the euphotic layer and thus fuel the phytoplankton productivity at depth (Gaube et al., 2014;McGillicuddy, 2016). Honda et al. (2018) have shown how cyclonic eddies may affect the POC flux in the oligotrophic region of the subtropical North Pacific Ocean through the shoaling of the pycnocline driven by the passage  [2010][2011][2012][2013] of cyclonic eddies over the observation area, resulting in increased Chl-a in the subsurface layer and enhanced POC flux to the deepsea region. This flux is positively related to the concentration of POC in the water column (Roca-Martı́et al., 2021). In contrast, anticyclonic eddies with clockwise geostrophic velocities (arrows) are convergence zones favoring the piling up of relatively warmer and lower density water masses at their core (Supplementary Figure 4), producing positive SSH ns as seen in Figure 3. Anticyclonic eddies are generally considered nutrientdepleted relative to cyclonic eddies (Biggs et al., 1988;Seki et al., 2001). However, the POC fluxes observed in anticyclones in the eastern tropical North Atlantic show the opposite behavior, with larger fluxes comparable to those observed in mesotrophic or eutrophic regions (Fiedler et al., 2016;Fischer et al., 2016). Here, stations at the CAE mostly showed lower POC average concentrations by~0.4 µmol L -1 (Figure 7) than at those in the BAE. In addition, near the core of the anticyclones there was lower density gradient with depth, lower fluorescence and a deeper MFD. A characteristic transect through an anticyclonic eddy shows a deepening (shallowing) of the pycnocline and MFD near the core (border) of these mesoscale features.
Previous work reported differences in cyclonic and anticyclonic eddies within the GoM. Lee-Sańchez et al. (2022) analyzed vertical nutrient profiles for two cruises, one in the winter (February-March) of 2013 and the other in the beginning of summer (June) of 2016. Their results showed that mesoscale eddies play a strong modulating role in the vertical distribution of nitrate and nitrite, particularly in the first 250 m. In the cores of recently detached Loop Current eddies, the nitracline reached maximum depths, and there were lower nutrient concentrations in the euphotic layer. Data from ARGO measurements in the GoM showed that the depth and the magnitude of the DFM or deep chlorophyll maximum are also strongly controlled by the mesoscale variability, showing consistently lower chlorophyll concentration in anticyclones (Pasqueron de Fommervault et al., 2017 andDamien et al., 2021). Based on a biogeochemical model, primary productivity in the Loop Current Eddies is higher than the average rates in the surrounding open waters of the GoM. This anomalous behavior is explained by the mixed layer response during winter convective mixing, which reaches deeper into nutrient-rich waters due to the very low-density gradient with depth in the the anticyclonic eddy. Although we didn't observe this process, since our data was always collected in the summer season, we coincide with their observation of the strong influence exerted by the mesoscale activity in the surface ocean chlorophyll concentration and more specifically within Loop Current Eddies (Damien et al., 2021).
All these observations show how mesoscale structures do exert some influence in the concentrations of POC in the upper layer. However, correlations of POC concentrations versus SSH ns , which in principle is related to the mesoscale structures, didn't show a significant relationship. To explore this apparent lack of correlation we calculated the Pearson correlation coefficient between satellite POC concentrations and SSH ns in every coordinate of the data mesh of the GoM, considering the average of 8 days for the three summer months of 2015, 2016 and 2017. The relationship is presented in Figure 9, and a clear negative relationship between the POC concentrations in the surface ocean and SSH ns was observed over a broad region in the central abyssal region, especially during 2015 and 2016. This is most likely related to the east-west trajectories of large, detached eddies from the LC. These negative correlations show lower POC concentrations associated with higher SSH ns , implying the importance of these mesoscale structures in controlling POC concentrations in the euphotic layer of the deepwater region of the GoM.
Although our observations show how the POC concentrations may be modulated by mesoscale dynamics in the southern deep-water region of the GoM, it is important to mention that other processes deserve further studies to explain the observed variability. For instance, primary production is associated with seasonal high freshwater inputs, including those from the Tamaulipas shelf in the west, the Mississippi-Atchafalaya River System in the north, and the outflows of the Grijalva-Usumacinta River System in the south (Figure 1). In addition, further studies should be focused on the possible effects of increasing the stability of the upper water column due to rising sea surface temperatures in the interior of the GoM, and how this warming affects mixing in mesoscale features, and consequently primary production and upper layer POC concentrations to better understand the biological carbon pump in the GoM. Spatial correlation maps for the summers of 2015, 2016, and 2017 between surface POC and SSH ns , both derived from satellite data. Pearson's correlation coefficients are shown in colors, with hot colors indicating positive correlation and cold colors negative correlation. The black contours represent a confidence level of 95%.

Data availability statement
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://zenodo.org/record/7302327.

Author contributions
YC: Formal analysis, methodology, visualization, data curation, investigation, conceptualization, writingoriginal draft preparation. SH: Project administration, funding acquisition, writing-review and editing, manuscript improvement. GV: Data analysis, manuscript improvement. JH: Conceptualization, supervision, investigation, methodology, data curation, resources, project administration, funding acquisition, writing-review and editing, manuscript improvement. All authors contributed to the article and approved the submitted version.

Funding
Research funded by the National Council of Science and Technology of Mexico -Mexican Ministry of Energy -Hydrocarbon Trust, project 201441. This is a contribution of the Gulf of Mexico Research Consortium (CIGoM). We acknowledge PEMEX's specific request to the Hydrocarbon Fund to address the environmental effects of oil spills in the Gulf of Mexico. Y.V. Contreras-Pacheco was supported by a CONACyT PhD fellowship.