Climate drivers of Arctic tundra variability and change using an indicators framework

This study applies an indicators framework to investigate climate drivers of tundra vegetation trends and variability over the 1982–2019 period. Previously known indicators relevant for tundra productivity (summer warmth index (SWI), coastal spring sea-ice (SI) area, coastal summer open-water (OW)) and three additional indicators (continentality, summer precipitation, and the Arctic Dipole (AD): second mode of sea level pressure variability) are analyzed with maximum annual Normalized Difference Vegetation Index (MaxNDVI) and the sum of summer bi-weekly (time-integrated) NDVI (TI-NDVI) from the Advanced Very High Resolution Radiometer time-series. Climatological mean, trends, and correlations between variables are presented. Changes in SI continue to drive variations in the other indicators. As spring SI has decreased, summer OW, summer warmth, MaxNDVI, and TI-NDVI have increased. However, the initial very strong upward trends in previous studies for MaxNDVI and TI-NDVI are weakening and becoming spatially and temporally more variable as the ice retreats from the coastal areas. TI-NDVI has declined over the last decade particularly over High Arctic regions and southwest Alaska. The continentality index (CI) (maximum minus minimum monthly temperatures) is decreasing across the tundra, more so over North America than Eurasia. The relationship has weakened between SI and SWI and TI-NDVI, as the maritime influence of OW has increased along with total precipitation. The winter AD is correlated in Eurasia with spring SI, summer OW, MaxNDVI, TI-NDVI, the CI and total summer precipitation. This winter connection to tundra emphasizes the role of SI in driving the summer indicators. The winter (DJF) AD drives SI variations which in turn shape summer OW, the atmospheric SWI and NDVI anomalies. The winter and spring indicators represent potential predictors of tundra vegetation productivity a season or two in advance of the growing season.


Introduction
Indicators provide a framework for understanding complex ideas using simplified metrics that can be standardized, used for long-term monitoring, and communicated effectively with the broader public. The climate indicators framework was developed in order to support climate assessments (Kenny et al 2016) and was successfully implemented in the U.S. Global Climate Research Program (www.globalchange.gov, NASEM 2017). As researchers are challenged to understand the rapidly changing Earth system, the indicators framework is increasingly being used to monitor and document change. One such example is the NOAA Arctic Report Card (e.g. Richter-Menge et al 2019) that develops and updates essays annually on the atmosphere, ocean, sea-ice (SI), vegetation, and other indicator variables.
The indicators framework has been applied to understand Arctic climate change and variability. Box et al (2019) analyzed key Arctic indicators in physical (e.g. SI area, temperatures, precipitation) and biophysical (e.g. tundra productivity, disturbances, riverine biogeochemistry) systems jointly to document their linkages. Overland et al (2019) combined a set of Arctic indicators into a composite Arctic Climate Change Index as a way to document the manifestation of global change in the Arctic. In this study, we apply the indicators framework to climate drivers of Arctic tundra vegetation variability and trends, enabling a diagnosis of relevant drivers including several that can serve as seasonal-scale predictors of NDVI variations.
The Normalized Difference Vegetation Index (NDVI) is commonly used to monitor vegetation greenness (Tucker 1979, Frost et al 2020 and has been applied to Arctic tundra from the first analyses indicating enhanced greening at high latitudes (Myneni et al 1997, Jia et al 2003 and continuing through annual monitoring reports (e.g. NOAA State of the Climate and Arctic Report Card). NDVI contrasts surface reflectance in red (R) and near infrared (NIR) wavelengths (figure 1(a)), calculated as (NIR − R)/(NIR + R). Most red light is absorbed by plants for photosynthesis, while NIR wavelengths are reflected by the complexities of the plant canopy. As vegetation cover increases particularly in less abundantly vegetated areas such as the Arctic, either spatially or temporally, additional red light is absorbed and more NIR radiation is reflected, increasing the NDVI ( figure 1(a)). Over the summer, Arctic NDVI increases rapidly after snowmelt in late May to early June, peaks in July and August, and decreases slowly until the first snowfall ( figure 1(b)).
The long-term satellite record indicates that 'greening' (increase in NDVI) is occurring in much of the Arctic tundra biome (non-glaciated land north of treeline), especially on Alaska's North Slope (Jia et al 2003, Frost et al 2020. NDVI trends, however, vary spatially and temporally and have recently shown a general stabilization or even decline especially in Arctic North America. At circumpolar and continental scales, tundra vegetation and productivity vary considerably across north-south climate gradients, coastal-inland continentality gradients, and with east-west floristic differences (CAVM Team 2003, Raynolds et al 2019. At regional scales, geology, elevation, and physiography strongly affect vegetation patterns, and at landscape and plot scales, small differences in microrelief, soil moisture, and disturbance strongly affect vegetation patterns . Field studies have documented the highly variable nature of the NDVI response to local temperature gradients, local precipitation, surface wetness, changes in the length of the snow-free period, growing season length, extreme weather events, and large disturbances such as fire, floods, and insect outbreaks (e.g. Trofaier et al 2013, Raynolds and Walker 2016, Bjerke et al 2017. This study applies an indicator framework to explore trends and variability in Arctic NDVI relative to a few key parameters, some used in previous work and some new. The Advanced Very High Resolution Radiometer (AVHRR) maximum annual NDVI (MaxNDVI) and summed growing season bi-weekly (i.e. time-integrated) NDVI (TI-NDVI)) time-series that now span from 1982 to 2019 form the core of this analysis, as they are the most direct observational indicators of tundra vegetation productivity for the full circumpolar Arctic using reasonably consistent sensor platforms with a long period of record. The other indicators that will be used in the study are chosen because of their established or hypothesized importance for NDVI. The amount of springtime SI in nearshore (100 km buffer) coastal waters influences the total summer open-water (OW) area and available summer warmth for vegetation productivity in the Arctic tundra, which is a maritime biome (Walker et al 2005). The decline in Arctic SI plays a key role in NDVI increases, since land surface temperatures become warmer once SI melts along the coast (Bhatt et al 2010, Dutrieux et al 2012). However, as the ice continues to decline, this relationship is weakening (Bhatt et al 2013) as other factors come into play. The total summer warmth depends on the timing of coastal SI melt and is described using the summer warmth index (SWI: sum of the degree months above freezing during April-September). Total summer precipitation (TSP) is a new indicator included in this context, because there is a growing body of evidence that the Arctic hydrological cycle is intensifying , will increase in the future (Bintanja and Selten 2014), and it influences tundra NDVI variations (Lara et al 2018, Campbell et al 2020. A continentality index (CI, the maximum minus minimum monthly temperature in a calendar year) is introduced to provide a measure of atmospheric variability across seasons, since winter warming has also been found to impact vegetation productivity (Bokhorst et al 2009). Additional indicators that are introduced in the context of tundra variability are the first two modes of sea level pressure variability in the Arctic, the Arctic Oscillation (AO) (Thompson and Wallace 1998) and the Arctic Dipole (AD) (Wang and Ikeda 2000) (explained below).
This study presents the pan-Arctic and continental scale relationships between key indicators and tundra vegetation. We discuss each indicator variable through its climatology (long-term average) and spatial trends as shown in circumpolar maps, time-series graphs, and through the correlations between each of the variables and vegetation for Arctic, Eurasian, and North American Arctic tundra. Our goal is to identify the important parameters affecting tundra vegetation productivity and to better map the climate drivers operating on Arctic tundra vegetation. The sequencing of the processes across seasons will point to the potential for predicting NDVI at ranges of a season or two.

Data
Spatially distributed time series were calculated for ten environmental indicators for the 1982-2019 period: two measures of vegetation, two measures of SI, two temperature data sets, continentality, precipitation, and two teleconnection indices. Remotelysensed 8 km resolution NASA GIMMS (Global Inventory Modeling and Monitoring System) biweekly composited maximum NDVI (MaxNDVI) data (Pinzon and Tucker 2014) from 1982 to 2019 were derived from AVHRR sensors on NOAA-7 through NOAA-19 satellites. The GIMMS NDVI3g V1.2 product corrected discontinuities in the GIMMS NDVI north of 72 • N and permitted the first comprehensive analysis of NDVI trends in the High Arctic (Bhatt et al 2010). The GIMMS dataset used Seaviewing Wide Field-of-view Sensor data for calibrating among the collection of used sensors (Pinzon and Tucker 2014). This study uses the 1/12 • resolution NDVI data to more closely match the grids used for SI and surface temperature. The MaxNDVI is the highest summer NDVI value, representing peak vegetation photosynthetic capacity, and serves as an indicator of peak tundra biomass (Tucker 1977, Myneni et al 1997. The TI-NDVI is the sum of biweekly values from May to September that exceed the threshold value of 0.05. TI-NDVI incorporates the length of the growing season and phenological variations, and represents net primary production better than MaxNDVI (Tucker and Sellers 1986).
Spring SI and OW were investigated using the Special Sensor Microwave Imager SI concentration data from 1982 to 2019 (Comiso and Nishio 2008). Spring SI concentration is the percentage of a pixel covered by ice during a 3 week averaged period centered on the week when that pixel has 50% climatological ice cover. The methodology to construct SI was chosen, because TI-NDVI is most strongly correlated with SI at 50% concentration (Bhatt et al 2010), and the timing of the 50% threshold varies across the Arctic. OW is the inverse, the portion of a pixel that was not covered by SI (i.e. 100%-SI concentration). Summer OW is the average of the weekly May through August OW percentage.
Both land surface and air temperatures were used to construct indicators in order to compare their variability and change. The AVHRR-derived landsurface-temperature data were corrected through effective cloud-masking techniques and calibrated using in situ surface-temperature data including temperatures from the Surface Heat Budget in the Arctic experiment conducted in the central Arctic from October 1997 through September 1998 and 2 m air temperatures from meteorological stations (Comiso et al 2003). The European Center Reanalysis version 5 (ERA5) (Hersbach et al 2020) was used for 2 m air temperature, and precipitation. The ERA5 shows improvements over ERA-Interim in the Arctic (Hersbach et al 2020), due to a reduction in background errors and differences in data quality control. The ERA5 temperature and precipitation have been shown to perform well in the eastern Arctic (Graham et al 2019) and Alaska (White et al 2020). Note that summer ground-surface temperatures are usually warmer than 2 m air temperature due to the absorption of solar radiation by the land surface (Raynolds et al 2008). SWI is calculated as the sum of average April to September monthly surface temperatures above freezing at each pixel, in units of • C months. SWI based on the AVHRR land-surface temperatures is identified as SWIs, and SWI based on 2 m ERA5 air temperature is called SWIa. The CI is the difference between the maximum monthly air temperature and the minimum monthly air temperature for each pixel in a calendar year, in units of • C using ERA5. There are various ways to calculate the CI (Vilček et al 2016), but we chose the simple index, which has been used successfully at high latitudes to compare the sensitivities of mass balance in maritime and continental glaciers (De Woul and Hock 2005). Summer precipitation (mm) is the sum over June, July and August for each pixel for each year. Finally, the AO and AD indices were calculated for winter (December-February) and summer (June-August) over the 1979-2019 period by applying empirical orthogonal analysis (Kutzbach 1967) on the seasonal ERA5 sea level pressure.

Analysis methods
Our study area was the Arctic non-alpine tundra region as defined in the Circumpolar Arctic Vegetation Map, with the southern boundary defined as treeline (CAVM Team 2003). Linear trends were calculated for time series averaged over oceanic regions within 100 km of the Arctic coastline and over the full tundra domains at elevations <300 m. Results are presented for the full Arctic, and divided into North America and Eurasia. Spatial trends in indicators are shown as maps of magnitude change over the 1982-2019 period, and the least-squares-fit method was used to determine the trends. Standard Pearson correlation coefficients and regression coefficients were calculated after all series were linearly detrended. The statistical significance of correlations and trends were assessed using the two-tailed Student's t-test with a 95% and 90% threshold for significance. Reduced degrees of freedom for trend significance followed Santer et al (2000).

Normalized difference vegetation index
The tundra MaxNDVI and TI-NDVI trends display increases over most of the Arctic but there are also regions with declines: southwestern Alaska, high Canadian Arctic, and scattered areas of western Eurasia (figures 2(a) and (b)). There are more pixels with negative trends for TI-NDVI than for MaxNDVI over the 1982-2019 period. MaxNDVI pixels with trends lower than −0.1 account for 7% of the total in both Eurasia and North America. In contrast, TI-NDVI pixels with trends below −0.4 account for 15% of the total in Eurasia and 23% in North America. The North Slope of Alaska is the region of the Arctic with the most homogenous positive trends in MaxNDVI and TI-NDVI trends over the 1982-2019 period (figures 2(a) and (b)). Arctic MaxNDVI and TI-NDVI time series display a steady increasing trend over the 1982-2019 period though TI-NDVI shows recent declines (figures 2(c) and (d)). TI-NDVI displays a declining trajectory over the last decade that results from recent record lows: the second-lowest value was in 2017 in Eurasia and the lowest value was in 2018 in North America. The standard deviation for MaxNDVI is similar for Eurasia and North America while that for TI-NDVI is slightly higher in Eurasia (0.16) than North America (0.13) (table 1, bottom). MaxNDVI and TI-NDVI are more strongly correlated to each other in North America (0.81) than Eurasia (0.59) (table 2).
It is generally agreed that, to the first-order, Arctic tundra plants are temperature limited (Bliss and Petersen 1992), but there is growing evidence that precipitation should also be considered in this context (e.g. Keuper et al 2012, van der Kolk et al 2016. In some regions of the Arctic, SWI trends are increasing (section 3.2) yet the corresponding NDVI trends are decreasing, suggesting additional drivers are in play. This motivates analysis of other potential climate factors, such as continentality and precipitation.

Spring SI and OW
Spring SI and summer OW both display their largest trends (decreasing for SI and increasing for OW) in the Beaufort, Chukchi, Laptev, and Kara/Barents Seas (figures 3(a) and (b)). Spring SI trends are increasing in the southern portions of the Bering Sea where North Pacific decadal variability resulted in above normal SI from 2006 to 2013 (Frey et al 2015). Spring SI has decreased and OW has increased over northern Baffin Bay along the northwest Greenland coast. The expansion of the ice edge poleward is greater in the western Arctic (Beaufort and Chukchi Seas) than in the Atlantic sector (Barents and Kara Seas) (figure 3(b)), based on comparing spring SI versus summer OW trends (figures 3(a) and (b)). However, the OW trends are largest in the Atlantic sector. The SI declines and OW increases are larger on the western side of Hudson Bay compared to its eastern. OW has decreased in northeast Greenland where ice is being transported out of the Arctic Basin.
Spring SI time series for a 100 km coastal zone in the Arctic, Eurasia, and North America indicate a decline from around 60% ice cover in the early 1980s to between 30% and 40% in recent years (figure 3(c)). The SI linear trend is −18.7%, −24.0%, and −17.6% for the Arctic, Eurasia, and North America, respectively (table 1, top). Spring SI in the 100 km coastal zone has leveled off for North America in the last decade but has continued to decline in Eurasia. OW time series in the 100 km coastal zone display an increase from 40%-47% in the early 1980s to 60%-65% in recent years (figure 3(d)). The OW trend is 16.0%, 20.9%, and 13.0% for the Arctic, Eurasia, and North America, respectively (table 1). SI and OW trends and interannual variability are larger in Eurasia than North America because of the perennial SI in the high Canadian Arctic.

Summer warmth index
Climatological SWI of surface temperature (SWIs) and 2 m air temperature (SWIa) (figures 4(a) and (d)) display similar overall patterns where terrestrial Arctic SWI is greater than that over the Arctic Ocean. SWIs trends for the polar region north of 55 • N are largest over the ocean areas adjacent to the Arctic land areas and display a band of decreasing trends between 60 • and 70 • N over Siberia and in the vicinity of the Mackenzie River ( figure 4(b)). SWIa trends are increasing over most of the study region and display weak trends over the Arctic SI ( figure 4(e)).
Since the focus of this study is Arctic tundra, subsequent discussion will focus on tundra areas shown in figures 4(c) and (f) for SWI trends though plots will show pan-Arctic values to facilitate the climate driver discussion. Regional linear trends in Eurasia are larger for SWIa than SWIs while trends in North America are larger for SWIs than SWIa (table 1, top). SWIs displays the largest increasing trends around Greenland, Baffin Island, and southwestern Alaska ( figure 4(c)). SWIa has large increasing trends along the coasts of the Laptev Sea, Chukotka and southern Canadian archipelago (figure 4(f)). SWIs and SWIa time series display increasing trends in the Arctic, Eurasia and North America (figures 4(g) and (h)). In Eurasia, mean tundra SWIs is 30.7 • C month and SWIa is 23.5 • C month, about a 7 • C month difference between the surface and 2 m temperatures. In contrast, in North America mean tundra SWIs and SWIa are 26.2 • C and 27 • C month, respectively (table 1, center). It is noteworthy that the SWIs mean is lower than the SWIa, in disagreement with previous comparisons (e.g. Bhatt et al 2017). The large 2 m air temperature-based SWI in North America arise from mean ERA5 SWIa over western Alaska tundra regions that are larger than mean SWIs in the same region (figures 4(a) and (d)). These SWI values if correct suggest large poleward heat transport in the Pacific sector but may also arise erroneously and warrant further investigation in the ERA5, which is beyond the scope of the current study. The interannual variability quantified by the standard deviation is slightly higher in Eurasia than for North America for both SWIs and SWIa (table 1, bottom). The SWIs over Eurasia and N. America decreased from 2000 to 2010, suggesting low-frequency climate variations present during the summer over Arctic tundra.

Continentality index
The Arctic tundra is a maritime biome (Walker et al 2005), and the trend in the CI (Tmax-Tmin of average monthly temperatures in a calendar year) reflects the increasing maritime influence caused by more nearshore OW. The climatology of the CI is characterized by small values (0 • C-20 • C) over the ocean and larger values (40 • C-65 • C) over land ( figure 5(a)). The climatological CI is the largest over Eastern Siberia where winters are extremely cold and summers are warm. The western portions of Eurasia and North America have a lower mean CI than regions farther east, consistent with the strength of maritime influence in generally eastward atmospheric flow. Regarding tundra, eastern Eurasia and the Canadian Arctic have the largest CI values. The CI has decreased over most of North America and has mixed trends over Eurasia ( figure 5(b)). Trends over tundra regions are generally decreasing.
Time series of the CI display decreasing trends over Eurasian and North American tundra regions ( figure 5(c)). The CI trend over North America is −3.07 • C/38 years (significant at the 95% or greater level) and over Eurasia is −1.83 • C/38 years (table 1, top). North American tundra has become more maritime over the 1982-2019 period than has Eurasian tundra. The mean CI index over tundra is larger for Eurasia than North America (table 1, middle) and the CI index standard deviation is slightly higher for Eurasia than North America (table 1, bottom).

Precipitation
Quantifying precipitation with confidence in the Arctic remains a challenge but it is an important indicator for vegetation and will increase in importance as the climate warms. This study employs precipitation from ERA5, which appears to be a substantial improvement over previous products in the Arctic and less prone to discontinuities in Alaska compared to station observations (White et al 2020). The June-August (TSP) ranges between 5 and 20 cm over the tundra ( figure 6(a)). Mean precipitation is higher on the western versus the eastern sides of Eurasia and North America, which is consistent with the mean continentality gradient ( figure 5(a)). Climatological precipitation is higher in those parts of the Arctic at higher elevations. Trends in precipitation over the 1982-2019 period display areas of increasing and decreasing trends in the Arctic but are generally increasing over the tundra regions, particularly in North America ( figure 6(b)). The trend in tundra precipitation over North America is 1.15 cm/38 years (significant at the 90% or greater level) (table 1). The trend over Eurasian tundra is weak and not significant (0.40 cm/38 years). Similar to other indicators, the standard deviation in Eurasia is higher than North America for TSP, 1.58 versus 1.24 cm (table 1).

AO and AD index
The discussion of indicators in the Arctic would be incomplete if indices describing variability of atmospheric circulation, such as the AO and the AD, were not explored. These teleconnection indices can provide insight into the large-scale atmospheric drivers (e.g. Macias-Fauria et al 2012) that force SI which then influences tundra vegetation. Both indices are linear decompositions of sea level pressure using empirical orthogonal functions (EOF) or principal component analysis (PCA). This method separates differ parts of a signal into a few patterns that maximize variability and is similar to ordination of vegetation using PCA (Anderson 1971). The AO (Thomson and Wallace 1998) is the first EOF of atmospheric sea level pressure and the variance explained is 31.6% in DJF and 30.5% in JJA. The AD is the second EOF and the variance explained is 13.8% in DJF and 12.3% in JJA, less than half that of the AO. The AO is correlated with many physical and biological metrics in the Arctic (e.g. Baltzer et al 2005, Rigor et al 2002) but was weakly correlated with tundra indicators in both seasons, so we only present the AD results. The AD (Wu et al 2006) explains less variance than the AO but is important for SI variations and particularly the transpolar movement of SI exiting the Arctic which coincides with the transition zone between the anomalous high and low (figures 7(a) and (b)). Note that in the northern hemisphere low pressure has counterclockwise circulation while high pressure has clockwise circulation. The time series in figure 7(c) represents how the patterns vary over time. The magnitude trend over the 38 years of the study period is −0.58 for DJF-AD and 0.95 JJA-AD. Cai et al (2018) showed that the AD is linked to observed temperature and precipitation variations using the ERA-Interim and evaluated these relationships in a suite of CMIP5 models.
One key finding of Cai et al (2018) was that the AD had a larger impact on precipitation variations than on temperature. The DJF-AD and the JJA-AD were significantly correlated with some of the tundra indicators and this is explored further in the discussion of covariability of indices.

Discussion of covariability among indicators
Based on linearly detrended correlations between the indicators, we find that SI continues to influence tundra vegetation. Above average spring SI corresponds to below average OW, SWIa and TI-NDVI (table 2), consistent with previous studies (e.g. Dutrieux et al 2012, Bhatt et al 2017. However, the influence of SI on tundra indicators is not as strong as it was before based on a comparison of correlations over 1982-2019 to 1982-2008 (table 3). In Eurasia, correlations between SI have weakened with SWIs and TI-NDVI (table 3). In North America, correlations are weaker between SI and SWIs but have increased slightly with TI-NDVI (table 3). Additionally, the correlation between SWIs and TI-NDVI has weakened substantially in both Eurasia and North America (table 3). This suggests a weakening influence of SI and land surface temperatures on tundra vegetation productively and the emergence of other forcing factors over the last decade. Atmospheric moisture is a likely factor that may be gaining importance in the tundra story. In North America, significant correlations demonstrate that anomalously high OW results in above normal MaxNDVI and greater precipitation. TSP is weakly correlated in Eurasia and North America with the other indicators, leading us to conclude that moisture as a climate driver for tundra variability has not yet emerged but will as tundra regions become more maritime as the OW increases and the CI decreases.
The CI index was significantly correlated with only Arctic MaxNDVI and TI-NDVI (table 2) at 0.31 and 0.30, respectively (significant at the 90%  or greater level). Correlations between the Arctic CI index and spatial TI-NDVI are increasing throughout most of the tundra region though weak (<0.3, not shown). The stronger correlation at the pan-Arctic scale than the continental scale highlights the notion that key driving processes vary with region and scale, which adds challenges to identifying drivers (Myers-Smith et al 2020). The CI does not seem to contribute substantially at present to tundra vegetation productivity variability. However, a continued decrease of continentality driven by SI decline could make CI a more important climate driver in the future (table 1). TSP is significantly correlated with SWIs (−0.29, ⩾90%) in Eurasia and SI (−0.35, ⩾95%) and OW (0.45, ⩾95%) in North America (table 2). This suggests that above normal precipitation corresponds to reduced SI and enhanced OW in North America. This interpretation is confirmed by regression coefficients of TSP on Eurasian and North American OW (figures 8(a) and (b)), where positive precipitation anomalies correspond to above normal OW in their respective regions. While the correlations of indicators are very weak, precipitation is positively correlated with Max-NDVI and TI-NDVI for the Arctic, Eurasia and North America, where the correlations are somewhat stronger in Eurasia than North America. We speculate that the importance of precipitation as a climate driver for tundra NDVI may increase in the future as the climate warms, evapotranspiration increases, and permafrost thaw alters near-surface soil hydrology.  A new finding of this study is that SI and OW are more strongly correlated with SWIa than SWIs in all of the regions (table 2), where above normal SI corresponds to below normal SWIa and SWIs. Although the datasets underlying the two SWI indicators come from very different sources, the much weaker correlations of SWIs relative to SWIa suggest that the strong influence of vegetation and seasonal snow on the Arctic surface energy budget has dampened this relationship between SWIs and the other indicators from previous analyses.
The DJF-AD was significantly (⩾90%) correlated with SI, OW, MaxNDVI, TI-NDVI, CI, and TSP, while the JJA-AD was correlated with OW and TSP. This appears surprising at first but is consistent with the following sequence: the DJF-AD forces spring SI anomalies that then drive summer temperature, precipitation, and tundra NDVI variations. Correlations (table 2) indicate that the positive phase of the DJF-AD is associated with greater spring SI in the Eurasian 100 km coastal zone, reduced summer OW, below average Max and TI-NDVI, and below normal precipitation. Regressions provide a spatial perspective: positive DJF-AD departures are associated in Eurasia with positive SI and negative MaxNDVI anomalies and in the Alaska seas with negative SI and positive MaxNDVI anomalies ( figure 9(a)). The DJF AD drives SI variations which then force SWIa and NDVI anomalies in the following summer. The DJF AD connection with summer NDVI motivates examining precipitation in preceding seasons in the future as another potential climate driver. The two-season lag relationships provide the basis for multiseasonal predictions for which skill measures (i.e. how well a forecast performs) will exceed those of random chance and of climatology.
The positive phase of the JJA-AD enhances offshore flow over Eurasia (figure 7(b)) to increase OW, which is consistent with the significant positive correlations between JJA-AD and OW (table 2). This is confirmed from the spatial perspective with regressions of OW on JJA-AD ( figure 10(a)). Correlations of MaxNDVI and JJA-AD are weak but when regressed on JJA-AD ( figure 10(a)) the patterns reveal a varied structure over Eurasia. Eurasia west of 120 • E shows that more OW is associated with MaxNDVI declines, which appears inconsistent with positive correlations between NDVI and OW found in table 2. This is explored further by constructing regressions of SWI and TSP on JJA-AD. A consistent explanation emerges when all the variables are examined in the context of the JJA-AD in western Eurasia: the above normal OW and negative MaxNDVI anomalies correspond to positive precipitation (figure 10(b)) and negative temperature (figure 10(c)) departures. In the Laptev Sea the regression on JJA-AD indicates that enhanced OW corresponds to positive MaxNDVI anomalies (figure 10(a)), above normal precipitation (figure 10(b)), and above normal SWIa (figure 10(c)). OW has a spatially varied but consistent covariability in Eurasia with adjacent SWIa, MaxNDVI, and TSP. This result highlights the value of indicators at the large-scale and the need for caution when trying to generalize relationships regionally. In particular, the spatial pattern of the JJA-AD is such that offshore flow, above-normal temperatures and reduced SI are favored under the positive JJA-AD in eastern Eurasia (Laptev, East Siberian Seas), while onshore flow and colder-than-normal temperatures are favored over eastern Eurasia (Barents, Kara Seas).

Summary and conclusions
In this study, climate indicators known to be relevant for tundra productivity (i.e. coastal SI, coastal summer OW, SWI), two additional indicators (continentality and summer precipitation), and Arctic teleconnection indices (AO and AD) are analyzed for trends and co-variability in the context of Max and TI-NDVI. Over the study period of 1982-2019, significant decreasing trends in spring SI and significant increasing trends in summer OW, SWI (SWIs, SWIa) and MaxNDVI present a consistent story over Eurasia and North America. Despite a decline in recent years, TI-NDVI still displays overall increasing trends, though they are only statistically significant in Eurasia. Spatial TI-NDVI trends show declines in the Canadian High Arctic and southwest Alaska. High Arctic TI-NDVI declines may be due to thawing permafrost, melting ground ice and subsidence arising from high sensitivity to increasing SWIs because vegetation is sparse and organic layers are thin (Farquharson et al 2019). The CI has decreasing trends over tundra in Eurasia and North America, with significance (⩾95%) only for the latter. CI's relevance for tundra may increase since it is expected to continue to decrease in the future as warming is projected to be greater in winter than summer. TSP is increasing over Eurasia and North America, with the latter having weakly significant trends (⩾90%).
Correlation analysis identifies large-scale covariability between the indicators. The most significant correlations are between TI-NDVI and the climate indicator variables: above normal TI-NDVI corresponds to below normal SI (a spring indicator), and above normal OW, SWIs, and SWIa. TI-NDVI better reflects total seasonal productivity rather than peak biomass and is also a better tool for exploring changes in the growing season associated with more OW and changes in the timing of greenup and senescence. Correlations between SI and summer warmth and TI-NDVI have weakened when analyzed over the 1982-2019 period compared to 1982(Bhatt et al 2010 suggesting other processes are operating to change the influence of SI on tundra vegetation. The atmospheric mode of variability named the AD is defined as the second EOF of sea level pressure and was more closely correlated to the tundra indicators than the AO. In particular, winter AD is significantly correlated in Eurasia with SI, OW, MaxNDVI, TI-NDVI, CI and TSP. At first a winter connection to tundra productivity is perplexing but the link is through SI. The DJF-AD drives SI variations which shape OW, SWIa and NDVI anomalies during the following summer. This connection with DJF-AD has potential for use in prediction of Arctic vegetation variability.
Spatial regressions on JJA-AD demonstrate how the large-scale circulation favors regional climate anomalies that have competing effects on tundra vegetation. Above normal OW in the Kara-Barents Seas is associated with below normal MaxNDVI on adjacent land, above normal summer precipitation, and below normal SWIa during the positive phase of the JJA-AD. In contrast, in the Laptev Sea, above normal OW is associated with above normal MaxNDVI, TSP, and SWIa during the positive phase of the JJA-AD. This hints at potential competing effects of OW on TSP and SWIa and supports examining the climate drivers of tundra vegetation further at the regional scale.
It is challenging to synthesize climate driver understanding across scales: from plot to landscape to region to continent to hemisphere. The importance of a given climate driver likely varies with scale, and further work to understand across-scale relationships (e.g. Assmann et al 2020) is needed. Regional scale studies using the indicators framework may be fruitful, and the inclusion of other parameters such as cloud and snow cover would be an important next step. Clouds and snow cover are likely changing in the Arctic and both are relevant for tundra vegetation productivity. Clouds, precipitation, and surface snow cover are among the most challenging parameters for earth system models (Kushner et al 2018, McIlhattan et al 2020 to simulate, let alone predict, yet these are key for vegetation productivity. Advancing our understanding of climate drivers of tundra vegetation is relevant for Arctic prediction on seasonal-to-decadal scales and needed to anticipate future water and carbon budgets. While the indicator framework presented in this study shows the potential for seasonal prediction, climate model simulations will be required to quantify the likelihood of future changes in Arctic vegetation over decadal and longer timescales.

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.