Vulnerability of Arctic Ocean microbial eukaryotes to sea ice loss

doi:10.21203/rs.3.rs-4448593/v1

Download PDF

Article

Vulnerability of Arctic Ocean microbial eukaryotes to sea ice loss

https://doi.org/10.21203/rs.3.rs-4448593/v1

This work is licensed under a CC BY 4.0 License

You are reading this latest preprint version

The Arctic Ocean (AO) is changing at an unprecedented rate, with ongoing sea ice loss, warming and freshening impacting the extent and duration of primary productivity over summer months. Surface microbial eukaryotes are vulnerable to such changes, but basic knowledge of the spatial variability of surface communities is not well understood. Here, we sampled microbial eukaryotes in surface waters of the Beaufort Sea from four contrasting environments: the Canada Basin (open ocean), the Mackenzie Trough (river-influenced), the Nuvuk region (coastal) and the under-ice system of the Canada Basin. Microbial community structure and composition varied significantly among the systems, with the most phylogenetically diverse communities being found in the more coastal systems. Further analysis of environmental factors showed potential vulnerability to change in the most specialised community, which was found in the samples taken in water immediately beneath the sea ice, and where the community was distinguished by rare species. In the context of ongoing sea ice loss, specialised ice-associated microbial assemblages may transition towards more generalist assemblages, with implications for the eventual loss of biodiversity and associated ecosystem function in the Arctic Ocean.

Biological sciences/Ecology/Microbial ecology

Biological sciences/Ecology/Climate change ecology

Biological sciences/Ecology/Molecular ecology

Earth and environmental sciences/Ecology/Biooceanography/Microbial biooceanography

Environmental perturbations induced by anthropogenic activities are impacting the Arctic Ocean (AO) at a faster rate than other regions of the planet ^1,2. Mean air temperatures in the Arctic have risen at least twice as much as the global average in the past 20 years ^3,4. This warming is causing sea surface temperatures to rise, leading directly to the loss of older sea ice through increased melting ⁵. The melting of thicker ice is transitioning the AO from a system dominated by multi-year ice (MYI) to one comprising primarily first-year ice (FYI) ⁶, with both sea ice extent and thickness decreasing ^7–9. These changes have been particularly apparent in the Pacific-influenced Western AO ¹⁰, where there have been marked reductions in sea ice cover ¹¹ and increased freshwater inflow from rivers due to a faster hydrological cycle and permafrost thawing in the High North ^12–14. The growing volume of terrestrial sources of freshwater flowing into the AO contributes to the ongoing freshening and heightened salinity stratification of the surface AO ¹⁵, and leads to greater suspended sediment loads in surface waters ¹⁶. These ongoing environmental perturbations are causing a regime shift towards a warmer, wetter and more unstable climate over the AO ¹⁷, which is already impacting its biogeochemistry and ecosystem productivity ^18–20.

Eukaryotic phytoplankton drive AO primary productivity ^21,22, and there is a pressing need to better understand how changing environmental conditions alter the biogeography of microbial eukaryotic communities, in particular at the air-sea interface. However, the phylogenetic composition and structure of surface microbial eukaryotic assemblages occurring across the range of AO ecosystems, including under sea-ice systems, remains understudied. The Beaufort Sea stretches from the coast along northern Alaska (USA) to Nunavut (Canada) and over the Canada Basin to nearly 80° N ²³. Environmental change here has been especially rapid ^24,25, making the Beaufort region a vanguard for understanding change elsewhere in the Arctic ²⁶. The Beaufort Sea comprises the Mackenzie Canyon, Beaufort Shelf, Amundsen Gulf and the Canada Basin (Fig. 1a) ^27,28. The Beaufort Gyre, one of the major AO features and reservoir of most of the AO’s freshwater ¹⁵, dominates the Canada Basin.

Here, we examine the taxonomic composition and phylogenetic structure of microbial eukaryotic communities recovered from AO surface waters spanning four distinct regions of the 310000 km2 Beaufort Sea. To determine whether the surface AO systems (coastal, river-influenced, open ocean and under-ice) harbour distinct microbial eukaryote communities, we sequenced the V4 region of 18S rRNA from samples collected across the Mackenzie Trough, the Utqiagvik (Barrow) region, the Canada Basin and samples collectively classified as under-ice. Canada Basin stations were predominantly open-ocean and influenced by the Beaufort Gyre and included stations with non-compacted (drift) sea ice. Molecular phylogenetics approaches can differentiate between cryptic species or taxa that are difficult to distinguish by microscopy or pigment-based approaches and can also identify the large non-photosynthetic components of the microbial eukaryotic community. Previous molecular taxonomic surveys of microbial assemblages in the Beaufort Sea have detected shifts in community composition associated with environmental perturbations and varied environmental conditions; however, these earlier studies covered small geographic areas [29, 30]. Using RNA as a template and targeting 18S rRNA to characterise microbial communities favours detecting the metabolically active, living community members in samples [31, 32], we show that microbial eukaryote communities from contrasting regions were phylogenetically distinct, and the under-ice community had a specific phylogenetic signature.

Sampling area

Samples were collected in the Beaufort Sea from the icebreaker CCGS Louis S. St-Laurent as part of the Joint Ocean Ice Studies (JOIS) program between 12 and 29/09/2018. Open water surface waters (≤ 5 m depth) were sampled with the ship’s underway pump, while under-ice surface waters were sampled from holes drilled through the ice. 5 L seawater were filtered first through a 50 µm mesh pre-filter, and the material was collected on a 0.4 µM polycarbonate filter. Filters were stored at 80°C with PowerSoil isolation buffer (Mo Bio Laboratories, Carlsbad, CA USA). Twenty-four stations, including three ice camps on the northern section of the Beaufort Sea, were sampled from four different systems: the Mackenzie Trough (n = 4), the Utqiagvik region (formally Barrow, n = 8), the Canada Basin open water (n = 9) and under-ice stations (n = 3) (Fig. 1a). Station coordinates and GSHHS shoreline data ²⁹ were used to map sampling stations in ArcMap. Daily sea ice concentration data during the sampling period were retrieved from the National Snow and Ice Data Center (NSIDC) ³⁰ and used to calculate the average sea ice concentration across the entire sampling period.

Nucleic acid extraction, amplification and sequencing

Total RNA was extracted from the stored polycarbonate filters with the RNeasy Mini Kit (QIAGEN, Hilden, Germany), DNase treated and reverse transcribed to cDNA using a SuperScript kit (Invitrogen, Waltham, MA USA). PCR amplification of the V4 region of the 18S rRNA was carried out using the E572F/E1009R V4-18S primers ³¹. Library preparation and Illumina MiSeq sequencing were carried out at Dalhousie University’s Integrated Microbiome Resource (Halifax, Nova Scotia, Canada). Sequences are available at the NCBI Sequence Read Archive (SRA) under the accession number PRJNA931954.

Environmental data

Ancillary data collected using a conductivity-temperature-depth (CTD) cast from each station were obtained from Woods Hole Oceanographic Institution (data are available at www2.whoi.edu/site/beaufortgyre/data/ctd-and-geochemistry/). Water samples drawn from Niskin bottles (at the same depth sampled for microbial diversity) were analysed on board for salinity, oxygen and nutrients. Unfiltered nutrients (nitrate, silicate and phosphate) samples were analysed within 12 hours of collection onboard using a three-channel SEAL Analytical nutrient AutoAnalyser 3, following the methods described by the manufacturer. Ammonium samples were collected along the BL and MK shelf lines, and ammonium concentrations were measured following the Holmes et al. protocol A (0 to 3 µm) ³², modified for 40 mL sample volume. Samples of 40.5 (± 0.5) mL of seawater were collected into 50 mL acid-washed Falcon tubes, in duplicate, from the Niskin bottles at each station. For chlorophyll and phaeopigments, water was filtered under low pressure onto 0.7 µm pore size GF/F 25 mm filters (Whatman GF/F, GE Healthcare) on board. If water could not be filtered immediately, samples were kept in the dark until filtered within 4 hours. Filters were then folded inside a 90 mm grade 1 filter paper (Whatman, GE Healthcare), wrapped in aluminium foil and stored at -80°C for later analysis onshore. Onshore, the pigment samples were extracted in glass scintillation vials with 90% acetone for 24 hours in the dark at -20°C. Samples were removed from the freezer and placed in the dark one hour before being read using a Turner Field Fluorometer. The environmental metadata (Table S1) were converted to z scores, and the resulting standardised metadata table (Table S2) was used to compute a Euclidean distance matrix with the vegdist function in Vegan v2.6-2 ³³. Hierarchical clustering was then performed on the distance matrix and used to construct a dendrogram. A heatmap of the standardised metadata table and a dendrogram was constructed using the R ComplexHeatmap package v2.12.1 ³⁴. For the under-ice samples and two BL stations, where ship CTD casts were not done, the environmental data were retrieved from the stations in the immediate vicinity (ICE1: station CB17; ICE2: station CB16; ICE3: station CB9; BL5: station BL4 and BL7: station BL8).

Sequence processing

Sequencing adapters were trimmed using TrimGalore v0.6.5 ³⁵ and primers were specifically removed using Cutadapt v2.10 ³⁶. The reads were then pre-processed using DADA2 v1.14 ³⁷; read quality profiles were visualised before and after filtering and trimming [truncLen = c(280,265), trimLeft = c(18,20), maxEE = c(2,2) and truncQ = 2]. Error rates for forward and reverse reads were then assessed using DADA2. Amplicon sequence variants (ASVs) were identified, and forward and reverse reads were merged. Chimeric ASVs were removed, and the PR² reference database v4.12.0 ³⁸ was used for taxonomic assignments with the RDP classifier ³⁹ (output available at 10.5281/zenodo.7612768). ASVs classified as multicellular taxa (Metazoa, Streptophyta) were discarded from this study.

Phylogenetic community analyses

ASVs were aligned to the 18S Silva reference alignment (Release 102) ⁴⁰ using MAFFT v7 ⁴¹ with the G-INS-i iterative refinement method. The ASV maximum-likelihood phylogenetic tree was reconstructed using IQ-TREE v2.0 ⁴² with 1000 non-parametric bootstrap replicates and using the TIM2 + F + R10 model, which was evaluated as the best-fit model by ModelFinder ⁴³ and BIC criterion. Bray-Curtis, unweighted and weighted UniFrac distance matrices were computed with Phyloseq v1.40.0 ⁴⁴. The homogeneity of variance (homoscedasticity) between the four communities was assessed using the betadisper function in Vegan and an ANOVA test. A PERMANOVA test (Vegan R package; adonis function) was used to assess if the community composition differed significantly between Beaufort Sea habitats. The distance matrices were then ordinated and plotted as a non-metric multidimensional scaling (NMDS) plot using the ordinate and plot_ordination functions in Phyloseq. The R hclust function was used for the hierarchical clustering of the distance matrices. Additionally, samples were clustered using the UPGMA method in QIIME 2 2023.9 ⁴⁵ (q2-diversity plugin) at a rarefaction depth of 20000 sequences with 1000 iterations.

Species richness measures were calculated on the ASV table with samples rarefied to an even depth of 25000 sequences. The mean number of ASVs (species richness; SR_25k), mean Faith’s phylogenetic diversity (PD_25k) and mean Simpson index were calculated from 100 rarefaction iterations over each sample, using the rarefy_even_depth and estimate_richness functions from Phyloseq, and the Picante R package v1.8.2 ⁴⁶ function for PD. We used Picante to calculate the standardised effect sizes of phylogenetic diversity (SES_PD) on rarefied and non-rarefied ASV tables to identify if rarefaction affected the alpha diversity estimates. Standardised effect size of mean pairwise distances (SES_MPD) ⁴⁷ were calculated per cluster with 1000 iterations using the non-rarefied ASV table (converted to relative abundances). Relative abundances per sample were randomised 1000 times and used to calculate the mean NRI (-1 x SES_MPD) ⁴⁸ of each community cluster.

Non-negative matrix factorisation

We used a non-negative matrix factorisation (NMF) approach to partition the samples using the NMF R package ⁴⁹. NMF is a dimension reduction method that decomposes the ASV abundance matrix in the product of 2 matrices. The coefficient matrix describes the overall structure of the abundance matrix with a small number of descriptors (the number of descriptors is the rank of the NMF decomposition). The basis matrix gives the weights of each original descriptor (ASV) on the new descriptors. NMF has been used to describe the spatial structure of microbial communities’ compositions (microbial species as descriptors ⁵⁰). NMF analysis was first computed with rank values ranging from 3 to 7, 100 runs, and various algorithms (nsNMF, Brunet, KL). The optimal results were obtained for the nsNMF algorithm, random seed of the factorised matrices, and a rank value of 4. We performed the final analysis with 200 runs, rank of 4, random seed and nsNMF algorithm.

The ASV indices were calculated by combining two methods described in Jiang et al. ⁵¹ and Kim et al. ⁵²: 1) the Spearman correlation coefficient between an ASV (i) and a descriptor (k) r_i,k was calculated using the abundance profile of the ASV (original abundance matrix) and SC (coefficients matrix) along all the samples; 2) The multidimensional projection between an ASV and a descriptor was calculated as the cosine of the angle between the vectors represented by an ASV abundance in the sample space and the vector represented by descriptors in the sample space. The abundance profiles of ASVs and descriptors were first normalized, and we calculated the multidimensional projection as:

$$\:{Cos\theta\:}_{i,k}=\sum\:_{j=1}^{n}{a}_{i,j}\times\:{s}_{k,j}$$

where Cosq_i,k is the multidimensional projection between the ASV i and the descriptor k, n is the number of samples, a_i,j is the normalized abundance of the ASV i in the sample j, and s_k,j is the normalized abundance of the descriptor k in the sample j. We then used the basis matrix to calculate the score of each ASV:

$$\:Score\left(i\right)=1+\frac{1}{{\text{log}}_{2}\left(q\right)}\sum\:_{k=1}^{q}p\left(i,k\right){\text{log}}_{2}\left(p\left(i,k\right)\right)$$

where i is the ASV, q is the number of descriptors (4 in our study), k is the descriptor number, p(i,k) is the probability of finding the ASV i in the descriptor k. The final index of ASV i on descriptor k (I_j,k) was obtained by multiplying the Spearman correlation coefficient, cos theta and the ASV score:

$$\:{I}_{i,k}={\rho\:}_{i,k}\times\:{Cos\theta\:}_{i,k}\times\:Score\left(i\right)$$

To select the ASVs with the greatest index values for each descriptor, the distributions of ASV indices for each descriptor were plotted (Fig. S1) and means and standard deviations calculated. The ASVs with index values greater than three standard deviations above the mean of each descriptor were classified as signature ASVs for their corresponding descriptor.

Functional community analyses

In addition to taxonomic and phylogenetic analyses, the data were also examined from a functional perspective. Using an 18S functional association database ⁵³ to obtain information on chloroplast presence/absence for each sequence identification number in the PR² database, ASVs were assigned as either ‘chloroplast-containing’ (including mixotrophs) or ‘heterotroph’ taxa. We carried out beta-diversity analyses on both subsets of ASVs.

Rare biosphere analyses

We defined rare ASVs as those with a total abundance of less than 0.01% across all samples ⁵⁴ and used UniFrac (weighted and unweighted) and Bray-Curtis distance methods to analyse the beta-diversity patterns of rare taxa within the dataset.

Distance-decay of community similarity

A distance matrix of geographic distances between sample locations was calculated in R using the geodist function from the geodist R package ⁵⁵. Geographic distances between samples were then plotted against community dissimilarities (weighted UniFrac distances) and the decay.model function in the betapart R package ⁵⁶ was used to fit a model of distance-decay to the data. A multiple regression on distance matrices (MRM) was carried out with the ecodist R package ⁵⁷ to assess the strength of the relationship between geographical distance and community dissimilarity.

Validation with operational taxonomic units

An operational taxonomic unit (OTU) approach was also used to verify the results of the ASV-based approach. Briefly, the sequences from the DADA2 sequence table were clustered at 98% similarity (TreeLine function) following alignment (AlignSeqs function) and distance matrix generation (DistanceMatrix function) using the DECIPHER v2.24.0 R package [40]. The resulting OTU table was filtered to remove Metazoa and Streptophyta, and an OTU alignment and phylogeny were generated using MAFFT v7 and IQ-TREE v2.0. The filtered OTU table was used to compute a weighted UniFrac distance matrix, which was ordinated and plotted as an NMDS using Phyloseq v1.40.0. Additionally, samples were clustered using the UPGMA method in QIIME 2 2023.9 ⁴⁵ (q2-diversity plugin) at a rarefaction depth of 20000 sequences with 1000 iterations.

Environmental context

Samples were collected from surface waters at 24 stations (Fig. 1a), representative of four habitats found across the Beaufort Sea: 1) open ocean Canada Basin (CB); 2) adjacent under-ice communities in the Canada Basin (ICE); 3) and 4) and two separate onshore to offshore river-influenced sites. Although the ICE stations were geographically close to several of the CB stations, the samples were taken from the water directly underneath the surface of the sea ice, whilst all other samples were retrieved from the ice-free surface water at the top of the polar mixed layer using the ship’s underway pump. The Mackenzie (MK) samples were coastal and within an area influenced by the Mackenzie River plume, adjacent to the Mackenzie Trough. The Barrow (BL) stations, situated off the coast of Nuvuk (formerly Point Barrow, Alaska, USA), were also coastal and near the Barrow Canyon. The BL stations were the most westerly stations of this study and the most influenced by Pacific Waters from the Bering Strait. Both canyon and trough systems are influenced one after the other by the east-flowing Beaufort shelf break jet as well as asymmetric upwelling of the canyons ⁵⁸. Average sea ice concentrations across the Beaufort Sea varied, with much greater concentrations of ice present in the eastern regions and more open water to the west (Fig. 1a).

The standardised environmental data heatmap (Fig. 1b) separated the most coastal, Pacific-influenced sites within the BL region (stations BL1, BL2 and BL3). A feature of these stations was their apparent O₂ depletion compared to the other BL stations. These waters also had higher concentrations of silicate, phosphate and available nitrogen, especially ammonium, with BL1 having the highest nitrate levels of all the stations in this study. BL1, BL2 and BL3 also had high chlorophyll a and phaeopigment concentrations. The other BL stations were further offshore and closer to the Barrow Canyon upwelling, with lower temperatures, salinities, phosphate and nitrogen concentrations than the inshore BL stations. Samples from the Mackenzie Trough region clustered together with stations CB2, CB3 and CB27 from the open water CB region. These CB stations were closer to the shelf break and had lower salinities than the more offshore CB stations. Stations CB28b and MK4 had lower salinities and higher temperatures, silicate, CDOM and phaeopigment concentrations than the other stations in the river-influenced cluster. Both the MK and CB regions were low in nitrate, in keeping with surface nutrient depletion. The northernmost CB stations (CB10, CB12, CB15 and CB16) had the highest phosphate concentrations outside of the BL system. The CB and ICE stations in the Canada Basin region formed an open ocean cluster characterised by greater O₂ levels than the BL stations and most of the river-influenced stations, as well as the lowest temperatures found in any of the clusters.

Phylogenetic biogeography and diversity

Taxonomic assignment and filtering of 18S rRNA ASVs resulted in a final sequence table with 1449 microbial eukaryote ASVs from ~ 1.4 million reads. Approximately 21% of ASVs could not have their taxonomy fully assigned (to species level). Variance within the four communities was homogeneous (betadisper and ANOVA; F = 0.56, p = 0.65), but the communities differed significantly in composition (PERMANOVA; R² = 0.713, F = 16.57, p < 0.001).

We assessed phylogenetic beta diversity using weighted UniFrac distances. NMDS ordination and UPGMA clustering were used to identify four distinct microbial communities (Fig. 2). The microbial assemblages collected from open-ocean waters grouped together, forming a distinct CB cluster, while samples taken from waters directly beneath the sea ice grouped apart, forming the ICE cluster. Samples from the MK transect grouped together, as did the BL samples. The OTU-based sample clustering yielded near identical results to the ASV approach in the weighted UniFrac NMDS and UPGMA sample clustering (Fig. S2). There was a positive trend between weighted UniFrac and geographic distances of the samples, as seen on the distance-decay plot (Fig. S3; p < 0.001). Furthermore, a multiple regression on matrices (MRM) carried out on the phylogenetic and geographic distance matrices also showed a significant relationship between the two matrices (MRM; R² = 0.107, F = 97.88, p < 0.05).

We analysed the community partitioning using unweighted UniFrac and Bray-Curtis distances (Fig. S4). Except for the BL samples, which clustered together, the samples did not form distinct clusters in the NMDS plot of unweighted UniFrac distances (Fig. S4a). Instead, the MK stations clustered with CB and ICE stations. In the Bray-Curtis NMDS plot (Fig. S4b), the sample clusters were more clearly defined, with the MK stations forming a cluster and most BL stations clustering together. To determine which community cluster (Fig. 2) was the most phylogenetically diverse, we calculated Faith’s phylogenetic diversity (PD_25k; Fig. 3a). Mean PD_25k varied significantly between clusters (Fig. 3a; ANOVA, p < 0.0001, F = 129). The MK community had the highest phylogenetic diversity compared with the other three communities (Tukey test, p < 0.0001) and the highest mean species richness (SR_25k; Fig. S5a). The lowest SR_25k was seen in the ICE community. The greatest range in PD_25k was observed in the CB cluster, but there was no significant difference in mean PD between CB and BL clusters. SES_PD values for each cluster were comparable between the rarefied and non-rarefied ASV tables (Fig. S6). We also calculated the net relatedness index (NRI) for each of the four clusters of samples in the weighted UniFrac plot (Fig. 3b). NRI is a measure of phylogenetic relatedness of species within a community ⁴⁸; NRI reflects the degree to which species (or here, ASVs) in a community are more closely related to each other than would be expected by chance. Positive NRI values indicate phylogenetic clustering and suggest that abiotic environmental filtering is a dominant process in shaping community composition ⁵⁹. Negative NRI values indicate phylogenetic overdispersion where either competitive exclusion acting on evolutionarily conserved traits, or environmental filtering acting on evolutionarily convergent traits plays an important role ⁵⁹. Mean NRI varied significantly between clusters (Fig. 3a; ANOVA, p < 0.001, F = 9.86). Both BL and CB clusters had positive mean NRI values, whereas the ICE and MK clusters had negative mean NRI values. There were significant differences in mean NRI between the BL and MK clusters (Tukey test, p < 0.01) and between the CB and MK clusters (Tukey test, p < 0.05). The ICE cluster also had significantly lower NRI values than those from BL and CB clusters (Tukey test, p < 0.01 and p < 0.001, respectively).

Functional and rare biosphere community analyses

To determine if proximate trophic modes influenced community composition, each ASV was labelled as either chloroplast-containing (inclusive of mixotrophs) (n = 227) or heterotrophic (n = 1222), based on its taxonomic classification. Additional phylogenetic analyses were carried out on these two categories. When considering only chloroplast-containing taxa, ICE1 is situated apart from all other samples in the weighted UniFrac and Bray-Curtis NMDS ordinations (Fig. 4). In ordinations, the BL and MK samples also formed distinct clusters, whereas the CB samples formed clusters with ICE2 and ICE3 samples. All community clusters (except BL) had negative mean NRI values when examining exclusively chloroplast-containing taxa. The heterotrophic taxa from MK samples clustered together in both weighted UniFrac and Bray-Curtis NMDS ordinations. All three ICE samples clustered together in the heterotroph-only weighted UniFrac NMDS, albeit with some overlap with some CB samples. Both ICE1 and BL8 were distinct from their respective clusters in the Bray-Curtis NMDS. The ICE and MK heterotrophic communities had negative mean NRI values in contrast with the BL and CB communities, the latter of which had the highest mean NRI of the heterotroph-only communities.

To test the hypothesis that rare taxa significantly influenced phylogenetic community structure, we carried out the same phylogenetic analyses on the rare biosphere, defined here as the ASVs with < 0.01% relative abundance across all samples (n = 1203) ⁵⁴. In the NMDS ordination of the weighted UniFrac distances, the MK and ICE clusters formed distinct and exclusive clusters. The CB samples were split across two clusters, however. Four CB samples clustered together in their own group, whereas the other five clustered with the BL samples. In the Bray-Curtis NMDS ordination, the BL and MK samples formed distinct clusters, and the CB samples clustered with ICE2 and ICE3, leaving ICE1 separate from the other clusters. The BL community had the highest mean NRI and both MK and ICE communities had negative mean NRI values. The mean NRI for the CB community was above zero but the range of values among samples went from negative to positive.

Signature sequence variants

We then applied a non-negative matrix factorisation (NMF) analysis to identify key taxa (signature sequence variants based on Fig. S1) in each community cluster (Fig. 5). This method groups samples based on the contribution of SCs descriptors in each sample, and we found that there was a strong agreement between the clusters identified using weighted UniFrac and the descriptors from the NMF analysis, with the only exception being that station CB16 aligned more closely with the ICE community than with the CB community, but also retained affinity to the CB community.

The BL community had 44 signature ASVs, more than any other community (Table S3; Fig. S5). Among the signature ASVs of the BL community, 17 originated from photosynthetic taxa. Micromonas clade B3 had the third greatest ASV index value and was unique to the BL signature community. Of the 33 signature ASVs in the CB community, the majority were dinoflagellates and included species in the genera Gymnodinium and Prorocentrum, with many other unclassified dinoflagellate ASVs (Table S4). The taxon with the greatest ASV index value for the CB community was Mantoniella, a chlorophyte alga belonging to the Mamiellophyceae ⁶⁰. Many of the 42 MK signature ASVs were heterotrophic species such as ciliates and choanoflagellates, as well as nominally parasitic dinoflagellates in the Syndiniales (Table S5). There were four diatoms among the MK signature ASVs, three of which were Chaetoceros spp., as well as three prymnesiophyte species. Among the 34 ICE community signature ASVs, the only chlorophyte was classified as Carteria sp. (Chlamydomonadales). The signature diatoms of the ICE community were an unclassified pennate diatom, Naviculales sp. and Bacillaria paxillifer. Most of the signature ASVs in the ICE community belonged to heterotrophic or mixotrophic taxa, almost all dinoflagellates and ciliates (Table S6).

Our analysis suggests that microbial communities under ice are particularly susceptible to ongoing environmental shifts, underscoring the vulnerability of these ecosystems in the face of change. By leveraging amplicon sequence data, we inferred species or ecotype presence, absence, and relative abundance, which allowed us to elucidate patterns of biodiversity and biogeography across the surface Arctic Ocean. Through sampling surface waters throughout the Beaufort Sea, we aimed to dissect the microbial community structure. Utilising a high-resolution phylogenetic technique, we were able to discern subtle diversity differences across regions, effectively distinguishing between under-ice stations and open-water microbial communities. Our results highlight the critical link between biodiversity and ecosystem resilience, suggesting that areas with less diverse assemblages, such as those under ice, are at greater risk of disruption due to environmental perturbations.

Challenges in microbial community studies

Using amplicon sequence data has proven to be a fundamental tool for understanding ecological dynamics across various environments ^61–63, and is of particular importance in under sampled polar regions to provide initial perspectives on community diversity, which may subsequently be validated by qPCR/dPCR ⁶⁴. Although the use of rRNA as a marker gene could skew diversity estimates due to copy number variation between taxa, the effect of multiple rRNA gene copies may in fact contribute little to the combined bias incurred from each methodological step in metabarcoding studies ⁶⁵. rRNA gene copy number tends to increase with phytoplankton cell size; here we focused on cells < 50 µm, which in theory harbour fewer rRNA gene copies than larger cells ^66,67.

In this study, we primarily use ASV over OTUs. Although OTUs may reduce the likelihood of discovering false biological diversity by combining similar amplicon sequences into taxonomic units, they may simultaneously overlook the fine-scale variation that exists between microbial strains, which can be distinguished with modern high-quality sequencing methods ^37,68. Intragenomic rRNA polymorphisms may also inflate the diversity detected in a dataset if different gene copies are classified as separate ASVs/OTUs. However, using even a strict cut-off of genetic distance of 0.01 to cluster OTUs should not create artificial diversity because the mean genetic distance between multiple 18S rRNA copies in the same genome is less than 0.003, meaning that multiple copies from the same genome should cluster as a single OTU ⁶⁹. We used rRNA instead of rDNA metabarcoding to examine the diversity of predominantly living and metabolically active community members ^70–72. Studies using both rRNA and rDNA have found general agreement and differences between the two methods ^72–75; some of the observed differences may have been overlooked if only rDNA had been studied ⁷⁰. Moreover, rRNA sequencing has been suggested to be superior to rDNA approaches ⁷⁶, for instance, in its sensitivity to environmental conditions ^77,78. This is an important factor to note in the context of a changing AO. Errors in RNA transcription could cause overestimations of biodiversity as these would likely not be detected and corrected by the DADA2 algorithm. However, transcriptional errors remain uncommon. In eukaryotes, most rDNA genes are transcribed by RNA polymerase I; a study on transcription errors in yeast found the error rate of RNA polymerase I to be 4.3×10^− 6 per base pair ⁷⁹.

Despite these challenges, our results remain robust when using OTUs clustered at 98% similarity, with the same overall patterns of diversity observed. Future microbial community studies of the AO may opt to utilise rDNA and rRNA in parallel to minimise the impact of sequence artefacts from dead cells while also allowing researchers to compare ‘total’ and ‘active’ microbial communities in contrasting environments and seasons ^70,75,80. There is, however, still a need for standardised methods in microbial ecology, from field sampling to laboratory and bioinformatic protocols ⁸¹.

Contrasting marine systems

Previous studies in the Canada Basin area focused on the vertical (different depths in the water column) diversity of microbial eukaryotic communities that vary depending on the water mass of origin and are relatively easy to identify in this region of the Arctic ^82,83. In contrast, our study focused on microbial eukaryotic communities' horizontal diversity across the Beaufort Sea. We identified three geographic regions using abiotic factors: 1) Canada Basin; 2) Mackenzie trough and Beaufort shelf stations; and 3) the coastal Barrow Canyon. The surface microbial communities only partially followed these divisions and could be separated as Barrow Canyon, Mackenzie-influenced, open-water and under-ice Canada Basin microbial communities. These results highlight the variability of microbial communities across the surface AO, even within the short 18-day period in which these samples were taken. Whilst our snapshot precludes predicting how these communities may adjust to ongoing ocean surface warming and freshening over an annual cycle, the study indicates community complexity and diversity associated with areas of the Beaufort Sea. The proportion of sequences (~ 21%) in our dataset whose taxonomy could not be assigned to species level by the reference 18S database has been noted previously in the extensive pan-Arctic sampling effort by Tara Oceans ⁸⁴, and suggests a high proportion of AO eukaryotes remains to be characterised. This highlights the need for ongoing culture and sequence validation as recently reported for deep branching predatory flagellates in the world ocean, including in the Arctic ⁸⁵. Culture-independent methods such as single-cell omics also play a vital role in characterising the extent of microbial eukaryote diversity, since these methods are not reliant on available reference genomes or cultured isolates ^86–89.

High Chl a and phaeopigments indicated an area of potential high productivity and potential phytoplankton turnover in coastal waters from the Alaskan Beaufort Sea (stations BL1, BL2 and BL3) ⁹⁰. The relatively higher nutrient concentrations in the coastal BL stations indicate an influence of nutrient-rich Pacific water flowing into the Beaufort Sea through the Bering Strait and Chukchi Sea ^91,92. However, our sampling was somewhat late in the season, and the high Chl a concentrations at the surfaces suggest ongoing nutrient input and local surface warming in the stratified water column. The Alaskan Beaufort Sea, where the BL samples were taken, experiences periodic upwellings driven by winds, particularly in the late autumn ^93,94. These upwelling events bring denser and more nitrate-rich Atlantic water to the surface of the Nuvuk region ⁹⁵. As surface nutrients usually are depleted during the productive spring and summer months, such upwelling events may locally replenish nutrients to the surface waters of the Nuvuk region in the late summer, allowing primary productivity to continue later than in other regions. This is consistent with our data, where saltier, nutrient-rich waters were seen at the coastal BL stations compared to other stations. The relatively fresher, cooler and more oxygenated waters of the offshore BL stations (BL4-BL8) suggest a region of upwelling close to the Barrow Canyon. The BL microbial community had the greatest number of photosynthetic signature ASVs, which, along with the greater concentrations of Chl a and phaeopigments, indicates an area of higher primary productivity and grazing and degradation ⁹⁶. The BL community also had more signature ASVs than any other communities, reflecting a high level of uniqueness of the Barrow taxa to these samples.

Despite the BL and MK communities both being coastal systems influenced by rivers and upwelling, they were distinct from each other. Although higher temperatures were also observed at MK4 and CB28b (the Mackenzie River offshore stations), lower salinity levels accompanied by high CDOM concentrations were consistent with an offshore upwelling of warmer, fresh water from the Mackenzie River plume, as previously suggested ⁹⁷.

The MK community had the second highest number of signature ASVs, mostly from heterotrophic taxa, in contrast with the BL community, dominated by chloroplastidic groups. In addition, the MK community had higher mean PD_25k and SR_25k values than the BL community but had a negative NRI, indicative of phylogenetic overdispersion (evenness) ^48,98,99. This may help explain why this community was more diverse than any others. Phylogenetic overdispersion commonly occurs when competition regulates community membership, and ecologically similar taxa cannot coexist, leading to a phylogeny where species are more distantly related than would be expected by chance (competitive exclusion) ^59,100. However, overdispersion may also result from distantly related species having convergent traits allowing them to survive in their environment, in which case habitat filtering may be driving MK community structure ^99,101. Another explanation could be that input of freshwater species from the Mackenzie River and the development of a brackish water community ¹⁶ led to phylogenetic overdispersion, as previously noted in Mackenzie Shelf bacterial communities ⁸⁰.

By the same measure, the positive NRI of the BL community indicated phylogenetic clustering, where species are more closely related than expected by chance; given the distinctive environmental conditions measured in the BL samples, habitat filtering was the more likely driver of community structure here. The BL community was also the only community that was distinct from the others in all three beta-diversity measures used, implying that its structure differed both qualitatively and quantitatively from the other communities ¹⁰². This may have been due to the influences of Pacific water and localised upwelling, which provided an influx of nutrients that sustained a community characterised by phytoplankton.

Functional perspectives on community structure

To relate phylogenetic community structure to ecological function, we carried out beta-diversity analyses on subsets of taxa based on their proximate trophic status. Chloroplast-containing (including mixotroph) taxa played essential roles in structuring the BL and MK communities, as these samples formed distinctive clusters in both weighted UniFrac and Bray-Curtis NMDS plots. ICE1 clustered apart from the other ICE samples in the chloroplast-only plots, indicating that its community composition was both qualitatively and quantitatively distinct from the communities recovered from other ice-covered surface waters. Carteria sp., a chlorophyte alga, appeared only in ICE1 and is therefore likely to have contributed to differentiating this sample from others. Carteria has previously been found to be abundant in Arctic summer sea ice melt ponds, which can become connected to seawater ¹⁰³. A possible explanation for the presence of Carteria in only one of the ICE communities is that the under-ice habitat here may have been receiving meltwater from the sea ice surface containing Carteria, suggesting connectivity between upper and lower sea ice surfaces. Another difference in chloroplast-containing taxa between ICE1 and ICE2/ICE3 was the presence of haptophyte species such as Phaeocystis spp. and Chrysochromulina spp. in ICE2/ICE3, neither of which were detected in ICE1. Sarcinochrysis sp. and two other unclassified pelagophytes were also found in ICE1 but not in ICE2 or ICE3. One of the unclassified pelagophytes in ICE1 (ASV_1365) was subsequently found to have > 99% sequence identity with the Arctic pelagophyte Plocamiomonas psychrophila (CCMP2097) ¹⁰⁴, which is a euryhaline species likely adapted to living in ice ¹⁰⁵. These findings are indicative of a unique assemblage in the ICE1 sample. Arctic genotypes of Phaeocystis and Chrysochromulina have been previously recorded ¹⁰⁶, and P. pouchetii dominated phytoplankton communities during the 2007 record ice melt ¹⁰⁷. Moreover, P. pouchetii abundances rapidly increased during sea ice retreat in the Barents Sea ¹⁰⁸, suggesting that this taxon may be associated with sea ice melt.

The predominance of heterotrophic and parasitic signature ASVs in the MK and ICE communities may result from environmental factors limiting phytoplankton growth, such as decreased light levels and nutrient availability. In coastal MK and under-ice ICE stations, the surface waters had low photosynthetically active radiation (PAR) measures. The sea-ice cover limited the light intensity for the ICE stations, whereas the turbidity introduced by the Mackenzie River plume likely limited light levels. This reduction in light constrains the growth of photosynthetic organisms, resulting in an ecological niche where heterotrophic and parasitic lifestyles become advantageous ^82,109. Such conditions underscore the need to consider the physical, chemical, and biological interplays, including nutrient influx and light availability, in shaping microbial community compositions.

Influence of the rare biosphere

To investigate if low abundance taxa were driving the overall diversity patterns, we also conducted beta-diversity analyses on the rare biosphere. The MK samples formed clusters in the weighted UniFrac and Bray-Curtis NMDS plots, clustering less clearly in the unweighted UniFrac NMDS, suggesting that the diversity of this community was driven more strongly by quantitative (abundance) than by qualitative (presence/absence) differences in rare species. The BL samples formed clear clusters in the unweighted UniFrac and Bray-Curtis ordinations but formed a cluster with some of the CB samples in the weighted UniFrac ordination, suggesting that there were similarities in rare taxa and taxa abundances within the CB community. Rare taxa did not drive the structure of the CB community, as the CB samples did not form distinct and exclusive clusters in any of the rare biosphere NMDS ordinations. Both ICE and MK rare biospheres exhibited phylogenetic overdispersion, in contrast to the BL rare biosphere, which was phylogenetically clustered. These data support the notion that the rare biosphere was a strong driver of phylogenetic community structure in the ICE and MK communities but much less so in the BL and CB communities.

Microbial communities are often dominated by low-abundance taxa ¹¹⁰. Taxa may be permanently or transiently rare, switch between abundance and rarity, or always be rare with periodic changes in abundance ^110,111. Permanently rare taxa are associated with distinct life history strategies and trade-offs and may occur only in specific ecological niches ¹¹¹. Previous research on the AO rare biosphere suggested that it is driven primarily by stochastic processes, in contrast to abundant communities, which were more influenced by deterministic selection ⁵⁴. Given that the AO under-ice habitat fluctuates between seasons and years and may require specific physiological adaptations to survive ¹¹², the dominance of the rare biosphere here is ecologically coherent. The influence of freshwater may also have played a role in the dominance of the rare biosphere, since the lowest salinities were observed in the Mackenzie River-influenced region, and the under-ice habitat may also have been fresher than surrounding open ocean surface waters. In the Arctic summer, meltwater from the sea ice bottom can accumulate and form a layer of fresher water immediately beneath the sea ice ¹¹³, and riverine inputs also intensify due to increased meltwater and rainfall ¹¹⁴. Given that both MK and ICE communities are likely to be subjected to both stochastic and seasonal environmental fluctuations and increased freshwater input due to environmental change, these results hint at possible disruption of these AO rare biosphere communities under climate change ¹¹⁵.

The role of sea ice in structuring microbial communities

The loss of multi-year ice (MYI) and decreasing duration of the ice-covered period across the Arctic impacts protist communities associated with sea ice ^31,116,117. The ice-water interface beneath first-year ice (FYI) has been found to contain relatively more sea ice protists (but lower protist diversity) than beneath MYI ¹¹⁸. This could result from increased light penetration to surface waters covered by FYI, which is less thick than MYI. The Western Arctic Ocean is predicted to become dominated by FYI, with an expected increase in primary productivity coinciding with a decrease in ice-associated algal diversity ¹¹⁹.

In this study, we sampled surface water at the ice-water interface immediately beneath MYI (> 1.5m thick). The influence of ice was evident by species recovered from samples immediately beneath sea ice, including typical ice-associated species such as Carteria ¹²⁰ and Polarella glacialis ¹¹⁶. Our phylogenetic approach detected the distinct under-ice community, keeping with sea ice's importance in structuring microbial assemblages. This unique ice community was also apparent when looking at heterotroph-only and rare biosphere data subsets, suggesting that these two groups were drivers of the observed structure. Although diatoms are typically the major primary producers found in sea ice habitats ^118,121, our under-ice samples contained a greater proportional abundance of ciliates and dinoflagellates than diatoms. This could reflect the time of year the samples were taken; there is little known about MYI communities, and light and nutrient limitation at the end of summer could also affect diatoms more severely. Overall the low surface nutrients would have created unfavourable conditions for primary productivity. The apparent dominance of ciliates and dinoflagellates in the under-ice community was also reflected in the signature ASVs, where most species were ciliates or dinoflagellates. Previous research on ciliates within and beneath sea ice has suggested that both ciliate diversity and abundance are greater in the water column immediately under the ice than within the ice itself ¹²². Alveolates, including ciliates and dinoflagellates, were previously found to be more dominant in the Canada Basin than in the Bering Strait region, where stramenopiles dominated ¹²³. There were three diatom species among the under-ice signature ASVs, all pennate. Pennate diatoms are more commonly associated with sea ice than centric diatoms ^117,124. Although the relative abundance of diatoms at the ice-water interface in our samples was low, their classification as signature taxa of the samples implies a high level of uniqueness to this community. Six of the signature ASVs of the under-ice community were assigned as Spirotrichea ciliates, which are obligate mixotrophs ¹²⁵. All ciliates are heterotrophic; however, many species can acquire and retain their prey (or chloroplasts from their prey) and therefore become mixotrophs ¹²². Spirotrichea ciliates have been commonly found in the Arctic ¹²⁶, including in melt ponds associated with sea ice ¹²⁷, so their presence in our samples is not unexpected. The top signature ASV (with the highest ASV index value) came from an unclassified ciliate in the clade Phyllopharyngea; species belonging to this group are generally parasitic or prey on other ciliates ¹²⁸. The species found in the under-ice community may have been a predator feeding on the abundant ciliates present in these samples.

Our results highlight not only the spatial heterogeneity of surface AO protist communities but also the contribution of the sea ice habitat to shaping microbial assemblages. Phylogenetic diversity varied by community; the greatest diversity was in the MK community, whilst the lowest diversity, indicative of strong selection pressure, was seen within the ICE community. The greatest mean NRI values were in the BL and CB communities and were indicative of phylogenetic overdispersion and competition-driven community assemblages ^59,100. By contrast, the ICE community had the lowest mean NRI, signifying that phylogenetic clustering and habitat filtering had occurred ⁵⁹. The ICE community structure was driven by heterotrophs and rare species, which were more closely related than expected by chance due to their similar ecologies and traits, which allowed them to occupy this environment ¹⁰⁰. This emphasises the degree of ecological specialisation in ice-associated microbial communities and the vulnerability of these communities to ongoing losses in MYI. With the continued and accelerated rates of sea ice loss in the Arctic and the transition from a MYI to a FYI-dominated system, shifts in microbial community assemblages may occur due to changing selection pressures. Phylogenetically selective ice-associated communities composed of rarer species could gradually be replaced with more generalist assemblages and common taxa as the ice habitat is lost, or entirely different assemblages may form. If this happens, AO microbial community structure and function will be altered ⁶², with potentially widespread impacts on the wider AO ecosystem.

Author Contribution

A.M. and C.L. designed the study; A.M. supervised the study; T.G. and D.F.L.O. conducted the sampling; E.S.K. and V.L.N.J. performed the sequencing preparation; V.L.N.J., T.G. and A.M. analysed the data; V.L.N.J. prepared the figures and tables with input from T.G.; V.L.N.J., C.L. and A.M. wrote the manuscript with inputs from all authors.

Acknowledgement

We thank the captain and crew of the CCGS Louis S. St-Laurent and JOIS 2018 chief scientist, Dr Sarah Zimmerman, and JOIS lead scientist, Dr Bill Williams (Fisheries and Oceans Canada). The ancillary data were collected and made available by the Beaufort Gyre Exploration Program (BGEP) based at the Woods Hole Oceanographic Institution (https://www2.whoi.edu/site/beaufortgyre/) in collaboration with researchers from Fisheries and Oceans Canada at the Institute of Ocean Sciences. CL received support from the program Sentinel North financed by the Canada First Research Excellence Fund (CFREF) and the Natural Sciences and Engineering Research Council of Canada (NSERC) discovery and northern supplement grant; additional infrastructure was supported by the fond de Recherche du Quebec Nature et Technologies (FQRNT) (Québec Océan) and the Network of Centres of Excellence ArcticNet. AM is funded by a Royal Society University Research Fellowship (UK).

Data Availability

The sequence data generated and analysed during the study are available at the NCBI Sequence Read Archive (SRA) repository under the accession number PRJNA931954. The ASV sequences, ASV phylogeny, ASV table and taxonomic classification data are available from the Zenodo repository (doi:10.5821/zenodo.7612768). The OTU (98% clustering) table and phylogeny are also available at this DOI.

Ardyna, M., Gosselin, M., Michel, C., Poulin, M. & Tremblay, J. Environmental forcing of phytoplankton community structure and function in the Canadian High Arctic: Contrasting oligotrophic and eutrophic regions. Mar. Ecol. Prog. Ser. 442, 37–57 (2011).
Screen, J. A. & Simmonds, I. The central role of diminishing sea ice in recent Arctic temperature amplification. Nature 464, 1334–1337 (2010).
Dai, A., Luo, D., Song, M. & Liu, J. Arctic amplification is caused by sea-ice loss under increasing CO2. Nature Communications 10, 121 (2019).
Meredith, M. et al. Polar Regions. in IPCC Special Report on the Ocean and Cryosphere in a Changing Climate (IPCC, 2019).
Stroeve, J., Holland, M. M., Meier, W., Scambos, T. & Serreze, M. Arctic sea ice decline: Faster than forecast. Geophysical Research Letters 34, (2007).
Lange, B. A. et al. Comparing springtime ice-algal chlorophyll a and physical properties of multi-year and first-year sea ice from the Lincoln Sea. PLOS ONE 10, e0122418 (2015).
Stroeve, J. et al. Arctic sea ice extent plummets in 2007. Eos, Transactions American Geophysical Union 89, 13–14 (2008).
Kwok, R. & Untersteiner, N. The thinning of Arctic sea ice. Phys. Today 64, 36–41 (2011).
Kwok, R. Arctic sea ice thickness, volume, and multiyear ice coverage: Losses and coupled variability (1958–2018). Environ. Res. Lett. 13, 105005 (2018).
Proshutinsky, A. et al. Analysis of the Beaufort Gyre freshwater content in 2003–2018. Journal of Geophysical Research: Oceans 124, 9658–9689 (2019).
Fujiwara, A., Hirawake, T., Suzuki, K., Imai, I. & Saitoh, S.-I. Timing of sea ice retreat can alter phytoplankton community structure in the western Arctic Ocean. Biogeosciences 11, 1705–1716 (2014).
Déry, S. J., Hernández-Henríquez, M. A., Burford, J. E. & Wood, E. F. Observational evidence of an intensifying hydrological cycle in northern Canada. Geophysical Research Letters 36, L13402 (2009).
Ahmed, R., Prowse, T., Dibike, Y., Bonsal, B. & O’Neil, H. Recent trends in freshwater influx to the Arctic Ocean from four major Arctic-draining rivers. Water 12, 1189 (2020).
Bintanja, R. The impact of Arctic warming on increased rainfall. Sci Rep 8, 16001 (2018).
Giles, K. A., Laxon, S. W., Ridout, A. L., Wingham, D. J. & Bacon, S. Western Arctic Ocean freshwater storage increased by wind-driven spin-up of the Beaufort Gyre. Nature Geoscience 5, 194–197 (2012).
Ardyna, M. et al. Shelf-basin gradients shape ecological phytoplankton niches and community composition in the coastal Arctic Ocean (Beaufort Sea). Limnology and Oceanography 62, 2113–2132 (2017).
AMAP. Snow, Water, Ice and Permafrost in the Arctic (SWIPA) 2017. (2017).
Brussaard, C. P. D. et al. Arctic microbial community dynamics influenced by elevated CO2 levels. Biogeosciences 10, 719–731 (2013).
Brown, K. A., Holding, J. M. & Carmack, E. C. Understanding regional and seasonal variability is key to gaining a pan-Arctic perspective on Arctic Ocean freshening. Front. Mar. Sci. 7, (2020).
Vincent, W. F. Microbial ecosystem responses to rapid climate change in the Arctic. The ISME Journal 4, 1087–1090 (2010).
Hoppe, C. J. M., Wolf, K. K. E., Schuback, N., Tortell, P. D. & Rost, B. Compensation of ocean acidification effects in Arctic phytoplankton assemblages. Nature Climate Change 8, 529–533 (2018).
Arrigo, K. R. & van Dijken, G. L. Continued increases in Arctic Ocean primary production. Progress in Oceanography 136, 60–70 (2015).
Ayles, G. B. & Snow, N. B. Canadian Beaufort Sea 2000: The environmental and social setting. ARCTIC 55, 4–17 (2002).
Bourgain, P. & Gascard, J. C. The Atlantic and summer Pacific waters variability in the Arctic Ocean from 1997 to 2008. Geophysical Research Letters 39, L05603 (2012).
Huntington, H. P. et al. Evidence suggests potential transformation of the Pacific Arctic ecosystem is underway. Nature Climate Change 10, 342–348 (2020).
Carmack, E. C. et al. Freshwater and its role in the Arctic marine system: Sources, disposition, storage, export, and physical and biogeochemical consequences in the Arctic and global oceans. Journal of Geophysical Research: Biogeosciences 121, 675–717 (2016).
Carmack, E. & Macdonald, R. Oceanography of the Canadian Shelf of the Beaufort Sea: A setting for marine life. Arctic 55, 29–45 (2002).
Carmack, E., Macdonald, R. & Jasper, S. Phytoplankton productivity on the Canadian Shelf of the Beaufort Sea. Mar. Ecol. Prog. Ser. 277, 37–50 (2004).
Wessel, P. & Smith, W. H. F. A global, self-consistent, hierarchical, high-resolution shoreline database. Journal of Geophysical Research: Solid Earth 101, 8741–8743 (1996).
Maslanik, J. & Stroeve, J. Near-Real-Time DMSP SSMIS Daily Polar Gridded Sea Ice Concentrations, Version 1. (12/Sep/2018-29/Sep/2018). (NASA National Snow and Ice Data Center Distributed Active Archive Center, Boulder, Colorado, USA, 1999).
Comeau, A. M., Li, W. K. W., Tremblay, J.-É., Carmack, E. C. & Lovejoy, C. Arctic Ocean Microbial Community Structure before and after the 2007 Record Sea Ice Minimum. PLoS One 6, e27492 (2011).
Holmes, R. M., Aminot, A., Kérouel, R., Hooker, B. A. & Peterson, B. J. A simple and precise method for measuring ammonium in marine and freshwater ecosystems. Can. J. Fish. Aquat. Sci. 56, 1801–1808 (1999).
Oksanen, J. et al. vegan: Community Ecology Package. (2022).
Gu, Z., Eils, R. & Schlesner, M. Complex heatmaps reveal patterns and correlations in multidimensional genomic data. Bioinformatics 32, 2847–2849 (2016).
Krueger, F. Trim Galore. (2015).
Martin, M. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet.journal 17, 10–12 (2011).
Callahan, B. J. et al. DADA2: High-resolution sample inference from Illumina amplicon data. Nature Methods 13, 581–583 (2016).
Guillou, L. et al. The Protist Ribosomal Reference database (PR2): A catalog of unicellular eukaryote small sub-unit rRNA sequences with curated taxonomy. Nucleic Acids Research 41, D597–D604 (2013).
Wang, Q., Garrity, G. M., Tiedje, J. M. & Cole, J. R. Naïve Bayesian classifier for rapid assignment of rRNA sequences into the new bacterial taxonomy. Appl Environ Microbiol 73, 5261–5267 (2007).
Schloss, P. D. Silva reference files. Mothur https://mothur.org/wiki/silva_reference_files/#release-102 (2019).
Katoh, K., Rozewicki, J. & Yamada, K. D. MAFFT online service: multiple sequence alignment, interactive sequence choice and visualization. Briefings in Bioinformatics 20, 1160–1166 (2019).
Nguyen, L.-T., Schmidt, H. A., von Haeseler, A. & Minh, B. Q. IQ-TREE: A fast and effective stochastic algorithm for estimating maximum-likelihood phylogenies. Molecular Biology and Evolution 32, 268–274 (2015).
Kalyaanamoorthy, S., Minh, B. Q., Wong, T. K. F., von Haeseler, A. & Jermiin, L. S. ModelFinder: Fast model selection for accurate phylogenetic estimates. Nature Methods 14, 587–589 (2017).
McMurdie, P. J. & Holmes, S. phyloseq: An R package for reproducible interactive analysis and graphics of microbiome census data. PLoS ONE 8, e61217 (2013).
Bolyen, E. et al. Reproducible, interactive, scalable and extensible microbiome data science using QIIME 2. Nat Biotechnol 37, 852–857 (2019).
Kembel, S. W. et al. Picante: R tools for integrating phylogenies and ecology. Bioinformatics 26, 1463–1464 (2010).
Kembel, S. W. Disentangling niche and neutral influences on community assembly: Assessing the performance of community phylogenetic structure tests. Ecology Letters 12, 949–960 (2009).
Webb, C. O. Exploring the phylogenetic structure of ecological communities: An example for rain forest trees. The American Naturalist 156, 145–155 (2000).
Gaujoux, R. & Seoighe, C. A flexible R package for nonnegative matrix factorization. BMC Bioinformatics 11, 367 (2010).
Cai, Y., Gu, H. & Kenney, T. Learning microbial community structures with supervised and unsupervised non-negative matrix factorization. Microbiome 5, 110 (2017).
Jiang, X. et al. Functional biogeography of ocean microbes revealed through non-negative matrix factorization. PLOS ONE 7, e43866 (2012).
Kim, H. & Park, H. Sparse non-negative matrix factorizations via alternating non-negativity-constrained least squares for microarray data analysis. Bioinformatics 23, 1495–1502 (2007).
de Vargas, C. et al. Eukaryotic plankton diversity in the sunlit ocean. Science 348, 1261605–1261605 (2015).
Galand, P. E., Casamayor, E. O., Kirchman, D. L. & Lovejoy, C. Ecology of the rare microbial biosphere of the Arctic Ocean. PNAS 106, 22427–22432 (2009).
Padgham, M., Sumner, M. D. & Karney, C. F. F. geodist. (2021).
Baselga, A. & Orme, C. D. L. betapart: an R package for the study of beta diversity. Methods in Ecology and Evolution 3, 808–812 (2012).
Goslee, S. C. & Urban, D. L. The ecodist package for dissimilarity-based analysis of ecological data. Journal of Statistical Software 22, 1–19 (2007).
Lin, P. et al. On the nature of wind-forced upwelling and downwelling in Mackenzie Canyon, Beaufort Sea. Progress in Oceanography 198, 102674 (2021).
Webb, C. O., Ackerly, D. D., McPeek, M. A. & Donoghue, M. J. Phylogenies and community ecology. Annual Review of Ecology and Systematics 33, 475–505 (2002).
Yau, S. et al. Mantoniella beaufortii and Mantoniella baffinensis sp. nov. (Mamiellales, Mamiellophyceae), two new green algal species from the High Arctic. Journal of Phycology 56, 37–51 (2020).
Freyria, N. J., Joli, N. & Lovejoy, C. A decadal perspective on north water microbial eukaryotes as Arctic Ocean sentinels. Sci Rep 11, 8413 (2021).
Jacquemot, L., Vigneron, A., Tremblay, J.-É. & Lovejoy, C. Contrasting sea ice conditions shape microbial food webs in Hudson Bay (Canadian Arctic). ISME COMMUN. 2, 1–10 (2022).
Lin, Y. et al. Decline in plankton diversity and carbon flux with reduced sea ice extent along the Western Antarctic Peninsula. Nat Commun 12, 4948 (2021).
Gold, Z. et al. Signal and noise in metabarcoding data. PLOS ONE 18, e0285674 (2023).
McLaren, M. R., Willis, A. D. & Callahan, B. J. Consistent and correctable bias in metagenomic sequencing experiments. eLife 8, e46923 (2019).
Zhu, F., Massana, R., Not, F., Marie, D. & Vaulot, D. Mapping of picoeucaryotes in marine ecosystems with quantitative PCR of the 18S rRNA gene. FEMS Microbiology Ecology 52, 79–92 (2005).
Godhe, A. et al. Quantification of Diatom and Dinoflagellate Biomasses in Coastal Marine Seawater Samples by Real-Time PCR. Applied and Environmental Microbiology 74, 7174–7182 (2008).
Callahan, B. J., McMurdie, P. J. & Holmes, S. P. Exact sequence variants should replace operational taxonomic units in marker-gene data analysis. The ISME Journal 11, 2639–2643 (2017).
Gong, W. & Marchetti, A. Estimation of 18S gene copy number in marine eukaryotic plankton using a next-generation sequencing approach. Frontiers in Marine Science 6, (2019).
Adamo, M., Voyron, S., Chialva, M., Marmeisse, R. & Girlanda, M. Metabarcoding on both environmental DNA and RNA highlights differences between fungal communities sampled in different habitats. PLOS ONE 15, e0244682 (2020).
Stoeck, T., Zuendorf, A., Breiner, H.-W. & Behnke, A. A molecular approach to identify active microbes in environmental eukaryote clone libraries. Microb Ecol 53, 328–339 (2007).
Wutkowska, M., Vader, A., Mundra, S., Cooper, E. J. & Eidesen, P. B. Dead or Alive; or Does It Really Matter? Level of Congruency Between Trophic Modes in Total and Active Fungal Communities in High Arctic Soil. Front. Microbiol. 9, (2019).
Dowle, E., Pochon, X., Keeley, N. & Wood, S. A. Assessing the effects of salmon farming seabed enrichment using bacterial community diversity and high-throughput sequencing. FEMS Microbiology Ecology 91, fiv089 (2015).
Cox, F., Newsham, K. K. & Robinson, C. H. Endemic and cosmopolitan fungal taxa exhibit differential abundances in total and active communities of Antarctic soils. Environmental Microbiology 21, 1586–1596 (2019).
Guardiola, M. et al. Spatio-temporal monitoring of deep-sea communities using metabarcoding of sediment DNA and RNA. PeerJ 4, e2807 (2016).
Li, R. et al. Comparison of DNA-, PMA-, and RNA-based 16S rRNA Illumina sequencing for detection of live bacteria in water. Sci Rep 7, 5752 (2017).
Kong, H. et al. RNA outperforms DNA-based metabarcoding in assessing the diversity and response of microeukaryotes to environmental variables in the Arctic Ocean. Science of The Total Environment 876, 162608 (2023).
Laroche, O. et al. First evaluation of foraminiferal metabarcoding for monitoring environmental impact from an offshore oil drilling site. Marine Environmental Research 120, 225–235 (2016).
Gout, J.-F. et al. The landscape of transcription errors in eukaryotic cells. Science Advances 3, e1701484 (2017).
Ortega-Retuerta, E., Joux, F., Jeffrey, W. H. & Ghiglione, J. F. Spatial variability of particle-attached and free-living bacterial diversity in surface waters from the Mackenzie River to the Beaufort Sea (Canadian Arctic). Biogeosciences 10, 2747–2759 (2013).
Yates, M. C., Derry, A. M. & Cristescu, M. E. Environmental RNA: A Revolution in Ecological Resolution? Trends in Ecology & Evolution 36, 601–609 (2021).
Monier, A. et al. Oceanographic structure drives the assembly processes of microbial eukaryotic communities. ISME J 9, 990–1002 (2015).
Lovejoy, C., Galand, P. E. & Kirchman, D. L. Picoplankton diversity in the Arctic Ocean and surrounding seas. Mar Biodiv 41, 5–12 (2011).
Dorrell, R. G. et al. Convergent evolution and horizontal gene transfer in Arctic Ocean microalgae. Life Science Alliance 6, (2023).
Tikhonenkov, D. V. et al. Microbial predators form a new supergroup of eukaryotes. Nature 612, 714–719 (2022).
Mauger, S., Monard, C., Thion, C. & Vandenkoornhuyse, P. Contribution of single-cell omics to microbial ecology. Trends in Ecology & Evolution 37, 67–78 (2022).
Kolisko, M., Boscaro, V., Burki, F., Lynn, D. H. & Keeling, P. J. Single-cell transcriptomics for microbial eukaryotes. Curr Biol 24, R1081-2 (2014).
Wideman, J. G. et al. Unexpected mitochondrial genome diversity revealed by targeted single-cell genomics of heterotrophic flagellated protists. Nat Microbiol 5, 154–165 (2020).
Pachiadaki, M. G. et al. Charting the Complexity of the Marine Microbiome through Single-Cell Genomics. Cell 179, 1623–1635.e11 (2019).
Terrado, R., Lovejoy, C., Massana, R. & Vincent, W. F. Microbial food web responses to light and nutrients beneath the coastal Arctic Ocean sea ice during the winter–spring transition. Journal of Marine Systems 74, 964–977 (2008).
Carmack, E. & Wassmann, P. Food webs and physical–biological coupling on pan-Arctic shelves: Unifying concepts and comprehensive perspectives. Progress in Oceanography 71, 446–477 (2006).
Wang, D., Henrichs, S. M. & Guo, L. Distributions of nutrients, dissolved organic carbon and carbohydrates in the western Arctic Ocean. Continental Shelf Research 26, 1654–1667 (2006).
Pickart, R. S. et al. Long-term trends of upwelling and impacts on primary productivity in the Alaskan Beaufort Sea. Deep Sea Research Part I: Oceanographic Research Papers 79, 106–121 (2013).
Pisareva, M. N., Pickart, R. S., Lin, P., Fratantoni, P. S. & Weingartner, T. J. On the nature of wind-forced upwelling in Barrow Canyon. Deep Sea Research Part II: Topical Studies in Oceanography 162, 63–78 (2019).
Beaird, N. L., Shroyer, E. L., Juranek, L. W., Hales, B. & Goñi, M. A. Nutrient-rich gravity current formed by upwelling in Barrow Canyon: High-resolution observations. Journal of Geophysical Research: Oceans 125, e2020JC016160 (2020).
Dobrzyn, P., Tatur, A. & Keck, A. Photosynthetic pigments as indicators of phytoplankton development during spring and summer in Adventfjorden (Spitsbergen). Oceanology 49, 368–376 (2009).
Garneau, M.-È., Vincent, W. F., Alonso-Sáez, L., Gratton, Y. & Lovejoy, C. Prokaryotic community structure and heterotrophic production in a river-influenced coastal Arctic ecosystem. Aquatic Microbial Ecology 42, 27–40 (2006).
Horner-Devine, M. C. & Bohannan, B. J. M. Phylogenetic clustering and overdispersion in bacterial communities. Ecology 87, S100–S108 (2006).
Kraft, N. J. B., Cornwell, W. K., Webb, C. O. & Ackerly, D. D. Trait evolution, community assembly, and the phylogenetic structure of ecological communities. The American Naturalist 170, 271–283 (2007).
Cooper, N., Rodríguez, J. & Purvis, A. A common tendency for phylogenetic overdispersion in mammalian assemblages. Proc Biol Sci 275, 2031–2037 (2008).
Cavender-Bares, J., Ackerly, D. D., Baum, D. A. & Bazzaz, F. A. Phylogenetic overdispersion in Floridian oak communities. The American Naturalist 163, 823–843 (2004).
Lozupone, C., Lladser, M. E., Knights, D., Stombaugh, J. & Knight, R. UniFrac: An effective distance metric for microbial community comparison. ISME Journal 5, 169–172 (2011).
Kilias, E. S., Peeken, I. & Metfies, K. Insight into protist diversity in Arctic sea ice and melt-pond aggregate obtained by pyrosequencing. Polar Research 33, 23466 (2014).
Daugbjerg, N., Lara, C., Gai, F. F. & Lovejoy, C. (in press) Plocamiomonas psychrophila gen. et sp. nov. (Pelagophyceae, Heterokontophyta), an Arctic marine nanoflagellate characterized by microscopy, pigments and molecular phylogeny. European Journal of Phycology (2024).
Freyria, N. J. et al. Salinity tolerance mechanisms of an Arctic Pelagophyte using comparative transcriptomic and gene expression analysis. Commun Biol 5, 500 (2022).
Terrado, R., Scarcella, K., Thaler, M., Vincent, W. F. & Lovejoy, C. Small phytoplankton in Arctic seas: Vulnerability to climate change. Biodiversity 14, 2–18 (2013).
Lasternas, S. & Agustí, S. Phytoplankton community structure during the record Arctic ice-melting of summer 2007. Polar Biol 33, 1709–1717 (2010).
Wassmann, P., Ratkova, T. & Reigstad, M. The contribution of single and colonial cells of Phaeocystis pouchetii to spring and summer blooms in the north-eastern North Atlantic. Harmful Algae 4, 823–840 (2005).
Kellogg, C. T. E., McClelland, J. W., Dunton, K. H. & Crump, B. C. Strong seasonality in Arctic estuarine microbial food webs. Frontiers in Microbiology 10, 2628 (2019).
Lynch, M. D. J. & Neufeld, J. D. Ecology and exploration of the rare biosphere. Nat Rev Microbiol 13, 217–229 (2015).
Jia, X., Dini-Andreote, F. & Falcão Salles, J. Community assembly processes of the microbial rare biosphere. Trends in Microbiology 26, 738–747 (2018).
Young, J. N. & Schmidt, K. It’s what’s inside that matters: physiological adaptations of high-latitude marine microalgae to environmental change. New Phytologist 227, 1307–1318 (2020).
Salganik, E. et al. Temporal evolution of under-ice meltwater layers and false bottoms and their impact on summer Arctic sea ice mass balance. Elementa: Science of the Anthropocene 11, 00035 (2023).
Beel, C. R. et al. Emerging dominance of summer rainfall driving High Arctic terrestrial-aquatic connectivity. Nat Commun 12, 1448 (2021).
Ardyna, M. & Arrigo, K. R. Phytoplankton dynamics in a changing Arctic Ocean. Nature Climate Change 10, 892–903 (2020).
Comeau, A. M. et al. Protists in Arctic drift and land-fast sea ice. Journal of Phycology 49, 229–240 (2013).
Piwosz, K., Wiktor, J. M., Niemi, A., Tatarek, A. & Michel, C. Mesoscale distribution and functional diversity of picoeukaryotes in the first-year sea ice of the Canadian Arctic. The ISME Journal 7, 1461–1471 (2013).
Hop, H. et al. Changes in sea-ice protist diversity with declining sea ice in the Arctic Ocean from the 1980s to 2010s. Frontiers in Marine Science 7, 243 (2020).
Lannuzel, D. et al. The future of Arctic sea-ice biogeochemistry and ice-associated ecosystems. Nat. Clim. Chang. 10, 983–992 (2020).
Kilias, E., Kattner, G., Wolf, C., Frickenhaus, S. & Metfies, K. A molecular survey of protist diversity through the central Arctic Ocean. Polar Biol 37, 1271–1287 (2014).
Leu, E., Søreide, J. E., Hessen, D. O., Falk-Petersen, S. & Berge, J. Consequences of changing sea-ice cover for primary and secondary producers in the European Arctic shelf seas: Timing, quantity, and quality. Progress in Oceanography 90, 18–32 (2011).
Majaneva, M., Rintala, J.-M. & Blomster, J. Taxonomically and functionally distinct Ciliophora assemblages inhabiting Baltic sea ice. Microb Ecol (2021) doi:10.1007/s00248-021-01915-4.
Liu, Q., Zhao, Q., McMinn, A., Yang, E. J. & Jiang, Y. Planktonic microbial eukaryotes in polar surface waters: Recent advances in high-throughput sequencing. Mar Life Sci Technol 3, 94–102 (2021).
Barber, D. G. et al. Selected physical, biological and biogeochemical implications of a rapidly changing Arctic Marginal Ice Zone. Progress in Oceanography 139, 122–150 (2015).
Onda, D. F. L. et al. Seasonal and interannual changes in ciliate and dinoflagellate species assemblages in the Arctic Ocean (Amundsen Gulf, Beaufort Sea, Canada). Frontiers in Marine Science 4, 16 (2017).
Lovejoy, C., Massana, R. & Pedrós-Alió, C. Diversity and distribution of marine microbial eukaryotes in the Arctic Ocean and adjacent seas. Applied and Environmental Microbiology 72, 3085–3095 (2006).
Xu, D. et al. Contrasting community composition of active microbial eukaryotes in melt ponds and sea water of the Arctic Ocean revealed by high throughput sequencing. Frontiers in Microbiology 11, (2020).
Terrado, R. et al. Protist community composition during spring in an Arctic flaw lead polynya. Polar Biol 34, 1901–1914 (2011).

No competing interests reported.

Download PDF

Editorial decision: Revision requested
31 Jul, 2024
Reviews received at journal
29 Jul, 2024
Reviews received at journal
22 Jul, 2024
Reviewers agreed at journal
12 Jul, 2024
Reviewers agreed at journal
10 Jul, 2024
Reviewers invited by journal
09 Jul, 2024
Editor assigned by journal
09 Jul, 2024
Editor invited by journal
09 Jul, 2024
Submission checks completed at journal
09 Jul, 2024
First submitted to journal
20 May, 2024

You are reading this latest preprint version

Vulnerability of Arctic Ocean microbial eukaryotes to sea ice loss

Status:

Version 1

Abstract

Figures

Introduction

Methods

Sampling area

Nucleic acid extraction, amplification and sequencing

Environmental data

Sequence processing

Phylogenetic community analyses

Non-negative matrix factorisation

Functional community analyses

Rare biosphere analyses

Distance-decay of community similarity

Validation with operational taxonomic units

Results

Environmental context

Phylogenetic biogeography and diversity

Functional and rare biosphere community analyses

Signature sequence variants

Discussion

Challenges in microbial community studies

Contrasting marine systems

Functional perspectives on community structure

Influence of the rare biosphere

The role of sea ice in structuring microbial communities

Declarations

Author Contribution

Acknowledgement

Data Availability

References

Additional Declarations

Supplementary Files

Status:

Version 1