Spring onsets of a young forest in interior Alaska determined based on time‑lapse camera and eddy covariance measurements

Spring phenology is essential in modeling the carbon balance of high‑latitude ecosystems and is possibly sensitive to climate change. In the present study, we evaluated the onset of the growing season for three species ( paper birch, bog blueberry, and bog Labrador tea ) in interior Alaska from 2012 to 2019 using photos taken using time‑lapse cameras. We also evaluated the onset of the growing season at the ecosystem scale from 2010 to 2019 on the basis of the CO 2 flux by the eddy covariance method at the site. On the basis of the growing degree‑day ( GDD ) model with the parameters estimated using the Bayesian approach, we found that the interannual variations in the spring onsets were explained by the model, and the thermal forcing requirement differed among the species. At the ecosystem scale, the spring onset was closely linked to the snow disappearance date. Under the possible future climate scenarios indicated by the representative concentration pathway 8.5 scenario, the spring onsets were predicted to be one to three weeks earlier than the present dates for the three species. The ecosystem‑scale onsets were also predicted to be five days to a little over a month earlier at the end of this century. The future spring onset is highly sensitive to the snow disappearance date for high‑latitude vegetation; thus, further understanding of climate change before snowmelting is required.


Introduction
Spring phenology is important in modeling the carbon balance of high-latitude ecosystems and is possibly sensitive to climate change Myneni et al., 2001;White et al., 1999 . Satellite remote sensing data indicate that the onset of the growing season became 0.32 days yr -1 earlier at 40 northward in North America from 1982 to 2005Zhang et al., 2007 . These changes are possibly caused by ongoing high-latitude warming Piao et al., 2020 . Despite its importance, the spring phenology of high-latitude vegetation is not fully understood Jarvis and Linder, 2000;Suni et al., 2003 , resulting in considerable uncertainties in predictions of the future carbon cycle Euskirchen et al., 2006;Holmberg et al., 2019 . Time-lapse cameras are useful for monitoring vegetation phenology at the species and ecosystem scales Ide and Oguma, 2010;Keenan et al., 2014;Kobayashi et al., 2016;Nagai et al., 2013;Xie et al., 2018 . Time-lapse cameras allow us to monitor phenology at smaller spatiotemporal scales than those monitored by satellite remote sensing and thus can be useful for filling the gap between satellite and field observations. The onset of the growing season as determined by camera observations was consistent with that detected by visual inspection at a site Keenan et al., 2014 , which was linearly related to the canopy leaf area index. Despite the wide availability of data, limited studies have used cameras for phenological studies in high-latitude vegetation Kobayashi et al., 2016;Nagai et al., 2013 . Previous studies revealed that the interannual variations in the spring onsets were controlled by spring temperatures at the ecosystem scale on the basis of analyses using camera images Ide and Oguma, 2010;Kobayashi et al., 2016 , but few studies analyzed such variations at the species scale Xie et al., 2018 . Eddy covariance measurements provide information on ecosystem-scale phenology. Previously, the CO 2 fluxes observed by the eddy covariance method were used to understand the mechanism of spring and autumn phenological transitions in 21 global forests Richardson et al., 2010and 29 forests Melaas et al., 2013. Melaas et al. 2013 found that the growing degree days GDD model could be useful to predict the spring onset of gross primary productivity GPP after considering the geographic variability in the photoperiod and thermal forcing accumulation required for spring onset. In boreal forests, Suni et al. 2003 evaluated the importance of five-day running air temperature for spring onsets on the basis of CO 2 fluxes using eddy covariance data. On the basis of eddy covariance measurements in three boreal forests, Barr et al. 2009 also found that soil thaw was an important driver of spring onset and that spring temperature primarily controlled annual production. In boreal forests, few studies have evaluated spring onset using both ecosystem-scale data from eddy covariance measurements and species-scale data from cameras. The use of both datasets could help us better understand each other s results. Boreal forests often have an open canopy structure, and the contributions of understories to CO 2 exchange are important Ikawa et al., 2015 . Thus, it is important to understand the spring onset of overstory and understory species as well as those at the ecosystem scale.
In this study, our objective is to quantify the temperature sensitivity of spring phenology at a young forest in interior Alaska. Toward the goal, we evaluated the spring onsets of individual species by camera images and ecosystem-scale CO 2 fluxes with eddy covariance measurements in a young forest in interior Alaska. Then, we explained the interannual variability in the onset for each species and ecosystem-scale CO 2 fluxes. The variability in the onset was evaluated on the basis of the GDD model with the parameters determined by the Bayesian approach. Finally, we predicted the spring onset at the end of this century. We discuss important environmental variables, such as air temperature and snow disappearance date, in future modeling activities for high-latitude regions.

Study area
The study site was a young forest after a wildfire in the Poker Flat Research Range in interior Alaska, USA 65 07'11"N, 147 25'44"W, altitude 491 m; hereinafter referred to as the US-Rpf site of AmeriFlux code; Ueyama et al., 2019 . The US-Rpf site was a burned black spruce forest after the wildfire continued from late June to early August 2004, where shrubs e.g., blueberry and Labrador tea and deciduous broad-leaved trees e.g., paper birch, aspen, and willow had regenerated by the time of observation. The deciduous trees grow in the upper layer, whereas evergreen trees black spruce , herbaceous plants, and mosses e.g., juniper hair cap grow in the lower layer. At the site from 2010 to 2019, the mean annual air temperature ranged from 2.6 C to 1.3 C, and the annual rainfall ranged from 231 mm yr -1 to 516 mm yr -1 . The snow disappearance date was generally observed from the end of April to the middle of May. Further details on the study site were available for CO 2 flux Iwata et al., 2011;Ueyama et al., 2019, energy flux Ueyama et al., 2020, and leaf photosynthetic characteristics Ueyama et al., 2018 .

Meteorological observations
We used meteorological and eddy covariance data from 2010 to 2019Ueyama et al., 2019 . We briefly describe the measurement system used in this study. Air temperature was measured at 1 m above the ground using a temperature and humidity sensor HMP45, Vaisala, Finland mounted in a ventilated radiation shield. Downward and upward shortwave radiations were measured using pyranometers PCM-01, CMP3, and CNR4, Kipp & Zonen, Netherlands . CO 2 flux was measured by the eddy covariance method using several types of sonic anemometers and gas analyzers. Intercomparisons between observation systems were conducted, and systematic differences were adjusted Ueyama et al., 2019 . The meteorological and turbulence data were collected using a data logger CR3000, Campbell Scientific Inc., USA .
Net ecosystem exchange NEE was calculated by accounting for a storage change determined from the concentration at the height of eddy covariance observation. The artificial negative CO 2 flux measured by the open-path gas analyzer LI-7500, Li-Cor, USA was corrected on the basis of Amiro 2010 . The correction was not applied for the new analyzer EC150, Campbell Scientific Inc., USA Helbig et al., 2016 . Further details for the observations and data processing are shown in Ueyama et al. 2019 .
We also used the precipitation from the closest weather station at Fairbanks International Airport, Alaska 40 km apart from the site from 2010 to 2019 by the National Climatic Data Center.

The time-lapse camera observation and data analysis
We used the photos taken by the GardenWatchCam time-lapse camera Brinno, Taiwan; resolution of 1280 1040 px from 2010 to July 2016 and the TLC2000 time-lapse camera Brinno, Taiwan; resolution of 1280 720 px from July 2016 to 2019. The former cameras were mounted horizontally at a height of 1 m. The latter cameras were mounted at heights from 2 m to 6 m facing down to the ground depending on the vegetation regrowth. Photos were taken every hour during the daytime using an automatic white balance. We extracted the region of interest from the photo where one species dominated, for each species Fig. 1 . The color digital values of red R , green G , and blue B in the region of interest were averaged and then used to calculate the green color index g cc : the relative proportion of G to the sum of R, G, and B Harazono et al., 2009;Keenan et al., 2014;Xie et al., 2018 . To reduce artificial fluctuations in g cc due to different brightness values in different weather conditions, the median g cc was extracted from a 3-day moving window for each day.
We used g cc data after the snow disappearance date because the onset of the growing season for each species occurred after the snow disappeared. The snow disappearance date was determined as the date when the measured albedo was continuously less than the typical growing-season albedo at this site 0.15 . We estimated the date of spring budburst or greening for the three species using time-lapse camera images. We used the definition of the spring onset on the basis of the so-called threshold method according to Xie et al. 2018 . The spring onset was determined as the first day when g cc normalized with the seasonal maximum value successively exceeded 0.05 in the period from the snow disappearance date to the end of July.

Study species
We investigated the spring onsets of paper birch Betula papyrifera , bog blueberry Vaccinium uliginosum , and bog Labrador tea Ledum groenlandicum on the basis of camera images. The paper birch is classified as deciduous broad-leaved trees, and bog blueberry is classified as a deciduous understory. Bog Labrador tea is an evergreen understory species but showed obvious greening in the spring and browning in the autumn.
Because of the data availability, the study period and sample size differed in the species. The camera images for the paper birch were available for one sample in 2012, 2013, and 2017, two samples in 2018, and three samples in 2019. Images of bog blueberries were available for one sample in 2013, 2014, 2015, and 2018. The images of bog Labrador tea were available for one sample in 2012, 2013, 2014, and 2015. Samples used for one species varied from year to year.

Data analysis for each species
We used the GDD model Cannell and Smith, 1983;Leinonen and Kramer, 2002;Melaas et al., 2013 to explain the onset of the growing season for each species. In the model, the onset was determined as the day when GDD exceeded the cumulative temperature required for spring onset GDD c ; hereinafter referred to as the thermal forcing requirement . The GDD on a given day t 1 was calculated by cumulating the mean air temperature above a degree-day base temperature T b after the day on which accumulation started t 0 . We set the starting date t 0 as the snow disappearance date. Although chilling requirement in winter is essential for spring onset in the temperate region Ettinger et al., 2020 , it could be easily met under long and cold winter in interior Alaska; hence, we did not consider the chilling requirement in this study.
We determined two parameters T b and GDD c for each species using the same data with the differential evolution method Das and Suganthan, 2010 in a preliminary experiment. Based on the experiment, the T b varied little among the three species, where the T b was estimated as 0.3 C for paper birch, 1.8 C for bog blueberry, and 0.3 C for bog Labrador tea, and the GDD c was estimated as 42.7 C days for paper birch, 15.9 C days for bog blueberry, and 25.5 C days for bog Labrador tea. Hence, we fixed T b for all species to the lowest T b , 0.3 C, and determined only GDD c that varied greatly among the three species to minimize equifinality in the parameterization. We estimated the probability density distributions of GDD c using a constraint of observed g cc -based spring onset on the basis of the Bayesian approach using the library PyMC3 Version 3.7 in Python. We determined the posteriori distributions of the parameter from the initial maximum posterior probability using the Metropolis -Hasting method. In the parameterization, we used 10,000 samples after 1,000 samples of warm-up and then evaluated the 95 Bayesian highest posterior density HPD intervals to identify uncertainties in the parameter.

Data analysis for the ecosystem scale
We determined the onset of the growing season for NEE or GPP to explain the onset at the ecosystem scale. We used the gap-filled NEE and GPP from 2010 to 2019, where the data gap was filled using a combined look-up-table method and nonlinear regression method Ueyama et al., 2019 . GPP was partitioned from NEE on the basis of a nighttime approach using the Flux Analysis Tool program version 2.0 Ueyama et al., 2012 . The onset of NEE was determined as the date when the seven-day moving mean NEE continuously showed negative values. The onset of GPP was the date when the seven-day moving mean GPP was continuously higher than 10 of its annual amplitude Melaas et al., 2013 . We applied a linear regression to explain the onsets of GPP and NEE using the snow disappearance date because the onsets were linearly correlated with the snow disappearance date shown in section 3.3 . We did not apply the GDD model to explain the ecosystem-scale onsets because 1 the onsets occurred earlier than the snow disappearance date but 2 the parameter t 0 was set as the snow disappearance date. We determined the linear regression based on the Bayesian approach with a no-U turn sampler. In the parameterization, we used 2,000 samples and then evaluated the 95 Bayesian HPD intervals.

Future predictions of spring onsets
We predicted the future spring onsets for each species, NEE, and GPP from 2010 to 2099 under the representative concentration pathway RCP 8.5 scenario. We chose the RCP8.5 scenario because our aim was to evaluate the possibility when global warming was the most severe. We used simulated climate data for the selected pixel for the study site, which included the predicted air temperature and precipitation, on the basis of five climate models HadGEM2-ES, IPSL-CM5A-LR, MIROC-ESM-CHEM, GFDL-ESM2M, and NorESM1-M produced in the Coupled Model Intercomparison Project 5 Hempel et al., 2013 . In the climate data, the mean air temperatures were 2.4 C-7.7 C higher in the 2090s than in the 2010s, whereas the annual precipitation was 55-157 mm yr -1 larger. We predicted the snow disappearance date using the snow model embedded in a terrestrial ecosystem model, VISIT Ito and Inatomi, 2012 , by inputting the climate forcing. The snow model estimates the snow water equivalent at the daily timestep using air temperature and precipitation. We found that the model accurately estimated the snow disappearance dates at the site from 2010 to 2018 Fig. 2 .

Spring onset based on g cc for each species
The g cc for the three species increased after snowmelt in the spring Fig. 3 , reflecting the onset of the growing season.  This result was caused by the budburst in paper birch and bog blueberry and by the greening of the leaves of bog Labrador tea. The background of the target plant, such as the forest floor, had little influence on the spring increase in g cc because we did not analyze plants whose background color changed significantly Fig. 1 . The timing of the onset due to the increasing g cc differed among the species. Bog Labrador tea started greening earliest, and the budburst of paper birch and bog blueberry followed. The different onsets represented the species-specific response to environmental conditions during the spring. Duplicated samples of the onset n = 2 or 3 were collected for only paper birch in individual years, showing that the difference in the onset among samples was up to three days within a year.

Parameterization of the GDD model for each species
We found that the interannual variations in the spring onset were well explained by the GDD model R 2 = 0.71; p < 0.05 for paper birch, R 2 = 0.65; p < 0.05 for bog blueberry, R 2 = 0.87; p < 0.05 for bog Labrador tea , and the thermal forcing requirements differed among the species. The mode and 95 Bayesian HPD intervals of the thermal forcing requirements were estimated to be 17.7 2.5-40.8 C day for bog Labrador tea, 70.0 51.5-80.4 C day for paper birch, and 156.0 116.8-193.0 C day for bog blueberry Fig. 4 . These intervals showed significant differences in the parameters among the species.
The GDD model reproduced the interannual variations in the observed onset for each species Fig. 5 . The earliest onset in 2016 was explained by the warmest air temperature in May 9.7 C , and the latest onset in 2013 was explained by the coldest temperature in May 4.2 C for all three species. The 95 Bayesian HPD intervals of the root-mean-square error were up to 3.9 days for bog Labrador tea, up to 2.7 days for paper birch, and up to 4.4 days for bog blueberry.
The GDD model generally reproduced the interannual variations in the observed onset for each species Fig. 5 . The respective standard deviations of the observed and calculated onsets were 5.4 and 5.2 days 5.0-5.5 days for paper birch, 7.5 and 8.8 days 8.4-9.3 days for bog blueberry, and 11.2 and 10.5 days 10.0-11.0 days for bog Labrador tea, where the values in the calculated onsets were the mean and the 95 Bayesian HPD interval. The interannual variability was consistent for paper birch and was slightly overestimated or underestimated for the other species. The interannual variability was smallest for paper birch, followed by bog blue berry and then bog Labrador tea, which were well reproduced by the model. Generally, all species basically showed synchronized interannual variations at spring onset. On the basis of a linear regression between mean air temperature from mid-April to May and the spring onset, we estimated the temperature sensitivity in the spring onset in terms of the interannual variations. The estimated sensitivity of the spring onset to air temperature from mid-April to May was 1.9 days C -1 for paper birch, 1.8 days C -1 for bog blueberry, and 2.9 days C -1 for Labrador tea.

Spring onset for the ecosystem scale
The onset of the growing season for NEE and GPP was positively correlated with the snow disappearance date Fig. 6 . The linear regression for the onset of NEE had a slope of 0.89 days day  linked to each other. The onset of the CO 2 fluxes occurred earlier than those for the deciduous species paper birch and bog blueberry Fig. 7 . The sensitivity of the ecosystem-scale onset to air temperature from mid-April to May was 3.5 days C -1 for NEE and 3.0 days C -1 for GPP. The similar relationships for NEE and GPP indicate small variability in ecosystem respiration in the spring. The contributions of GPP before budburst of birch to the annual GPP were negatively correlated with the mean air temperature in May Fig. 8 . We calculated the contributions as a ratio of GPP before the budburst of birch Fig. 5a to the annual GPP. The contributions varied each year, ranging from 38.6 in the cold year 2013 to 7.8 in the warm year 2019 , showing that a cooler spring induced higher contributions of understory in the spring, but a warmer spring caused lower contributions 5.2 C -1 ; R 2 = 0.91 .

Future predictions of the spring onset
Under the RCP8.5 scenario, the spring onset for the three species was predicted to be one to three weeks earlier at the end of this century Fig. 9a . The advanced onset in the future was predicted to be 17.8 days 8.0-24.1 days for paper birch, 18.1 days 8.1-22.3 days for bog blueberry, and 18.1 days 7.9-27.3 days for bog Labrador tea; the mean and the 95 confidence interval were estimated by considering the uncertainties in the determined parameter and climate change simulations. Because of the large variabilities in the climate by the different models, the differences in the advanced onset were not clear among the species.
The advanced spring onset for the three species was predicted to occur due to earlier snow disappearance because the snow disappearance date was estimated to be 7.8-33.7 days earlier based on the snow model Fig. 10 . By contrast, the days required for accumulating thermal forcing did not substantially change under the scenarios. The mean and 95 confidence interval of the days for accumulating thermal forcing were predicted to be 1.4 days 6.8-2.2 days for paper birch, 0.9 days 7.9-3.7 days for bog blueberry, and 0.9 days 4.2-1.4 days for bog Labrador tea. The spring onset at the ecosystem scale was predicted to  be five days to a little over a month earlier at the end of this century under all RCP8.5 scenarios Fig. 9b . Considering the uncertainty in the climate simulations, the mean and 95 confidence interval of the spring onset were predicted to be 18.7 days 5.6-38.6 days for NEE and 16.4 days 4.9-34.3 days for GPP. This earlier onset was predicted because of the earlier snow disappearance because the onset at the ecosystem scale was explained through the linear regression using the snow disappearance date Fig. 6 . The advancements for the ecosystem-scale onset were similar to those for each species Fig. 9a .

Discussion
The combination of a time-lapse camera with g cc was a useful approach for determining the species-specific onset of the growing season Fig. 3 . The g cc was previously used to evaluate ecosystem-scale phenology Ide and Oguma, 2010;Keenan et al., 2014;Kobayashi et al., 2016 andspecies-specific phenology Xie et al., 2018 and to estimate CO 2 flux Ahrends et al., 2009;Harazono et al., 2009;Ide et al., 2011;Ueyama et al., 2013 . Since the CO 2 balance in boreal forests is often contributed by not only overstory vegetation but also understory vegetation Ikawa et al., 2015 , the species-specific phenological onset for both overstory and understory species could be useful for understanding the CO 2 balance in boreal forests.  The GDD model successfully explained the spring onset for the boreal species Fig. 5 . We confirmed that the species-scale onsets determined by the GDD model were better reproduced than those determined by the regression using the snow disappearance date. The onsets for the three species were well regulated by air temperature after snowmelt Fig. 5 , meaning that warming could similarly influence the onsets of the three species Fig. 9a . Photoperiod could have little effect on the onsets at this site because photoperiod around the spring onset at this site is longer than previously reported requirements 9.6-10.5 h for boreal deciduous forests Melaas et al., 2013 . Despite the simple structure, the GDD model was used to successfully evaluate spring onset in various plants Vitasse et al., 2011 and eddy covariance-based fluxes Melaas et al., 2013 . The good performance of the model was reported for high-latitude ecosystems on the basis of a satellite-based regional-scale study Fu et al., 2014a . The ecosystem-scale spring onset was closely linked to the snow disappearance date Fig. 6 and was earlier than that for the deciduous species Fig. 7 . This difference was possibly because evergreen shrubs including black spruce seedlings and mosses could start photosynthesis during and immediately after snowmelt at the study site. After snowmelt, surface soils ~10 cm depth immediately thawed and reached above-zero temperatures. The onset of plant productivity is known to be linked to liquid water availability, which is related to the presence of above-zero soil temperatures Jarvis and Linder, 2000 . Hence, the ecosystem-scale onset was sensitive to snowmelt, which determines energy loading into soils and light availability for shrubs and mosses. This result indicated that the ecosystem-scale onset was determined by the meteorological conditions in late winter. Conversely, the onset for deciduous species was regulated by cumulating air temperature after snowmelt. The species-scale onset was determined by the snow disappearance date and air temperature after snowmelt. Since disappearing snow increases net radiation and regional air temperatures Ueyama et al., 2020 , the warming rate before and after snowmelt could differ even under similar forcing changes. This result indicated that the prediction of species-scale onsets could be more complex than that of ecosystem-scale onsets. Hence, the different sensitivities of the species-and ecosystem-scale onsets could be essential in the boreal CO 2 balance.
The warming spring temperatures could regulate photosynthesis in the understory due to low light availability during the early spring. In the period between snowmelt and budburst in the future, daily incoming solar radiation will decrease and day length will shorten because the period will move forward from summer solstices. After the expected earlier budburst of overstory vegetation, light availability for the understory will be reduced for the rest of the growing season. This phenomenon could result in decreased contributions of understory productivity into the warmer spring Fig. 8 .
The advanced future onset of the growing season was consistent with those predicted by land surface models and field-based estimates. Holmberg et al. 2019 predicted the earlier onset of boreal regions at the end of the 21st century by approximately three weeks using the model. Euskirchen et al. 2006 predicted 43.4 days longer growing seasons in boreal and Arctic regions from 1976 to 2100 on the basis of the model, mostly due to the earlier spring onset. Based on our study, the predicted future sensitivity 28.4 days advancement under a 1.6-7.2 C increase; hence, up to 3.9 days per degree Celsius was consistent with those estimated on the basis of trans-Alaska camera analysis 2.5 to 3.9 days C -1 ; Kobayashi et al., 2016 . On the basis of the sensitivity by Kobayashi et al. 2016 with the RCP8.5 climate scenario, the spring onset was roughly assumed to be a week to a month earlier in the boreal region at the end of this century, which is similar to our estimates. We found that snow disappearance is an essential process for determining future spring onset for boreal forests. The earlier onset could accelerate GPP and CO 2 uptake at the ecosystem scale because of the prolonged growing season Euskirchen et al., 2006 , although earlier onset possibly results in earlier autumn senescence Fu et al., 2014b;Keenan and Richardson, 2015 . The earlier snow disappearance leads to an earlier decrease in albedo in the spring, resulting in further warming at the regional scale with changes in surface energy exchange Euskirchen et al., 2007;Ueyama et al., 2014Ueyama et al., , 2020 .

Conclusions
We evaluated the species-specific onset of the growing season for three species in interior Alaska using time-lapse cameras and the ecosystem-scale onset of CO 2 fluxes. We found that the species-specific onset was explained by the GDD model with a species-specific thermal forcing requirement. The species-specific phenology indicated that careful considerations of species-specific parameters could be necessary for estimating the consequences of spring warming in each plant, although the predicted future change in the onset was unclear because of high uncertainties in climate scenarios. The ecosystem-scale onset was solely explained by the snow disappearance date, which suggested that the snow disappearance date and the onset of the CO 2 fluxes were closely linked to each other. This result was possibly because thawing of surface soils immediately after snow disappeared could be important for photosynthesis. The understory contributions before the budburst of the birch decreased with warming spring temperatures, suggesting that the earlier budburst with a shorter day length could increase the regulation of understory production. Under the RCP8.5 scenario, the onset was predicted to be one to three weeks earlier for the three species and five days to a little over 1 month earlier at the ecosystem scale at the end of this century. The future onset is highly sensitive to the snow disappearance date, and hence, we recommend further understanding of climate change during the snowmelting season.