Coastal upwelling systems as dynamic mosaics of bacterioplankton functional specialization

Coastal upwelling areas are extraordinarily productive environments where prokaryotic communities, the principal remineralizers of dissolved organic matter (DOM), rapidly respond to phytoplankton bloom and decay dynamics. Nevertheless, the extent of variability of key microbial functions in such dynamic waters remains largely unconstrained. Our metatranscriptomics analyses of 162 marker genes encoding ecologically relevant prokaryotic functions showed distinct spatial-temporal patterns in the NW Iberian Peninsula upwelling area. Short-term (daily) changes in speci ﬁ c bacterial functions associated with changes in biotic and abiotic factors were superimposed on seasonal variability. Taxonomic and functional specialization of prokaryotic communities, based mostly on different resource acquisition strategies, was observed. Our results uncovered the potential in ﬂ uence of prokaryotic functioning on phytoplankton bloom composition and development (e.g., Cellvibrionales and Flavobacteriales increased relative gene expression related to vitamin B12 and siderophore metabolisms during Chaetoceros and Dinophyceae summer blooms). Notably, bacterial adjustments to C-or N-limitation and DMSP availability during summer phytoplankton blooms and different spatial-temporal patterns of variability in the expression of genes with different phosphate af ﬁ nity indicated a complex role of resource availability in structuring bacterial communities in this upwelling system. Also, a crucial role of Cellvibrionales in the degradation of DOM (carbohydrate metabolism, TCA cycle, proteorhodopsin, ammonium, and phosphate uptake genes) during the summer phytoplankton bloom was found. Overall, this dataset revealed an intertwined mosaic of microbial interactions and nutrient utilization patterns along a spatial-temporal gradient that needs to be considered if we aim to understand the biogeochemical processes in some of the most productive ecosystems in the world´s oceans.


Introduction
Coastal upwelling areas are complex and dynamic environments where winds displace surface water off the coast, allowing nutrientrich deep water to reach the surface and driving remarkable productivity at all trophic levels (ultimately sustaining major portions of global fisheries) (Kämpf and Chapman, 2016;Stock et al., 2017) [see Kämpf, J. and Chapman, P 2016 and references therein].Bacterioplankton is a crucial component of the ocean biome that has long been recognized for its principal role in regulating the biogeochemical cycling of carbon (Rivkin and Legendre, 2001), but also the remineralization and utilization of elements like N and P that are vital for ocean productivity (Caron, 1994).Moreover, recent findings emphasize the role of bacteria and their interactions with phytoplankton in providing vitamins and growth factors necessary for organisms throughout the food web (Shilova et al., 2014;Amin et al., 2015).Still, the specific adjustments in bacterial metabolism that regulate such processes in upwelling areas, and the taxa involved, remain largely unexplored.
The usage of 'meta-omics' techniques has allowed the scientific community to gain comprehensive understanding on the succession of specific prokaryotic taxa and their metabolisms operating during the evolution of phytoplanton blooms.Thus, a recent study based on mesocosm experiments carried out in summer in the NW Iberian Peninsula upwelling system, showed important daily changes in the transcription of glycoside hydrolases, peptidases, and transporters, supporting the hypothesis that resource partitioning is affected by temporal changes in the availability of dissolved organic matter (DOM) (Pontiller et al., 2022).Similarly, recurrent patterns at the functional level, particularly regarding substrate-induced responses, have been observed throughout bacterioplankton succession during phytoplankton blooms in the southern North Sea (Teeling et al., 2012;Teeling et al., 2016) and in the Monterey upwelling system (Nowinski and Moran, 2021).Although geographically distant, there are interesting similarities in the involvement of bacterial taxa in these studies, with pronounced contributions of Flavobacteriaceae (Bacteroidota; recognized for their polymer degrading capacities) along with Alphaproteobacteria (e.g.Roseobacter) and Gammaproteobacteria (e.g.Alteromonadales and Cellvibrionales) that differ in preference for different quantity and quality of resources.Nevertheless, and surprisingly, little is known about the temporal shifts in the functional profiles of different bacteria associated with changes in nutrient availability, hydrodynamic conditions, or phytoplankton community composition and activity along the coast to offshore spatial gradient in upwelling systems.
The environmental conditions for life in upwelling areas vary substantially along both spatial and temporal dimensions (Kämpf and Chapman, 2016).A key axis in the spatial structuring of upwelling ecosystems is the transition from the coastal zone, where the upwelling occurs, to offshore areas, where surface waters are displaced over time.This axis is strongly influenced by different coastal features like interactions between winds and the coastal topography (Figueiras et al., 2002), potential riverine inputs (Teixeira et al., 2018), and tidal forcing (Souto et al., 2001), generating spatial heterogeneity at various scales.In this regard, Satinsky et al. (2017) suggested that changes in microbial element cycling along the coastal system of the Amazon River Plume are largely associated with shifts in specific metabolisms.Accordingly, important spatial and temporal variability in the structure and function of microbial communities in large and heterogeneous upwelling areas can be expected (Zdanowski and Figueiras, 1997;Arıśtegui et al., 2009;Cury et al., 2011;Bertrand and Allen, 2012;Joglar et al., 2021a).
Short-term variability in wind regimes superimposed on seasonality induces patchiness in the prevalence of upwelling and downwelling events (typically 2-5 days (Pisareva et al., 2019)), particularly in the coastal zone (Blanton et al., 1984).Measurements of microbial bulk abundance and activity (e.g., production, respiration, enzymatic activities) show that such short-term temporal dynamics are associated with phytoplankton blooms and affect microbial community biomass, activity, and composition (Souto et al., 2001).The spatial and temporal scales of variability of phytoplankton growth and microbial food web functioning in upwelling systems are, therefore, inherently complex (Kerkhof et al., 1999;Barbosa et al., 2001;Bergen et al., 2015;Wear et al., 2015;Hernańdez-Ruiz et al., 2018), and constraining this complexity remains an important challenge.
The coastal area of the NW Iberian Peninsula is a highly dynamic and productive ecosystem, characterized by a light-limited winter and a nutrient-limited summer period, and by diatom and dinoflagellate blooms during spring and autumn, respectively (Figueiras and Rıós, 1993).On top of these seasonal patterns, this region is affected by the intermittent upwelling of cold and inorganic nutrient-rich Eastern North Atlantic Central Water (Nogueira et al., 1997).Although upwelling-favorable northerly winds prevail from March to September and downwelling-favorable southerly winds the rest of the year, out-of-season upwelling or downwelling events have been frequently recorded.During downwelling periods, surface, nutrientpoor seawater moves from the ocean to the coast.The extensive knowledge of the physical, chemical, and biological oceanography in the area provides a solid background against which to assess aspects of the microbial ecology (Souto et al., 2001;Figueiras et al., 2002;Teira et al., 2015;Teixeira et al., 2018;Joglar et al., 2021a).Thus, we approached this ecosystem with the objective of refining our understanding of microbial and biogeochemical processes in upwelling areas.In the current work we aimed to capture the range of variability in frequencies and amplitudes of key microbial functions.Our comparative metatranscriptome analysis comprised a collection of 162 genes previously used as informative markers for nutrient (e.g., C, N, P, S, Fe) fluxes and microbial metabolism (e.g., response to stimuli and stress, phototrophy, and synthesis of vitamins) (Satinsky et al., 2014;Robidart et al., 2019) and covered different oceanographic and biological settings.We hypothesized that short-term oceanographic features as wells as differences in upwelling intensity between coastal and offshore areas will have distinct effects on the diverse set of key microbial functions studied.

Study site and sampling
Three 10-day cruises were conducted on board B/O Ramoń Margalef in February (winter), April (spring), and August (summer) 2016 to cover a wide range of hydrographic and ecological conditions.Samples were taken at two locations in the eastern Atlantic Ocean, in the upwelling system near the Rıá de Vigo: the Coastal station (88 m depth, 42.14°N, 8.95°W, sampled between 8:00 and 8:20 am) and the Offshore station (260 m depth, 42.10°N, 9.60°W, sample between 6:30 and 6:45pm) (Figure 1A).It is important to note that the different time of sampling may have biased the metatranscriptome comparisons between both sampling stations, as it has been previously shown that prokaryotic transcriptomes from natural samples may exhibit diel periodicity, particularly those of Cyanobacteria in open ocean waters (Ottesen et al., 2013;Ottesen et al., 2014).
Detailed information on the experimental setup and specific methods below are given in Supplementary Material.

Hydrographic survey
Vertical profiles of temperature (°C), salinity, turbidity (NTU), total chlorophyll fluorescence, and photosynthetically active radiation (PAR) were obtained using a Seabird CTD rosette.Inorganic-nutrient determinations of collected water were done with a Bran + Luebbe segmented flow analyzer (Hansen and Grasshoff, 1983).

Chlorophyll a and prokaryotic abundance
Chlorophyll a (Chl-a) concentration was measured by the nonacidification technique as a phytoplankton biomass proxy (Welschmeyer, 1994).Fluorescence was determined with a TD-700 Turner Designs fluorometer (absorption coefficient as 87.7 at 663 nm) (Lorenzen and Newton Downs, 1986).Prokaryotic abundance in seawater samples was calculated as described by Gasol and Del Giorgio (2000) using a Becton Dickinson FACSCalibur flow cytometer (488 nm light produced by argon-ion laser).

Microbial community composition
Community composition was assessed from 24 samples by sequencing the 16S rRNA gene for prokaryotes and the 18S rRNA gene for eukaryotes (0.22 -3μm size-fraction and > 3μm size-fraction, respectively).DNA was amplified using the universal -and NH₄⁺divided by PO₄³⁻), and silicate (SiO 2 ) at surface (5m) during winter (blue), spring (red) and summer (green) at surface (5m) at days 1, 3, 5 and 7 at the offshore and coastal stations.
primers 515F and 926R for prokaryotes (Parada et al., 2016) and TAReuk454FWD1 and TAReukREV3 for eukaryotes (Logares et al., 2014).Amplified regions were sequenced in an Illumina Miseq platform.The sequencing depth for amplicon 16S and 18S was 10,000 reads per sample The sequences obtained were analyzed with the software DADA2 for amplicon sequence variants (ASVs) (Callahan et al., 2016) using the SILVA reference database (v138.1)for taxonomic assignment of 16S (Quast et al., 2013), and the databases PR2 (Guillou et al., 2012) and the marine protist from the BioMarKs project (Massana et al., 2015) for taxonomic assignment of 18S ASVs.

Metatranscriptomic analyses: prokaryotic community gene expression
RNA samples (< 3μm size-fraction) were taken during the experimental period at the surface in the coastal and oceanic stations (a total of 24 samples).Sampling for gene expression analyses was carried out on days 1, 3, 5, and 7 of each cruise.mRNA samples were extracted and sequenced using Illumina Miseq sequencing technology at the Science for Life Laboratory in Stockholm (SciLifeLab; www.scilifelab.se).The sequencing depth was 70,000 sequences.Illumina HiSeq 2500 raw paired-end reads (2x125bp) were quality checked with FastQC and MultiQC (Ewels et al., 2016), primer sequences were removed with cutadapt (Marcel, 2011), reads were trimmed with Sickle (Joshi and Fass, 2011) and filtered with ERNE (Del Fabbro et al., 2013) against a "contamination reference" database to filter reads designated as artifacts.Subsequently, forward and reverse reads were merged with PEAR (Zhang et al., 2014).All Merged reads were aligned with DIAMOND (Buchfink et al., 2014) against the NCBI RefSeq protein database (O'Leary et al., 2016).Afterward, all the sequences were compared to the 162 protein families (Pfams and TIGRfams) by running MEGAN (Huson et al., 2016) to annotate taxonomy and all possible transcripts to select marker genes relevant to a variety of ecological processes.The 162 genes found out of the 198 target protein profiles were grouped into nine categories and 44 subcategories (Table S1).To make comparisons between samples, considering that the transcripts have different lengths, and the libraries differ in size, we used a Transcript Per Million (TPM) normalization.

Statistical analyses
Non-parametric Mann-Whitney U Test analyses were performed to compare the salinity, nutrient concentrations, bacterial abundance and chlorophyll a concentration during the different seasons at both stations.T-tests analyses were performed to compare the contribution of each taxonomical group to the 16S dataset and to the metatranscriptome.To test the null hypothesis that data came from a normally distributed population, a Kolmogorov-Smirnov & Lillierors test was used.
Analysis of Non-parametric MultiDimensional Scaling (nMDS) and PERmutational Multivariate ANalysis Of VAriance (PERMANOVA) were used to test the grouping of metatranscriptomes by seasons and stations.ReDundancy Analyses (RDA) were applied to show the variability in the metatranscriptomes that is explained by environmental data and by the composition of phytoplankton and bacteria communities and PERMANOVA was used to test the significance of explanatory variables.
A Benjamini-Hochberg correction was applied to the Spearman rank correlation analyses.

Physicochemical and biological conditions
The results of the hydrographic surveys showed that downwelling conditions prevailed during winter and spring (Figures 1B; S1-S2).During this period, the important vertical mixing and reduced light availability prevented phytoplankton blooms and relatively high concentrations of inorganic nutrients were measured in surface waters (NO 3 -, NO 2 -and SiO 2 concentrations and the N:P ratio during winter and spring at both stations were significantly higher than those measured in summer, Mann-Whitney U Test, p-value < 0.05; Table S2; Figure S3).At the coastal station, a persistent surface halocline associated with rain and terrestrial runoff was observed during winter and spring (during winter and spring surface salinity was significantly lower than during summer, Mann-Whitney U Test, p-value < 0.05) (Figures S1, S3; Table S2).By contrast, during summer, the upwelling brought deep, salty, and nutrient-rich water to the surface.In summer, a massive phytoplankton bloom was observed at the coast and spread to the offshore station: at the surface significantly higher Chl-a concentrations accompanied by significantly higher bacterial abundance were measured compared to spring and winter (Mann-Whitney U Test, p-value < 0.05) (see satellite images, Figure S4).The observed very low nutrients concentrations during summer in surface waters at the coastal station are probably due to the massive consumption of nutrients by the phytoplankton bloom (Figures 1B; S1-S3; Table S2).At the coastal station, pronounced changes in Chl-a concentration were observed at the end of the winter cruise (Figures 1B; S1, S5) and even more so, at the end of the summer cruise (likely due to the rapid wind relaxation that interrupted the upwelling and disrupted the bloom (Figure S1), see Joglar et al., 2020(Joglar et al., 2020).
Collectively, there was more daily variation in physicochemical and biological oceanographic conditions (e.g.salinity, nutrients, chlorophyll a and bacterial biomass) in the coastal compared to the offshore and in summer compared to winter (Figures S1-S3, S5; Table S2).An interesting exception to this was temperature, which registered a similar large daily variation in the offshore and coastal stations in summer due to the intense upwelling of cold water that precisely promoted the changes in the rest of variables.
The eukaryotic community at both stations was dominated by dinoflagellates and diatoms (Dinophyceae, Thalassiosira, and Chaetoceros, particularly in summer) (Figure S6A) (Joglar et al., 2020).Among the bacteria, Pelagibacterales and Synechococcales showed higher relative proportions at the offshore station, particularly in winter, while other Alphaproteobacteria and Flavobacteria reached higher relative abundances in spring and summer on the coast.At both stations, Oceanospirillales were relatively more abundant in winter, whereas Cellvibrionales increased in summer (Figure S6B) (Joglar et al., 2020).

Relationships between metatranscriptomes and environment variables
Our metatranscriptomics analysis showed a dominance by Alphaproteobacteria (mainly Pelagibacterales and Rhodobacterales), Gammaproteobacteria (including Alteromonadales and Cellvibrionales), and Bacteroidota (mainly Flavobacteriales) in both stations and in the three sampled seasons (Figure 2A).The taxonomic distribution of transcripts for functional genes was similar in winter and spring, except for higher levels of Cyanobacteria at the offshore station in winter.During the summer phytoplankton bloom at the coastal station there was an important increase in the relative transcription levels of Cellvibrionales (Figure 2A).A relatively higher representation in the metatranscriptomics than in the DNA-based 16S rRNA gene amplicon data was observed for both Pelagibacterales (offshore station) and Other Gammaproteobacteria and Cellvibrionales (coastal station) during the three different seasons (t-test, p-value < 0.05) (Figure S7).

Analysis of Non-parametric MultiDimensional Scaling (nMDS)
and PERmutational Multivariate ANalysis Of VAriance (PERMANOVA) of bacterioplankton metatranscriptomes revealed a grouping of samples based on seasons and sampling stations (PERMANOVA, p-value < 0.05) (Figure S8).Redundancy Analysis (RDA) showed that functional profiles of winter and spring bacterial communities from both stations were significantly (PERMANOVA, p-value < 0.001) associated with high concentrations of specific nutrients (NO 3 -and SiO 2 ) and with fairly high relative abundances (estimated from DNA-based 16S rRNA gene amplicon data presented in Figure S6) of the eukaryotic phytoplankton Ciliophora, Cryptophyceae, and Dictyophiceae (PERMANOVA, p-value < 0.005) and the bacteria Oceanospirillales, Pelagibacterales, and Alphaproteobacteria (PERMANOVA, p-value < 0.001) (Figures 2B-D).During the summer bloom, the bacterial functional profiles showed higher variability and appeared to be linked to increased bacterial abundance, Chl-a concentration, PO³, and NH 4 + concentration (PERMANOVA, p-value < 0.001) as well as with Chaetoceros, Dinophyceae, and Fungi at both stations (PERMANOVA, p-value < 0.005) (Figures 2B-D).At the offshore station, bacterial functioning during summer was associated with high relative abundances of Amylibacter, Flavobacteriales, and Rhodobacterales (PERMANOVA, p-value < 0.001), while at the coastal station, it was

Spatial-temporal variability in expression of ecological marker genes
We next assessed the microbial functions related to nutrient fluxes (C, N, P, S, and Fe) and to vitamin synthesis, response to stress and stimuli, and phototrophy in this upwelling system.Different seasonal (Figures 3; S9-S17), daily (Figures 4; S18-S23) and spatial patterns were observed in the different functions studied.

Carbon, nitrogen and phosphorous metabolisms
In the carbon metabolism category, high relative transcription of carbohydrate metabolism and TCA cycle genes, particularly by Cellvibrionales and some Alphaproteobacteria, was observed in both stations as the seasons progressed, with higher values in summer during the phytoplankton bloom (Figures 3; S9, S18).Isocitrate dehydrogenase (icd) and succinate dehydrogenase (sdhA) were the TCA cycle genes with highest relative transcription (Figures S9,  S18).While icd, (which produces 2-oxoglutarate, also a precursor in amino acid metabolism, see below) was mainly transcribed by Gammaproteobacteria (Cellvibrionales), sdhA was transcribed by several taxa, including Bacteroidetes (Flavobacteriales) and Gammaproteobacteria at the coast and Alphaproteobacteria (Pelagibacterales) at the offshore (Figures S9, S18).
In the nitrogen metabolism category, ammonium uptake and amino acid metabolism were the subcategories with the highest proportion of transcripts, while urea utilization and nitrogen fixation remained low (Figures 3; S10, S18).There were pronounced seasonal changes in the relative expression of ammonium transporter (amt) genes (increasing from winter to summer, when low ammonium concentrations were measured during the phytoplankton bloom in the coast) and different taxa contributed to the transcription of amt-1 (mainly Cellvibrionales) and amt-2 (Pelagibacterales and Rhizobiales) (Figures S10, S18).The dominant genes in amino acid metabolism were also involved Relative transcript abundance (TPM) of selected key biological process: carbohydrate metabolism, TCA cycle, amino acid metabolism, ammonium uptake, phosphate metabolism, DMSP assimilation, siderophore metabolism, iron storage, vitamins, motility, photosynthesis, and proteorhodopsin.Color bars represent the relative contribution of main taxonomic groups to transcription at the costal and offshore station during winter (WIN), spring (SPG), and summer (SUM).The mean of four replicates (days 1, 3, 5 and 7) is shown.Note that different number of genes are included in each category (see supplementary material) and that different scales are used on the y-axis.
in nitrogen assimilation, such as the glnA gene, also importantly transcribed by Cellvibrionales and encoding glutamine synthetase that is the key enzyme that catalyzes the incorporation of ammonium or amino groups from amino acids into glutamine, and the global nitrogen-metabolism regulatory protein P-II (glnB) that modulates the expression of both glutamine synthetase and ammonium transporter genes (Figures S10, S18).
The phosphorus metabolism category was dominated by the expression of genes for phosphate metabolism, with higher relative values in spring and summer than in winter, especially at the coast (Figures 3; S11, S18).Pelagibacterales and Rhodobacterales (offshore and coastal station, respectively) and Gammaproteobacteria (particularly Cellvibrionales in summer at the coast) were the main groups accounting for this expression (Figures 3; S11, S18).Thus, phosphate import ATP-binding protein (pstB) was mainly transcribed by Pelagibacterales at the offshore station and by Rhodobacterales and Cellvibrionales at the coast and the transcription (mainly by Cellvibrionales) of the high affinity phosphate transporter pstS relatively increased in the coast in summer, while the low-affinity phosphate transporter (pit), which peaked in spring at the coast, was mainly transcribed by Gammaproteobacteria at both stations (Figures S11, S18).
Several of the C, N and P metabolism genes at both stations showed important daily changes in relative transcription levels over the course of the 7-day cruises (particularly during summer at both stations, Figures S21-S23) coinciding with short-term oceanographic variability (Figures 5; S5).As an example, individual phosphate acquisition genes showed different dynamics during the last part of the summer cruise at the coastal station (i.e., pstB increased while pstS decreased) and during winter and spring at the coastal station (i.e., pit, pstB, pstS) (Figures 4; S21), coinciding with changes in phosphate availability in the field (Figures 1; S5).Interestingly, daily changes (calculated as the ratio of the daily value divided by the mean value during the four days of the specific season and station) in PO 4 -3 concentration observed in the field were inversely correlated with the day-to-day changes in dmdA1, amt1, amt2, glnA, pstS and pstB (p < 0.05) (Figure 5).A similar (and counterintuitive) pattern was observed between the daily changes in these genes and those in the N/P ratio (Figure 5).It is important to note here that daily changes in the N/P ratio were mostly driven by changes in NO 3 -concentrations and not in PO 4 -3 availability (Figure 1B).On the other hand, daily changes in bacterial biomass were inversely correlated with the day-to-day changes in icd, zwf and amt2 and directly correlated with those in sdhA, pit, dmdA1 and fbpA (p < 0.05) (Figure 5).Finally, daily changes in temperature were directly correlated with the day-to-day changes in icd, zwf, amt1, amt2, glnA, glnB and pstS and inversely correlated with those in pit (p < 0.05) (Figure 5).

Sulfur, iron and vitamin metabolisms
Dimethylsulfoniopropionate (DMSP) assimilation, especially Pelagibacterales dmdA1 and dmdA2 transcription, was the category with the highest relative expression levels of all categories studied (Figures 3, 4; S12).A large peak in DMSP assimilation was observed in summer at the offshore station coinciding with a relative increase in the contribution of Dinophyceae to the eukaryotic community (Figure S6).Large short-term changes were observed in the transcription of DMSP assimilation genes during the final phase of the summer cruise at both the coastal and offshore stations (Figures 4; S19, S22).Relative expression of genes for sulfite and sulfate assimilation also reached fairly high levels, accounted for by Pelagibacterales and Gammaproteobacteria for sulfite (e.g., adenylsulfate reductase genes aprAB) and by a variety of taxa for sulfate (e.g., cysN) (Figure S12).Synechococcales and Cellvibrionales increased their relative contribution to cysN transcription in winter at the offshore station and during the summer bloom at the coastal station, respectively (Figure S12).
Vitamin synthesis marker genes selected were on one hand, thiL and thiC as marker genes of B1 synthesis metabolism and, on the other hand, cbiX (as marker gene of the anaerobic pathway) and, cobG, cobN and cobT (as marker genes of the aerobic pathway) for B12 synthesis metabolism.Vitamin synthesis genes generally showed higher relative expression levels in summer at the offshore station and during spring and summer at the coastal station (Figure 3).Vitamin synthesis genes were mainly related to Alphaproteobacteria (mostly Pelagibacterales at the offshore station), Gammaproteobacteria (particularly Cellvibrionales at the coastal station) and Synechococcales and Flavobacteriales at the offshore and at the coastal station, respectively (Figure 3).Specific genes related to vitamin synthesis presented different spatialtemporal patterns and were transcribed by different taxonomic groups (Figures 4; S14).Still, there were marked spatial-temporal patterns for different vitamins and bacterial taxa.For example, vitamin B1 gene thiL and B12 gene cobG were mainly transcribed by Cellvibrionales and Flavobacteria, with increased relative transcription in summer, particularly at the coastal station (Figures 4; S14, S19).Pelagibacterales instead dominated the transcription of vitamin B12 gene cobT, particularly at the offshore station during summer (Figures S14, S19).By contrast, vitamin B12 gene cbiX was transcribed by multiple taxa with no clear patterns (Figures S14, S19).Highest daily variability in the relative transcription of genes related to vitamins was observed during winter at the offshore station (Figure S22).

Motility, phototrophy and response to stress
Transcription of stress response genes was distributed among Pelagibacterales (offshore) and Cellvibrionales (coast) and increased at the coastal station in summer and at the offshore station in spring (Figures S15, S20).The expression of stimuli-related genes (e.g., motility genes like fliC) was mostly associated with Cellvibrionales, Alphaproteobacteria, and Verrucomicrobia.(Figures S16, S20).Highest daily variability in the relative transcription of genes related to response to stress and stimuli were observed during summer at both stations (Figure S22).
In the phototrophy category, the relative transcription of our selected photosynthesis genes (up to ~5500 TPM; Figure 3) was dominated by Synechococcales and was generally higher at the offshore station where it peaked in winter (Figures 3; S17, S20).This included genes for light harvesting antenna (e.g., cpcAB for phycocyanin), the photosystems (psa and psb) and carbon fixation (rubisco subunits).In contrast to the variability in photosynthesis genes, the proteorhodopsin (prd) gene expression was relatively stable and high (~4000) across all samples (Figure 3).The prd expression was generally dominated by Pelagibacterales with fair contributions also by other taxa -especially Cellvibrionales at the coastal station in summer (Figures 3; S17, S20).Bacteriochlorophyll gene (puf) expression was below 5 TPM (Figures S17, S20).Highest daily variability in the relative transcription of genes related to phototrophy was observed during winter at the offshore station (Figure S23).

Relationships between bacterial functioning and phytoplankton community composition
Since eukaryotic phytoplankton are key components of upwelling ecosystemspotentially shaping the conditions for bacterial activitieswe carried out Spearman rank correlation analyses (applying a Benjamini-Hochberg correction) between the relative abundance of the major eukaryotic taxa in the 18S rRNA gene amplicon data and the relative transcription levels of the Spearman rank correlation analysis between transcription of selected prokaryotic gene subcategories and (A) the relative abundance of the major eukaryotic taxa (18S) and; (B) environmental variables.The metatranscriptome data was clr-transformed before the analysis and a Benjamini-Hochberg method was used to adjust the p-values.Dots indicate p < 0.05.
studied prokaryotic functional categories (Figure 6A) and environmental variables 6B).Notably, bacterial temperature, and salinity displayed significantly positive correlations (p-value < 0.05 marked with asterisks in Figure 6) with categories related to bacterial uptake and utilization of resources (e.g., metabolism of carbohydrates, amino acids, phosphate and sulfate, and ammonium uptake) while the opposite pattern was observed for nutrients concentration (pvalue < 0.05 Figure 6B).Significant positive correlations were also found between the categories related to bacterial uptake and utilization of resources and the relative abundances of Stramenopiles (particularly Chaetoceros), Dinophyceae, marine stramenopiles (MAST), and Fungi (p-value < 0.05) (Figure 6A).A positive relationship was found between the transcription of DMSP assimilation genes by the prokaryotic community and the abundance of Dinophyceae (p-value = 0.052) and marine alveolate (MALV) (p-value < 0.05) (Figure 6A).In contrast, other metabolic traits, such as sulfite assimilation, iron storage, and photosynthesis genes, showed a significantly positive correlation (p-value < 0.05) with Chlorophyta (particularly Ostreococcus and Micromonas) and Cryptophyceae (Figure 6A) and with low chlorophyll a values (pvalue < 0.05) (Figure 6B).

Discussion
Coastal upwelling systems: dynamic grids of bacterioplankton functional specialization Our analyses provided novel information about the potential mechanisms governing microbial interactions, nutrient utilization, and energy and matter fluxes in a representative upwelling system.We reason that it is essential to identify and characterize the large diversity of microbial processes and patterns for two important reasons.First, it partly results from two-way interactions with the inherently diverse distribution patterns of phytoplankton and bacteria within the microbial communities.Second, because of the complex interdependency between different physicochemical conditions and the myriad of ecological processes driven by different microbes in the ocean.Thus, focusing on the transcriptional investment by any one taxon in any particular metabolic pathway or a specific season of choice, would unavoidably have resulted in underestimating the magnitude of dynamics and distribution patterns characteristic of coastal upwelling systems.In the following we will discuss these findings with the intention to understand potential linkages between physicochemical and biotic processes.

Resource availability as a powerful driver in structuring microbial functioning
The present work is, to our knowledge, the first studying high resolution spatial and temporal changes in bacterial functioning by means of metatranscriptomics in an upwelling system.The coastal area of the NW Iberian Peninsula is a highly dynamic ecosystem, affected by the intermittent upwelling of cold and inorganic nutrient-rich water that promotes phytoplankton and bacterial growth (Figueiras and Rıós, 1993;Nogueira et al., 1997;A ́lvarez-Salgado et al., 2002).Our dataset suggests that the coast to offshore gradient and the temporal changes in environmental and biological conditions promoted by the upwelling are key drivers of the observed pronounced taxonomic and functional bacterial diversification in this system.Several examples extracted from our results may be given.
The present dataset exemplifies some specific functional mechanisms behind the important role of nitrogen (N) and phosphorus (P) in structuring the functioning of microbial communities previously suggested to exist in this ecosystem (Martıńez-Garcıá et al., 2010;Teira et al., 2016).Regarding nitrogen, the observed high expression of the amt, glnAB, and icd genes by Gammaproteobacteria (particularly Cellvibrionales) during summer at the coastal station was remarkable.The enzymes encoded by these genes are linked in the process of nitrogen assimilation from ammonium in that amt encodes an ammonium transporter, glnAB encode glutamine synthase that synthesizes glutamine from ammonium and glutamate (glutamine can next be used for amino acid synthesis), and icd encodes the TCA cycle enzyme isocitrate dihydrogenase that produces 2-oxoglutarate, which is the precursor of glutamate.This may suggest that the competitive ability of Cellvibrionales during the phytoplankton bloom in part lies in their ability to gear their metabolism towards utilizing N from ammonium as a means to utilize the C that they process.This is in line with the emerging view that the dialogue between the metabolic N and C cycles is largely mediated at the intercept between the TCA-cycle and the synthesis of amino acids (Forchhammer et al., 2022).At the offshore station, on the other hand, Pelagibacterales showed a relatively higher transcription of ammonium transporters compared to other taxa, suggesting this group is more capable of assimilating ammonium at the low concentrations typically found in oceanic areas.We think analysis of adjustments of central metabolism between biosynthetic as compared to energy-yielding pathways (such as respiration) has a large potential to inform on how marine bacteria meet challenges in elemental stoichiometry of their resources (as imposed through e.g.C-or N-limitation).Especially if advances are made in characterizing gene systems like amino acid transporters and extracellular proteases involved in the utilization of dissolved organic nitrogen.
The high expression of the gene that encodes the transport system for orthophosphate (pstB) during summer, when phosphate concentration was low, is coherent with the suggested role of the Pst system (pstABCS) in inorganic phosphate (Pi) uptake only under Pi-limiting conditions (Harke and Gobler, 2013;Jimenez-Infante et al., 2017).On the other hand, the transcription of the pit gene, a low-affinity phosphate membrane transporter typically expressed under high phosphorus concentrations and post-bloom conditions (Alonso-Saéz et al., 2020), peaked at the coast in spring, coinciding with the highest ammonium and phosphate concentrations registered.Note here that, even though the phosphate concentration was similar in winter and spring at the coastal site, dissolved inorganic N, particularly ammonium, was considerably higher in suggesting that higher of the pit gene in spring than in winter might be related to the higher ammonium to phosphate ratio.Overall, the different spatial-temporal patterns of variability observed in genes with different phosphate affinity in this upwelling system, regardless of the relatively limited variability in phosphate ambient concentration, reinforces the hypothesis of a complex role of phosphate in structuring bacterial functioning (Harke and Gobler, 2013;Satinsky et al., 2014).
Our results support the idea of a taxonomic diversification in iron (Fe) provision strategies in this upwelling system.Thus, while Cellvibrionales and Flavobacteriales are suggested to mostly use siderophore-associated Fe, Pelagibacterales probably use ferric ions and Synechococcales accumulate Fe as ferritin in this ecosystem (Ahlgren et al., 2020).On the other hand, a widespread genetic capacity for the uptake and utilization of sulfate in surface waters off NW Iberian Peninsula can be deduced from our data, as opposed to previous knowledge limiting sulfate use during phytoplankton blooms to a few bacterial groups (Zhou et al., 2020).Interestingly, our results also highlight a probable dominant role of Pelagibacterales in sulfite assimilation in the system, extending to upwelling areas previous findings in offshore epipelagic layers (Needham et al., 2017).The present dataset suggests for the first time the relevant role of organic sulfur (S) compounds as sources of reduced S for marine prokaryotes in this upwelling system, confirming previous studies in other areas (Kiene et al., 2000;Moran et al., 2003).Interestingly, the high relative abundance of dmdA transcripts at the offshore station during the summer upwelling period coincided with a dominance of the phytoplankton community by Dinophyceae.Dinoflagellates are known to have high dimethylsulfoniopropionate (DMSP) production rates (Zhao et al., 2021), and previous work has suggested that their DMSP-producing ability could affect the structure of their associated bacterial community (Lin et al., 2021).In this regard, the observed dominant role of Pelagibacterales in the transcription of DMSP utilization genes expands the work by Varaljay et al. (2012) in the North Pacific subtropical gyre to this upwelling system.Furthermore, DMSP has been previously reported to have a photoprotective function (Archer et al., 2009), which may be related to the high radiation levels and the strong thermal stratification registered in the offshore station during summer.
Finally, the present dataset suggests that different strategies for using light as a source of energy by bacteria occur in this upwelling system.Thus, the relative transcription of photosynthesis related genes by Synechococcales was more important in the offshore station during winter, where the abundance of this group was relatively higher (Joglar et al., 2020).It is important to note here that, as explained in the Methods section, the coastal and offshore stations were sampled at different time points in the daily cycle.Thus, the high relative transcription of photosynthesis-related genes by Synechococcales during winter in the offshore station may be not only related to the broad differences in physicochemical and biological conditions (higher relative abundance of Synechococcales during winter in the offshore station, see Figure S6) but also partly to an increase in the investment in light capture due to the sampling time.In this regard, previous studies have shown that prokaryotic transcriptomes from natural samples may exhibit diel periodicity, particularly those of Cyanobacteria in open ocean waters (Ottesen et al., 2013;Ottesen et al., 2014).Interestingly, the data on the transcription of the prd gene suggest that photrophy mediated though proteorhodopsin was relevant and widespread among different taxa, including Pelagibacterales, in accordance with the results by Olson et al. (2018) in oligotrophic waters and Gammaproteobacteria during the coastal phytoplankton bloom.

Short-term dynamics in microbial functioning
Previous studies have shown recurring daily patterns in the functioning of heterotrophic bacteria from relatively stable oligotrophic environments (Ottesen et al., 2014;Muratore et al., 2022).It is also well known that the metabolism of naturally occurring bacteria experience inter-and intra-seasonal changes associated to phytoplankton blooms (Teeling et al., 2012;Teeling et al., 2016).However, this is to our knowledge the first time that systematic changes in specific microbial functions (identified from gene expression analyses) associated to short-term oceanographic variability are shown.Our results evidence important short-term dynamics in bacterioplankton functioning superimposed on seasonal variability that coincided with changes in hydrography, nutrient availability and biotic variables in this system.Interestingly, this day-to-day changes in transcripts for specific functions related to nutrient fluxes and response to stimuli were even more pronounced when daily variability in oceanographic conditions (e.g.temperature, phytoplankton biomass, resource availability) was highest, that is: close to the coast during summer.
In summer, a wind relaxation event interrupted the ongoing upwelling, which modified the vertical density gradient at the coastal station, altered the hydrographic conditions and nutrient availability, and disrupted the phytoplankton bloom (see (Barbosa et al., 2001) and (Joglar et al., 2021b)).Concomitantly, a sharp decrease in the transcription of functions related with bacterial utilization of resources derived from the phytoplankton bloom was observed, including amino acid metabolism (glnAB), TCA cycle (icd), and ammonium uptake (amt).This suggests that bacteria quickly responded to the changes in resource supply and potentially also in the C/N ratio of resources, as previously suggested by others (see review by Forchhammer et al. ( 2022)).Also, the relative importance of high-affinity phosphate uptake genes (e.g., pstS and pstB, which are considered an indicator of microbial P limitation and that varies seasonally in relation to P availability (Dyhrman et al., 2007)) decreased quickly after bloom disruption (on day 5), this probably due to the cease of the competition with phytoplankton.In contrast, the transcription of low-affinity phosphate uptake gene pit followed the opposite pattern and relatively increased at the end of the summer sampling at the coast when phosphorus availability increased.Similarly, bacteria shifted from a strategy based on iron storage to an increase in the relative importance of iron uptake (as suggested by the high relative transcription of ftnA and fbpAB, respectively) when the phytoplankton bloom was interrupted.Interestingly, the important short-term changes observed in the transcription of DMSP assimilation genes coincided with day-to-day peaks in the abundance of dinoflagellates the final phase of the summer cruise at both stations.Overall, these results suggested a close coupling between the physicochemical conditions and the bacterial and algal component of microbial communities from this upwelling system.

Unraveling functional interactions between bacteria and coastal phytoplankton
Overall, our results demonstrate that specific functional processes of the prokaryotic community identified from gene expression analyses) are associated with different phytoplankton groups in the upwelling system off the NW Iberian Peninsula.The presence of different phytoplankton and bacterioplankton populations during periods with different environmental conditions has been suggested in previous studies both in the area (Joglar et al., 2020;Costas-Selas et al., 2023) and in other marine ecosystems (Pinhassi et al., 2004;Needham et al., 2017;Bunse et al., 2019).Recent work in the NW Iberian Peninsula suggests a key role of biotic interactions as structuring factors of the eukaryotic community, mostly driven by positive associations between phytoplankton and bacteria (Costas-Selas et al., 2023).These findings may suggest some coevolution between the algae and the associated bacterial community that extends beyond the ability to transform phytoplankton-derived macromolecules by bacteria.The present dataset unravel the specific metabolic processes and bacterial and phytoplankton taxa that may be behind these associations through space and time.It is important to note here again that the interactions between phytoplankton and bacteria may change throughout the light cycle as it is known to affect the functioning of photosynthetic organisms (Ottesen et al., 2013;Ottesen et al., 2014).Therefore, the sampling at different times in both stations might influence to some extent the results discussed below.
Two interesting examples of possible phytoplankton-bacteria linkages observed during the productive season were associated to the bacterial transcription of genes related to the metabolism of DMSP (see above) and to the synthesis of B-vitamins.The observed increase in the relative abundance of transcripts of vitamin-related genes during summer is likely related to the high phytoplankton and bacterial requirements characteristic of blooming populations (Croft et al., 2005;Koch, 2012;Sañudo-Wilhelmy et al., 2012;Koch and Trimborn, 2019).Interestingly, genes belonging to the aerobic pathway of B12 synthesis (e.g.cobG and cobT), which is mostly associated with phytoplankton blooms (Zhou et al., 2020), increased during summer, while transcription of cbiX, which belongs to the anaerobic pathway of B12 synthesis that has been previously suggested to be independent of phytoplankton activity (Grossman, 2016), did not follow a clear seasonal pattern.The relative increase in thiL and cobG transcription by Cellvibrionales and Flavobacteria during summer suggests that these groups provided B1 and B12 vitamins to blooming eukaryotic phytoplankton.Interestingly, this hypothesis is also supported by microcosms experiments stimulating phytoplankton growth and including B1 and B12 vitamins additions carried out in summer in the NW Iberian Peninsula upwelling system that showed relatively important contributions of Cellvibrionales and Flavobacteria to the transcription of genes related to vitamin metabolism (Joglar et al., 2021a).
Finally, our results suggest strong functional linkages between bacteria and phytoplankton during non-blooming winter conditions.In this regard, the association between transcripts related to bacterial metabolic traits such as sulfite assimilation, iron uptake and storage, and photosynthesis, and Clorophyta and Cryptophyceae may reflect the close association between these phytoplankton groups and bacteria from the orders Pelagibacterales and Synechococcales, which played prominent roles in the relative transcription of genes associated with these metabolic functions.The close connection between these phytoplankton groups and Synechococcus sp., also described by Costas-Selas et al. (2023) in the same study area, were hypothesized to be linked with B12 provision.Our results point to an additional interaction mechanism between these groups, which could be related to the high capacity of Synechococcus to store Fe.As Clorophyta and Cryptophyceae may have mixotrophic behavior (Stoecker et al., 2017), predation on Synechococcus might be a source of Fe for them.A previous experimental study observed a sharp decrease in genes related to Fe uptake in a mixotrophic Haptophyta fed with bacteria, suggesting that bacterivory may provide Fe to the algae (Liu et al., 2015).
The close correlation between Chlorophyta (particularly Ostreococcus and Micromonas), Cryptophyceae, and Pelagibacterales suggests a link mediated by the sulfur metabolites, given the strong correlation between these phytoplankton groups and the relative abundance of sulfite assimilation genes.In this regard, a recent study suggests that adenylylsulfate reductase genes (aprAB) in Pelagibacterales clades could be related to taurine metabolism, oxidizing the sulfite derived from the conversion of taurine into acetyl-CoA (De Corte et al., 2021;Ruiz-Perez et al., 2021), a result consistent with the demonstrated ability of Ostreococcus and Micromonas to produce high quantities of taurine (Durham et al., 2019).

Cellvibronales: an unexpected key player in upwelling systems
Knowledge of the ecophysiology and ecology of Cellvibrionales is still scarce, but they appear to reach highest relative abundances in surface waters and coastal areas (Cheńard et al., 2019;Pajares et al., 2020;Reji et al., 2020).Although Cellvibrionales is not the most abundant bacterial group in this system, it appears to play a relevant and unexpected role during the summer phytoplankton bloom formation and development.
First, our field data suggest a prominent role of Cellvibrionales in degradation of phytoplankton-derived organic matter.Thus, during the Dinophyceae and Chaetoceros summer bloom, Flavobacteriales and Cellvibrionales dominated the relative expression of genes related to resource utilization and remineralization.It is well known that the primary role of Flavobacteriales in the degradation of polymeric organic matter (Teeling et al., 2012;Fernańdez-Goḿez et al., 2013) is important during diatom in upwelling ecosystems (Klindworth et al., 2014).However, the high transcription by Cellvibrionales of carbohydrate and TCA cycle-related genes as well as ammonium utilization and amino acid synthesis genes during the phytoplankton bloom presented here is remarkable.This behavior may be related to efficient energy production (Wang et al., 2022) and to adjustments to C or N-limitation by this group when resources are available during summer phytoplankton blooms (see above).Thus, the increase in the transcription of C metabolism-related genes by Cellvibrionales is consistent with recent findings in concurrent mesocosms experiments (performed in the same sampling area during the present study) (Pontiller et al., 2022) and in other regions (Poretsky et al., 2010;Liu et al., 2020), suggesting a quick response of Cellvibrionales to phytoplanktonderived organic matter.In surface waters of coastal systems like the NW Iberian Peninsula, the Gullmar Fjord and the German Bight, Cellvibrionales have been shown to importantly invest energy in the expression of different polymer degrading enzyme systems including for example laminarases that allow bacteria to decompose the algal storage glucan laminarin (Teeling et al., 2012;Pontiller et al., 2021;Pontiller et al., 2022) emphasizing their significant role in the turnover of phytoplankton-derived organic matter in these systems.Similarly, the relatively abundant Cellvibrionales transcription of motility-related genes during the summer bloom is in accordance with their description as a copiotroph group, known to quickly respond to the increase of DOM surrounding the microenvironment of phytoplankton cells in productive environments (Smriga et al., 2016).In addition, the high transcription of the proteorhodopsin gene (i.e.prd) by Cellvibrionales during summer may indicate that this group uses proton pumps as a supplemental energy source to uptake recently released complex carbon compounds during phytoplankton blooms (Pinhassi et al., 2016).
Our data suggest that the association between blooming phytoplankton and Cellvibrionales may be more complex than just the ability to transform phytoplankton derived macromolecules.We found higher levels of expression of genes for vitamins B1 and B12 synthesis by these bacteria associated with phytoplankton groups dominant during the blooms (i.e., Chaetoceros and Dinophyceae).This finding suggests that Cellvibrionales have an important role as suppliers of micronutrients needed for phytoplankton growth and bloom development.Similarly, the increase in the expression of siderophore-uptake related genes by Cellvibrionales associated with bloom-forming phytoplankton may indicate that these bacteria compete with Chaetoceros and Dinophyceae algae for organically-bound iron during the bloom season (Amin et al., 2009).Therefore, our data suggest that Cellvibrionalles develop a mutualistic relationship with blooming phytoplankton by producing growth-promoting compounds like vitamin B12 and facilitating phytoplankton functioning.They also show, as explained above, that Cellvibrionalles seem to act as competitors of blooming phytoplankton for iron provision in this productive upwelling system.
In conclusion, the present work highlights the mosaic of microbial processes governing energy and matter transformations in this representative upwelling system.Uncovering the linkages between bacterial activity and spatial-temporal variability in environmental conditions represents a relevant goal for microbial oceanography The systematic patterns that emerge from this work are encouraging in the context of future modeling attempts aiming to predict biogeochemical processes in such highly productive ecosystems.
FIGURE 2 (A) Contribution of the major taxonomic groups to the bacterial metatranscriptome at surface (5m) during winter (WIN), spring (SPG) and summer (SUM) at the coastal and offshore stations.The mean of four replicates (days 1, 3, 5 and 7) per station and season is shown.Redundancy analysis (RDA) of bacterioplankton metatranscriptome data constrained by: (B) environmental variables, (C) relative abundances of main eukaryotic groups (18S), and (D) prokaryotic groups (16S).The metatranscriptome data was clr-transformed to calculate Euclidean distances and to perform the RDA.A PERMANOVA analysis (p-values provided in each plot) was conducted to determine the significance of the analysis.Black arrows indicate the variables that significantly constrain each model.The percentages representing variability are displayed on the axes.Circles: coast, triangles: offshore.Winter: blue, spring: red and, summer: green.

FIGURE 4
FIGURE 4Daily variability of selected genes at the costal and offshore stations during winter (WIN), spring (SPG), and summer (SUM).Relative transcript abundance (TPM) of selected key biological process is provided.Color bars represent the relative contribution of main taxonomic groups to transcription.The value for each sampling day (days 1, 3, 5 and 7) is shown.Note that different scales are used on the y-axis.
FIGURE 5 (A) Daily variability (estimated as the ratio of the daily value of each variable divided by the mean value of that variable during the four days of the specific season and station) of phosphate (PO₄³⁻) and of pstB transcripts measured at the costal and offshore stations during winter (WIN), spring (SPG), and summer (SUM).The color and the vertical position of the bubble (y-axis) represent the value of the ratio (B) Correlation between values presented in panel (A), Spearman coefficient and p-value of the correlation are provided.(C) Spearman rank correlation analysis between daily variability of the different environmental variables and that of selected genes.A Benjamini-Hochberg method was used to adjust the p-values.Dots indicate p < 0.05.NP stands for N:P ratio.