Field sites

A total of 40 red snow samples was collected from four well-separated Arctic locations on 16 glaciers and permanent snow fields: Svalbard (SVA, n=12), Northern Sweden (TAR, n=24), Greenland (MIT, n=2) and Iceland (ICE, n=2; Fig. 1; Supplementary Table 1). These localities were chosen as they represent different geographical settings including low (67.9°N) versus high (78.9 N) latitude, low (150–400 m) versus high (∼1,200–1,400 m) elevation, and maritime versus continental settings. Vestre Brøggerbreen, Midtre Lovénbreen, Austre Brøggerbreen, Pedersenbreen, Austre Lovénbreen and Feiringbreen in Svalbard were sampled in July and August 2013. Samples from Storglaciären, Rabot, Liljetopsrännan, SE-Kasskasatjåkkå, Björling and nearby permanent snow fields in Northern Sweden were collected in July 2013 and July 2014. Mittivakkat glacier in Greenland and the glacier Drangajökull and permanent snow field Laugafell in Iceland were sampled in July 2012. Red snow samples were collected late in the melt season as those are the typical snow algal blooms that will have the largest effect on albedo.

Field sampling and measurements

All sampling, field measurements and most analyses have previously been described in full detail11,33. Here we summarize previously employed methods and give full details of new methods. At each sampling site, we measured pH, conductivity and temperature with a daily calibrated metre (Hanna instruments, HI 98129) before sampling. Photosynthetic active radiation, ultraviolet radiation and surface albedo (400–700 nm range) were measured using a radiometer (SolarLight, PMA2100) with specific photosynthetic active radiation (PMA2132), ultraviolet-A (PMA2110) and ultraviolet-B (PMA2106) sensors. Albedo was calculated by taking the ratio of reflected to incident radiation (400–700 nm range) and measuring the values always in the same position to the sun. The reading of the sensor was not affected by shading by the observer. Measurements were carried out with the sensors held at 30 cm above the snow surface (field of view 160°). At first, the sensor was pointed upwards (incident radiation) and then towards the snow surface (reflected radiation). Five measurements for incident and reflected radiation were acquired each, and the average was taken to avoid measuring bias. The standard deviations for each measurement set was below 10%. Data of the relative contribution of pigments, other light impurities (that is, mineral dust, black carbon) and snow metamorphism is lacking42. However, based on qualitative microscopic observations in the field, particularly mineral impurities (most often light coloured quartz and feldspars) were less important in changing albedo measurements in red snow surface samples in comparison to the pigment distribution. Moreover, Aoki et al.35 and Painter et al.34 showed that red snow has much higher light absorption below 600 nm because of the algal pigments in comparison to mineral dust or black carbon. Samples were collected in sterile centrifuge tubes or sterile Whirl-Pak bags and in pre-ashed glass jars (450 °C, >4 h) for organic analyses. All samples were slowly melted at room temperature within ∼6 h, and processed and preserved (for example, filtered, acidified) within ∼8 h after collection. Samples for DNA and organic analyses were flash-frozen in liquid nitrogen and stored at −80 °C until analysed. Inorganic samples were stored cold (4 °C) and in the dark. All analyses were carried out in the Cohen Laboratories at the University of Leeds unless stated otherwise.

Aqueous analyses

Aqueous analyses were carried out by Ion Chromatography (Dionex, anions), by inductively coupled plasma mass spectrometry (Agilent, cations, at the University of Sheffield), on a total organic carbon analyser (Shimadzu, for DOC contents, DOC, at the Plymouth University), and by segmented flow-injection analyses (AutoAnalyser3, Seal Analytical, dissolved phosphate).

Particulate analyses

Particulates in the samples were analysed for δ15N and δ13C by a Vario Pyro Cube elemental analyser (Elementar Inc) coupled to an Isoprime mass spectrometer. Samples were combusted in tin capsules at 1,150 °C, and gases were separated using temperature-controlled adsorption/desorption columns. Carbon analyses were calibrated with in-house C4-sucrose and urea standards assigned values of −11.93‰ and −46.83‰, respectively via calibration with the international standards LSVEC (−46.479‰), CH7 (−31.83‰), CH6 (−10.45‰) and CO-1 (+2.48‰). Nitrogen isotope values were calibrated using the international standards USGS-25 and USGS-26 with assigned values of −30.4‰ and +53.7‰, respectively. Total carbon, total nitrogen and total sulphur were derived from the thermal conductivity detector in the elemental analyser and calibrated using a sulphanilamide standard. Particulate phosphorus was extracted by ashing of the samples at 550 °C for 2 h and incubating in 1 M HCl for 16 h according to extraction step V in Ruttenberg et al.44.

Algal biomass

Algal cells were imaged on a Leica DM750 microscope equipped with a × 63objective and counted with a haemocytometer in triplicate. For cell size analyses, 100 cell diameters per sample were measured in ImageJ. Cell volumes were calculated assuming a perfect spherical shape (V=4/3*π*r3). Total algal biomass was calculated using the average cell volume and cell abundance.

Pigment analysis

Carotenoid and chlorophyll contents in the samples were quantified by high-pressure liquid chromatography (HPLC) and a modified carotenoid/chlorophyll-specific extraction protocol45. Cells were disrupted by shock freezing in liquid nitrogen for 10 min followed by grinding with a Teflon mortar and pestle. The resulting powder was re-suspended in 1 ml of dimethylformamide and 1.0 mm glass beads and horizontally shaken on a laboratory shaker (MoBio Vortex Genie 2) at maximum speed (3,000 r.p.m.) for 10 min, followed by centrifugation for 5 min at 10,000 r.p.m. The supernatant was separated from the debris by filtering through a 0.45-μm Teflon filter and the filtrate was mixed with methanol (25 vol%).

Extracted samples were analysed immediately on an Agilent Technologies 1200 Infinity HPLC instrument with a gradient pump, an autosampler, a variable wavelength detector and ODS Hypersil column (250 × 4.6 mm2; 5 μm particle size). Two solvents were used: solvent A consisted of a mixture of acetonitrile/water/methanol/hexane/tris buffer at ratios of 80:7:3:1:1, whereas solvent B was a mix of methanol and hexane at a ratio of 5:1. The HPLC was run at a flow rate of 1 ml min−1, with an injection volume of 25 μl. Spectra were recorded from 200 to 800 nm. Chromatograms were quantified at 450 nm for carotenoids and 660 nm for chlorophyll a and b. Run time was 60 min. The protocol required a 15-min run with 100% of solvent A followed by a linear gradient from 100% solvent A, to 100% solvent B between 32 and 45 min, and finally with 15 min of column re-equilibration through a 5-min linear gradient from solvent B back to 100% solvent to A, followed by a further column conditioning with 100% solvent A for 10 min. The following commercially available standards were used for peak identification: chlorophyll a, chlorophyll b (Sigma), violaxanthin, neoxanthin, antheraxanthin, lutein, β-carotene, trans-astaxanthin and cis-astaxanthin (Carotenature).

Fatty acids analysis

Fatty acids were extracted according to the method described by Wacker and Martin-Creuzberg46. Briefly, 20 ng of internal standard (tricosanoic acid methyl ester) were added to each sample, followed by ultrasonic extraction using dichloromethane:methanol (2:1 (v:v)), and centrifugation to remove particulates and evaporation of solvent from the supernatant. Fatty acids were transesterified by adding methanolic HCl to the dried extract and heating at 60 ° C for 20 min. After cooling, fatty acid methyl esters were extracted in isohexane, the solvent was removed under nitrogen and the sample resuspended in isohexane for analysis.

Analysis of fatty acid methyl esters was carried out using a Trace 1300 gas chromatograph with flame ionization detector (Thermo Scientific), equipped with a non-polar-fused silica capillary column (CPSil-5CB, 50 m × 0.32 mm × 0.12 mm, Agilent Technologies). Samples (1 μl) were injected in splitless mode, with the injector maintained at 200 °C. Carrier gas was helium, with a constant flow rate of 1.5 ml min−1. The following temperature programme was used: initial temperature 40 °C, rising to 140 °C at 20 °C min−1, then rising to 240 °C at 4 °C min−1, holding at 240 °C for 5 min. Fatty acid methyl esters were identified by comparison of retention time with those of reference compounds (Supelco) and by gas chromatography mass spectrometry (GC–MS). GC–MS was carried out using the gas chromatograph and column previously described, with identical operating conditions, coupled to an ISQ mass spectrometer (Thermo Scientific). The transfer line and the ion source were maintained at 300 °C. The emission current was set to 50 mA and the electron energy to 70 eV. The analyser was set to scan at m/z 50–650 with a scan cycle time of 0.6 s.

DNA sequencing

Total DNA was extracted from pelleted biomass using the PowerSoil DNA Isolation kit (MoBio Laboratories). 16S rRNA genes were amplified using bacterial primers 27F (5′-AGAGTTTGATCMTGGCTCAG-3′) and 357R (5′-CTGCTGCCTYCCGTA-3′; tagged with the Ion Torrent adapter sequences and MID barcode) spanning the V1-V2 hypervariable regions. 18S rRNA genes were amplified using the eukaryotic primers 528F (5′-GCGGTAATTCCAGCTCCAA-3′) and 706R (5′-AATCCRAGAATTTCACCTCT-3′; Cheung et al., 2010 (ref. 48); tagged with the Ion Torrent adapter sequences and MID barcode) spanning the V4-V5 hypervariable region. PCRs were performed using Platinum PCR SuperMix High Fidelity according to the manufacturer’s protocols. Initial denaturation at 95 °C for 5 min was followed by 30 cycles of denaturation at 95 °C for 30 s, annealing at 60 °C for 30 s and elongation at 72 °C for 30 s. Final elongation was at 72 °C for 7 min. Archaeal 16S rRNA genes were amplified following a nested PCR approach. The first PCR reaction was carried out using primers 20F (5′-TCCGGTTGATCCYGCCRG-3′) and 915R (5′-GTGCTCCCCCGCCAATTCCT-3′). Initial denaturation at 95 °C for 5 min was followed by 35 cycles of denaturation at 95 °C for 30 s, annealing at 62 °C for 30 s and elongation at 72 °C for 180 s. Final elongation was at 72 °C for 10 min. The PCR product was used as template for the second PCR reaction with primers 21F (5′-TCCGGTTGATCCYGCCGG-3′) and 519R (5′-GWATTACCGCGGCKGCTG-3′; tagged with the Ion Torrent adapter sequences and MID barcode) spanning the V1-V2 hypervariable region. Initial denaturation at 95 °C for 5 min was followed by 30 cycles of denaturation at 95 °C for 30 s, annealing at 60 °C for 30 s and elongation at 72 °C for 30 s. Final elongation was at 72 °C for 7 min. Detailed information on the sequencing primers can be found in the Supplementary Information. All PCRs were carried out in triplicates to reduce amplification bias and in reaction volumes of 1 × 25 μl and 2 × 12.5 μl. All pre-amplification steps were done in a laminar flow hood with DNA-free certified plasticware and filter tips. The pooled amplicons were purified with AMPure XP beads (Agencourt) with a bead–to-DNA ratio of 0.6 to remove nucleotides, salts and primers. Quality, size and concentration were determined on the Agilent 2100 Bioanalyser (Agilent Technologies) with the High-Sensitivity DNA kit (Agilent Technologies). Sequencing was performed on an Ion Torrent Personal Genome Machine using the Ion Xpress Template Kit and the Ion 314 or Ion 316 chips following the manufacturer’s protocols. All PCR amplifications and sequencing were carried out at the Aberystwyth University and the University of Bristol. The raw sequence data were processed in QIIME47. Barcodes and adapter sequences were removed from each sequence. Filtering of sequences was performed using an average cutoff of Q20 over the full sequence length (350 bp). Reads shorter than 200 bp were removed. Operational taxonomic units (OTUs) were picked de novo using a threshold of 97% identity. Taxonomic identities were assigned for representative sequences of each OTU using the reference databases Greengenes for bacteria and archaea. The Silva database (ref. 49; extended with additional 223 sequences of cryophilic algae kindly provided by Dr Thomas Leya from the CCCryo—Culture Collection of Cryophilic Algae, Fraunhofer IZI-BB) was used for eukaryotes. Data were aligned using PyNAST and a 0.80 confidence threshold. Singletons were excluded from the analysis. For bacterial sequence matching, plant plastids were removed from the data set before further analysis. For eukaryotic sequence matching, Chloroplastida were pulled out of the data set and stored in a separate OTU table. In order to focus upon algal diversity, sequences matching Embryophyta (for example, moss, fern) were removed from the data set. For archaea, sequences matching bacteria were removed. Finally, for further analyses, samples were rarefied to the minimum library size and Shannon indices were calculated in QIIME. All analyses were conducted at the 97% OTU level. A matrix of each OTU table representing relative abundance (raw data) was imported into PAST v3.06 (ref. 50) for multivariate statistical analyses (principal component analysis, canonical correspondence analysis) and Pearson correlations. One-way analysis of variance test was done in SPSS v19 (IBM).

Sequencing primers

Primers targeting the 18S rRNA gene were chosen because there are more sequences in the databases for green algae (that is, Chlorophyta, Charophyta) than for rbcL or internal transcribed spacer (ITS). Before sequencing, we carried out an in-silico investigation including 18S rNRA sequences from 218 snow and permafrost algae in order to make sure that the chosen primers are suitable for green algae and that there is enough variability in the chosen region (v4-v5) to distinguish between species.

Previous studies51 have found that one primer pair is not sufficient to recover all eukaryotic groups in one sample. However, we chose our primer pair based on one group we were specifically targeting, that is, the green algae. We do not claim to have equally recovered all other groups among the micro-eukaryotes such as fungi or the ‘SAR’-group. Furthermore, they found that libraries derived from different primer pairs grouped together for individual samples with no significant differences. Based on our in silico test of 218 snow and permafrost algae and the rarefaction curves (Supplementary Fig. 6), we are fairly confident that the choice of our primer pair has resulted in a good coverage of the algal diversity. However, we acknowledge that PCR-based approaches will always introduce a certain amount of bias.

This is similar for the archaea, which show no biogeographical patterns in our samples. The primers used are specific for archaea and since they are not the focus of this study and only the associated microbiome, we did not explore other primer possibilities. However, the results match what other studies have found before in cryo-environments28,52.

Overall sampling design

All samples for DNA and aqueous analyses were analysed in exactly the same way for all samples from Greenland, Iceland, Svalbard and Sweden. Pigment and fatty acid data are only shown for the samples from Svalbard and Sweden because for Greenland and Iceland these data have previously been published11,33. The samples from Greenland have been excluded from the comparison here because the pigment and fatty acid data have been collected and quantified in a different way. The pigments were normalized to chlorophyll a, whereas in all other study area they were quantified with the appropriate pigment standards. The fatty acid components were analysed by GC–MS, whereas all samples from Svalbard and Sweden were also quantified by flame ionization detection. The pigment and fatty acid data from Iceland were also excluded from the comparison, as in all samples large amounts of moss (identified by microscopy and DNA sequences) that could not be separated from the algae before pigment and fatty acid extraction were present. This moss contribution would strongly ‘skew’ the data and thus these were excluded.

In addition, only selected samples in Greenland11 and Iceland33 were included in the comparison. This is because at both sites samples were collected at different stages in the melt season. The study in Greenland was conducted at the onset of melting and over a 3-week period when snow algae just started to bloom and a decrease in relative chlorophyll content and increase in carotenoid content could be observed. This led to our conclusion that there is a great heterogeneity in pigment composition both in space and time11. However, the few samples collected at the end of the study showed similar carotenoid contents. This end-of-season homogeneity in the red snow samples was the reason why in the current study we focused solely on samples from late in the melt season, which is the dominant red snow stage with the largest impact on albedo. Thus, we only included two DNA samples from Greenland. Similarly in Iceland33, most of the samples were collected earlier in the year (June—July) and those samples were described as less ‘typical’ of red snow patches33. Off all samples from Iceland again only the two samples that were collected late in the melt season and therefore matched the conditions of the samples in the current study were used for comparisons.

Integrated albedo change

Using our mean, minimum and maximum measured albedo values for wet clean snow and red snow and literature data22 for clean dry snow, we used a simple one-dimensional moving boundary approach that allows us to predict the effect of red pigmented snow algae on albedo. This approach is valid under the assumption that the snow and ice surfaces melt downwards relative to a fixed depth, and that at the same time such a change is accompanied by changes in albedo38. The parameters, equations and boundary conditions used are as follows:

Table 1 shows measured minimum, maximum and average albedo values for dry clean snow before the onset of melting11, wet clean snow (no visual presence of algae) at the onset of melting and red snow (full red pigmented snow algal bloom). We used these values to derive linear regressions for albedo changes over a 100-day melt season (Supplementary Fig. 4). A conservative 100-day scenario was chosen, as this encompasses all our albedo measurements in the current and previous studies (June—August)3,11,33. In addition, this also corresponds to the number of days with mean air temperatures above 0 °C in the same period (Ny Alesund: 116 days in 2013 and 105 days in 2014, kindly provided by Dr Marion Maturilli and Siegrid Debatin, AWI; Storglaciären: 132 days in 2013 and 108 days in 2014; kindly provided by Dr Peter Jansson, Stockholm University; data are also publicly available at http://bolin.su.se/data/tarfala/). We compare a benchmark case of purely physical-driven albedo change (that is, changes in snow crystal sizes and shapes and increasing water content, scenario 1) with albedo change due to red pigmented algal growth (scenario 2).

Scenario 1 considers the transition over 25 days from clean snow to a wet melting surface without algal growth and with an albedo of 0.80 (a minimum value), which with continued melting results in an albedo of 0.75 (an average value) after 50 days and 0.70 (a maximum value) after 100 days (Table 1). Our benchmark case (scenario 1) shows albedo (α) changes with time and fits the equation:

Scenario 2 considers the transition from clean snow to a surface where the growth of algae after 25 days produces light red snow with an albedo of 0.77 (a minimum value), and continued melting produces darker red snow with an albedo of 0.65 (an average value) after 50 days and 0.53 (a maximum value) after 100 days (Table 1). The albedo changes with time for this scenario fit the equation:

These two equations can be integrated to obtain the cumulative effects of albedo (α) changes with time to give:

Subtracting equations (3) and (4) gives

which represents the albedo changes attributable to algae growth alone. For a melt season of 100 days Δα=1.85+11.3=13.15≈13.

In order to assess the error of our analysis, we carried out a sensitivity analysis using the data below (see also Supplementary Fig. 5 for details):

Comparing minimum and average albedo values:

Subtracting gives Δα=0.0409t–0.0038t1.5

So when t=100, Δα=4.09–0.0038 × 1,000=4.09–3.8=0.29

Comparing average and maximum values:

Subtracting gives Δα=0.0408t–0.0038t1.5

So when t=100, Δα=4.08–0.0038 × 1,000=4.08–3.8=0.28

So our sensitivity test is giving a crude range of ∼0.3 about the mean.

Data availability

DNA sequences have been deposited to the European Nucleotide Archive (ENA) under accession number PRJEB10548. All other data are available in the Supplementary Information.