Modeling seasonal growth of phototrophs on bare ice on the Qaanaaq Ice Cap, northwestern Greenland

Abstract

Glacier phototroph blooms on the surfaces of ice sheets and glaciers cause albedo reduction, leading to increased melting rates. We observed seasonal changes in the abundance of phototrophs on the Qaanaaq Ice Cap in northwestern Greenland from June to August 2014, and reproduced these changes using numerical and empirical models. The phototroph community on the ice surface mainly consisted of the glacier alga Ancylonema nordenskioldii and the cyanobacterium Phormidesmis priestleyi. The glacier alga appeared on the ice surface in late June, after which its abundance increased exponentially throughout the melting period. A logistic growth model designed for snow algal growth reproduced the measured exponential increases, suggesting that growth could be explained using the model as a function of the ice melting duration. Cyanobacteria appeared and their abundance increased in late July but did not change exponentially thereafter. The abundance of cyanobacteria was explained with an empirical model expressed as a function of the amount of mineral dust on the bare ice surface. Our numerical and empirical models for reproducing glacier algae and cyanobacteria could be useful for quantifying the albedo reduction caused by their growth and the melt rates of the Greenland ice sheet and glaciers in the future.

1. Introduction

The Greenland ice sheet (GrIS) is the second-largest continuous body of ice worldwide; however, it has lost mass because of a decrease in the surface mass balance and an increase in ice discharge over the last two decades (Shepherd and others, 2020; IPCC, 2022). The decrease in the surface mass balance is mainly attributed to increased air temperature and decreased surface albedo (Riihelä and others, 2019; Hanna and others, 2021). In the last decade, ice surface darkening, which involves decreased surface albedo caused by biological activity on the bare ice surface, has been observed in the southwestern GrIS, and has contributed to the accelerated surface mass loss of the ice sheet based on field and satellite observations (Wientjes and others, 2011; Yallop and others, 2012; Shimada and others, 2016; Wang and others, 2018, 2020; Cook and others, 2020).

Such darkening phenomena are mainly caused by an exponential increase in the population (phototroph bloom) of cold-adapted autotrophic microbes on the snow and ice surfaces of ice sheets and glaciers. For example, Sanguina nivaloides, a surface snow alga, is commonly found in the upper snow area of glaciers. Blooming of this species results in red surface snow because of its secondary carotenoid pigments (Cook and others, 2017; Ganey and others, 2017; Procházková and others, 2019; Onuma and others, 2020). In contrast, more diverse phototrophs have been observed on the ice surface at lower elevations (Lutz and others, 2014; Tanaka and others, 2016; Uetake and others, 2016; Anesio and others, 2017; Takeuchi and others, 2019). Algae found on the bare ice surface are distinct from those found on snow surfaces and are referred to as glacier algae; blooming populations can be identified as purple or brown areas on ice surfaces (Cook and others, 2020; Di Mauro and others, 2020; Williamson and others, 2020; McCutcheon and others, 2021). Furthermore, filamentous cyanobacteria enable the entanglement of mineral and organic particles and form dark-colored aggregates, known as cryoconite granules, on the ice surface (Hodson and others, 2010; Langford and others, 2010; Takeuchi and others, 2010; Wientjes and others, 2011; Uetake and others, 2019).

Blooms of glacier phototrophs (both glacier algae and cyanobacteria), including Ancylonema nordenskioldii, Ancylonema alaskana (recently renamed from Mesotaenium berggrenii by Procházková and others, 2021), Cylindrocystis (Cyl.) brebissonii and Phormidesmis priestleyi commonly occur in ice sheets and glaciers (Uetake and others, 2010; Takeuchi and others, 2014; Stibal and others, 2017; Lutz and others, 2018; Williamson and others, 2018). Blooms of glacier algae, which have dark-reddish pigments in their cells, have been observed on the ice surface during the summer. These pigments are likely to play an adaptive role in extracting the liquid water required for the growth of glacier algae through ice melting (Dial and others, 2018). Additionally, cryoconite granules formed by filamentous cyanobacteria on ice surfaces can sink into glacial ice and form cryoconite holes by effectively absorbing solar radiation on the ice surface because of their dark color, which affects the ice surface albedo on the GrIS and glaciers (Chandler and others, 2015; Cook and others, 2016; Takeuchi and others, 2018).

The habitats of glacier phototrophs on ablating glacier surfaces are topologically heterogeneous and can be classified into three main types: (1) bare ice surfaces, (2) cryoconite holes and (3) meltwater streams (Irvine-Fynn and others, 2014; Holland and others, 2019; Tedstone and others, 2020). Among these habitats, the bare ice surface is typically the most extensive in the GrIS and glaciers (Chandler and others, 2015; Ryan and others, 2018; Tedstone and others, 2020). On bare ice surfaces, glacier algal cells typically appear to be present in the thin layer of liquid meltwater, which increases during the melting season, and are washed out by heavy rainfall (Stibal and others, 2017; Williamson and others, 2018). Glacier algae are concentrated in rough local depressions on the ice surface and in porous ice (weathering crust), suggesting that ice surface roughness is an important factor that affects the abundance of glacier algae (Tedstone and others, 2020). In contrast, a previous study suggested that the abundance of cyanobacteria is closely related to the amount of mineral dust rather than to nutrient conditions on the bare ice surface (Uetake and others, 2016).

To predict the growth of microbes and their effects on surface albedo and melt rates of glacier surfaces using a climate model for the GrIS (Langen and others, 2015; Fettweis and others, 2017; Niwano and others, 2018; Noël and others, 2018), a numerical model that accounts for the effects of the growth of glacier phototrophs is essential. Several models have been proposed to reproduce temporal changes in microbial abundance on snow and ice surfaces. In mountain snowpacks or upper accumulation areas of glaciers, the cell abundance of green algae on the snow surface exponentially increases with snow melting, as reported by field observations worldwide, such as in Greenland, Alaska and Japan (Takeuchi, 2013; Onuma and others, 2016, 2018); the authors proposed that temporal changes in the abundance of snow algal cells on surface snow could be expressed using a simple differential logistic growth equation known as the ‘snow algae model’. The latest version of the snow algae model is used to reproduce global spatiotemporal changes in snow algae abundance (Onuma and others, 2022b). This model can reproduce the exponential growth of snow algae during snow melting using their initial cell concentration, growth rate and carrying capacity, which limits algal growth. Field observations in the ablation area of the GrIS have reported that glacier algal abundance exponentially increases on the ice surface during summer (Stibal and others, 2017; Williamson and others, 2018). Their studies suggested that the growth rate of glacier algae can be represented using the doubling time, which is the inverse of the cell division rate of algae. In addition, two models for glacier algal blooms have been proposed on the basis of field observations in the southwestern GrIS: a model to reproduce the increase in glacier algal abundance in response to the recession of the snow line (Williamson and others, 2018) and a model to reproduce glacier algal abundance using solar radiation, air temperature and snow depth on the bare ice (Williamson and others, 2020). However, studies on numerical modeling of the growth of glacier phototrophs, especially glacier cyanobacteria, are still limited. Furthermore, a numerical model for the growth of glacier algae is necessary based on a different concept from previous studies in order to evaluate its performance by inter-comparison using several models in the future. Therefore, we applied the snow algae model to glacier phototroph growth on the ice surface.

In this study, biological, chemical and meteorological observations were conducted on the ablating ice surface of the Qaanaaq Ice Cap in northwestern Greenland to establish a numerical model to reproduce the temporal changes in the abundance of glacier phototrophs throughout the melting season. First, we describe temporal changes in algal abundance at four locations in the ablation area of the glacier from June to August 2014 and evaluate the seasonal growth patterns of the glacier phototrophs in Section 3.1. In Section 3.2, we describe the numerical model established for glacier algal growth based on periodic observations, and model evaluation using glacier algal abundance observed for the three seasons. In Section 3.3, we describe and evaluate the empirical model established for glacier cyanobacterial growth, similar to the previous section.

2. Study sites and methods

2.1. Site location and its biological features

The study was conducted at the Qaanaaq Ice Cap in northwestern Greenland (Figs 1a, b) from June to August of 2012, 2013 and 2014. This ice cap, which is located on a small peninsula in northwestern Greenland, covers an area of 286 km², with an elevation range of 0–1200 m a.s.l. (Sugiyama and others, 2014, 2017). We selected five study sites at different elevations (S1, S2, S3, S4 and S5) on an outlet glacier of the Qaanaaq Ice Cap that is easily accessible by foot from the Qaanaaq village. These sites are located in the middle of the glacier at elevations of 247, 441, 672, 772 and 944 m a.s.l. from sites S1 to S5, respectively. The equilibrium line altitude of the glacier ranged between 862 and 1001 m a.s.l. for the three seasons (Tsutaki and others, 2017); thus, the ice surfaces at these sites were exposed during the melting season.

Fig. 1. Map of Greenland (a) and Qaanaaq Ice Cap in northwestern Greenland (b) and typical bare ice surface for sample collection in this study (c). Sampling sites are shown in (b). The snow lines for 2012, 2013 and 2014 reached elevations above S5, above S2, and between S4 and S5 during the observation seasons, respectively.

Previous studies have revealed five major taxa of phototrophs on the glacier, including three algae (A. nordenskioldii, A. alaskana and S. nivaloides) and two cyanobacteria (P. priestleyi and Chroococcaceae (Chr.) cyanobacterium) (Uetake and others, 2010, 2016; Takeuchi and others, 2014). The community structure varied with elevation: the glacier alga A. nordenskioldii was dominant at the lower sites (S1 and S2), whereas both A. nordenskioldii and P. priestleyi were present in the middle parts of the glacier (Uetake and others, 2016). The albedo of the ice surface was significantly lower at the midstream site at S2, S3 and S4 than at the clear ice surface at S1 because of the abundance of dark impurities, including glacier phototrophs (Takeuchi and others, 2014).

2.2. Surface ice collections and analysis methods

The samples were collected from bare ice surfaces at the study site (Fig. 1c). Field observations suggest that the total area of the bare ice surface accounts for more than 80% of the ablation zone of the Qaanaaq Ice Cap and the southwestern GrIS (Chandler and others, 2015; Ryan and others, 2018; Takeuchi and others, 2018; Tedstone and others, 2020). Therefore, we focused on glacier phototrophs only on the bare ice surface, which likely has the greatest effect on the surface albedo reduction in the bare ice area of the GrIS.

Surface ice collection was conducted periodically from day 176 (25 June) to day 215 (3 August) in 2014. During this period, ice was collected three times at all sites from day 176 to day 197 (16 July), and five times at all sites from day 198. Details of sample collection frequencies are listed in Supplementray Table 1. The surface ice of the glacier was covered with winter snow until late June. It was first exposed at the lower sites (S1 and S2) in late June, and then at the higher sites (S3 and S4) in mid-July. As the ice surface at the highest site S5 was not exposed until the end of the season (Onuma and others, 2018, 2020), no sample was collected from this site.

Surface ice was collected in late July 2012 and 2013 (days 200 and 206, respectively). These samples were collected five times at sites S1–S5 in 2012, but only at S1 and S2 in 2013 because of the presence of snow cover at site S3 in late July. Phototroph abundances and mineral weights quantified from samples collected at the study sites in 2012 have been previously published (Takeuchi and others, 2014; Uetake and others, 2016). In this study, we compare these observational results with those of late July 2014 (days 201 and 203) to focus on the annual changes in the glacier phototroph community.

Samples were collected from 3 to 5 randomly selected surfaces (depth of the top 1 cm) using a stainless-steel scoop at each study site. The sampling area ranged from 100 to 400 cm² and was recorded for each collection. All ice samples were preserved in plastic bags (Whirl-Pak® bags; Nasco, Fort Atkinson, WI, USA). The electrical conductivity (EC) and pH of the collected samples were measured using a portable pH-conductivity meter (F-54, HORIBA, Kyoto, Japan) after the samples were melted in the Qaanaaq village. The samples used for mineral abundance and algal cell analyses were collected separately. These samples were melted at room temperature and immediately preserved in 3% formalin in 30 mL clean polyethylene bottles before being transported to Chiba University, Japan, for analysis.

The abundance of glacier phototrophs (bio-volume) was quantified using the direct cell counting method and represented as the bio-volume per unit ice surface area, according to previous studies (Tanaka and others, 2016; Onuma and others, 2018). The sample water (2–100 μL) was filtered through a hydrophilized PTFE membrane filter (pore size 0.45 μm; Omnipore JHWP, Merck Millipore, Japan). The number of algal cells on the filter was counted using an optical microscope (BX51; Olympus, Tokyo, Japan). Filamentous cyanobacteria were counted in units of length (1.4 μm). Counting was conducted three times for each sample, and the cell concentration (cells mL⁻¹) in the sample was obtained from the mean of the cell counts and filtered sample volume. Cell numbers per unit area (cells m⁻²) were calculated using the cell concentration and area of sample collection. To determine the bio-volume, the mean cell volumes were geometrically estimated by measuring the size of 5–80 cells for each species under a microscope. The total bio-volume per unit area (μL m⁻²) for each taxon was obtained by multiplying the cell number and cell volume.

The abundance of mineral dust on the bare ice surface was quantified using the combustion method (Takeuchi and Li, 2008; Onuma and others, 2018). The samples were melted at room temperature in the Qaanaaq village and their masses were measured on a weight scale. The dust precipitated in the bag was preserved in 30 mL clean polyethylene bottles and transported to Chiba University, Japan for analysis. The samples were dried (60°C, 24 h) in pre-weighed crucibles and combusted at 500°C for 3 h in an electric furnace to remove organic matter. The abundance of mineral dust per unit area (g m⁻²) was obtained from the combusted sample weight and collected sample area, as only mineral dust remained after combustion.

The statistical values for the 2012, 2013 and 2014 seasons regarding the analyzed EC, pH, bio-volume and mineral dust abundance are summarized in the supplemental information (Supplementray Tables 1 and 2).

2.3. Meteorological observation

Meteorological data used in the study were collected using an automatic weather station (AWS) at the SIGMA-B site (77°31′N, 69°04′W; 944 m a.s.l.), which is the same position as site S5 in this study, installed on day 201 (July 19) in 2012 (Aoki and others, 2014). The surface air temperature and solar radiation data were collected hourly during the summer season (June to August) from 2012 to 2014 using the AWS. The temperature sensor and pyranometer of the AWS were placed at heights of 3.0 and 2.5 m above the snow surface. Aoki and others (2014) have described the AWS in detail. Seasonal changes in the measured meteorological conditions for the three seasons are shown in Figure 2 (solid black lines). Meteorological observations indicated that 2012 was the warmest season and 2013 was the coldest season among the three years.

Fig. 2. Temporal changes in meteorological conditions for 2012–2014 seasons at S5 (944 m a.s.l.) in Qaanaaq Ice Cap. Blue and black solid lines indicate atmospheric reanalysis and observational data, respectively. Blue solid line and shading indicate average values and minimum or maximum values derived from the three reanalysis data sets (WFDEI, GSWP3-FD and CRUJRA), respectively. R and RMSE in the figures mean Pearson correlation coefficient and root mean square error between reanalysis and automatic weather station data, respectively.

2.4. Atmospheric datasets and ice melting simulation

As we did not measure precipitation levels with AWS observations throughout the 2012–2014 seasons, we additionally used reanalysis datasets obtained from atmospheric simulations to ensure consistency with meteorological observations, which created gridded historical datasets at global or regional scales. We used three different reanalysis datasets, WFDEI (Weedon and others, 2014), GSWP3-FD (van den Hurk and others, 2016; Kim and others, 2017) and CRUJRA (Harris, 2019), which were derived from global reanalysis data near the land surface. The datasets include three- or six-hourly information on surface air temperature, surface air pressure, downward radiation (shortwave and longwave), humidity, wind speed and precipitation rate. As the datasets have a low horizontal resolution (0.5° × 0.5° globally), they were corrected with elevation information at the study sites following the protocol of Onuma and others (2022b). The temperature lapse rate for correction was assumed to be −7.80 × 10⁻³ K m⁻¹, which was estimated from field observations of the Qaanaaq Ice Cap in late July 2012 (Sugiyama and others, 2014). Rainfall, snowfall and ice surface temperature at the study sites were derived from numerical simulations using the land surface model described in the following paragraph. Seasonal changes in atmospheric conditions at S5 and S1 for the three seasons derived from the reanalysis datasets are shown in Figure 2 and Supplementary Figure 1 (blue solid lines), respectively.

To estimate the duration of ice melting at the study sites, we conducted numerical simulations using a land surface model with reanalysis datasets. The land surface model used in this study was the minimal advanced treatments of surface interaction and runoff (MATSIRO; Takata and others, 2003; Nitta and others, 2014, 2020). This model is developed to simulate land-based physical processes in a general circulation model, and version six, MATSIRO6, involving up to three snow layers, six soil layers (14 m in total) and a single canopy layer, is used for interdisciplinary research on climate (MIROC6; Tatebe and others, 2019).

MATSIRO6 can simulate spatiotemporal changes in snow water equivalent, snow cover fraction and snow or ice surface temperatures. The snow water equivalent is simulated based on the water balance and is derived from the snowfall rate, snow sublimation, snowmelt and refreezing of rainfall and snowmelt. To express the presence or absence of snow cover at each study site, the areal fraction of snow cover (a value from 0 to 1) in a spatial grid of the simulation is determined from the calculated snow water equivalent. The snow and ice surface temperatures are simulated using the thermal conductivity equation. The detailed method used to calculate these physical properties is described by Takata and others (2003) and Nitta and others (2014). We used ice surface temperatures calculated using MATSIRO6 to quantify the ice-melting period at each site elevation.

To simulate the ice surface temperature, rainfall and snowfall amounts at the study sites, we conducted 1-D numerical simulations using MATSIRO6, with the reanalysis datasets described in this section as input data. One-dimensional simulations were conducted at all study sites (ablation areas) using the reanalysis datasets (WFDEI, GSWP3-FD and CRUJRA) described in the first paragraph of this section. The calculation period for the simulations was from 1 January 2011, to 31 December 2014. Simulations in 2011 were conducted to obtain the initial conditions of snow depth (in water equivalent) in 2012. To derive the rainfall and snowfall amounts at each site elevation, these amounts were corrected from the total precipitation data using the ratio of rain to snow, which was obtained using the surface air temperature, air pressure and specific humidity provided in MATSIRO6.

2.5. Model description and experimental design

2.5.1. Modification of the snow algae model (glacier algae model)

To reproduce the growth of glacier phototrophs on the bare ice surface, we applied a numerical model, referred to as the ‘snow algae model’, which utilizes a general differential equation of microbial growth to the observed algal growth curve (Onuma and others, 2018, 2020; Onuma, 2022). A previous study suggested that the exponential growth of S. nivaloides (snow algae) on snowpacks on the Qaanaaq Ice Cap can be expressed using a differential logistic growth equation as follows:

(1)$$\displaystyle{{{\rm d}X} \over {{\rm d}t}} = \left\{{\matrix{ {\mu X\left({1-\displaystyle{X \over K}} \right)\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{\rm if}\;T_{{\rm snow}} \ge 0\;[^{{\circ}}C ] } \cr {0\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{\rm if}\;T_{{\rm snow}} < 0\;[^{{\circ}}C ] } \cr } } \right., \;$$

where X represents the population density of microbes at growth period GP and μ indicates the hourly microbe growth rate. K denotes the ice surface alga carrying capacity. In this model, X increases when the snow surface temperature T _snow is above 0°C (as algal growth only occurs on the melting snow surface). When algae first appeared on the snow surface, X represents the initial population density, X ₀.

To adapt the snow algae model to the growth of glacier algae (A. nordenskioldii) on the bare ice surface in this study, we modified the snow surface temperature, T _snow, in Eqn (1) to the ice surface temperature, T _ice.

We refer to the model for A. nordenskioldii on a bare ice surface as the glacier algae model. In this model, the population density increases when T _ice is above 0°C (as glacier algal growth likely occurs on the melting ice surface). The cumulative hours exceeding T _ice of 0°C are defined as the algal growth period (GP) to quantify the duration required for the specific algal cell concentration. In the land model MATSIRO used in this study, GP can increase even with snow cover, because T _ice may reach 0°C when the snow layer on the ice surface is thin. In this study, the initial value of X at GP = 0 is represented as the initial bio-volume, X ₀, and it is assumed to be the surviving cells on the ice surface under a snowpack since the previous summer.

2.5.2. Biological parameters for the glacier algae model

None of the biological parameters, including initial bio-volume (X ₀), growth rate (μ) and carrying capacity (K), are known for glacier algae, and thus must be estimated prior to numerical simulations using the glacier algae model. X ₀ and μ were obtained by fitting Eqn (1) to the bio-volume of A. nordenskioldii observed in 2014 using the least square method after logarithmic conversion of the observed bio-volumes. The regression method has been applied to obtain biological parameters for snow algae in northwestern Greenland (Onuma and others, 2018). Regression was conducted using the bio-volume observed at S1 and S2 because the growth curve of the algae was observed at these sites (Section 3.1). In the regression, the maximum bio-volume observed for 3 years was assumed to be K. GP in 2014 was derived from T _ice, which was the simulation result from 1 January 2014, to 31 December 2014, using MATSIRO6 as described in Section 2.4. The reanalyzed datasets (WFDEI, GSWP3-FD and CRUJRA) were used to simulate T _ice as the model input data rather than the observed meteorological data. The three GPs were obtained using T _ice derived from the numerical simulations using each reanalysis dataset. Subsequently, biological parameters were obtained by fitting the observed bio-volume to the glacier algae model using each GP.

2.5.3. Model simulation setting for three seasons

In this study, the bio-volumes of A. nordenskioldii in 2014 were simulated and compared with the observational data at S1, S2, S3 and S4. The biological parameters and GP used in this simulation were obtained using the regression and MATSIRO6 simulations described in Section 2.5.2, respectively.

To further evaluate the glacier algae model obtained in this study, the bio-volumes of A. nordenskioldii in 2012 and 2013 were simulated and compared with the observational data. The simulation was conducted at S1 and S2 only, as there were no observational data at S3 and S4 in 2013. The biological parameters and GP were obtained in the same manner as those in the 2014 simulation. We assumed that the glacier algal abundance increased with the duration of ice melting until the end of the summer. To simulate algal growth in 2012 and 2013, we reset GP and bio-volume to zero and X ₀, respectively, when the snow cover fraction in the simulation using MATSIRO6 reached a maximum value of 1.0, indicating that snow completely covered the glacier surface.

3. Results and discussion

3.1. Differences in seasonal growth pattern between glacier algae and cyanobacteria

Microscopic analysis revealed that A. nordenskioldii was the dominant taxon across the bare ice surface of the Qaanaaq Ice Cap, similar to the southwestern GrIS where significant darkening has recently been reported (Yallop and others, 2012; Stibal and others, 2017; Lutz and others, 2018; Wang and others, 2018, 2020; Williamson and others, 2018; Cook and others, 2020). Microscopic observations identified three major taxa of algae (A. nordenskioldii, A. alaskana and S. nivaloides), and two major taxa of cyanobacteria (P. priestleyi and Chr. cyanobacterium) at the Qaanaaq Ice Cap study site (Table 1). These taxa were identified according to the shape, size and pigmentation of the taxon previously observed on this glacier in 2007 and 2012 (Uetake and others, 2010, 2016; Takeuchi and others, 2014). Seasonal changes in bio-volume showed that the proportion of A. nordenskioldii to the total bio-volume was more than 50% at all study sites throughout the growth season in 2014 (Fig. 3). The dominance of glacier algae (A. nordenskioldii or A. alaskana) on the Qaanaaq Ice Cap is consistent with the results of a previous report in 2007 (Uetake and others, 2010) and our observations in 2012 and 2013 (Fig. 4, top bar graphs). The cell length and width of A. nordenskioldii observed in the Qaanaaq Ice Cap were smaller than those observed in the southwestern GrIS, Patagonia and Switzerland-Italy (Takeuchi and Kohshima, 2004; Di Mauro and others, 2020; Williamson and others, 2020) and are similar to those observed in Svalbard and Siberia (Remias and others, 2012; Tanaka and others, 2016).

Fig. 3. Temporal changes in algal biomass, EC, pH and mineral abundance at each site from the top to the bottom. Legend for the top figures indicates each taxon observed at study sites. The dashed lines indicate that snow line reached the lower elevation sites (S1 and S2) on day 176 and the higher elevation sites (S3 and S4) on day 197.

Fig. 4. Inter-annual changes in the algal community, total algal biomass, electrical conductivity, pH and mineral abundance across the four study sites on the glacier. Legend in the top figure is the same as that in the top figure of Figure 3.

Table 1. Descriptions of glacier phototrophs observed in the Qaanaaq Ice Cap

The exponential increase in A. nordenskioldii abundance on the bare ice surface suggests that a differential logistic growth equation can be applied to reproduce seasonal bio-volume changes. The bio-volume of A. nordenskioldii exponentially increased with the duration of ice melting until the end of the observation period at S1 and S2 (Fig. 3). An exponential increase in algal abundance of A. nordenskioldii has also been observed previously on the bare ice surface of the GrIS (Stibal and others, 2017; Williamson and others, 2018). Meteorological conditions during the study period indicated that the ice surfaces continuously melted without snowfall events at S1, allowing for the growth of phototrophs (Supplementary Figure 1, top and bottom right panels). EC and pH often affect the growth of snow and glacier algae (Remias and others, 2012; Hoham and Remias, 2020); however, EC and pH in meltwater on the bare ice surface did not show seasonal trends that explained the change in algal growth in this study (Fig. 3, Supplementary Table 1). These observational results suggest that a major factor affecting the growth of A. nordenskioldii is the duration of ice melting. Continuous ice melting may supply nutrients for glacial algal growth. The algal abundance of S. nivaloides (snow algae) can be reproduced as a function of the snow melting duration using the differential logistic growth equation described in Section 2.5. Although glacier algae are a different species compared to snow algae, an increase in their population by cell division can generally be assumed to follow this equation during the growth stage. Therefore, this equation can likely be used to evaluate seasonal changes in A. nordenskioldii on a bare ice surface.

In contrast to the lower sites, no significant increase in A. nordenskioldii was observed at the upper sites S3 and S4 on the glacier, although the bare ice surface continuously melted after snow disappearance (Fig. 3), suggesting that the A. nordenskioldii abundance had already reached the carrying capacity at the time of snow disappearance at the sites. The durations of ice melting from day 197 (16 July) to day 215 (3 August), which were simulated using the land surface model with reanalysis data, were 392 ± 23 and 372 ± 47 h at S3 and S4, respectively. This result indicates the continuous melting of ice at the sites. However, there was no significant difference in the bio-volume at the sites during the observation period (Fig. 2). As previous studies suggest that resource limitation for glacier algal growth exists on the bare ice surface (Williamson and others, 2018), the bio-volumes at S3 and S4 might have already reached the carrying capacity on day 197, which is the first date for our sample collection at the sites. Glacier algae are likely to be able to grow on ice surfaces even under thin snow cover because their dark-colored cells absorb blue light and heat effectively (Dial and others, 2018). Model simulations showed that the duration of ice melting from day 182 (July 1) to day 215 was 663 ± 71 and 593 ± 76 h at S3 and S4, respectively, suggesting that the liquid water required for glacier algal growth had been supplied on the ice surface since day 182. Considering that the presence of liquid water is one of the limiting factors for snow algal growth (Ganey and others, 2017), glacier algae may be able to grow on the ice surface under favor of liquid water caused by the ice melting before the disappearance of snow cover.

In contrast to the growth pattern of A. nordenskioldii, that of cyanobacteria (P. priestleyi) was not exponential, suggesting that their seasonal change is distinct from that of glacier algae. At S1 and S2, the bio-volume of P. priestleyi increased only slightly from day 176 (25 June) to day 214 (2 August), although the bare ice surface melted continuously during this period (Fig. 3). In contrast, the bio-volumes of P. priestleyi at S3 and S4 drastically increased from day 197 to day 203 in 2014 and then showed no significant change until day 215 (Fig. 3), suggesting that the abundance of P. priestleyi cannot be simply determined by the duration of ice melting, as observed for A. nordenskioldii. Therefore, a different model from the logistic growth equation may be necessary to reproduce the seasonal growth of P. priestleyi.

When the bio-volume of P. priestleyi drastically increased at S3 and S4 after snow disappearance, the abundance of mineral dust also increased, possibly because of the redistribution of cryoconite granules on the bare ice surface induced by the collapse of cryoconite holes. Uetake and others (2016) suggested that the growth of P. priestleyi, which actively entangles mineral and organic particles in cryoconite holes, is not associated with nutrient conditions but with an abundance of mineral dust. Figure 3 shows the synchronization of the variation in the bio-volume of P. priestleyi and the abundance of mineral dust from day 197 (July 16) to day 203 (July 22) at S3 and S4. The cryoconite holes collapsed under cloudy conditions on days 206 and 208 in 2012 in the middle part of the ice cap near sites S3 and S4 (Takeuchi and others, 2018). As the reanalysis data showed similar meteorological conditions in late July 2014 (Fig. 2), the collapse of cryoconite holes likely occurred during this period, coinciding with the timing of the increase in P. priestleyi and mineral dust particles. These results indicate that the drastic increase in filamentous cyanobacterial abundance can be reproduced by projecting the collapse of cryoconite holes. A detailed discussion regarding numerical modeling of filamentous cyanobacterial growth is described in Section 3.3.

3.2. Numerical simulation of glacier algal growth on bare ice surface

3.2.1. Biological parameters for glacier algal growth

There was no significant difference in the initial bio-volume X ₀ between S1 and S2, as determined by regression analysis. The X ₀ obtained by fitting the glacier algae model was 9.7 × 10⁻³ and 2.7 × 10⁻¹ μL m⁻² at S1 and S2, respectively (Table 2). As the X ₀ at the sites are assumed to be the initial bio-volume under snow cover before ice melting, the comparison of the X ₀ with the bio-volume observed at the sites and reported by previous studies is difficult in this study. However, the bio-volumes observed at S1 and S2 immediately after snow disappearance (1.8 ± 3.1 × 10⁻² and 1.2 ± 1.8 × 10⁻¹ μL m⁻², respectively) are of similar order to the X ₀ at the sites. Further observations are necessary to quantify the initial algal abundance on bare ice under winter snow cover on the Qaanaaq Ice Cap and other glaciers.

Table 2. Biological parameters of the glacier algae model

Initial bio-volume and growth rate were obtained by the fitting observed bio-volume of Ancylonema nordenskioldii to the glacier algae model using growth periods. The maximum bio-volume of the algae at each site was assumed to be the carrying capacity. The growth periods were derived from T _ice simulated using a land surface model with three reanalysis data sets. The parameters derived from observed cell concentrations per volume (cells mL⁻¹) are shown in Supplementary Table 3.

The growth rate obtained by regression also did not differ significantly among the sites, suggesting that the growth rate of A. nordenskioldii was constant across the ice cap. The mean growth rate (μ) of A. nordenskioldii obtained using regression analysis was 0.010 h⁻¹ at both sites (Table 2). Previous studies suggested that the doubling time for their growth is 132.0 ± 40.8 h⁻¹ (Stibal and others, 2017) and ranges from 2.6 to 172.3 h⁻¹ in southwestern GrIS (Williamson and others, 2018). Doubling time DT converted from the growth rate μ in our study via DT = ln(2)/μ for S1 and S2 was 69.4 ± 3.7 and 82.4 ± 5.8 h⁻¹, respectively, which are consistent with the results of their study. However, further in situ and cell culture studies on the growth rate of this alga are required.

Although we assumed the maximum bio-volume observed at each site to be the carrying capacity in this study, further investigation is needed to reveal the major factor(s) affecting this parameter. A previous study suggested that the doubling time depends on nutrient limitations (Williamson and others, 2018). In a logistic growth equation, the effect of such a limitation on the growth rate is implicitly defined as carrying capacity (Onuma and others, 2018). As an increase in phosphorus, which is supplied from mineral dust on the ice surface, is likely to induce the growth of glacier algae in GrIS (McCutcheon and others, 2021), the concentration of phosphorus on the ice surface might be a key factor controlling the carrying capacity and maximum algal abundance during the growth season. The maximum abundance of glacier algae varies significantly between glaciers and ice sheets worldwide. The maximum cell abundance is 1.8 × 10⁵ cells mL⁻¹ in southwestern GrIS (Stibal and others, 2017), 2.4 × 10⁵ cells mL⁻¹ in Svalbard (Remias and others, 2012), 2.5 × 10⁵ cells mL⁻¹ in Siberia (Tanaka and others, 2016) and 4.7 × 10⁵ cells mL⁻¹ in northwestern GrIS (for this study, see Supplementary Table 3), respectively. In addition, the maximum abundance is likely to vary among elevations, even for the same glacier, as shown in Tables 2 and Supplementary Table 3. Further investigation to quantify the relationship between the maximum abundance of glacier algae and nutrient concentration, especially phosphorus, is necessary in the future.

3.2.2. Evaluation of the glacier algae model

The glacier algal simulations using biological parameters indicate that the model can reproduce the seasonal changes in bio-volume on the bare ice surface at the lower sites of the Qaanaaq Ice Cap in 2014. Figure 5 shows the seasonal changes in the observed and simulated bio-volumes of A. nordenskioldii at the study sites in 2014, using X ₀ and μ at S1, as shown in Table 2. The coefficients of determination in the regressions (R ²) with the periodical observation in 2014 were 0.93 (P < 0.05), 0.95 (P < 0.05), 0.71 (P > 0.05) and 0.66 (P > 0.05) at S1, S2, S3 and S4, respectively. The numerical simulations using X ₀ and μ at S2 in Table 2 also showed similar coefficients at the sites (Supplementary Figure 2), indicating that the glacier algae model can reasonably reproduce the growth of A. nordenskioldii at the lower sites (S1 and S2). Although heavy rainfall washed out algal cells on the bare ice surface in the southwestern GrIS during the 2014 season, as reported previously (Stibal and others, 2017), no heavy rainfall occurred on the ice cap (northwestern GrIS) during the melting period of this study (Supplementary Figure 1). As extreme rainfall events in the northwestern GrIS during the summer of the past decade have been obviously fewer than those in the southwestern GrIS, as suggested by Niwano and others (2021), the effect of heavy rainfall on the abundance of algal cells in the Qaanaaq Ice Cap may be less than that in the southwestern site under current climate conditions. Therefore, the loss of algal cells due to rainfall may have been insignificant in the present study.

Fig. 5. Seasonal changes in the algal abundance of Ancylonema nordenskioldii during summer in 2014 simulated using the glacier algae model. (a) S1, (b) S2, (c) S3, (d) S4. Each solid line and shade indicate averaged and maximum or minimum bio-volume simulated with the glacier algae model using the different meteorological conditions, respectively. The initial bio-volume and growth rate at each site are values obtained by the fitting the glacier algae model to the observation at S1. The maximum bio-volumes at each site for three seasons are assumed to be the carrying capacity for (a), (b), (c) and (d).

The numerical simulations using the glacier algae model suggested that the initial abundances of A. nordenskioldii at the higher sites were greater than those at the lower sites of the Qaanaaq Ice Cap. The bio-volume of A. nordenskioldii simulated using X ₀ at S2 in Table 2 agreed with that observed at S3 and S4 rather than that simulated using X ₀ at S1 in Table 2 (Supplementary Figure 2 and Fig. 5, respectively). However, the simulated bio-volumes shown in Supplementary Figure 2 still underestimated the observed values. These results suggest that the initial abundance of A. nordenskioldii on the ice surface at the higher sites was greater than that at the lower sites (S1 and S2). At higher elevation sites, glacier algal cells on the ice surface may be less washed out by meltwater during the snow-melting season because the ice surfaces are less likely to melt. Indeed, the bio-volume observed at S2 (441 m a.s.l.) immediately after snow disappearance was greater than that observed at S1 (247 m a.s.l.). In addition, meltwater, heavy rainfall and collapse of cryoconite holes may redistribute cells on the ice surface or supply those found downstream of the glacier. Although it is difficult to quantitatively discuss the horizontal movement of glacier algal cells in this study, further field observations to quantify the spatial distribution of glacier algal abundance are necessary. In addition, the effect of horizontal movement on algal abundance should be incorporated into a numerical model for glacier algal growth in the future.

The bio-volume in late July 2012 in the simulation agreed well with the observations at the study sites, suggesting that the glacier algae model can reasonably reproduce the bio-volume over multiple seasons on the Qaanaaq Ice Cap. Model simulation using the three reanalysis datasets of meteorological conditions showed that the algal bio-volume at S1 was 19.8–74.5 μL m⁻² (mean: 46.0 μL m⁻²) on day 204 (22 July) of 2012 and was close to the observed values (Fig. 6). Similarly, the bio-volumes at S2 in 2012 in the simulation agreed with the observations. However, the model simulation significantly overestimated the bio-volume in 2013, particularly at site S1. The reanalysis datasets showed that snowfall frequently occurred from days 180 to 210 in 2013 (Supplementary Figure 1). The shade in Figure 6 shows that the uncertainty of the bio-volume derived from the model simulation was larger in 2013 than in 2012, suggesting that frequent snowfall events during summer caused an increase in the uncertainty of the simulated algal abundance. In contrast to 2013, minimal snowfall occurred after the ice surface was exposed in 2012. The meteorological conditions were similar to those of 2014. To evaluate the model performance, satellite observations, which can estimate glacier algal abundance based on spectral light absorption by them (Di Mauro and others, 2020; Bohn and others, 2022), are necessary, as well as further periodical field observations over multiple seasons. However, our results suggest that this model can reproduce A. nordenskioldii blooms during seasons with little snowfall, which is not in every season.

Fig. 6. Seasonal changes in snow water equivalent (top), growth period of Ancylonema nordenskioldii (middle) and bio-volume of the algae (bottom) simulated using the glacier algae model at S1 (left) and S2 (right) during summer in 2012 and 2013 seasons. Each solid line and shade indicate averaged and maximum or minimum values simulated with the land surface model or glacier algae model using the different meteorological conditions, respectively. Black marks indicate observed bio-volume.

Numerical simulations using the glacier algae model suggest that the difference in the abundance of A. nordenskioldii among the three seasons can be explained by the duration of ice melting. The bio-volume of A. nordenskioldii in 2012 was the greatest among the three seasons at all study sites (Fig. 4). As discussed in Section 3.1, the abundance of A. nordenskioldii likely increased exponentially on the bare ice surface during the melting season. The GP in the simulation indicated that the duration of ice melting was significantly longer in 2012 than in 2013 (Fig. 6). The duration of ice melting each year may indicate the growth of glacier algae during the season. As the duration can be simulated using global or regional climate models, interannual trends in glacier algal abundance can be estimated from the past to the future.

Although the numerical simulations agreed well with the observations, further investigation of the biological and physical processes in the glacier algae model is required. We assumed that the simulated algal bio-volume was reduced to its initial abundance when the summer ended each year. Field observations support this assumption, as the bio-volumes observed at S1 and S2 in 2013 were significantly lower than those in 2012 (Fig. 6). Most A. nordenskioldii cells were in the vegetative form, which is intolerant to freezing, and thus cannot survive on the glacier during winter. Although they typically form zygotes (pre-akinete-type cells) under conditions inadequate for growth, the zygote of A. nordenskioldii was rarely observed on the ice surface (Remias and others, 2012). Therefore, the abundance of algae at the end of the growth season may have been inherited in the next season. The formation of pre-akinete type cells at the end of summer may be important for determining the initial cell abundance in the next melting season. However, further investigation of surviving glacier algal cells during winter is necessary to quantify their initial abundance, because there is little information on the cells during winter. In addition to such biological processes, glacier algae cells may move on bare ice surfaces or in gaps of ice in weathering crust with meltwater. However, the model does not include the horizontal displacement of algal cells on the glacier surface. Furthermore, hydrological and topographical processes may change with time and space, as described by the concept of bio-cryomorphology (Cook and others, 2015; Irvine-Fynn and others, 2021). Thus, further biological and physical processes should be incorporated into glacier algae models.

3.3. Empirical modeling of cyanobacterial abundance on bare ice surface

3.3.1. Glacier filamentous cyanobacteria model

As the increase in P. priestleyi was not exponential, the abundance of cyanobacteria appeared to be reproduced empirically from the abundance of mineral dust, rather than using the logistic growth equation for glacier algae. As described in Section 3.1, the abundance of P. priestleyi was not associated with the duration of ice melting, suggesting that the abundance of cyanobacteria was limited by growth conditions; that is, it had already reached the carrying capacity. The relationship between the abundance of mineral dust and P. priestleyi at all study sites in 2014 is shown in Figure 7a and exhibited a significant positive correlation with the polynomial approximation (R ² = 0.88). Based on this relationship, the following regression curve was obtained:

(2)$$X_{\rm P} = 0.017{\rm M}{\rm D}^2 + 0.068{\rm MD}, \;$$

where X_P and MD are the bio-volumes of P. priestleyi (μL m⁻²) and mineral dust weight (g m⁻²) on the bare ice surface, respectively. This empirical model (glacier filamentous cyanobacteria model) reproduces the abundance of cyanobacteria from mineral dust abundance on the bare ice surface.

Fig. 7. Relationship between mineral dust weight and Phormidesmis priestleyi bio-volume. (a) Fitting result with observed bio-volume in 2014. (b) Comparison of the bio-volume simulated using the glacier filamentous cyanobacteria model (Eqn (2)) with that observed in 2012 and 2013.

3.3.2. Evaluation of the glacier filamentous cyanobacteria model

The glacier filamentous cyanobacteria model can estimate yearly trends in cyanobacterial abundance. In this section, we evaluate the applicability of our empirical model (Eqn (2)) to the 2012 and 2013 seasons using observed mineral dust data. Figure 7b shows that the simulated cyanobacterial abundance underestimated the observed value, particularly for higher mineral dust weights. However, the calculated cyanobacterial abundance was consistent with that observed in 2012 and 2013 (R ² = 1.0), indicating an obvious relationship between cyanobacterial abundance and mineral dust weight. A previous study suggested that small mineral particles with diameters of 30–249 μm are bases or footholds for the growth of filamentous cyanobacteria because they provide a source of phosphate and are generally less influenced by meltwater (Uetake and others, 2016). In contrast, filamentous cyanobacteria may also play a role in avoiding the washing of mineral particles from the ice surface. Thus, the relationship between the abundance of cyanobacteria and mineral particles was interactive. Although the quantitative relationship between cyanobacteria and mineral dust is not well understood, the mineral components in cryoconite granules may also affect the growth of P. priestleyi.

The transportation of mineral dust through the dynamics of cryoconite holes onto a bare ice surface may be a key process in reproducing seasonal changes in the abundance of P. priestleyi. Previous studies suggested that the mineral dust on GrIS was mainly supplied from the local moraine via the atmosphere, and its deposition was significantly correlated with the snow cover duration on soil in the coastal area of Greenland (Nagatsuka and others, 2014, 2016, 2021). Furthermore, mineral dust can also be derived from outcropping englacial ice with its melting (Wientjes and Oerlemans, 2010; Goelles and Bøggild, 2017). In contrast, meltwater running on the ice surface can wash mineral particles out of the glacier (Goelles and Bøggild, 2017). Therefore, the mineral dust abundance on the ice surface depends on the balance between supply and loss on the ice. According to these studies, the abundance of mineral dust on the bare ice surface is likely to increase gradually during the ice-melting season after the disappearance of snow in the coastal glacier of Greenland (Qaanaaq Ice Cap). However, observations revealed that dust abundance drastically changed over time on the bare ice surface (Fig. 3), likely because of the collapse and formation of cryoconite holes (Takeuchi and others, 2018). Therefore, to simulate seasonal changes in P. priestleyi on a bare ice surface, it is necessary to quantify seasonal changes in the abundance of mineral dust while considering the dynamics of cryoconite holes.

4. Conclusions

We investigated seasonal and inter-annual changes in the glacier phototroph community over 3 years (2012–2014) on an ice cap in northwestern Greenland. Field observations revealed that the abundance of A. nordenskioldii, a dominant algal species on this glacier, exponentially increased with the duration of ice melting. The exponential increase in snow algae (S. nivaloides) can be reproduced using a logistic growth equation (referred to as the snow algae model) proposed based on snowfield observations in the same ice cap (Onuma and others, 2018). Hence, we applied the same logistic growth equation to the growth of glacier algae on a bare ice surface using the duration of ice melting as the input data for algal simulation. The numerical glacier algae model showed that the simulated algal abundance agreed well with observations from 2012 and 2014, when minimal snowfall occurred. The results indicate that the logistic growth equation can reproduce algal growth on bare ice and snow surfaces. In the case of P. priestleyi (filamentous cyanobacteria), the increase in cyanobacterial abundance likely synchronizes with increases in the mineral dust weight on the bare ice surface, suggesting that cyanobacterial abundance can be reproduced using an empirical model based on mineral dust abundance on the ice surface. Recently, many studies have reported that albedo reduction of the ice surface on GrIS and other Arctic glaciers is caused by blooms of glacier phototrophs (Wientjes and others, 2011; Shimada and others, 2016; Cook and others, 2020; Wang and others, 2020). As albedo reduction highly depends on phototroph abundance (Cook and others, 2020; Williamson and others, 2020), projecting spatiotemporal changes in glacier phototrophs abundance is important for evaluating the mass balance of the GrIS and glaciers. Our study provides a foundation for evaluating the mass balance using numerical models for the growth of glacier algae and cyanobacteria.