Multi-Source Remote Sensing Based Modeling of Vegetation Productivity in the Boreal: Issues & Opportunities

Abstract

Understanding the processes driving terrestrial vegetation productivity dynamics in boreal ecosystems is critical for accurate assessments of carbon dynamics. Monitoring these dynamics typically requires a fusion of broad-scale remote sensing observations, climate information and other geospatial data inputs, which often have unknown errors, are difficult to obtain, or limit spatial and temporal resolutions of productivity estimates. The past decade has seen notable advances in technologies and the diversity of observed wavelengths from remote sensing instruments, offering new insights on vegetation carbon dynamics. In this communication, we review key current approaches for modeling terrestrial vegetation productivity, followed by a discussion on new remote sensing instruments and derived products including Sentinel-3 Land Surface Temperature, freeze & thaw state from the Soil Moisture & Ocean Salinity (SMOS) mission, and soil moisture from the Soil Moisture Active Passive (SMAP) mission. We outline how these products can improve the spatial detail and temporal representation of boreal productivity estimates driven entirely by a fusion of remote sensing observations. We conclude with a demonstration of how these different elements can be integrated across key land cover types in the Hudson plains, an extensive wetland-dominated region of the Canadian boreal, and provide recommendations for future model development.

RÉSUMÉ

Il est essentiel de comprendre les processus régissant la dynamique de la productivité de la végétation terrestre des écosystèmes boréaux pour pouvoir évaluer avec précision les flux de carbone. Le suivi de cette dynamique nécessite la fusion d‘observations de télédétection à grande échelle, d‘informations climatiques, et d‘autres données géospatiales, lesquelles comportent des erreurs inconnues, sont difficiles à obtenir, et limitent les résolutions spatiales et temporelles des estimations de la productivité. Au cours de la dernière décennie, des progrès technologiques notables ont diversifié les longueurs d‘onde observées par les instruments de télédétection, ce qui a permis de mieux comprendre la dynamique du carbone liée à la végétation. Dans cette communication, nous passons en revue les principales approches pour la modélisation de la productivité de la végétation terrestre, suivies d‘une discussion sur les nouveaux instruments de télédétection et leurs produits dérivés, notamment la température de la surface terrestre extraite de Sentinel-3, l‘état de gel ou de dégel du sol dérivé de la mission SMOS, et l’humidité du sol dérivée de la mission SMAP. Nous décrivons comment ces produits peuvent améliorer le détail spatial et la représentation temporelle des estimations de la productivité boréale entièrement basées sur une fusion d‘observations de télédétection. Nous concluons par une démonstration de la façon dont ces différents éléments peuvent être intégrés dans les principaux types de couverture terrestre des plaines hudsoniennes, une vaste région de la forêt boréale canadienne dominée par des zones humides, et nous formulons des recommandations pour de futures améliorations du modèle.

Introduction

Forests cover approximately 31% of the world’s terrestrial area (Keenan et al. 2015), sequestering ∼30% of annual anthropogenic carbon emissions from the atmosphere (Canadell and Raupach 2008; Xu et al. 2019). Between the tropical, subtropical, temperate, and boreal biomes, the boreal is the most northern and largest terrestrial biome in the world (Brandt et al. 2013; Douglas et al. 2014). Canada’s boreal zone accounts for 28%, or 552 Mha, of the global boreal ecosystem, and is characterized by a variety of forest stand ages, complex disturbance histories, a small number of dominant species, and high-latitude continental climatic conditions (Brandt et al. 2013; Kurz et al. 2013; Price et al. 2013). With regional variation in species distribution, their productivity, and underlying soils, the total boreal ecosystem is estimated to sequester an average of 8850 TgC yr⁻¹ (1E + 9 Kg of Carbon per year)^, 62% of which is during the growing season (Virkkala et al. 2021). Due to these complex heterogenous patterns and drivers of land cover and land cover change, quantifying the rate of carbon sequestration in these boreal environments is challenging. Overcoming these challenges is imperative for the continued improvement of robust boreal vegetation productivity models and the interannual estimates of regional carbon sequestration they provide. These improved productivity estimates and their impact as a key policy driver forms the motivation of this communication.

The boreal is a patchwork of fine scale vegetation structure and composition reflective of the heterogeneous climate, moisture regimes, permafrost and natural & anthropogenic disturbances (Bonan and Shugart 1989; Chapin et al. 2010). Boreal ecosystem function and structure is primarily driven by climate and weather through temperature, precipitation, and disturbances such as wind, fire, and insect & disease infestation (Brandt et al. 2013; Friend et al. 2014). In the boreal forests of North America, the majority of the characteristically modest precipitation falls in the brief summer period, whilst the cold & long winters are typically dry (Bonan and Shugart 1989; Brandt et al. 2013). The boreal forest is typically underlain by low relief terrain, acidic and nutrient-poor soils, and discontinuous permafrost, with continuous permafrost present north of the treeline (Smith et al. 2010). Permafrost regions may be defined as “substrate remaining below 0 °C for more than 2 consecutive years” (Douglas et al. 2014) and account for approximately 388 Mha, or 70%, of the total Canadian boreal zone (Kurz et al. 2013). Notably, permafrost constrains the movement of water and nutrients to the active layer, which freezes and thaws annually, strongly impacting soil and vegetation processes linked to the carbon cycle (Douglas et al. 2014; Petrone et al. 2006). Changing seasonal dynamics and increased permafrost thaw under warming climate scenarios may thus lead to alterations in the distribution, function and structure of Canada’s northern forests (Mekonnen, Grant, and Schwalm 2018; Price et al. 2013; Friend et al. 2014).

Despite the limited productivity of boreal vegetation, warming climates will likely significantly impact the distribution of key land cover types for carbon storage such as wetlands. Boreal wetlands are recognized as a key carbon sink, storing approximately 20-30% of the world’s carbon (Lal 2008; Kendall et al. 2021; Coletti et al. 2013). The Canadian boreal is home to some of the world’s largest wetlands, and includes treed swamps, shrub swamps, treed fens, and bogs (Kendall et al. 2021). Capturing the heterogeneity of these environments in modeling efforts, particularly as components of the terrestrial carbon cycle, is difficult. Using traditional field-based methods and in-situ data has been limited given the vast and remote boreal expanse, which presents many challenges for geospatial modeling.

Under rapidly changing climatic conditions, the interplay between permafrost thaw, wetland expansion and the northward movement of the treeline will translate into dynamic carbon productivity conditions. In order to monitor the carbon cycle of Canadian boreal forested ecosystems, it is imperative that we understand the underlying environmental dynamics and drivers of change, and develop robust terrestrial carbon estimates at appropriate temporal and spatial scales. Impacts of changing climate have complex impacts on these landscapes, and failure to capture changing conditions and resultant vegetation and landscape patterns may translate into considerable uncertainty in modeling estimates (Srinet et al. 2022; Lees et al. 2021; Xie and Li 2020).

Our aim is to review the key drivers of boreal forest productivity, assess the models used to estimate carbon dynamics and land cover change, and highlight current gaps in model development which could be improved by considering a number of recently available satellite datasets and products. To do so, we first review and summarize key drivers of land cover change in the Canadian boreal, and discuss methods for modeling boreal vegetation productivity in section 2. In section 3 we then discuss several novel and thematic remote sensing datasets that can be leveraged to address key uncertainties identified in current modeling approaches. Specifically, new remote sensing datasets, and remote sensing derived products, are highlighted to address key knowledge gaps in current modeling frameworks.

In the 4th section of the review, we demonstrate a conventional remote sensing-based productivity model across the Hudson Plains Ecozone, characterized by boreal forests and wetlands, and the most southern continuous permafrost in Canada. Following this, we propose a new modeling framework which integrates a number of these new satellite datasets for improved model performance in boreal environments.

Drivers of change in the Canadian boreal

In order to capitalize on newly available remote sensing products, and produce a robust approach for modeling vegetation carbon dynamics at appropriate spatial and temporal resolutions, it is critical to identify the key drivers of change in boreal ecosystems. Anthropogenic and natural disturbances drive mortality, succession, and carbon dynamics in boreal ecosystems (Maynard et al. 2014; Price et al. 2013). The most prominent of these disturbances in the boreal terrestrial ecosystem are fire, insects, and disease, which are integral to biogeochemical cycling and productivity maintenance (Maynard et al. 2014). The impacts of these types of disturbances include single tree mortality and reduced growth rates for low severity disturbances, but in severe scenarios, may also lead to stand replacement (Kurz et al. 2013). The degree of disturbance impact subsequently results in the redistribution of carbon in the ecosystem (Kurz et al. 2013). Stand replacing disturbances turn the carbon sinks associated to a landscape into a carbon source, and have complex periods of recovery time, determined by factors such as drought, soil type, disturbance frequency and intensity, and competing vegetation (Douglas et al. 2014). Disturbance regimes and species adaptations are thus key in understanding carbon dynamics, the most prominent of which are discussed in more detail below.

Fire

Coupled with warming climates and increased periods of drought, fires are a critical driver behind boreal species composition (Price et al. 2013). For example, well-adapted deciduous species such as Aspen and Birch are unlikely to experience difficulty in post-fire recovery given their ability to directly regrow their vegetation. Conversely conifer species may be outcompeted, as increasingly hot and dry conditions, along with increasing disturbance frequency, limit new seedling establishment (Price et al. 2013, 2001). In addition to the implications for forest species distribution, there is an increasing trend of fires in peatlands, which not only contribute to increased permafrost degradation, but also the release of considerable amounts of stored carbon, influencing the global carbon cycle (Price et al. 2013). During the 1990-2008 period, boreal fires in Canada are reported to have released an average of 23 TgC yr⁻¹ (Stinson et al. 2011). These fires are the source of 9% of global fire emissions, and are the primary driver behind stand renewal in the Canadian boreal (van der Werf et al. 2010; Brandt et al. 2013).

Insects and disease

Insects and pathogens regulate primary productivity through their impact on the structure and composition of the forest, in turn, they are themselves heavily influenced by fire disturbance and climate (Fleming and Volney 1995; Price et al. 2013). The steadily warming climate across the Canadian boreal is projected to intensify current disease dynamics involving forest pathogens and insect species (Price et al. 2013), particularly through the expansion of insect species range (Kurz et al. 2013). Warming climates are expected to increase the rate of insect species development and reproduction, whilst also reducing winter mortality, which will greatly impact insect distribution as they expand their ranges northwards (Harrington, Fleming, and Woiwod 2001; Price et al. 2013). Insect outbreaks and disease create complex mosaics of impacted forest stands across the landscape (Maynard et al. 2014), many impacting at the individual tree or small patch scale. Consequently, we must carefully consider the spatial scale of modeling efforts moving forward, and critically assess whether they are able to capture the spatial patterns of prevalent disturbances.

Climate, temperature, and drought

Moisture is known as a key control factor governing growth, mortality, and biomass accumulation (Hogg and Wein 2005; Price et al. 2013). Soil moisture in particular, is a primary driver behind vegetative processes like water use, photosynthesis and the decomposition of organic materials, all of which are key components to carbon flux (Jones et al. 2017; He et al. 2017). The drying climate is impacting forests directly through temperature-induced drought stress and decreasing regeneration, and indirectly through drought-related increases in mortality related to disturbances like fire and disease (Price et al. 2013; Hogg and Wein 2005). The Canadian boreal zone is likely to experience increased frequency and severity of drought in the 21^st century, with drought related mortality expected to increase significantly (Hely et al. 2001; Kurz et al. 2013; Michaelian et al. 2011; van der Molen et al. 2011). Under an intermediate future emissions scenario (scenario ‘RCP 4.5′, Taylor et al 2012), extreme drought events are projected to increase by a factor of 3.1 whereby the contribution of drought on global GPP loss is increased from 28% during 1850-1999 to 49% by 2075-2099 (Xu et al. 2019). Future increases in drought are predicted to result in significant alterations in species distribution (Coops and Waring 2011), as well as an overall decrease in productivity of surviving species in the boreal forest, particularly in western North America (van Mantgem et al. 2009; Michaelian et al. 2011).

Recent warming trends have also led to observed increases in productivity, however, largely due to increases in growing season length and available atmospheric carbon in boreal forests (Briffa et al. 2008; Greenwood and Jump 2014; Magnani et al. 2007; Mekonnen, Grant, and Schwalm 2018). The impacts of climate change on boreal forest composition will be more pronounced along the ecotones than the center of a given species’ survival range (Stratton, Price, and Gajewski 2011), and it is expected that climate-sensitive species such as boreal conifers will gradually disappear from the southern boundaries of their current extent by 2100 (Price et al. 2013). With the aforementioned warming climate and changes in permafrost thaw, Canada’s forests are expected to move northwards in latitude and upward in elevation (Price et al. 2013; Greenwood and Jump 2014). Consequently, modeling changes in terrestrial productivity in these transition zones is an area of considerable uncertainty. Remote sensing data products that can be leveraged to address these uncertainties will be revisited later in this communication.

Remote sensing approaches to support modelling of boreal vegetation carbon dynamics

Conventional GPP modelling

Measures of forest vegetation growth and carbon accumulation have conventionally been performed using field inventory plots (Asner et al. 2013). In response to the challenges and labor associated with plot-based measurements, eddy covariance data from tower networks became a prominent method by which to measure temporal dynamics of net CO2 flux and derived GPP. Similarly to field plot data, however, remote environments and the significant installation and maintenance cost of these towers serve to limit the number of eddy covariance sites (Srinet et al. 2022; Sims et al. 2008; Asner et al. 2013; Lees et al. 2021). In addition, regional (e.g., AmeriFlux) and global (i.e., FLUXNET) flux tower networks, via which these measurements are standardized and archived, under-sample the global range of vegetation, climate, and disturbances, leading to considerable bias in measurement-based estimates and model verification (Jones et al. 2017; Schimel et al. 2015). Additionally, GPP estimates from flux towers are themselves a model output with inherent bias and error, and represent the mean flux of the environment surrounding the sensors and tower (Lees et al. 2019), posing challenges as a data source for highly heterogeneous landscapes.

Basic remote sensing approaches for modelling vegetation carbon dynamics

Due to the high spatial and temporal variability of GPP, remote sensing is a critical component toward building regional and global GPP datasets (Verma et al. 2015; Xie and Li 2020). Many models designed to estimate GPP include the fusion of remote sensing data with other geospatial environmental layers and climate data. Some of these datasets, such as Leaf Area Index (LAI) and gridded products derived from ground-based meteorological observations, are difficult to obtain and include substantial errors, making the models challenging to apply accurately over large spatial extents (Flack-Prain et al. 2019; Liu et al. 1997; Z. Wang et al. 2021). In response, remote sensing has become increasingly central to carbon estimation efforts (Asner et al. 2013; Lees et al. 2019; Ciais et al. 2014). Remote sensing methods represent an approach which can directly assess changes in carbon dynamics and land cover transitions without reliance on ground-based inputs, and offer an opportunity for improving our projections of future changes in these environments (Sims et al. 2008; Verma et al. 2015; Schimel et al. 2015; Lees et al. 2021).

Some of the simplest model formulations for GPP prediction are based on the correlation between GPP and a remotely sensed greenness index (Zhang and Ye 2021). Greenness indices from remote sensing observations commonly make use of Near Infrared (NIR) and Red wavelengths of the electromagnetic spectrum. Changes in these indices reflect vegetation growth and vigor, and are key in GPP predictions (Tang et al. 2021). The Enhanced Vegetation Index (EVI) for example, has become increasingly prominent in estimations of GPP. It is commonly applied in models such as the Vegetation Photosynthesis Model (VPM), Temperature and Greenness (TG) model, Greenness and Radiation (GR) model, and the Vegetation Index (VI) model (Jinwei Dong et al. 2015). Research has demonstrated that these models typically have good agreement with eddy-covariance data, yet are commonly noted for their inaccurate performance under drought conditions as a result of lacking sensitivity to water stress related variables (Jinwei Dong et al. 2015).

More recently, advances in microwave remote sensing technologies have presented an alternative to the use of the direct relationship between greenness indices and GPP through the introduction of Vegetation Optical Depth (VOD) products. VOD captures the radiation absorbed by vegetation as observed by microwave observations, which is primarily driven by water content, biomass, vegetation type, and vegetation density (Wild et al. 2022). Notably, microwave observations are not affected by the impact of cloud cover as optical remote sensing products like greenness indices are (Kumar et al. 2018). Due to this key advantage, VOD has become a commonly used proxy for above-ground biomass, and has been used alongside flux tower observations to model GPP (Wild et al. 2022; Scholze et al. 2019).

In addition to greenness indices and microwave observations of VOD, land surface temperature (LST) is a remotely sensed variable which has proven valuable in GPP estimation. LST is an essential climate variable, driven by land surface type and emissivity properties, and with its detection impacted by illumination, and atmospheric conditions amongst other considerations (Yang et al. 2020). Because of the complex relationship between LST and these factors, it is challenging to collect large in-situ LST datasets capable of accurately representing large areas, highlighting the potential for satellite-based observations (Yang et al. 2020).

An example of a relatively simple approach to estimating productivity which integrates the drivers of greenness and temperature is presented by Sims et al. (2008), who developed the temperature-greenness (TG) model for GPP prediction, based solely on remotely sensed inputs. The TG model uses two key drivers: vegetation greenness and LST, which are scaled using reference GPP information. In their 2008 demonstration of this model framework, Sims et al. leveraged EVI and LST from MODIS along with reference GPP from eddy covariance flux tower data to model GPP across their 11 sites. To account for the inactive winter period, Sims et al. set GPP to 0 when LST falls bellow 0 °C. To account for drought stress, the maximum tolerable LST is set to 30 °C, and GPP is set to 0 if temperatures exceed 50 °C. Parameterizing the TG model then relies on estimating the slope between the combined EVI & LST information and the reference GPP dataset (Sims et al. 2008).

The TG model was found to be successful at predicting GPP across a wide range of vegetation types in North America, reporting mean coefficients of determination (r²) of 0.83 for evergreen needleleaf forest sites, and 0.91 for broadleaf forest sites relative to flux tower GPP estimates from 2003-2005 across four and three sites, respectively (Sims et al. 2008). The TG model, like many relatively simple spectral index based approaches, has difficulty in predicting productivity in environments which are heavily influenced by moisture dynamics, principally those affected by drought or high moisture wetland environments (Lees et al. 2019; Sims et al. 2008; Srinet et al. 2022; Tang et al. 2021). The inclusion of additional remote sensing variables to improve model performance and address such limitations has therefore become a key area of exploration as remote sensing technologies and data products continue to develop. For example, in their assessment of simple Light Use Efficiency (LUE) remote sensing-based models for productivity estimation, Wang et al. (2021) found that, compared to eight other model formats, their two new proposed models including correction factors for moisture out-performed the other models they investigated. They implemented moisture scalars via the use of evaporative fraction and the normalized difference water index (NDWI), respectively, reporting r² of 0.86 and 0.82 in relation to a global network of 151 flux towers (Z. Wang et al. 2021).

Comprehensive approaches to modelling vegetation carbon dynamics

In addition to the development of remote-sensing based models, the past decade has seen development of comprehensive terrestrial remote-sensing based carbon dynamics models designed for carbon cycle estimation (see Li et al., 2018 for a comprehensive review of several key satellite-based models). One such modeling approach is the Carnegie Ames Stanford Approach (CASA), which estimates net primary productivity (NPP) through the use of absorbed photosynthetically active radiation (APAR) and LUE (Y. Wang et al. 2019). A limitation of CASA is that it relies on APAR as a function of NDVI, which may result in underestimation in densely vegetated areas due to saturation issues with optical spectral indices (Fang et al. 2021). In addition, the original CASA model has notable issues regarding its soil moisture parameter, and fails to distinguish between vegetation types in its use of LUE (Y. Wang et al. 2019). In contrast, boreal biome-specific models such as the Boreal Ecosystem Productivity (BEPS), explicitly incorporate soil water holding capacity, vegetation and soil information to model productivity in boreal environments (Liu et al. 2002, 1997). Similar to other models, however, BEPS relies on flux data, and is limited by the spatial resolution of meteorological information.

Data from the MODerate-resolution Imaging Spectroradiometer (MODIS) onboard the Terra and Aqua platforms has also been utilized to derive GPP estimates via the use of the MOD17 GPP algorithm (Running et al. 2004). The MOD17 is based on LUE models, and uses APAR and land cover information with constants representative of physiological attributes such as biome specific scaling for LUE (Schubert et al. 2012). Evaluations of the MOD17 GPP product have noted uncertainties related to both the APAR and input meteorological data, leading to updated versions of the product being released. Although it is regarded as a useful measure of GPP, and often used for model comparison, MOD17 is limited by the uncertainty and spatial resolution of its inputs. Most notably, the input MODIS 1-km landcover product used in MOD17 is spatially coarse and therefore not well suited for accurately capturing the heterogenous landscapes of Canadian boreal environments.

Although these aforementioned models focus directly on the estimation of GPP, there have also been models dedicated to monitoring carbon sinks and sources that incorporate estimates of GPP as an intermediate variable. Incorporating microwave-derived VOD introduced in the previous section, recent works such as Scholze et al. (2019) highlight the potential of remotely sensed atmospheric carbon and biomass estimates in quantifying carbon sinks through an adaptation of a diagnostic biosphere model from Kaminski et al. (2017). This biosphere model calculates carbon uptake through the LUE approach, incorporating a vegetation specific light-use efficiency parameter, a water stress factor, and incident solar radiation (Kaminski et al., 2017). Assimilating in-site atmospheric carbon into this model framework along with soil moisture and VOD, Scholze et al. (2019) demonstrate the potential of new remotely sensed estimates of carbon dynamics drivers to produce estimates of carbon sinks that are consistent with more traditional modeling approaches. In particular, they found that the inclusion of microwave derived soil moisture and VOD can reduce flux uncertainty by an additional 15-25% relative to models that only assimilate CO₂.

In addition to the continued development of productivity models, new remote sensing missions have been developed to provide new insights on soil moisture, landscape freeze/thaw state, boreal carbon dynamics and their drivers, as well as deepening long-term data archives on land cover and land cover change dynamics. These missions present opportunities for new model formulations which utilize a wider range of remote sensing input, allowing finer spatial and temporal time scales to be modeled and with a reduced set of input variables thereby increasing its applicability and ease of use over larger areas. Whilst remote sensing-only model formulations have biases of their own, they may offer the most viable means for assessing long-term carbon dynamics assuming desired accuracies and appropriate lengths of observations can be achieved (Schimel et al. 2015).

New remote sensing datasets for boreal productivity modelling

Whilst the basic model formulation using remote sensing observations of landscape greenness and temperature is common to many productivity modeling frameworks, many studies have investigated the need for additional input variables (Lees et al. 2019; Sims et al. 2008; Xie and Li 2020) to ensure the models produce accurate results over a variety of land cover types (Jiaqi Dong et al. 2017; Lees et al. 2021). These additional variables are intended to represent key drivers of vegetation productivity either directly, or through proxy as scalar correction factors. Capturing these drivers of carbon flux at sufficiently fine spatial scales from remote sensing observations has remained a challenge, whilst robust, fine scale predictions of a variety of these additional model drivers are needed. Such predictions may be pursued using several new and exciting satellite missions; a number of which are described more fully below.

Remote sensing of water availability status

Small changes in soil moisture and water stress have large effects on carbon dynamics, and are required to ensure accurate GPP modeling, especially in boreal forest and wetland environments (Ju et al. 2006; Chen et al. 2007; Jones et al. 2017; Flack-Prain et al. 2019; Tang et al. 2021; Jinwei Dong et al. 2015). A viable method to capture seasonal variability in GPP is through the introduction of soil moisture as a proxy for drought; a well-documented control on GPP not captured by LST (Stocker et al. 2019; Verma et al. 2015; Xu et al. 2019). Ground-based soil moisture observations have insufficient spatial density to accurately represent the heterogeneity of soil moisture at spatially extensive scales, and remote sensing measurements of soil moisture are thus a promising alternative (Kumar et al. 2018).

SMOS was launched in November 2009 and was the first satellite mission to provide global measurements of L-band brightness temperatures, from which soil moisture and ocean salinity data products were derived. SMOS is equipped with the Microwave Imaging Radiometer using Aperture Synthesis (MIRAS), which is a passive microwave 2-D interferometric radiometer operating at 1.413 GHz (Mecklenburg et al. 2016). L-band microwave observations are highly sensitivity to soil moisture (Jackson and Schmugge 1991), which along with near-transparency of the atmosphere in the L-band spectral range, allows for both accurate and consistent global observations of emitted radiation from land and ocean surfaces (Mecklenburg et al. 2016).

SMOS measures the brightness temperature emitted from the earth at L-band across a swatch of approximately 1000-km, with a spatial resolution of 35-50-km, and a revisit time of 3 days at the equator (Mecklenburg et al. 2016). Soil Moisture products are available in several products with varying processing levels, with gridded spatial resolutions ranging from 15-25km.

The Soil Moisture Active Passive (SMAP) mission, launched in January 2015, commenced operations in April 2015. Notably, the SMAP radar instrument (providing 1-km resolution measurements) failed after three months of science data collection, necessitating sole use of the radiometer measurements (at 9 to 36-km resolution). SMAP retrieves soil moisture from measurements of naturally emitted L-band (1-2 GHz) microwave radiation, the observed radiance of which is highly sensitive to soil moisture in the top 5 cm of the soil column (Reichle et al. 2019). Notably, these microwave wavelengths are not limited by traditional optical remote sensing concerns such as solar illumination and cloud cover (Kumar et al. 2018). SMAP provides global coverage at a temporal resolution of 2 to 3 days, and can provide insights into both soil moisture status as well as landscape freeze/thaw state. Global grid-projected estimates of soil moisture are available as L3 products and consist of daily global composites of landscape soil moisture, derived from the radiometer brightness temperature measurements at 36-km. Level 4 products are also available, which model soil moisture through the use of SMAP observations and the NASA Catchment land surface model. Beyond an estimate of surface (0-5cm) soil moisture, the Level 4 algorithm also produces an estimate of root-zone (0-100cm) soil moisture at 9-km resolution (Reichle et al. 2019). Although coarse resolution remote sensing information is able to capture regional scale soil moisture anomalies, stand-level information on the interplay between drought, disturbances, and productivity can be difficult to estimate.

Remote sensing of near-Surface soil state and growing season length

In addition to providing information on water availability and drought, another key use of L-band measurements is for monitoring landscape freeze/thaw state. Seasonal freeze/thaw (F/T) dynamics are key drivers of vegetation productivity and soil respiration in boreal and arctic regions, where vegetation growth and photosynthetic activity is heavily constrained by the short growing season, or non-frozen period (Kim et al. 2012; Kalantari et al. 2015). Soil freeze is notably not a primary cause for decreasing plant productivity, yet serves as a suitable proxy for delimiting the end of the vegetation active period - defined as a dormant state for photosynthetic activity and carbon assimilation during winter conditions (Böttcher et al. 2018), and is the primary driver of growing season length in boreal and arctic regions (Kalantari et al. 2015; Derksen et al. 2017). The surface F/T state is monitored by the combined use of active and passive microwave instruments, each with their own unique advantages and limitations (Rautiainen et al. 2014; Derksen et al. 2017). Active RADAR is able to detect freeze/thaw state at high spatial resolution, but is limited by low temporal resolutions (Podest, McDonald, and Kimball 2014). Passive microwave instruments conversely have high temporal resolution, but are limited by their coarse spatial resolution.

Derived from the aforementioned SMOS L-band radiometer observations, the SMOS Level 3 Freeze/Thaw product is gridded at a spatial resolution of 25-km at 3-day time-steps (Böttcher et al. 2018; Kim et al. 2011; Rautiainen et al. 2014). Freezing of moisture in the first 5 cm of ground surface and in vegetation increases the measurements of brightness temperature captured by the L-band due to the high contrast in permittivity between frozen and thawed states (Rautiainen et al. 2014; Roy et al. 2020; Mecklenburg et al. 2016). Capitalizing on this relationship, NASA’s SMAP mission also provides users with a Level 3 Freeze/Thaw product, derived from 36-km spatial resolution L-band observations, and gridded with a spatial resolution of 9 km. Derksen et al. (2017) note several key uncertainties for L-band retrievals of F/T state such as an insensitivity to the soil state when wet snow is present (Derksen et al. 2017).

Fine scale land cover and land cover change datasets

The relationships between greenness, temperature and GPP differs across land cover types due to the variable conditions under which the vegetation is photosynthetically active. In order to understand the impacts of heterogeneous land cover on these relationships, fine-scale land cover information is critical. Canada-wide land cover information is available annually at 30-m from 1984 to present via the National Terrestrial Ecosystem Monitoring System (NTEMS). Landsat annual image composites, produced following the Composite2Change (C2C) approach (Hermosilla et al. 2016) utilize the free and open Landsat archive (White et al. 2014) to produce annual, gap-free, surface-reflectance image composites with a 30-m spatial resolution as well as the detection and characterization of forest disturbance information (Hermosilla et al. 2015). The Landsat archive allows pixel scoring functions to be developed to create best-available pixel image (BAP) composites (White et al. 2014). The pixel scoring function ranks and selects the optimum annual pixel observation to be used in an image composite based on sensor type, day of acquisition, distance to cloud or cloud shadow, and atmospheric opacity (White et al. 2014). The Virtual Land Cover Engine (VLCE) is a methodological framework that utilizes time series of Landsat surface reflectance to enable change-informed annual land cover mapping across Canada’s forested ecosystems (Hermosilla et al. 2018). The VLCE 2.0 extends on the original framework through the refinement of calibration data, the integration of new predictor variables and the optimization of regional model development and implementation (Hermosilla et al. 2022). Annual land cover maps for 1984 to 2019 were produced following the VLCE framework to map 12 land cover classes across Canada at a 30 m spatial resolution (Hermosilla et al. 2022). Out of the 12 classes, there are four vegetated treed classes: deciduous, coniferous, mixedwood and wetland-treed, as well as a wetland class. This land cover information can be used to assess complex and heterogenous landscape patterns and trends of change.

Sentinel 3-Derived land surface temperature

In addition to these new data products, refined and modern estimates of previously explored variables have also been developed over the past decade. Notably, the European Space Agency has introduced new estimates of LST via the Sentinel-3 program, with Sentinel-3A being launched in February 2016, and Sentinel-3B in April 2018. Sentinel-3 carries the Sea and Land Surface Temperature Radiometer (SLSTR) sensor, which is an infrared radiometer with six channels covering visible, near infrared, and shortwave infrared wavelengths, and three channels covering thermal infrared wavelengths (Yang et al. 2020). The Level-2 LST product is created using the SLSTR sensor (channels S7, S8 and S9), and is derived from a split window algorithm using information on surface biome, water vapor content, and vegetation fraction cover (Pérez-Planells et al. 2021; Yang et al. 2020). With a sun-synchronous orbit and high observation frequency, Sentinel-3 provides daily global coverage of LST at 1-km spatial resolution. Given the importance of LST as an indicator of environmental factors such as hydrology and drought, and the documented difficulty of collecting in-situ LST information, Sentinel-3 LST provides modern, high-quality insight into the drivers behind terrestrial carbon sequestration.

Knowledge gaps in predicting vegetation carbon dynamics in the boreal

Given existing modeling approaches, the advent of new remote sensing technologies and products, and the drivers of change in boreal environments, several key knowledge gaps stand out from the literature. Despite the aforementioned efforts for modeling vegetation carbon dynamics, global carbon observation efforts have been described as a “sparse, exploratory framework” (Ciais et al. 2014), with the various carbon monitoring models and products often diverging on estimates (Schimel et al. 2015). In the Canadian boreal, key sources of uncertainty are: the impact of the warming climate on the northern treeline, carbon estimates in wetland landscapes, and the implications of permafrost thaw.

The Canadian treeline

Changes in the timing of seasonal freeze and thaw transitions has a major impact on the growing season length and carbon flux in boreal regions (Kalantari et al. 2015; Derksen et al. 2017; Madani et al. 2021; Pulliainen et al. 2017). With warmer temperatures and longer growing seasons, we anticipate more rapidly regenerating forests, and northward movement of forest vegetation along the treeline (Beck and Goetz 2011; Nemani et al. 2003; Cazzolla Gatti et al. 2019). Although the rate at which species will migrate northward as climate warms is highly variable (Greenwood and Jump 2014), average northward shifts of climate envelopes in North America are predicted to reach as much as 700-km by 2071-2100 (Price et al. 2013; D. W. McKenney et al. 2007).

In the boreal zone, this migration rate is heavily constrained by seed dispersal distance, soil suitability, water, and the suitability of soil conditions (Asner et al. 2009; Price et al. 2013; D. McKenney, Pedlar, and O’Neill 2009; Higgins et al. 2003). In the Alaskan boreal for example, GPP has been increasing by 3.2 +/− 1.9 Tg C yr⁻¹ over the past 20 years, an increase attributable largely in part due to increased productivity early in the growing season (Madani et al. 2021). Although not directly indicative of treeline migration, this increase in GPP does offer some indication of increasingly suitable conditions for vegetation growth. Overall, the advancement of the treeline and the resultant impact on carbon dynamics will be slow and small per unit area (Limpens et al. 2021), yet can collectively magnify carbon accumulation at national levels (Greenwood and Jump 2014).

Wetlands

Despite their widespread distribution in Canada’s boreal (Wulder et al. 2018), the carbon flux related to boreal wetlands is poorly documented, and remains largely unknown (Magnan et al. 2020; Junttila et al. 2021; Lees et al. 2019). Forested wetland ecosystems across Canada accumulated an estimated 70.5 Tg C yr⁻¹ during the 20^th century (Ju et al., 2006) through the slow accumulation of carbon in wetlands due to decomposition of organic matter under anaerobic conditions. Sequestration rates between wetland types are strongly influenced by moisture availability, drainage, vegetation composition, and nutrient availability. Ju et al. (2006) note that well-drained regions accumulate an average of 10.8 kg C m⁻²yr⁻¹, whilst poorly drained areas accumulate an average of 175 kg C m⁻²yr⁻¹. Even within one wetland type, however, notable variability can exist. A study of carbon sequestration in forested peatland in Quebec by Magnan et al. (2020) suggests wetlands store between 62 to 172 kg C m⁻² in their organic layers, with an average accumulation rate across six sites of 15.9 g C m⁻²yr⁻¹.

Wetlands in particular are highly sensitive to disturbance, and large amounts of carbon can be released into the atmosphere through permafrost thaw and severe forest fires (Frolking and Roulet 2007; Magnan et al. 2020; Oris et al. 2014; Turetsky et al. 2011). Without accounting for disturbances, studies suggest that carbon sequestration still exceeds carbon release in forested wetlands by vast margins, despite changing climate (Kendall et al. 2021). Considering the Canadian boreal contains approximately 127 Mha of wetlands (Ju et al. 2006), addressing knowledge gaps in GPP modeling for this specific landcover type is a particularly important consideration for Canada-wide boreal modeling efforts.

Permafrost thaw

Carbon stored in boreal forests contain an estimated 10-20% of the global soil carbon pool (Alster, Allison, and Treseder 2020), yet these estimates do not account for the significant amounts of carbon stored in permafrost. Using the Northern Circumpolar Soil Carbon Database, estimates by Tarnocai et al. (2009) indicate an organic soil carbon pool of 1,672 billion metric tons attributable to the northern permafrost region, representing 50% of global belowground organic carbon. Permafrost temperatures in the northern circumpolar permafrost region have increased over the last several decades (Smith et al. 2010; Biskaborn et al. 2019). Carbon storage in these landscapes is highly sensitive to warming climates (O’Donnell et al. 2012), and may decrease due to rising microbial and enzyme activity, or growing litter inputs and root production, both resultant from climatic warming (Alster, Allison, and Treseder 2020).

Thaw-induced waterlogging is a key driver of forest loss in permafrost transition zones, and one of the dominant drivers behind landcover change in the subarctic boreal, and must be considered in tandem with the aforementioned northward treeline migration in order to capture the full context of landcover change in Canada’s northern boreal (Carpino et al. 2018). Changes in fire regimes and temperature induced thaw are likely to result in significant changes to sequestration ability and hydrology, which may result in increased carbon release from permafrost landscapes back into the environment (Douglas et al. 2014; O’Donnell et al. 2012; Williams et al. 2020; Helbig et al. 2017). As permafrost thaw changes under climate warming, many of the aforementioned conditions for tree migration are expected to change, and thaw dynamics alongside growing season changes may be key indicators of how Canada’s treeline is responding to changing climatic conditions.

Downscaling

In order to address these knowledge gaps using the new remote sensing datasets identified earlier in this communication, the issue of spatial-temporal scale disparities must be addressed. Due to the orbit characteristics of satellite-based sensors, coarse spatial resolution data products are typically delivered with high temporal revisit, as is the case with the Sentinel 3 LST product and the SMOS F/T products. Conversely, remote sensing vegetation greenness indices such as the NDVI used in the TG modeling framework can be derived from Landsat and Sentinel 2, resulting in much finer spatial resolution, but with significantly longer revisit times. As such, optimizing the use of thematic data products representing complex drivers of vegetation productivity requires the implementation of data integration and downscaling approaches. Downscaling studies typically involve the use of multispectral imagery or topographic information with high spatial resolution to model downscaled interpretations of medium or coarse-scale data products. While not the focus of this review, challenges of downscaling approaches are considered in the proposed model framework presented later in this communication.

Demonstration of a simple model for estimating boreal productivity

In order to demonstrate the identified issues and knowledge gaps in modeling estimates of GPP in boreal environments, we apply the existing TG model (Sims et al. 2008) described in section 2.2 over the Hudson Plains Ecozone. Through this demonstration, we can identify the inherent model biases and limitations of a simple GPP model prior to discussing additional model inputs and improvements.

The Hudson Plains Ecozone is a 36,408,957-ha region in Canada containing the largest wetland complex in the country, and the most southern continuous permafrost in North America (Abraham and McKinnon 2011). Encompassing sub-arctic ecoregions in the North and mid-boreal ecoregions in the South, the climate of the Hudson Plains is typical of the Canadian boreal, with cold and dry winters, and cool short summers when most precipitation occurs (Brandt et al. 2013; Magoun et al. 2005). Due to characteristically wet soil conditions, forests in the Hudson Plains are typically found in small clusters across the landscape, interspersed with ponds, lakes, and forested & non-forested wetlands, as shown in Figure 1 (Hermosilla et al. 2022, 2018). Denser forests are present in southern parts of the Ecozone, where soil drainage is higher (Abraham and McKinnon 2011).

Figure 1. National Terrestrial Ecosystem Monitoring System (NTEMS) 30-m land cover map representing the Hudson Plains Ecozone, 2018 (Hermosilla et al. 2022, 2018).

Data

We utilized daily Sentinel 3 SLSTR LST & NDVI data for 3 years spanning Jan 1^st 2017 through Dec 31^st 2019. For model evaluation and parameterization, we utilized the SMAP L4C data product. Delivered as part of the SMAP suite of products, the L4C product provides daily model estimates of GPP based on a complex integrated model with inputs from the Goddard Earth Observing System Model (GEOS) surface meteorology, MODIS fractional photosynthetically active radiation (fPAR) climatology, and SMAP-derived soil moisture (Jones et al. 2017). The SMAP L4C product delivers daily estimates of GPP (g C m⁻²day⁻¹) at 9-km grid spacing. In their assessment of SMAP L4 carbon products, Jones et al. (2017) reported a mean RMSE between L4C Ops products and 26 flux tower validation sites of 1.04 g C m⁻²day⁻¹. In addition, they note that, relative to other products, SMAP may be of particular use in areas impacted by drought due to the inclusion of SMAP soil moisture.

To enable the comparison, we resampled the 1-km Sentinel 3 SLSTR LST & NDVI to 9-km spatial resolution through bilinear interpolation using the resample function from the raster package in R (Hijmans and van Etten 2012). In order to ensure GPP is compared across relatively homogenous land cover types we extract NTEMS land cover information from 2018 using the ‘exact_extract’ function from the exactextractr package in R (Baston 2022), restricting analysis to 9-km cells which had a clear majority (>70%) of a single land cover type. Sentinel 3 SLSTR LST and NDVI were resampled and aggregated from their native resolution and used to predict 16-day mean GPP at 9-km resolution for 2018 & 2019 with reference to SMAP L4C, similar to the approach used for model parameterization per the concepts outlined by Sims et al. (2008). For an overview of input data sources, see Table 1.

Table 1. Input data sources Sentinel 3 SLSTR LST & NDVI are resampled using bilinear interpolation to match the spatial resolution and grid of SMAP L4C GPP.

Product	Description	Spatial Resolution	Temporal Resolution
Sentinel 3 SLSTR LST	Land Surface Temperature	1-km	Daily
Sentinel 3 SLSTR NDVI	Normalized Difference Vegetation Index (Ancillary)	1-km	Daily
SMAP L4C GPP	Gross Primary Productivity	9-km	Daily
NTEMS	Land Cover	30-m	Yearly

Landcover composition from NTEMS is extracted for each 9 km SMAP grid pixel.

TG model & adaptation

Following the TG concept, we utilized Sentinel 3 SLSTR LST and NDVI to create the product of (greenness * temperature) which was then scaled using the slope between this product and the commensurate SMAP L4C GPP (g C m⁻²day⁻¹), stratified by pixels containing homogenous (>70%) land cover types. Scalars were estimated based on the relationship between the combined NDVI & LST and SMAP GPP observations from 2017, and then set as constant for the 2018 and 2019 GPP simulations. A conceptual diagram of the TG model for the Hudson Plains Ecozone is shown in Figure 2. We then calculated Pearson’s correlation coefficient ‘r’ between model and SMAP GPP for all 16-day mean estimates, across all pixels per stratified area. As shown in Figure 3, correlations were generally high (>0.8), although notable patterns of seasonal over and under estimation are present, particularly in wetland and forested wetland environments.

Figure 2. 16-Day mean (2017/04/06 to 2017/04/22) LST (degrees Kelvin) (top-left), NDVI (bottom-left) for the Hudson Plains Ecozone (unstratified), combined following the TG concept into an initial estimate of TG GPP (g C m⁻²day⁻¹) (middle-right).

Figure 3. Estimated 16-day mean Gross Primary Productivity (GPP) per stratified land cover type across the Hudson Plains (2018 & 2019), where A = Coniferous > 70% (n = 598), B = Mixed-Wood > 70% (n = 6,440), C = Wetland-Treed > 70% (n = 84,318), and D = Wetlands > 70% (n = 9,016).

Amongst the four dominant land cover classes identified in the Hudson Plains, correlation coefficients (r) were higher for coniferous (Figure 3A, r = 0.90) and mixed-wood forests (Figure 3B, r = 0.93) than forested wetlands (Figure 3C, r = 0.85) and wetlands (Figure 3D, r = 0.84). The model shows higher overall prediction errors during spring and summer (Figure 3), and consistently overestimates GPP during winter and fall (Figure 4). These trends are summarized in Figure 5.

Figure 4. Temporal variation of 16 day mean Gross Primary Productivity (GPP), averaged per stratified land cover type across the Hudson Plains (2018 & 2019) predicted by SMAP L4C or the adapted TG model. A = Coniferous > 70%, B = Mixed-Wood > 70%, C = Wetland-Treed > 70%, and D = Wetlands > 70%.

Figure 5. General trends of over and under-prediction of seasonal GPP per the adapted TG model, relative to SMAP L4C GPP based on model performance shown in Figure 3.

A Proposed new model framework

Overcoming conventional model limitations

To improve model estimates of GPP in boreal environments, the knowledge gaps and conventional model limitations discussed in this communication will need to be meaningfully integrated into a robust new framework. Given the significant difference between TG and SMAP L4C model complexity, the results shown in Figures 3 and 4 highlight the strong performance of the core TG modeling concept as the basis for this framework. The TG model framework demonstrated it was able to capture key trends present in the L4C product across various land cover types in the Canadian boreal, using significantly less input data.

Results reaffirmed that the relationship between greenness & temperature and GPP is not as strong in land cover types that are heavily influenced by moisture dynamics, suggesting there is significant variation in GPP not captured by LST and NDVI, and echoing the notion that the addition of land cover specific moisture scalars and parameterization is critical for optimizing TG model performance and constraining model variance (Jiaqi Dong et al. 2017; Sims et al. 2008). In addition, results shown in Figures 3 and 4 suggest the TG model yields estimates of GPP during winter periods not shown in the SMAP estimates which suggests marking the winter inactive period by 0 °C may not be appropriate in boreal environments. Instead, constraining winter inactive periods through the use of F/T data (section 3.2) is expected to more accurately represent constraints on vegetation productivity in these northern landscapes.

Over and under-prediction of GPP seasonally, and by land cover type, is evident in Figure 3 and likely due to non-linear relationships between GPP and (scaled NDVI * scaled LST), particularly in the wetland-treed and wetland dominated pixels (Figure 3C and 3D), and during the spring and summer seasons. The non-linearity of this relationship is apparent given the inequality of variance visible across seasons in Figure 3. If GPP were predicted across the entire Hudson Plains ecoregion without land cover specific stratification and model parameterization, there would be significant over and under prediction across the various land cover types (Figure 5), given the land cover-specific differences in TG model parameters. Improved model performance will rely both on the implementation of appropriate model parameterizations to account for land cover type, as well as improved spatial resolution in model estimates of GPP through the application of downscaling approaches (Section 4.4).

Data & model overview

To address spatial and temporal heterogeneity, identify responses to disturbance, and integrate moisture information for improved GPP estimates, we propose a new modeling approach based on Sentinel-3 LST and greenness estimates combined with SMAP and SMOS derived soil moisture and freeze/thaw data. This approach capitalizes on the proven capabilities of a remote sensing-driven model using temperature and greenness, but augments this model with the required correction factors to address known uncertainties, and constrain estimates of GPP in boreal environments.

Soil moisture data will be downscaled using topographic information from a Digital Elevation Model (DEM) and used to compute a soil moisture modifier, whilst the wetland-oriented scalar NDWI will be derived from Sentinel 2, and introduced to capture the spatial variability in the heterogenous wetland landscapes. F/T data products will be downscaled via a Hidden Markov Model (HMM) algorithm to model downscaled transitions between F/T on a yearly basis to compute a growing season length modifier. Both SMAP and SMOS Soil Moisture and F/T products will be used to compute these model correction factors, and will be assessed for their ability to improve model performance.

Downscaling considerations

As outlined in the downscaling section earlier in this communication, the disparity in resolution between optical, topographic, temperature, and microwave remote sensing sensors and products translates into a need for the use of downscaling methods. In order to adequately capture drivers of GPP at appropriate spatial scales in heterogenous environments like wetlands and treeline transition zones, the resolution of key thematic datasets such as F/T have to meaningfully translated into fine-scale model inputs. At its core, downscaling methods revolve the identification of the relationship between a course-resolution dataset and finer-scale datasets that capture proxy indicators of the variable of interest. The downscaling approach then looks to minimize the error and maximize correlation and determination coefficients between the downscaled and reference dataset in order to provide a finer-scale product of the variable of interest. Recently, random forest regression models have become a popular method by which to achieve this. Wu & Li. (2019) used a random forest regression model to test downscaling of LST through the use of surface reflectance variables, as well as spectral indices, terrain factors, and land cover types. Similarly, Wei et al. (2018) used random forest regression to downscale SMAP soil moisture products with 1 km resolution MODIS optical and thermal information. Vergopolan et al. (2020) were able to improve the resolution of SMAP soil moisture through downscaling even further by applying a merged framework incorporating a land surface model, radiative transfer model, and Bayesian scheme to produce 30 m downscaled SMAP L3/4 soil moisture. Although outside of the scope of this communication, determining optimal downscaling approaches and output product resolution will have to carefully consider works such as those mentioned here, and how they impact model robustness, performance, and error propagation.

CAN-TG model of GPP

Once derived, these downscaled Soil Moisture and F/T products will be leveraged alongside Sentinel LST and Greenness in the newly developed CAN-TG GPP model, which is expected to result in improvements of interannual and intra-annual estimates of GPP at finer spatial scales than currently available (Li et al. 2018). Model parameterization will be achieved by first testing on homogenous areas to determine land cover specific slopes between the combined SM, FT, LST, & Greenness and the reference GPP. In principle, model testing in such areas would not require extensive downscaling, and could be used to build model confidence. Key downscaling methods can then be tested to bring the modeled GPP estimates down to appropriate spatial resolutions in key areas of interest, and assess changes in GPP along the northern Canadian tree line (Böttcher et al. 2018). For a proposed dataset overview, see Table 2; for an overview of proposed model structure, see Figure 6.

Figure 6. Proposed Model Framework Structure for CAN-TG. Red boxed indicate primary and secondary data sources. Green box indicates National Terrestrial Ecosystem Monitoring System (NTEMS) land cover attribution information.

Table 2. Primary data sources for proposed modeling framework.

Satellite Dataset	Analysis Level	Spatial Resolution	Temporal Resolution	Available Since	Data
SMOS	⅔	25-km	Daily	2010	Daily Freeze/Thaw & Soil Moisture
SMAP	¾	9-km	3 hours	2015	Daily Freeze/Thaw & Soil Moisture
Sentinel 3 SLSTR	2	1-km	Daily	2016	LST & Greenness
Sentinel 2	2	10 – 60-m	10 days	2015	NDWI

When adding variables to a model, their correlation is important to consider (Sims et al. 2008), therefore, the degrees of independence for the input variables will be investigated in order to present the most robust final modeling approach possible. Reducing model complexity in turn reduces the propagation of errors and uncertainties from the input data (Z. Wang et al. 2021), and will lead to improved GPP estimations in comparison to more complex models. In order to assess whether these simplified approaches can provide improved estimates of GPP across the heterogenous landscapes of the Canadian boreal, we recognize the outlined approach will need to be assessed against a suite of reference GPP datasets that accurately represent key boreal landscapes such as wetlands and forested wetlands. This assessment should use both collections of in-situ datasets and existing modeled estimates, and will demonstrate the potential advantages of CAN-TG over competing model frameworks such as SMAP L4 GPP.

Once verified, the new model will deliver spatial estimates of carbon and moisture dynamics from 2016-2020, as well as spatial information on drivers of northern treeline dynamics and wetland GPP accumulation. Due to the fine spatial and temporal resolution of the suggested model, this spatial information will also be appropriate for investigating the heterogenous impacts of disturbances such as wildfires and insect outbreaks on vegetation productivity. In addition, the model will be used to project future likely land cover changes by using existing climate change scenarios.

Conclusion

The Canadian boreal is a prominent component of the global terrestrial carbon cycle, sequestering significant amounts of carbon during its growing season, and accumulating vast stores of carbon in its wetlands. GPP, the rate at which vegetation sequesters carbon, is a key indicator of the complex carbon cycle, and may be modeled through various approaches. In this communication, we outline a variety of conventional and current frameworks for GPP modeling, and identify multi-source remote sensing frameworks based on the TG model ((scaled greenness * scaled land surface temperature) * GPP scalar) as a promising direction for continued GPP model development.

Estimates of GPP in the Canadian boreal are uniquely difficult to formulate due to the considerable uncertainty revolving wetland landscapes, a rapidly warming climate, and changing freeze/thaw dynamics. To address these uncertainties, we will utilize SMOS freeze/thaw, SMAP soil moisture, and Sentinel-3 greenness & land surface temperature information for the development of a remote sensing-driven way forward in GPP modeling.

Using the newly available greenness and land surface temperature information from Sentinel 3, we constructed this demonstration to provide initial GPP information for model parameterization and comparison, with reference to SMAP L4C. Although a highly complex model with its own uncertainties and limitations, these reference estimates provide an early indicator of how well trends in GPP may be captured via the approach we present. Model parameterization through the use of spatially explicit remote sensing datasets such as SMAP GPP provide a distinct advantage over traditional reference GPP data from flux towers, as they provide daily and spatially explicit reference data. Although both SMAP and flux-derived GPP are based on complex models, the spatial coverage of the remotely sensed SMAP GPP allows us to investigate drivers of carbon dynamics in unique Canadian boreal environments under changing climates on a pixel-by-pixel basis. This allows us to investigate the parametrization of robust model structures such as the one proposed in this communication with use of reference GPP which is not captured adequately in the flux station network.

Via the continued development of the adapted TG model through the proposed framework, we expect to capture greater levels of the variability present in SMAP L4C at fine spatial and temporal resolutions, leading to a robust, remote sensing driven model for accurate estimation of GPP in boreal environments. The inclusion of growing season length scalars is expected to constrain observed seasonal over-prediction of GPP in winter and fall, whilst the soil moisture modifier is expected to constrain underestimation of spring GPP, and overall variation of estimates during summer (Figure 5). The addition of these scalars is also expected to significantly improve estimates of GPP in environments dominated by wetlands and treed wetlands, particularly during spring and summer.

We reaffirm the need for continued advancement of exclusively remote sensing-driven GPP modeling efforts. Per the identified drivers of carbon dynamics, present knowledge gaps, and proposed model framework, we expect that estimates of GPP through the continued development of the TG model concept will result in robust, high-resolution information on boreal carbon dynamics.

Acknowledgements

The authors acknowledge the help of Dr. Sara Knox for her input on the preparation of the manuscript.

Disclosure statement

No conflict of interest was reported by the authors(s).

Additional information

Funding

This work was supported by the Canadian Space Agency under Grant 21SUESECDL; and in part by NSERC under Grant RGPIN-2018-03851.