Seawater samples were collected in May 2010 from two ice-covered stations, Stations 1 and 2, both located north of Spitsbergen (Figure 1). Station 1 in this paper was sampled between St. 2 and St. 3 of Johnsen et al. (submitted). After having drifted out of the sea ice, as described by Johnsen et al. (submitted), the ship re-located further NE into the sea ice, where our Station 2 was sampled. Sampling depths were selected based on hydrography and fluorescence profiles aiming to sample potentially different habitats for microbial protists. Water samples from five depths (0, 10, 45, 65, and 200 m) were collected at Station 1 and from four depths (0, 15, 120, and 200 m) at Station 2. From each depth, seawater was sequentially filtered through 3 μm filters (Isopore membrane filter, Millipore, USA) onto 0.22 μm filters (Durapore membrane hydrophilic PVDF filter, Millipore, USA) using a peristaltic pump at 40 rpm. One L was filtered for DNA extraction and 2 L were filtered for RNA extraction. Immediately after filtration, the RNA filters were submerged in 1 mL of Trizol (Invitrogen), the tube was shaken vigorously and kept at room temperature for 5 min. Filters were snap-frozen in liquid nitrogen and kept at −80°C until extraction. To measure total and >3 μm Chl a biomass, 250 mL seawater was filtered in triplicate onto GFF filters (Whatman, England) and 3 μm Isopore membrane polycarbonate filters (Millipore, USA), respectively. The filters were extracted in methanol for 12–24 h at 4°C in the dark and analyzed using a 10-AU-005-CE fluorometer (Turner, USA). Seawater for nutrient analyses was collected in acid-washed and sample-rinsed bottles and kept frozen until analysis. Nitrate + nitrite (NO₃+NO₂), phosphate (PO₄), and silicate (Si(OH)₄) were analyzed in triplicate by standard seawater methods using a Flow Solution IV analyzer from O.I. Analytical, USA. Photosynthetically active radiation (PAR) from the surface down to 40 m was measured using a cabled, flat Li-Cor sensor (LI-1400). The temperature, salinity and density measurements taken from the onboard SeaBird CTD was plotted in a TS diagram, and the water masses were characterized according to Aagaard and Carmark (1989); Manley (1995), and Cokelet et al. (2008). The thickness of the sea ice and the snow depths on top of the sea ice were measured in 6–10 drilled holes per station.

DNA was extracted using a modified CTAB protocol according to Vader et al. (2014) and each filter was cut into two halves and extracted individually. After thawing the filters sampled for RNA, the Trizol was transferred to a new tube and the filter was washed with 1 mL of fresh Trizol. In both tubes, the cells were lysed by beating two times 1 min with 200 μm zirconium beads (Molecular grade, pre-filled tubes, OPS Diagnostics) in a bead beater at 1/22s. Subsequently total RNA was extracted according to the manufacturers' protocol, except that 10 μg glycogen (RNA grade, Fermentas) was added to each tube as a carrier during precipitation of the RNA. The RNA was dissolved in DEPC-treated water, and DNAse-treated with Turbo DNAse (Ambion) according to the manufacturers protocol. RNA was reverse transcribed using the Retro Script kit (Ambion) with random decamers at 44°C for 1 h, after an initial denaturation step of 3 min at 75°C.

The hypervariable V4 18S rDNA region was amplified from the rRNA gene and from cDNA from reverse transcribed rRNA (referred to DNA and RNA, respectively) using eukaryote-specific 454 primers (Comeau et al., 2011). Samples were amplified in triplicate 25 μl reactions containing 1xHF-buffer, 0.2 μM of each dNTP, 0.2 μM of each primer, 0.5 U Phusion Hot start HF II polymerase and 1 μl undiluted DNA or RNA. All products were cleaned using AMPure XP beads (Agencourt, Beckman Coulter) according to the manufacturers instruction using a 4:5 ratio of beads to PCR product. Subsequently all PCR products were pooled at equimolar amounts and sequenced unidirectionally on a Roche 454 GS-FLX Titanium platform at the IBIS/Plateforme d'Analyses Génomiques de l'Université Laval.

Sequence analyses were done in QIIME v1.5.0 and v1.8.0 (Caporaso et al., 2010). Low quality reads were removed using default parameters except for min quality score 30, min and max sequence lengths 330 and 550 bp, max homopolymer length 7, and no ambiguous bases. Reads that passed the quality control were denoised using DeNoiser v0.851 (Reeder and Knight, 2010) as implemented in QIIME v1.5.0. Chimeras were removed in MOTHUR v1.32.1 (Schloss et al., 2009) with uchime (Edgar et al., 2011) using self as reference and default settings. Demultiplexed reads were clustered into operational taxonomic units (OTUs) at 98% sequence similarity (as often used in protist metabarcoding; cf. Comeau et al., 2011) using UCLUST (Edgar, 2010) with default parameters and the taxonomy of the most abundant representatives from each OTU were assigned using BlastN (e-value 0.00001) against the curated protist databases of Lovejoy et al. (2016) and PR2 (gb203; Guillou et al., 2013). Global singletons and reads assigned to Metazoa were removed. The number of sequence reads per sample was rarefied by random sampling to the sample having the lowest read number (6066). Diversity analyses were performed using the Vegan (version 2.0–10) package in R (Version 0.98.994–RStudio, Inc.). The Shannon diversity index and Pielou's index of evenness were calculated to compare alpha-diversity between the samples. Rarefaction analyses were conducted and computed every 100 sequences using 5,000 subsampling iterations.

Beta diversities were analyzed by clustering samples using UPGMA (unweighted pairwise grouping method with averaging) based on Bray-Curtis dissimilarities. To identify a possible differentiation of the communities under the constraint of environmental factors, canonical correspondence analysis (CCA) was done using the PAST software (v3.04, Hammer et al., 2009). CCA was performed on a Bray-Curtis dissimilarity matrix of the samples. A permutation test was performed to calculate the significance of constraints on the data for the environmental parameters (fluorescence, density, depth, temperature, salinity, PAR and concentrations of nitrate/nitrite, phosphate, and silicate). The environmental parameters were tested for statistical dependence using the Spearman rank correlation coefficient (data not shown). For two highly co-varying parameters (salinity co-varied with density and silicate co-varied with phosphate) only one parameter was shown.

To visualize the differences among the libraries, a heat map of relative abundances was generated (in R) for the abundant OTUs (>1% of total reads) and the five most abundant OTUs in each library, and overlaid the UPGMA dendrogram.

At the time of sampling both ice stations were snow-covered (30 ± 9.7 cm and 23 ± 6.3 cm at Station 1 and Station 2, respectively) and had thick overlying sea ice (121 ± 11.7 cm and 160 ± 32.3 cm at Station 1 and Station 2, respectively). The water masses were identified based on their potential temperature—salinity characteristics; Atlantic Water (AW), Polar Water (PW), and Arctic Surface Water (ASW; Aagaard and Carmark, 1989; Cokelet et al., 2008; Figures 2, 3). The hydrographic profiles of the two stations were similar, although the PW was deeper at Station 2 (90 m) compared to Station 1 (45 m; Figure 3). The upper samples from both stations were all collected within the PW (S < 34.4; θ < 0°C), and designated PW1 (0, 10, and 45 m from Station 1) and PW2 (0 m and 15 m from Station 2; Table 1). The samples from 200 m at both locations displayed Atlantic Water properties. The samples from depths 65 m at Station 1 and 120 m at Station 2 were collected beneath the halocline and thermocline at the two stations, and represented Arctic Surface Water (ASW1 and ASW2 from stations 1 and 2, respectively). The two stations both consisted of ice floes with leads in between, and the euphotic zone (1% of surface PAR) reached 27 m at Station 1 and 34 m at Station 2 (Figure 3).

A distinct algal bloom had developed at Station 1, with Chl a concentrations of 10.7–16.0 μg L⁻¹ in the PW (Table 1). Station 2 saw an order of magnitude lower Chl a biomass at 1.0–1.3 μg L⁻¹ in the PW. The size fractionated Chl a-values showed that cells >3 μm dominated the photosynthetic biomass at both stations; at Station 1 the autotrophic biomass of the picoplanktonic fraction was low or below detection, and at Station 2 it was very low (0.06–0.08 μg L⁻¹).

The depth distribution of the measured nutrients (NO₃ + NO₂, PO₄, Si(OH)₄) was similar at the two stations, with comparatively low nutrient concentrations in the PW, minimum concentrations at intermediate depths (120 m at Station 2) and considerably higher concentrations in the AW (Table 1).

A total of 286,979 reads were obtained, which were reduced to 224,272 after quality filtering and chimera checking. Clustering the reads at 98% similarity produced 1,547 OTUs, which after removing OTUs classified as Metazoa, global singletons and subsampling to the smallest sample size (6,066 reads) were further reduced to 1,288 OTUs. The rarefaction curves (Figure S1) showed that all libraries had been sequenced to saturation. There was a general trend of increasing diversity with depth both in OTU richness and evenness, also reflected in the Shannon-Wiener diversity index (Table 1). A higher eukaryotic diversity was seen in the DNA samples compared to the RNA samples.

Hierarchical cluster analyses as well as the CCA showed a distinct clustering according to a combination of water mass and local environmental characteristics (Figures 4, 5) where the two first axes in the CCA described 54.2% of the variation in the data. The two samples collected from ASW were clearly associated with the communities of either AW (the ASW2 sample from 120 m depth; Figures 4, 5) or PW (the ASW1 sample clustered with the PW1 samples, in particular with the 45 m PW1 sample; Figure 5). The permutation test performed on the CCA confirmed the importance of water mass as a structuring factor, as significant interactions were identified between the community composition of individual samples and density (R² = 0.13, p = 0.01). In addition, the ANOVA identified fluorescence (R² = 0.19, p = 0.001) as the other significant factor structuring the communities. While density mostly separated along the first CCA axis (which explained 34.8% of the variability in the data), fluorescence separated the communities along the second CCA axis (which explained 19.4% of the variability in the data) thus separating the surface samples of the two stations. Other parameters explaining the variability seen in the community composition (although not statistically significant) were temperature, depth, nutrient concentrations and PAR (Figure 4).

Alveolates were abundant in all samples, recruiting 30–80% of reads in individual libraries (Figure 6). While Ciliophora dominated in PW (7–73% of library reads), Dinophyceae and MALV were more common in AW (including the 120 m ASW sample from Station 2). The Dinophyceae reads were dominated by one abundant OTU (assigned to Dino_clone_NIF_3A1; Table 2, Figure 5) that accounted for 3–44% of the total reads in individual libraries and was present in all libraries. While the maximum number of reads was detected in the 0 m DNA sample from Station 1 (1D0 library), RNA abundances were highest in the deep samples (all AW and 65 m from Station 1, Figure 5). Reads assigned to Woloszynskia halophila was also highly abundant in the 1D0 library, but did not have a similar read abundance in the RNA library (1R0; Figure 5). Reads assigned to MALV were more abundant in the DNA libraries than in the RNA libraries, while the dominant Ciliophora and Dinophyceae OTUs displayed the opposite pattern, recruiting a larger proportion of reads within the RNA libraries (Figure 5).

Three Ciliophora OTUs were among the abundant taxa; Choreotrichia-1 was abundant in the 45 and 65 m libraries from Station 1, whereas Strombidiidae and Oligotrichia were prevalent in PW, with especially high read numbers of the latter in RNA from Station 1 (Figure 5). Four additional Ciliophora OTUs were among the five most abundant in single libraries; Strombidium sp. did not show a very distinct distribution pattern, but was most abundant in PW, Monodinium sp. was particularly abundant in the 0 m RNA sample of Station 1, Mesodiniidae was abundant in the RNA libraries from AW, and Protostomatea-1 was found in the PW (including at 65 m) of Station 1.

The haptophytes were dominated by an OTU assigned to Phaeocystis sp. (Table 2), which was the most abundant OTU in the total dataset (16% of all reads) and was especially prevalent in the PW (and ASW) at Station 1 (6–49% of library reads). The same OTU also recruited 5–21% of all reads in the deeper AW sample from the same station. Another dominant group was the pico-sized green algal phylum Mamiellophyceae which recruited 35–45% of reads from the PW samples of Station 2 (Figure 5). Nearly all of these sequences were assigned to the Arctic ecotype of Micromonas (Arctic Micromonas; Table 2).

A Bacillariophyceae OTU assigned to Thalassiosira antarctica was found at both stations, but showed highest relative abundance in PW at 45 m of Station 1. OTUs assigned to marine stramenopiles were found throughout the water column at both stations, with the most abundant OTU (MAST 7; ANT12_19) along with an OTU assigned to Chrysophyceae-Synurophyceae Clade-H predominating at 45 m and 65 m in Station 1 (Figures 5, 6). Another stramenopile OTU, assigned to Labyrinthulaceae, was only abundant in AW. An OTU assigned to the cryptophyte Teleaulax gracilis was abundant in PW at both stations, but had higher RNA read numbers at Station 1 (Figure 5). Picozoa were found at all depths at both stations, but were most abundant in the libraries from intermediate and deeper depths (Figure 6). An OTU assigned to the rhizarian Cryomonadida (OTU 20) was highly abundant in the surface DNA library from Station 1. Other abundant taxa in the deep AW included two OTUs assigned to the radiolarian Taxopodida B (Table 2, Figures 5, 6).

The results of the UPGMA clustering indicated that among the deeper samples (120 and 200 m), libraries made from the same molecule type (DNA vs. RNA) were more similar than libraries from the same location or depth (Figure 5). In particular, reads assigned to the two Taxopodida B and the five abundant MALV OTUs were more abundant in the DNA libraries compared to the RNA libraries. This effect was, however, not apparent in the CCA analysis where libraries prepared from the same depth and station were most similar (Figure 4). A higher number of rare OTUs were identified from the DNA libraries compared to the RNA libraries, whereas the numbers of abundant OTUs were similar in the two molecule types (data not shown).

The two stations showed limited sympagic growth underneath the melting sea ice (personal observation), and the depth of the salinity-reduced surface water (PW; 45 m at Station 1 and 90 m at Station 2; Figure 3) suggest that sea ice melting had been ongoing for a while, although typical sea ice algae as well as the dinoflagellate Polarella glacialis with resting cysts were identified in the bottom 3 cm of ice cores taken in the vicinity (Johnsen et al., submitted). The high Chl a concentrations along with reduced nutrient concentrations and the algal community composition suggest that Station 1 was in a peak bloom phase (cf. Kristiansen et al., 2000; Johnsen et al., submitted). Station 2, on the other hand, had an order of magnitude lower Chl a concentrations while nutrient levels were similar to those of Station 1 (Table 1). Although a permanent stratification with surface mixed layer nitrate values around 5 μM may exist in the region in the pre bloom phase (Randelhoff et al., 2016), the low Chl a and nutrient concentrations at Station 2 suggest that this area had already developed into a post-bloom stage. The more shallow euphotic zone at Station 1 was probably due to shading from the ongoing bloom (Table 1, Figure 3).

The pelagic bloom in Station 1 was dominated by large and typical spring-bloom species of diatoms and had probably advected into the area (Johnsen et al., submitted). Common spring bloom species of Arctic waters such as Phaeocystis sp. and Thalassiosira antarctica (von Quillfeldt, 1997; Degerlund and Eilertsen, 2010) were among the most abundant OTUs in our dataset (Figures 5, 6), and were also abundant in phytoplankton samples from the area studied by microscopy (Johnsen et al., submitted). An advection of Phaeocystis sp. into ice-covered polar water in the Fram Strait was suggested also by Metfies et al. (2016), although a recent study report an under-ice bloom dominated by Phaeocystis sp. developing in situ (Assmy et al., 2017). Under-ice blooms may become more common in a scenario of thinner sea ice and a changing under-ice light climate (see also Arrigo et al., 2012). The clear dominance of photoautotrophs >3 μm in the advected bloom (Station1, Table 1) is similar to previous reports from Svalbard waters (Hodal et al., 2011; Marquardt et al., 2016). The OTU assigned to the cold water dinoflagellate Woloszynskia halophila, a bloomforming species of the Baltic Sea (Kremp et al., 2005), was highly abundant in the DNA library from the surface (0 m), while very few reads were found in the corresponding RNA library. As 18S sequence divergence is known to not always distinguish species of Dinophyceae (Logares et al., 2007), it is possible that the identified OTU represents another Suessiales species. Interestingly, Suessiales cysts (identified to Polarella glacialis; Johnsen et al., submitted) were found in sea ice cores as well as in the pelagic, and this OTU may represent cysts of a Suessiales taxon which is not metabolically active and thus not well represented in the RNA library.

Other potentially photosynthesizing plankton among the abundant OTUs include oligotrich ciliates (Table 2; Figure 5) where several representatives have sequestered chloroplasts that are retained and function as kleptoplastids, enabling photosynthesis by the host (reviewed by Esteban et al., 2010). In support, the very abundant Oligotrichia OTU had particularly high RNA read abundances at 0 and 10 m depth at Station 1, suggesting high metabolic activity of this taxon within the photic zone at this station. However, other Strombidiidae are heterotrophic, and we cannot exclude the possibility that the abundant oligotrichs are bacterivorous as suggested for the abundant Strombidinium in the Amundsen Gulf (Terrado et al., 2011).

The community composition of Station 2 was distinctly different from the advected bloom community identified in Station 1 (Figures 4, 6), and although most of the abundant OTUs were present at both stations their relative abundances differed. The arctic clade of Micromonas sp. (Arctic Micromonas) recruited over 30% of the reads in both RNA and DNA libraries from PW in Station 2, in line with this cold-and shade adapted ecotype of Micromonas (Lovejoy et al., 2007) being important in blooms developing in colder water (see also Metfies et al., 2016). In a recent temporal study from Isfjorden, West Spitsbergen, Arctic Micromonas had especially high relative read abundances during pre- and post-bloom stages (Marquardt et al., 2016), in agreement with Station 2 being in a post-bloom (alternatively pre-bloom) stage. The dominance of Arctic Micromonas over Phaeocystis sp. in the colder station 2 suggests that development in this area, located further to the northeast, was less influenced by the large pelagic bloom advected under the ice at Station 1.

An OTU assigned to the cryptomonad Teleaulax gracilis was abundant in surface waters of both stations, supporting the developed stages of the blooms as cryptophytes often appear late in the bloom (Leu et al., 2006). RNA read counts suggested that the cryptomonad was especially active at Station 1, in agreement with the ongoing bloom. Fluorescence and PAR were the two environmental parameters identified to separate the two surface communities (Figure 4), suggesting that succession and light conditions were the main drivers of the community differentiation in the surface.

The samples collected around the halocline (45 and 65 m in Station 1) had increased abundances of OTUs assigned to the autotrophs Phaeocystis sp. and T. antarctica probably reflecting sinking of the huge bloom (Figure 5). The same pattern was seen for an OTU assigned to another bloom species, Chaetoceros socialis (data not shown), which also was abundant in microscopy counts from the vicinity (Johnsen et al., submitted). Abundant heterotrophic flagellates most likely feeding on the sinking bloom included the ciliates Choreotrichia-1 and Protostomatea, Picozoa, and MAST 7 (ANT12_10; Figure 5). MAST 7 is more often found in and below the subsurface chlorophyll maximum layer compared to the upper mixed layer in the Arctic Ocean (Monier et al., 2013).

As commonly found in studies of picoplankton diversity based on filtration approaches, sequences from larger protists and metazoans were recovered, probably due to cell breakage and deformation of flexible walled cells allowing their DNA and RNA to pass through the 3 μm filters (Terrado et al., 2011; Sørensen et al., 2013). Another technical issue using 18S sequence data to compare protist communities is the potential for overestimation of certain OTUs due to high rDNA copy numbers or artifacts of the sequencing procedure (e.g., Zhu et al., 2005; Potvin and Lovejoy, 2009).

The increased diversity and evenness found in the AW microbial communities compared to the PW surface samples (Table 1), suggest that these deep communities are both diverse and metabolically active. Also, the increased representation of rare taxa in AW libraries may represent a seeding stock that could function as a biological buffer to environmental change (Sogin et al., 2006; Pedros-Alio, 2012). The deep, dark and nutrient-rich Atlantic current extends from the Norwegian Coastal Current, flowing in a 200–500 m deep layer. This deep AW layer may access the sunlit part by upwelling processes, as was shown to occur in January 2012 onto the shelf of northern Spitsbergen (Falk-Petersen et al., 2015). The AW harbors a distinct community of microbial protists brought northwards with the North Atlantic Current from more southern waters as well as Arctic species introduced through mixing with local water. The increased temperature identified in the waters of outer Isfjorden over the last decade (Pavlov et al., 2013) suggests that the organisms transported northwards with the West Spitsbergen Current may have increased survival, both in the deep and if they are transported to the surface. As argued by Pedros-Alio (2012), the rare biosphere may be composed of species whose individual requirements do not fit the current environment, but that may be able to bloom under different conditions.

There were distinct differences in the abundances of several OTUs occurring in the deeper AW between the DNA and RNA libraries (Figure 5). Most of the OTUs predominating in the AW (MALV, RAD, Labyrinthulaceae, Mesodiniidae) have heterotrophic or parasitic lifestyles, except the Mesodiinidae which have representatives at different trophic levels. Radiolarians are known to be associated with a large diversity of marine alveolates of which they are often hosts (Guillou et al., 2008; Bråte et al., 2012). The radiolarian OTU RAD B (Taxopodida; Suzuki and Not, 2015) and the five most abundant marine alveolate OTUs (Figure 5) were all abundant in the AW DNA libraries, but were not active (i.e., not abundant in the RNA libraries), as has previously been shown in other studies (Not et al., 2009; Hu et al., 2016; Marquardt et al., 2016). The limited activities shown by the parasitic protists (RAB B and MALV) in AW suggests that these taxa occur in less active life cycle forms.

Radiolaria are known to host a wide diversity of eukaryote symbionts, and were shown to harbor both the MALV I and MALV II lineages (Bråte et al., 2012). Whether the association between RAD B and MALV lineages in the deep AW is due to a parasite-host relationship is unknown, but a similar correlation of MALV sequence abundances with radiolarians was found also in samples from the Western Antarctic Peninsula with increased abundances in deep waters (Cleary and Durbin, 2016).

Other OTUs displayed depth or station related patterns in their RNA to DNA ratios, suggesting that their activity levels were modulated by changes in the environment. One such case was the very abundant dinoflagellate Dino_Clone_NIF_3A1 that only had high activity levels in the deep, suggesting it to be a heterotrophic taxon (Figure 5).

All of the abundant phototrophic OTUs were detected in the AW, albeit at reduced read numbers. Interestingly Phaeocystis sp. was also abundant in the AW at Station 1, both in the DNA and RNA libraries, showing that the detected cells were alive. This finding is in agreement with previous reports on viable Phaeocystis sp. in the deep (Vader et al., 2014), suggesting that these haptophytes could have a mixotrophic lifestyle as was shown for Arctic Micromonas (McKie Krisberg and Sanders, 2014).

In conclusion, the microbial eukaryote communities of the two under-ice stations differed due to a combination of their water mass history and local processes, the homogenizing effect of water currents being less important in surface waters. Furthermore, the combined use of RNA and DNA libraries allow interpretation on the activity levels of the identified taxa. Thus, an extended use of RNA libraries in similar studies may improve our understanding of the active part of these communities.

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.