Glaciers and ice sheets currently cover ~10% of the surface of continents and contain >30 million km³ of ice. They contain distinct ecosystems that harbor diverse microbial communities and are places of significant biological activity (Hodson et al., 2008). The largest supraglacial (glacier surface) ecosystem on Earth is found on the Greenland ice sheet (GrIS). This is due to its high surface melting that may reach nearly 100% of the surface during extreme events (Nghiem et al., 2012). Supraglacial ecosystems in Greenland and elsewhere are changing rapidly due to climate warming, which causes retreat of the margins of glaciers and ice sheets and a potential inward expansion of the biologically active ablation areas.

Microorganisms in supraglacial ecosystems cycle carbon and nutrients (Hodson et al., 2008; Stibal et al., 2012a), and may also affect the physical behavior of glaciers by increasing melt via lowering the surface reflectivity of ice (Takeuchi et al., 2015; Lutz et al., 2016; Musilova et al., 2016). Supraglacial ecosystems may significantly impact neighboring terrestrial and marine ecosystems via meltwater export of live cells, organic carbon (OC), and other nutrients (Lawson et al., 2014; Hawkings et al., 2015; Cameron et al., 2017).

Insights into the supraglacial ecosystem have been driven mostly by empirical approaches relying on field sampling and laboratory measurements. However, these are often constrained by patchy spatial and temporal coverage, limiting their appropriateness for upscaling. Mathematical models are ideally suited to studying the supraglacial ecosystem since they can disentangle and quantify the interplay of various biogeochemical and physical processes, bridge the vast spatial scale of ice sheets, make scenario-based predictions, identify gaps in the current understanding, and thus help design future sampling and laboratory analyses, as reviewed by Bradley et al. (2016a). Despite increasing data on sources, sinks, and transformations of carbon and nutrients in the supraglacial ecosystem, few attempts of linking them together into an ecosystem model have been made to date (e.g., Hodson et al., 2010; Cook et al., 2012). As a result, estimates of microbial activity and associated carbon and nutrient transformations on a large scale are highly uncertain, and predictions of future ecosystem change are virtually impossible.

The aim of this perspective paper is to provide a theoretical framework of the supraglacial ecosystem in order to facilitate ecological modeling as a tool for understanding present day and future ecosystem dynamics. To do so, we focus on the largest supraglacial ecosystem on Earth—the GrIS surface—as an important example and present three conceptual models of the GrIS microbial ecosystem. However, these models are not Greenland-specific and can be applied to other supraglacial ecosystems. The models focus on organic carbon transformations but differ in complexity by the detail that biological processes are represented. We discuss the strengths and weaknesses of each modeling approach, their role in improving the understanding of carbon cycling and ecosystem dynamics, and the necessary integration of these modeling approaches with field data to parameterize and validate numerical output. This perspective is intended to guide future supraglacial ecosystem model development and field data collection, and encourage cross-disciplinary collaboration between modelers and experimentalists.

The supraglacial ecosystem comprises the top layer of the ice which is in contact with the atmosphere and receives solar radiation and atmospheric deposition of dust, aerosols, and microbial inocula. Three distinct but inter-connected habitats host living organisms which contribute to biogeochemical activity in the system: melting snow, bare ice, and surface debris (cryoconite). Meltwater, a key constraint on biological activity on glacier and ice sheet surfaces, is present in all three habitats during the melt season and also provides a transient habitat through which viable microorganisms can be transported.

Melting snow is an ephemeral habitat colonized by diverse microbial communities, including photoautotrophic snow algae (Uetake et al., 2010; Cameron et al., 2015). Microbial cells deposited with snow may remain on the ice surface after the snow has melted and become part of the surface ice community, or they can be flushed from the ice sheet surface with meltwater (Cameron et al., 2015). Bare ice is exposed seasonally around the margins of the ice sheet and hosts a high abundance of algae (Uetake et al., 2010; Yallop et al., 2012) and other microorganisms (Stibal et al., 2015a; Cameron et al., 2016). The algae are important primary producers in the ecosystem and may contribute to surface melting via darkening the ice due to pigment production (Cook et al., 2012; Yallop et al., 2012). Abundance of other, mostly heterotrophic, microbes in surface ice has been found to correlate with dust concentration (Stibal et al., 2015a). Cryoconite, usually concentrated in cryoconite holes (Hodson et al., 2010), hosts diverse and highly active microbial communities that consist of photoautotrophic cyanobacteria and a range of heterotrophic bacteria (Cameron et al., 2012, 2016; Stibal et al., 2012b, 2015b; Edwards et al., 2014; Uetake et al., 2016). Cryoconite holes provide a stable and nutrient-rich habitat within the supraglacial ecosystem, and are considered hotspots of microbial activity on glacier surfaces (Anesio et al., 2009; Cook et al., 2016a). Microbial activity in cryoconite may also contribute to surface darkening (Takeuchi et al., 2015; Musilova et al., 2016). The interconnected nature of these three habitats allows for transportation of biomass and nutrients between them over the course of an annual cycle (Cameron et al., 2016).

There is a wealth of data on the principal carbon cycling processes (primary production, secondary production, respiration) from the supraglacial ecosystems of the GrIS and other glaciers and ice sheets (e.g., Bagshaw et al., 2007, 2016; Foreman et al., 2007; Hodson et al., 2007, 2010; Stibal et al., 2008, 2012b; Anesio et al., 2009, 2010; Telling et al., 2010; Cook et al., 2012, 2016b; Bellas et al., 2013; Chandler et al., 2015; Rassner et al., 2016; Smith et al., 2016). However, data on rates of microbial exudation (biotic release of dissolved organic carbon (DOC) from living cells), decomposition of particulate organic carbon (POC), and cell mortality on glacier surfaces are currently lacking.

Numerous environmental controls on supraglacial microbial activity have been identified (as reviewed in Hodson et al., 2008; Anesio and Laybourn-Parry, 2012; Stibal et al., 2012a). First, liquid water, a key prerequisite for biological activity, is generated during the ablation season by snowmelt and ice melt. Data on meltwater quantity can be obtained from existing surface mass balance (SMB) models (e.g., Fettweis et al., 2013; Langen et al., 2017, for the GrIS). Second, sunlight provides energy for phototrophic growth, and melts snow and ice producing the water needed to support biological activity. Wintertime biological activity has not been measured in supraglacial ecosystems to date, however due to almost complete absence of light and liquid water during polar night, most biological processes are likely to stop or decline considerably. Third, the availability of organic carbon and dissolved inorganic nutrients is known to limit supraglacial microbial activity. Nitrogen is supplied to the ice surface in significant amounts by precipitation (e.g., Fischer et al., 1998; Geng et al., 2014). However, it may be limiting in glacier and ice sheet marginal areas, where microbial nitrogen fixation has been measured (Telling et al., 2011, 2012). Phosphorus, a critical macronutrient derived mainly from rock weathering, is likely to be a key factor limiting biological activity in the supraglacial ecosystem (Stibal et al., 2009).

In order to provide a theoretical framework of the supraglacial ecosystem of the GrIS and other glaciers and ice sheets, we present three process-based conceptual models of increasing complexity (Figure 1). The ultimate choice of model complexity will depend on the nature of the research question(s), the scale that is being resolved, the existing knowledge and data, and the potential of new knowledge gained. Each unique scientific question may have several best-possible approaches integrating models and data (Bradley et al., 2016a). Insight can be gained by exploring multiple approaches and developing models to varied levels of detail (Sierra and Müller, 2015). Here, we focus on a process-based approach and present the conceptual models as suggestions, intended to illustrate the varying levels of complexity that may be incorporated into each model formulation. Other modeling approaches are discussed in Bradley et al. (2016a). First, we discuss the major OC fluxes individually (see Table 1), before describing how they are applied and combined in each model formulation (Figure 1).

The foci of the presented models are the sources, sinks, and transformations of OC and biomass (Table 1). OC reservoirs and fluxes form the basis of these models, since most biological processes are in essence OC transformations, and much of the data available from field studies which may be used to inform and validate model predictions are presented as such (e.g., Hodson et al., 2010; Stibal et al., 2012b; Yallop et al., 2012; Chandler et al., 2015). Transformations of other elements such as nitrogen and phosphorus may be incorporated into model formulations using fixed or variable stoichiometries (e.g., Soetaert and Herman, 2009; Bradley et al., 2015), provided that appropriate knowledge and empirical/experimental data exist to inform and validate predictions. Conventional units in which data from the supraglacial ecosystem are reported are usually related to surface area (cells cm⁻², g C m⁻²; e.g., Cook et al., 2012; Yallop et al., 2012). Data reported as chemical concentrations or weight per volume of snow/ice (e.g., Chandler et al., 2015) can be converted to area by integrating over the “active” surface layer, while data from cryoconite, usually expressed per gram of debris, can be converted to area using measurements of surface cryoconite coverage (Bøggild et al., 2010; Hodson et al., 2010; Stibal et al., 2012b; Takeuchi et al., 2014).

There are three abiotic OC fluxes resolved within the suggested framework. Atmospheric deposition (W) of POC and microbial biomass is an important OC source for the GrIS supraglacial ecosystem (Edwards et al., 2014; Stibal et al., 2015b; Cameron et al., 2016). Direct measurements of W are challenging; however, first successful attempts from an Arctic glacier showed a cell flux of 10⁷ cells per m² per day (Irvine-Fynn et al., 2012). Melting of POC out of old ice (M) is a well-known phenomenon from the GrIS surface (Wientjes et al., 2012), the flux of which is dependent on the OC concentration in the ice, which can be directly measured, and the amount of ice melt, which can be obtained from SMB models (Fettweis et al., 2013; Langen et al., 2017). Meltwater flushing (F) is also dependent on surface melt rate, which can be obtained from SMB models, and on the concentration of OC in the meltwater leaving the ice sheet surface, which has been measured (Lawson et al., 2014).

The following biotic fluxes are resolved within the suggested framework. Net primary production (G_a), or growth of autotrophic microbes (bio_a), is an important source of OC on glacier surfaces, and has been shown to be in the range of 10⁻⁴–10⁻² g C m⁻² d⁻¹ for supraglacial habitats (Hodson et al., 2007; Anesio et al., 2009; Telling et al., 2010; Stibal et al., 2012b; Yallop et al., 2012). G_a is dependent on the abundance of BIO_a, the incident photosynthetically active radiation, and the concentration of limiting nutrients (Yallop et al., 2012; Bagshaw et al., 2016), and can be numerically formulated using e.g. Monod kinetics (Soetaert and Herman, 2009; Bradley et al., 2016b). Secondary production, or heterotrophic growth (G_h), represents the consumer component of the ecosystem (BIO_h). BIO_H consists of a range of heterotrophic microorganisms in all three environments. Due to the nature of the supraglacial environment, these are likely to be mostly aerobic chemoheterotrophs, i.e., oxidizing OC compounds by atmospheric oxygen to gain energy while also using OC as a carbon source. Typical rates of secondary production in supraglacial environments are between 10⁻⁷ and 10⁻⁴ g C m⁻² d⁻¹ (Hodson et al., 2007; Anesio et al., 2010; Bellas et al., 2013). Monod kinetics may also be appropriate for formulating G_H, since heterotrophic growth on GrIS is likely dependent on the concentration of DOC and limiting nutrients (Anesio et al., 2010). Respiration (R) is an important biotic sink of OC for the supraglacial ecosystem of the GrIS, with field measurements in the range of 10⁻⁴ and 10⁻² g C m⁻² d⁻¹ (Anesio et al., 2009; Telling et al., 2010; Stibal et al., 2012b). Estimates of bacterial growth efficiency (e.g., Anesio et al., 2010) may be used to explicitly relate bacterial production to respiration (e.g., Bradley et al., 2015). Similarly, exudation (E) can be modeled as a fraction of microbial growth (Bradley et al., 2015). Decomposition of POC (L) can be modeled using Monod kinetics as dependent on the proportion of BIO_H capable of decomposing POC, and POC concentration. Cell death rates (D) are problematic to define experimentally (Toal et al., 2000), and may be simplified as a density-dependent fraction of BIO_A and BIO_H (Bradley et al., 2015).

The simplest model (M1; Figure 1) adopts a simplistic approach and only resolves dissolved OC (DOC) and particulate OC (POC). POC comprises microbial active and inactive biomass and various types of organic residue. Deposition of wind-borne organic materials (W_poc), melting of POC out of old ice (M), autotrophic growth (G_a), and heterotrophic growth (G_h) contribute to the POC pool. G_A comprises net primary production (i.e., CO₂ fixed from the atmosphere during photosynthesis less CO₂ respired back) of snow algae, ice algae, and cryoconite cyanobacteria, whereas G_H serves as the conversion of DOC to biomass (POC) during heterotrophic growth. POC is depleted by decomposition (L), respiration (R), bacterial exudation (E), and meltwater flushing (F_poc). The DOC pool is depleted by heterotrophic uptake (G_H) and meltwater flushing (F_doc), while decomposition of POC (L) and exudation (E) provide a source of DOC. This approach may be useful for testing hypotheses concerning the GrIS supraglacial ecosystem as net carbon sink or source using data derived from field measurements of carbon fluxes, similar to prior upscaling attempts (Hodson et al., 2010; Cook et al., 2012; Chandler et al., 2015). If fluxes are formulated such that they are sensitive to perturbations in environmental conditions, the model could address such questions as: How will future climate warming affect the OC reservoir in the GrIS surface ice? However, biological processes are implied, rather than explicitly simulated, and insights into the biotic drivers of ecosystem dynamics are limited.

The next model (M2; Figure 1) adopts a more complex process-focussed approach, and distinguishes biomass (BIO) from non-living organic matter (POC). BIO captures living microbes in all three supraglacial environments (melting snow, bare ice, and surface debris). Abiotic fluxes are adjusted from M1 as follows: wind deposition of OC is split between POC and BIO (W_poc, W_bio), meltwater flushes material from all OC pools (F_doc, F_poc, F_bio), and no live biomass is assumed to melt out of old ice. For the biotic component, G_a, G_h, R, and E apply to the BIO pool. Cell death (D) results in the flux of OC from BIO to POC. This level of complexity may be suitable for answering questions related to the overall quantity and potential activity of microbes within the system. For example, will microbial biomass at the surface of the GrIS increase with enhanced surface melting? However, microbial functional diversity is not represented in sufficient detail to resolve the different controls of the principal carbon cycling processes, since all microbes are lumped into a single state-variable (BIO).

Finally, M3 (Figure 1) adopts a process-focussed approach, but differs from M2 by explicitly resolving the main functional groups of biomass, depending on principal carbon cycling pathways (autotrophic BIO_a and heterotrophic BIO_h). For abiotic fluxes, wind deposition of OC is split between POC, BIO_A, and BIO_H (W_poc, W_a, W_h, respectively), and all pools lose OC to meltwater flushing (F_doc, F_poc, F_a, F_h). For biotic fluxes, G_h and R only apply to the BIO_H pool, while G_a applies only to BIO_A. D and E are separate for autotrophs (BIO_A) and heterotrophs (BIO_H). A model of this level of complexity contains the basic components of the carbon cycle at the GrIS surface and may be useful in testing hypotheses concerning the impacts of various factors on carbon cycling processes such as primary production and respiration. For example, this model may provide answers to questions such as: May increased surface melting on the GrIS cause the system to switch between net autotrophy (carbon sink) and net heterotrophy (carbon source)? However, the model will require more robust and comprehensive datasets to inform and support predictions to an acceptable level of confidence.

The conceptual models presented in Figure 1 merge all the supraglacial habitats described above (snow, surface ice, and cryoconite). If hypotheses about specific components are to be tested, each component can be resolved individually, possibly at the expense of increased model complexity.

Field measurements and controlled experiments provide useful data that can inform model values, including initial conditions and parameters. Environmental controls on OC fluxes can then be included mathematically in a number of ways, as reviewed by Bradley et al. (2016a). If insufficient experimental data exist to directly inform parameters, their value may be estimated based on calibration to data from similar process-based models, or measurements from comparable low-temperature ecosystems such as Arctic soils or lakes (e.g., Jones et al., 2000; Mindl et al., 2007; Hollesen et al., 2011; Bradley et al., 2016b). Sensitivity analysis may then be used to determine the relative importance of obtaining an accurate value for each parameter, thereby informing future experimental work.

Forcing data can be obtained from existing monitoring networks and regional climate models. For example, the Programme for Monitoring the Greenland Ice Sheet (PROMICE) provides surface energy budget closure via automated weather stations distributed across the ice sheet ablation area (Ahlstrøm et al., 2008; van As et al., 2016). Models such as MAR (Fettweis et al., 2013) and HIRHAM5 (Langen et al., 2017) are useful to provide parameters influencing biological activity, such as surface temperature and precipitation (and so the availability of liquid water), under different climate scenarios. Thus, they may provide estimates of the spatial and temporal extensions of the biologically active area on the ice sheet and so enable predictions of the fate of the GrIS supraglacial ecosystem under a number of scenarios.

Spatial variability in carbon and nutrient concentrations and microbial abundance and activity has been reported from the GrIS (e.g., Hodson et al., 2010; Stibal et al., 2010, 2012b, 2015a; Telling et al., 2012; Yallop et al., 2012). However, direct measurements exist from very few discreet locations, and upscaling from point measurements to the whole ice sheet may introduce large errors. Therefore, it is crucial to identify processes and parameters for which spatial variability is a concern. For example, a small-scale zero-dimensional model simulating a unit surface area (e.g., a SMB model grid point) will be relatively easy to constrain in terms of forcing and validation data, thereby increasing confidence in model simulations. Under the assumption of homogeneity across the entire GrIS, estimates of whole ice sheet scale fluxes could be derived from a single 0D model simulation. Such a model will be tractable but would require substantial generalizations, thus possibly introducing errors. Alternatively, whole ice sheet scale fluxes can be derived from independently run model simulations with spatially discreet initial conditions, forcings, and parameter values, or by including transport terms between model grid points informed by existing ice flow models (e.g., Goelles et al., 2015).

Model validation is crucial for determining that model dynamics accurately represent the conceptual description and specifications, and presents an opportunity for collaboration between modelers and experimentalists. Where appropriate, model results may be validated using data from a short/small scale observations, upon which meaningful upscaling might provide robust predictions. Such data for the GrIS comprise measurements of microbial abundance and activity, and OC and nutrient concentrations, amongst others. Future efforts to provide useful field data to inform and validate the next generation of supraglacial ecosystem models should be a priority.

Modeling the supraglacial ecosystem will provide insights into present day dynamics going beyond simple extrapolation of in situ measurements, and enable the prediction of future ecosystem change in a warming climate. Iterative model development will also generate new hypotheses and improve sampling and experimental designs for in situ work. Exploring models of varying complexity and spatial and temporal scales, as presented here, will help answer specific questions of interest to the scientific community, including accurate quantification of OC fluxes, and the future fate of carbon and biomass on the GrIS surface and elsewhere. Close collaboration between modelers and experimentalists will be key to designing appropriate models and collecting useful data. Thus, ecological modeling forms a valuable contribution to the mosaic of cryosphere research currently underway.

All authors contributed to writing of the manuscript.