A series of three incubation experiments were conducted during a mission in the southern Gulf of St. Lawrence on the RV Coriolis II in August 2019 using a flow-through incubator system. Within the 6 flow-through incubators (Coleman^® coolers) were 18 enclosures. The enclosures were composed of PTFE lined buckets equipped with a PVC sampling pipe ( Figure 1 ). The pipe supported a HOBO^® temperature and light logger (Onset, Bourne, MA) and aeration was provided by an aquarium pump attached to a silicone tube that passed through the sampling pipe to an aquarium air stone resting at the bottom of the enclosure. Pumps were adjusted manually to provide gentle bubbling at approximately even rates across all enclosures. Bubbling provided some mixing energy and ensured that the enclosures did not stratify. A submerged pump on the vessel was connected to the incubators that were joined in 2 series of 3 to enable surface water to flow through the incubators at ~ 1-2 L min-1 to moderate the temperature in the enclosures.

Each experiment was designed to include a series of incubation periods and dilutions over 45 h, with sampling occurring at 8 time points ( Figure 1 ). For two experiments, samples were collected following known surface currents (Exp1 and Exp2) to mimic how a surface spill may be diluted ( Supplementary Figure 1 ). The third experiment aimed to follow the Orphan Bank Gyre, but likely sampled against the main flow of the Gaspe Current (Exp3). For each experiment and dilution, 120 L of surface water (2 m) was collected using a CTD-rosette system. In situ temperature, salinity and chlorophyll a (chl a) concentrations estimated from fluorescence were recorded for each sampling. Water from multiple 12 L Niskin bottles were pooled and distributed into 9 enclosures. At the beginning of the experiment (T0), one of three treatments was randomly assigned to each enclosure, with each treatment occurring once per incubator. For the oil treatments, either Bunker C or Heidrun was added to the surface of the seawater and the enclosure sampled immediately following the addition. Bunker C has been a commonly used fuel for decades, while restrictions on sulfur content in exhaust are driving a shift from bunker fuels to alternative, many vessels still use this fuel. Heidrun is a medium crude oil that represents the types of crude oils that are imported to Canada for refining. Both could be released from vessels transiting the Gulf of St. Lawrence. The control enclosures received no additions. Sampling was carried out from enclosure 1 to 9 each time. To facilitate addition of oil, small vials of Bunker C and Heidrun were warmed to 50°C and poured into each enclosure. Heating decreased the viscosity of the oil, but may increase volatilization when the vial was opened and the oil poured. To determine the amount of oil added, each vial was weighed before and after addition of the oil.

Samples were collected by placing a silicone tube through the sampling pipe and pumping water into sample vessels with a peristaltic pump from below the surface. The sampling pipe enabled the tubing to be lowered to the bottom of the enclosure without passing through any surface slick. After the initial addition of oil and the first sampling, the enclosures were incubated for 9 h before being sampled again (Incubation period 1 = I1, Figure 1 ). At 12 h, a second CTD cast collected surface water from a location down current from the first site ( Supplementary Figure 1 ) and water was distributed into the second set of 9 enclosures. Using the peristaltic pump at a flow rate of 1.9 L min^-1, 5 L of water from the first enclosure was pumped into an enclosure containing 7 L of freshly collected surface water, resulting in a 2.4-fold dilution of the original incubation. The diluted enclosure was immediately sampled (Dilution 1 = D1) and then the remaining incubations were diluted and sampled in order. The new incubations were sampled after an additional 9 h (21 h post oiling = I2) and then diluted 3 h later (24 h post oiling = D2). This was repeated at 33 h post oiling (I3), 36 h (D3) with a final sample collected at 45 h post oiling (I4). The final sample (45 h) was not collected for Exp2 due to rough weather and rolling of the ship that caused waves to wash over the enclosures.

At each time point, several different types of samples were collected, with <2.5 L collected at each sampling time. To analyze monocyclic aromatic compounds [benzene, toluene, ethylbenzene and xylenes (BTEX)], 40 mL of sample were collected into purge and trap vials containing 40 µL of HCl and stored at 4°C prior to analysis. Samples for total petroleum hydrocarbons (TPH) and quantification of specific polycyclic aromatic hydrocarbons (PAHs) were collected separately by pumping approximately 800 mL of water into two 1 L amber bottles. Dichloromethane (DCM) was added and each bottle was shaken well before being stored and transported to the Bedford Institute of Oceanography (BIO), Dartmouth, Nova Scotia, at 4°C. Approximately 450 mL of water was then collected in sterilized media bottles and subsampled in the wet lab on the vessel.

To determine the abundances of prokaryotes, duplicate 980 µL subsamples were pipetted into 2 mL cryovials containing 25 µL of 25% EM-grade glutaraldehyde (Cobanli et al., 2022). Duplicate 3955 µL subsamples were pipetted into 5 mL cryovials containing 40 µL of 25% EM-grade glutaraldehyde and 4 µL of 10% Pluronic F68 to determine the abundance of small phytoplankton (Marie et al., 2014). After 10 min on ice, the vials were flash frozen in liquid nitrogen and stored at -80°C. To measure chl a concentrations, water was filtered through a 25 mm GF/f filter (nominal pore size ~ 0.7 µm). The volume was recorded to calculate the sample concentrations. Each filter was placed in a glass scintillation vial and stored at -80°C. Approximately 30 mL of sample was poured into acid-cleaned Nalgene bottles for analysis of total inorganic (TIC) and organic carbon (TOC). Bottles were stored at -20°. Duplicate subsamples were poured into 10 mL vials and frozen at -20° to measure the concentrations of inorganic nutrients. To characterize the microbial diversity using amplicon sequencing, 100 mL were filtered through a 47 mm, 0.22 µm Millipore Durapore^® PVDF filter. The filter was placed in a sterile 50 mL polycarbonate tube and frozen at -80°C. A final subsample of water was collected and analyzed using the Microtox^® M500 system (Modern Water, New Castle, DE) on the vessel. All frozen samples were transported on dry ice to the BIO for processing.

Measured parameters were analyzed to identify if there were significant treatment effects using generalized additive mixed models (GAMM) in the mgcv 1.9.0 package (Wood, 2004; Wood, 2011; Wood et al., 2017) in R Statistical Software (v4.3.0; R CoreTeam, 2023). TPH, ΣPAHs and Σalkylated PAHs were log transformed prior to analysis after adding the reporting limit (TPH + 0.05, ΣPAHs/Σalkylated PAHs + 32.92) to compensate for the zero values. As the main question focused on the effects of oil addition and not on the specific patterns within each experiment, the data was pooled, with 27 mesocosms (9 from each experiment) over 3 experiments (n_{total samples}=207). The gam function was used with mesocosm nested in experiment and fit as random factors. Four separate models ( Supplementary Table 1 ) were fit for each parameter using Restricted Maximum Likelihood and the models compared using ΔAIC (AIC_model-AIC_min). The simplest model with a ΔAIC <2 was selected. All models included the random terms as well as a smooth term for time as a continuous variable. Additional models included a smooth term for time depending on treatment (control, Bunker C or Heidrun), treatment as a fixed parametric term or both treatment terms (Alguero-Muniz et al., 2017; Calderó-Pascual et al., 2021). The reference treatment was control. Models for the abundance of prokaryotes and the three phytoplankton populations as well as the concentration of chl a were analyzed using a Gamma distribution with log as the link function. The other models used the default Gaussian distribution with identity as the link function. Pearson correlations between parameters were calculated using the Rcmdr v2.8.0 package (Fox et al., 2023) with the Holm’s p-value correction using only data from the Heidrun and Bunker C enclosures.

Changes in prokaryote and microeukaryote community structure due to treatment were detected using Permutational Multivariate Analysis of Variance (PERMANOVA). The avgdist function in vegan (Oksanen et al., 2022) was run with 100 interactions to resample 20,000 or 15,000 reads from the 16S and 18S rRNA amplicon tables, respectively. The Bray-Curtis distance was calculated for each iteration and final distance matrices were calculated based on the average distances. Prior to PERMANOVA, the betadisper function was used to compare the communities at T0 for the three experiments. While distances between the centroids within an experiment were low, larger values between experiments suggested that the different T0 locations resulted in different initial communities. To increase the likelihood of detecting treatment effects, PERMANOVA was conducted separately for each experiment. The average distance matrices were passed to the adonis2 function where fixed factors included treatment, time and the interaction term and permutations were restricted based on the mesocosm (Simpson, 2022). To visualize how community structure varied with time and experiment, the 10 most abundant OTUs from each experiment were identified for both prokaryotes and microeukaryote communities. The average relative abundances at the end of each incubation period were plotted over time for each treatment and experiment. Overlap in abundant OTUs from the three experiments meant that there were 17 and 20 unique OTUs visualized for prokaryotes and microeukaryotes, respectively.

Two different oils were used in these experiments ( Table 1 ). Bunker C represents a heavy fuel oil compared to Heidrun, which is a medium crude oil from Norway. Both products have low percentages of BTEX (<0.5% of total content). n-alkanes were over 2-fold higher in Heidrun compared to Bunker C, which had higher amounts of parent and alkylated PAHs on a per g basis. Most PAHs in Bunker C were ≥4-ring, while Heidrun contained mostly 2-ring PAHs. Densities of both oils indicate they should float on seawater. Heating of oils facilitated pouring into enclosures. This approach allowed effective addition of Heidrun in all experiments (4.23 ± 0.05 g, n=9), with more variability in the amount of the more viscous Bunker C added (3.52 ± 0.91 g, n=9).

The experiments were all carried out using surface water collected from the southern Gulf of St. Lawrence, but the in situ conditions at the different sample locations varied (for details see the Supplementary Materials ). Exp1 was influenced by the Resitgouche River flowing into Chaleur Bay, while Exp2 and Exp3 were collected from the open Gulf ( Supplementary Figure 1 ). This had impacts on the salinities and chl a concentrations in these samples ( Supplementary Figure 2 ), although the total salinity range was small. Chl a was generally highest at the beginning of Exp1. Temperature in the mesocosms was slightly higher than measured in situ and varied across the samples ( Supplementary Figure 3 ). There was no significant impact of treatment on temperature or light for any of the experiments ( Supplementary Table 2 ). Experiments were started in the evening, with I2 and I4 occurring during the day ( Supplementary Figure 3B ).

Location of sampling also had an impact on the concentrations of NH₄ ⁺, silicate, PO₄ ^-3, TOC or TIC. Dilution appeared to have a larger impact on concentrations compared to addition of oil ( Supplementary Table 2 , Supplementary Figure 4 ). Concentrations of NO₃ ^- were found to differ significantly due to treatment (Model 2, Supplementary Table 2 ), with overall concentrations in Heidrun enclosures higher compared to the control or Bunker C enclosures, although this was not as obvious in Exp1 as in Exp2 and Exp3 ( Supplementary Figure 5 ). In all experiments and treatments, NO₃ ^- concentrations were low, with NH₄ ⁺ dominating the nitrogen pool. Nutrients were not limiting during the experiments.

Only 7 samples (n=206) had BTEX concentrations above the detection limit. None of these samples were from Exp2, and only 1 (Bunker C, 9 h) was from Exp3. From Exp1, samples with detectable BTEX included the three enclosures in 1 incubator at 9 and 12 h, suggesting possible contamination of the enclosures or the samples. Because of the few samples with detectable BTEX, no further analysis could be carried out, but low concentrations likely reflect the low initial concentrations of BTEX in the oil and potential volatilization following heating.

Water column hydrocarbons increased in the Bunker C and Heidrun enclosures, with no corresponding increase in the control enclosures ( Supplementary Materials ). TPH, ΣPAHs and Σalkylated PAHs showed similar patterns in all three experiments ( Figures 2A–C ). The concentrations increased during I1, decreased with dilution and then decreased further during I2. GAMM analyses indicated that the patterns in TPH, ΣPAHs and Σalkylated PAHs over time and absolute concentrations were affected by treatment (Model 4, Supplementary Table 2 ). TPH was highest in Heidrun treatments initially before becoming similar to Bunker C after I2. The higher TPH measured for Heidrun samples compared to Bunker C reflects the higher contribution of n-alkanes in the fresh Heidrun. There were only minor differences in the patterns of ΣPAHs between the two oils over time, with Heidrun slightly higher at the beginning of the experiments compared to Bunker C. While there was no difference between Heidrun and control concentrations after D2, ΣPAHs with Bunker C remained elevated relative to control enclosures up to I3. Similar patterns were observed for Σalkylated PAHs.

Both ΣPAHs and Σalkylated PAHs were dominated by 2-ring PAHs at I1 and D2 when concentrations were the highest ( Supplementary Table 3 ). The ≥4-ring PAHs that make up most of the PAHs of the Bunker C ( Table 1 ) have reduced solubility in water (Abdel-Shafy and Mansour, 2016), resulting in most of the PAHs detected in the Bunker C enclosures being naphthalene and its alkylated homologs and similar concentration between the two oils. Thus, while the overall concentrations of PAHs in oil can be useful in predicting water column concentrations after a spill, the composition of the PAHs and their solubility in the water column is an important consideration (Stoyanovich et al., 2022).

Although the concentration of hydrocarbons peaked at the end of I1, it is likely that more dissolution and dispersion would have occurred if the incubation period was longer. Following D1, no surface hydrocarbon slick was observed in any of the enclosures and no additional hydrocarbons could enter the water. During previous spills, concentrations of hydrocarbons in the water column have been observed to decrease over time. Concentrations of PAHs in the water column under the oil spilled from the Exxon Valdez reached 12 µg L^-1 7 days after the spill, but most measurements outside of boomed areas averaged 0.06 ± 0.4 µg L^-1 throughout the following month (Boehm et al., 2007). Elevated concentrations were observed under oil slicks, indicating continual inputs from surface oil. In the month following the sinking of the Prestige (December 2022), hydrocarbon concentrations were determined to range from 0.19 to 28.8 µg L^-1 of oil equivalents (González et al., 2006). By February 2003 concentrations had decreased to <0.15 to 17.4 µg L^-1 of oil equivalents. Unfortunately, during a spill samples are rarely collected during the initial hours or days, so the dilution rates and peak hydrocarbon concentrations are unknown. Two controlled experimental spills conducted in the North Sea determined that maximum concentrations of hydrocarbons in the water below the slick occurred between 10 and 24 h, with decreasing concentrations in the surface layer after this time (Brussaard et al., 2016).

Light inhibition measured with the Microtox^® method, represents the potential toxicity or bioavailability of hydrocarbons to prokaryotes. Values were similar in Heidrun and Bunker C enclosures ( Figure 2D ). GAMM analysis determined that the best model included treatment as a factor affecting the absolute amount of light inhibition and its pattern over time (Model 4, Supplementary Table 2 ). From the GAMM analysis, inhibition in Bunker C was lower than in Heidrun over the first part of the experiments, but the overall difference between the oils likely had a minimal ecological effect. Both oils had higher potential toxicity than controls until the end of I2. Significant correlations were identified between light inhibition and Σalkylated PAHs (ρ = 0.88, p<0.0001, n=207) and ΣPAHs (ρ = 0.90, p<0.0001, n=207), with a weaker correlation with TPH (ρ = 0.71, p<0.0001, n=207).

Increasing potential toxicity with increasing PAH concentration has previously been observed in other oil exposure experiments (Mueller et al., 2003; Brussaard et al., 2016; Ortmann et al., 2020), with the specific PAHs being important in determining potential toxic effects (Di Toro et al., 2007; Ball and Truskewycz, 2013). The strong correlation between light inhibition and PAHs is driven by the dominance of 2-ring PAHs in the water column ( Supplementary Table 3 ), which are known to have higher toxicity than PAHs with more rings (Di Toro et al., 2007). While the maximum light inhibitions measured in these experiments were too low to estimate an EC₅₀, total PAH concentrations at I1 averaged 39.1 ± 9.5 µg L^-1 in Heidrun enclosures and 31.8 ± 11.7 µg L^-1 with Bunker C, exceeding concentrations expected to have non-lethal impacts on a range of organisms (Almeda et al., 2013; Barron et al., 2018). However, with high rates of dilution, exposure times to these high concentrations would be relatively short, likely limiting impacts. As each dilution represented a 2.4-fold dilution, the original seawater was diluted 13.8 times after D3. Based on the potential toxicity measured at D3 (Control = 2.82 ± 3.49%, Bunker C = 1.35 ± 2.05%, Heidrun = 1.06 ± 2.12%, light inhibition, n=9 per treatment), the conditions in the oil treatments reflected background conditions within the Gulf. The average concentration of all measured PAHs was similar in the oil treatments compared to the control treatments (Bunker C = 0.81 ± 0.24 µg L^-1, Heidrun = 0.64 ± 0.34 µg L^-1, Control = 0.65 ± 0.38 µg L^-1), below concentrations that would be expected to be toxic.

Abundances of prokaryotes were not affected by treatment ( Figure 3 ; Supplementary Table 2 ). Generally, prokaryote abundances increased during I1, I3 and I4, but decreased through I2 and with most dilutions. The short time frame of the incubation period may have prevented the abundance of prokaryotes with oil from increasing relative to the controls. In many oil exposure experiments, a lag between the peak hydrocarbon concentrations and the increase in prokaryotes has been observed (Gutierrez et al., 2013; Chronopoulou et al., 2015; Ryther et al., 2022). During this lag, hydrocarbon degrading organisms increase, eventually outcompeting the non-hydrocarbon degraders present in the initial community. Previous studies have identified delays in the response of prokaryote communities to hydrocarbons, with lag times of 1 – 4 days (Jung et al., 2010; Lekunberri et al., 2010; Ortmann et al., 2012; Ryther et al., 2022), compared to the 9 h incubations in this study.

Prokaryote community diversity was monitored using 16S rRNA gene amplicons, representing both Archaea and Bacteria, although the communities were dominated by Bacteria. Richness, measured using the Chao1 estimate, did not show any significant treatment effect, with the average Chao1 richness being 311.55 ± 45.13 (n=206) ( Supplementary Table 2 ; Figure 4A ); however, there was a significant effect of treatment on the pattern of the Shannon diversity for the prokaryote community over time (Model 3, Supplementary Table 2 ; Figure 4C ). Diversity of prokaryotes in Heidrun enclosures tended to be lower than that of the controls initially (Heidrun_T0-I2 = 3.59 ± 0.13, Control_T0-I2 = 3.61 ± 0.15), and higher than controls after D2 (Heidrun_D2-I4 = 3.57 ± 0.11, Control_D2-I4 = 3.49 ± 0.12), although the effect was small. Richness and diversity tended to increase or decrease after each dilution due to differences between the communities in the incubations and those in the seawater in situ.

The initial prokaryote community structures in the three experiments were determined to be quite different from each other, with Exp2 and Exp 3 being more similar to each other than to Exp1 ( Supplementary Table 4 ). As differences between experiments were higher than within experiments, each experiment was analyzed separately to identify treatment effects. There were significant effects of both treatment and time on the structure of the prokaryote communities based on PERMANOVA analysis, but the interaction term was not significant for any of the analyses ( Supplementary Table 5 ). Time had the largest impact on the community structure, with 24 – 52% of the variation attributed to time and <7% of the variation due to treatment.

An examination of the most abundant OTUs across all three experiments, highlights the similarities in the communities across the Gulf and among treatments, even with differences in the overall community structures. Several of the abundant genera include known or putative hydrocarbon degraders. The prokaryote communities were always dominated by SAR11 Clade 1a (OTU-0), with Synechococcus (OTU-01), Planktomarina (OTU-02) and Actinomarina (OTU-03), occurring within the 5 most abundant OTUs in all experiments and treatments ( Figure 5 ). Planktomarina has been shown to have genes associated with hydrocarbon degradation, although degradation activity has not been experimentally shown (Voget et al., 2015; Peeb et al., 2022). Several abundant OTUs were identified as members of the Flavobacteriaceae (OTU-04, OTU-07, OTU-14 and OTU-18). Overall, the relative abundances of these OTUs were constant across treatments. Several members of the Flavobacteriaceae have been shown to degrade hydrocarbons, with many associated with degradation of high molecular weight organic matter (Vila et al., 2010). Three Alteromondaceae OTUs were observed in the experiments, with one increasing in abundance in all three treatments, especially in Exp3 (OTU-05), another identified as Alteromonas (OTU-15), reaching higher average relative abundances with Heidrun and Bunker C and Glaciecola (OTU-12) increasing in all treatments, with the highest average relative abundances with Bunker C. Alteromonas, along with other members of the Alteromonadaceae, are known to be important degraders of PAHs (Jin et al., 2012; Martinez-Varela et al., 2022). While Glaciecola is frequently observed to increase in oil exposed seawater, its role in hydrocarbon degradation has not been confirmed (Sauret et al., 2014; Martinez-Varela et al., 2022; Peeb et al., 2022). An OTU identified as Marinomonas (OTU-30) also increased in all treatments over time in Exp1, but reached higher average relative abundances with oil and may represent a facultative PAH degrading organism (Melcher et al., 2002). Although many of the OTUs that showed increases in relative abundances with oil addition can degrade hydrocarbons, none are known to be obligate hydrocarbon degraders.

Several of the other, poorly described clades (e.g., SAR86, SAR116 and OM43) may also have the potential to degrade hydrocarbons, but are unlikely to be obligate hydrocarbon degraders (Doyle et al., 2018; Sieradzki et al., 2021; Vigneron et al., 2023). The small amount of variability within the communities attributed to treatment may be a result of an initial community that already contains facultative hydrocarbon degraders, thus no wide spread change in the community structure was observed. In the absence of a radical change in the prokaryote community, decreases in hydrocarbons during I2 could be due to biodegradation being carried out by organisms already at high abundances.

Oil addition had a significant impact on some populations of phytoplankton based on GAMM analyses ( Supplementary Table 2 ). The temporal trends in the chl a concentrations were influenced by treatment (Model 3, Supplementary Table 2 ), but the absolute concentrations were not ( Figure 6A ). While chl a concentrations decreased with oil compared to controls during I1, I2 saw a larger decrease, with the biggest decrease in Bunker C. After D2, abundances increased in oil treatments, eventually reaching the same concentrations as control enclosures by the end of the experiments. Eight different populations could be differentiated by flow cytometry, which were combined into three groups for analysis ( Supplementary Figure 6 ). Five populations dominated by chl a were combined (eukaryote phytoplankton) and showed similar patterns to the chl a concentrations with the temporal patterns affected by treatment (Model 3, Supplementary Table 2 ; Figure 6B ). This suggests that many of the small eukaryotic phytoplankton were well represented in the total chl a estimate, which is supported by the strong correlations between the abundance of the eukaryote phytoplankton and chl a concentrations (ρ = 0.8025, p<0.0001, n = 205). Phytoplankton enumerated by flow cytometry were dominated by Synechococcus-type cyanobacteria, which were likely not captured by the chl a analysis (ρ = -0.05, p = 0.47, n = 205) ( Figure 6C ). The abundance of Synechococcus varied over time, but there was no significant effect of treatment on abundances ( Supplementary Table 2 ). The patterns observed here differ from those observed during 2 controlled spill experiments, where natural dilution occurred (Brussaard et al., 2016). In those spills, no decrease in to phytoplankton abundances were observed; however, decreases in chl a were detected under the slick 1 day after the spill compared to samples collected from outside the slick.

Two of the phytoplankton populations were impacted by nanodroplets (~ 1 – 5 µm) of oil, especially in the Heidrun treatments ( Supplementary Figure 6 ). These populations, when pooled, were strongly correlated with TPH (ρ = 0.7980, <0.0001, n = 205), Σalkylated PAHs (ρ = 0.6570, <0.0001, n = 205) and ΣPAHs (ρ = 0.5271, <0.0001, n = 205), and showed similar patterns to the hydrocarbon concentrations ( Figures 2 , 6D ). A weak correlation between these populations and chl a (ρ = 0.3363, p<0.0001) indicates that some of these particles are phytoplankton cells; however, cells could not be differentiated from the microdroplets. As for the hydrocarbon concentrations, there was a significant effect of treatment on the trends in abundance of these populations over time (Model 3, Supplementary Table 2 ), with significantly higher abundances in Heidrun from the beginning of the experiment through to D2, after which abundances were lower than those observed in control and Bunker C enclosures. The mechanism, fate and impacts of these oil nanodroplets remains to be determined. They are smaller than the droplets implicated in marine snow formation (Passow, 2016), but could coalesce over time, increasing their diameter and either surfacing or contributing to marine snow.

A decrease in richness and diversity in oil treatments at I2 can be observed for the microeukaryote community, but this was not significantly different from the patterns observed in the control enclosures ( Supplementary Table 2 ; Figures 4B, D ). This was likely due to the relatively large variation across replicates. PERMANOVA identified significant differences in microeukaryote community structure for all three experiments, with main effects but no interactions ( Supplementary Table 5 ). As for the prokaryotes (24 – 52%), variation attributed to time was high (19 – 46%) with <3% of the variation attributed to treatment for the microeukaryote community.

The largest impact of oil on the abundance, richness and diversity of the microeukaryote community was at I2. Photooxidation may have contributed to the decrease in hydrocarbons during I2 ( Figure 2 ) and may have produced toxic byproducts that reduced the abundances of phytoplankton and the diversity of the microeukaryotes in oil treatments (Maki et al., 2001; Prince et al., 2003; Aeppli et al., 2012). The experiments were started in the evening, with I1 occurring during the night. In the dark, with high hydrocarbon concentrations, there was little impact on the abundance or diversity of the microeukaryote community. I2 occurred under full sunlight. While not as dramatic as the decrease observed for chl a and small eukaryotic phytoplankton, there was a small decrease in the abundance of prokaryotes ( Figure 3 ) over this period, which could have been due to an increase in toxic photoproducts.

Treatment had a low overall impact on the relative abundances of the most abundant microeukaryotic OTUs ( Figure 7 ), with the largest differences in the relative abundances of some taxa at I2 with oil compared to within the control enclosures. The largest differences were observed for Bunker C. The most abundant sequence across all experiments was a dinoflagellate (OTU-01), with the exception of the I2 oil samples from Exp1 when Chaetoceros (OTU-11) became most abundant in both Heidrun and Bunker C enclosures. Bathycoccus (OTU-03), a green algae, was abundant, with relatively consistent abundances across time and treatment in all experiments. In the North Sea experimental spills, Chlorophyceae were detected 2 h after the spill at 1.5 m, but were not detected after 24 h when PAH concentrations peaked (Brussaard et al., 2016). Concentrations were higher in the study during I1 compared to the North Sea spills and it is possible that Bathycoccus (OTU-03) is adapted to higher PAHs and may contribute to PAH degradation in conjunction with co-occurring bacteria (Borde et al., 2003). Another group that was shown to be highly sensitive to oil in the North Sea experiments were the Cryptophyceae (Brussaard et al., 2016). In this study, Teleaulax (OTU-10), affiliated with the Cryptophyceae, decreased in all oil treatments at I2, further supporting the sensitivity of this phyla to oil spills.

There were two abundant OTUs affiliated with the parasitic Syndiniales groups. One, affiliated with group 1 (OTU-02), decreased at I2 with Bunker C relative to controls in all experiments, with smaller decreases in response to Heidrun. The other OTU, affiliated with group II (OTU-08), was more abundant in Exp2 and Exp3 compared to Exp1, but decreased in all experiments in the enclosures with oil. Decreases in the abundance of parasitic could reflect decreases in their hosts in response to hydrocarbons or high susceptibility of their zoospores to oil.

Few abundant OTUs were observed to increase in abundance with oil. Leucocryptos (OTU-13) was elevated in the oiled enclosures in all experiments at I2 compared to the control enclosures, although its relative abundance only reached ~10% of the total population. This organism is a phagotrophic flagellate that was also observed to increase in response to diluted bitumen exposure (Cobanli et al., 2022). This organism may be responsible for grazing on the prokaryote community, potentially preventing an increase in the abundance of prokaryotes in response to oil. This increase in a heterotrophic flagellate is in contrast to that observed in the North Sea experimental spills, where both nanoflagellates and ciliates decreased significantly in response to oil (Brussaard et al., 2016). Ciliates can also contribute to high grazing rates on both prokaryotes and phytoplankton. In this study, Pseduotontonia (OTU-06) showed different responses to oil across the experiments. While abundances remained relatively constant in Exp1, they decreased in Exp2 and Exp3, with slightly larger decreases with oil. Cryptocaron (OTU-14) showed similar patterns, except for Exp1 where relative abundances were often below detection. The reduction in ciliates, would suggest that any grazing that may have limited the abundances of prokaryotes or phytoplankton was likely due to flagellates or dinoflagellates. Three dinoflagellate OTUs were abundant throughout these experiments, including the most abundant OTU, which has the potential to exert significant grazing pressure on the microbial communities. Polykrikos (OTU-04) also increased with Bunker C in Exp2 and may have contributed to the observed decrease in phytoplankton through grazing (Matsuyama et al., 1999; Lee et al., 2015). Biecheleria (OTU-31), the third dinoflagellate, was relatively constant over the experiments, with a slight decrease at I2, but as a mixotroph, may also have exerted grazing pressure on the prokaryote and phytoplankton communities (You et al., 2023); however, Biecheleria is also known prey item for Polykrikos (Yoo et al, 2013), highlighting the complexity of interactions in the marine microbial communities.

In Exp1, Chaetoceros, a diatom, had higher relative abundances in both the Bunker C and the Heidrun enclosures at I2 compared to the control enclosures. Chaetoceros was rarely detected in Exp2 and Exp3 and may have been restricted in abundance outside of Chaleur Bay. While the diatom Pseudo-nitzschia (OTU-09) increased in the control and Heidrun incubations in Exp2 at I2. Chaetoceros and Pseudo-nitzschia are both diatoms and grew well in some of the oil enclosures as has previously been observed (Ozhan and Bargu, 2014; Melliti Ben Garali et al., 2021). It has been suggested that diatoms may be somewhat resistant to oiling (Quigg et al., 2021), although Pseudo-nitzschia appeared to be inhibited by Bunker C during Exp2 and Thalassiosira (OTU-24) grew poorly with Bunker C in Exp1. These patterns may be a result of competition between diatom species for resources. Bacteria associated with the phycosphere of some diatoms have been shown to be capable of degrading PAHs, which could contribute to their success in oiled conditions (Gutierrez et al., 2014), while more susceptible diatoms may lack hydrocarbon degraders (Severin and Erdner, 2019). Diatoms did not appear to be well supported by the environmental conditions in Exp3, highlighting that variations in environmental conditions can select for different organisms in response to oil spills.