How temperature regimes near the equinox synchronize spring biological events

Many biological processes, including plant leafout and flowering, occur once cumulative temperatures reach a threshold (the thermal-sum model). In this way, temperatures are thought to coordinate the timing of biological events. But growing evidence suggests that as climates warm, both the advancement of spring has slowed (declining sensitivity) and the variance in the timing of spring events has increased (declining synchrony)—raising questions about the resilience of temperature-based coordination to anthropogenic climate change. To answer these questions, researchers have complicated the thermal-sum model, introducing additional factors and mechanisms. We consider whether such complexity is necessary. Using results from the theory of stopped random walks, we show that sensitivity and synchrony are exactly as predicted by the basic thermal-sum model. The theory suggests a nonlinear relationship between temperatures and both the timing and synchrony of biological events. In particular, it predicts that as temperatures increase and springtime events shift from the equinox toward the solstice, the events themselves become less coordinated and more variable. We verify these predictions using experimental and real-world data, including 10,000 observations of common lilacs (United States, 1956-2025). We conclude that the theory provides a powerful tool for understanding the thermal-sum model—particularly when considering additional complexity.

Introduction

Many biological events are triggered not by the temperature on a single day, but by the temperature accumulated over many days (Larcher, 1980). A canonical example is the law of the flowering plants, which states that a bud blooms in the spring once cumulative temperatures exceed a plant-specific threshold. This and similar mechanisms are thought to synchronize a wide variety of seasonal events in an ecosystem—for instance, ensuring that pollinators emerge when flowers are in bloom. They also support an array of forecasting tools in which scientists track cumulative temperatures or ‘growing degree days’ to predict harvest dates, manage pest emergence, and assess how shifting temperature patterns reshape ecological communities with climate change (Schwartz, 2024).

Shifts in the timing of biological events have been the clearest fingerprint of anthropogenic climate change across both natural and managed systems (Pörtner et al., 2022). But these shifts have proven difficult to fully predict or explain. The timing of phenological events have occurred earlier (i.e., advanced) as climates warm, however, the rate of advancement has slowed in recent years (Fu et al., 2015). At the same time, the variance in the timing of phenological events has in some cases declined and in others increased, without a clear pattern (Vitasse et al., 2018; Stemkovski et al., 2023). Higher variance has also manifested as previously unknown or rare phenological events becoming more common (Chuine et al., 2025). While research has provided myriad explanations—often with contrasting predictions and logic (Fu et al., 2015; Wolkovich et al., 2021; Gao et al., 2024)—no work has examined the extent to which these observations are consistent with the basic thermal-sum model in which a biological event is triggered once cumulative temperatures reach a threshold.

Scientists have used the basic thermal-sum model for nearly two centuries (Quetelet, 1849), but have rarely if ever integrated the relevant mathematical tools developed to characterize its behavior—specifically results from renewal theory designed for analyzing stopped random walks. Here we draw on these results to predict how events triggered by a cumulative temperature sum will shift with warming. Our model explains both the observed deceleration in the advancement of springtime events as well as the shifts in the variance.

Throughout we emphasize the critical role of the temperature regime under which the biological events take place. For example, the variance in the timing of biological events can decrease or increase with warming depending on whether temperatures during the period of accumulation are stationary (as in the winter and summer near the solstice) or increasing (as in the spring near the equinox). We find that climate change can (i) move spring events earlier at a faster or slower rate, and (ii) increase the variance of their timings, producing less synchronized events as they arrive earlier, a phenomenon we refer to as a ‘scattered spring.’

Results & Discussion

Many springtime events occur when cumulative temperatures first reach or ‘hit’ a plant-specific temperature threshold. We use results from renewal theory (Siegmund, 1967; Woodroofe, 1982; Gut, 2009) to approximate the distribution of event times (hitting times) and thus predict the behavior of temperature-triggered biological events.

Consider one such event, the bloom date of a bud. We model the effective temperature $X_{i}$—the daily temperature experienced by the bud on day $i$—in two temperature regimes, visualized in Figure 1. The first is the stationary temperature regime in which the average daily temperature $\alpha$ > 0 is constant. The stationary regime is commonly seen in the winter months following the solstice, and we refer to the temperatures in this regime as ‘winter temperatures.’ The second regime is the increasing temperature regime in which the average daily temperature increases at constant rate $\beta$ > 0. The increasing temperature regime is commonly seen in the spring months closer to the equinox, and we refer to it as ‘spring warming.’ In both regimes, the effective temperature is the daily average temperature plus noise $\epsilon_{i}$. That is, $X_{i}=\alpha+\beta\,i+\epsilon_{i}$ where $\beta=0$ in the winter temperatures regime and $\beta>0$ in the spring warming regime.

Figure 1: We model the daily temperature (solid line) as trend (dashed line) plus noise. The timing of springtime events in temperate climates depends on whether temperature is accumulated during one of two regimes, 1. the stationary temperature regime (relatively flat near the winter solstice), which we refer to as ‘winter temperatures,’ or 2. the increasing temperature regime (relatively linear near the spring equinox), which we refer to as ‘spring warming.’

Noise in this formulation could represent many biological factors depending on the scale and our developing understanding of what triggers springtime events. At the bud-level, it represents idiosyncratic variation at that bud, which could be driven by microclimate variation (Schwartz et al., 2014), bud color variation (Peaucelle et al., 2022), or internal differences that vary when the bud starts accumulating temperatures (van der Schoot et al., 2014; Pan et al., 2023). We assume the $\epsilon_{i}$ are independent with mean $0$ and variance $\sigma^{2}$, although this assumption can be relaxed.

The bloom date $\nu(\tau)$ is the day $n$ that the cumulative temperatures $Z_{n}\;=\;\sum_{i=1}^{n}X_{i}$ first reach or ‘hit’ a plant-specific temperature-sum threshold $\tau>0$, i.e. $\nu(\tau)\;=\;\min\{m:\,Z_{m}>\tau\}$. This formulation is common in growing degree day models, and more generally, $Z_{n}$ is the ‘forcing’ part of most process-based models of plant phenology (Schwartz, 2024). Standard results for stopped random walks yield asymptotic normal approximations for the bloom date $\nu(\tau)$ when the thermal-sum threshold $\tau$ is sufficiently large that accumulation takes place over many days.

These results predict that the mean and variance of bloom dates depend on the underlying temperature regime (i.e., values of $\alpha$ and $\beta$) during which the plant accumulates temperature. In the winter temperatures regime, average temperatures are relatively flat during accumulation ($\beta=0$), and cumulative temperature grows at an approximately linear rate. It follows that expected bloom date advances with climate change in step with increasing winter temperatures ($\mathbb{E}[\nu(\tau)]\approx\tau/\alpha$) as is often described in the literature (Vitasse et al., 2022). In the spring warming regime, average temperatures increase during accumulation ($\beta>0$), and cumulative temperatures grow quadratically. In this regime, bloom date advances as a function of the rate of spring warming ($\mathbb{E}[\nu(\tau)]\approx\sqrt{2\tau/\beta}$).

Within either regime, the relationship between bloom dates and temperature is non-linear as observed in recent papers (Fu et al., 2015; Wolkovich et al., 2021). Increasing temperatures advance bloom dates at a diminishing rate as warming (measured by the winter temperatures $\alpha$ or spring warming $\beta$) increases (for the winter temperatures regime: $\frac{\partial}{\partial\alpha}\,\mathbb{E}[\nu(\tau)]\approx-\tau/\alpha^{2}$; for the spring warming regime: $\frac{\partial}{\partial\beta}\,\mathbb{E}[\nu(\tau)]\approx-\frac{1}{2}\sqrt{2\tau/\beta^{3}}$). But when climate change shifts the period when plants accumulate temperatures from the spring warming regime to the winter temperatures regime, advances may actually accelerate, particularly if increases in winter temperatures $\alpha$ outpace increases in spring warming $\beta$.

Importantly, the two regimes make very different predictions regarding the synchrony of events—that is the variance in event timings within the same year and location. Both predict decreased variance with warming (for the winter temperatures regime: $\mathrm{Var}(\nu(\tau))\approx\sigma^{2}\,\tau/\alpha^{3}$; for the spring warming regime: $\mathrm{Var}(\nu(\tau))\approx\sigma^{2}/\sqrt{2\beta^{3}\tau}$), but the influence of noise is strikingly different as reflected by whether the threshold $\tau$ is in the numerator or the denominator.

In the spring warming regime, the threshold $\tau$ is in the denominator, and thus the variance is small when the threshold is large, reflecting the fact that noise is quickly overwhelmed by the quadratic increase in cumulative temperatures that occurs during spring warming. This effectively acts to synchronize plants (and buds within plants) with the same threshold—making idiosyncratic differences such as microclimate negligible. In contrast, under the winter temperatures regime, the threshold $\tau$ is in the numerator, and thus the variance is large when the threshold is large. The interpretation is that noise can remain relatively influential so that idiosyncratic differences translate into larger differences in threshold-crossing times.

The fact that the variance is predicted to be small in the spring warming case and larger in the winter temperatures case suggests that the observed variance will increase when anthropogenic warming shifts accumulation from the spring warming regime to the winter temperatures regime. The shift in regimes is observable in data as a ‘scattered spring’ in which springtime events simultaneously advance and become less coordinated.

Not all warming scatters spring. Whether warming leads to a more synchronized or more scattered spring depends on the local joint evolution of winter temperatures and spring warming ($\alpha$ and $\beta$) and on the threshold ($\tau$) relevant for a given species and life-cycle event. In some locations, climate change may manifest primarily as higher winter temperatures (an increase in $\alpha$) with relatively little change in the rate at which temperatures rise during spring (a modest change in spring warming $\beta$). In others, winter temperatures may change less while spring transitions steepen (a larger increase in spring warming $\beta$).

The stopped-random-walk framework offers a simple way to translate this spatially varying climate change into spatially varying predictions of both the advancement of spring and its variability. To demonstrate the framework, we model bloom dates from the lilac–honeysuckle phenology network compiled by Rosemartin et al. (2015), collected across the continental United States from 1956–2014. We limit our analysis to the purple common lilac (Syringa vulgaris) and the phenophase full bloom, and we update the dataset using USA National Phenology Network records from the same sites and species until 2025 (Switzer et al., 2025).

Figure 2: Mean and log standard deviation of lilac bloom dates (n = 10,613 from Rosemartin et al. 2015, updated to 2025 using USA National Phenology Network data) by winter temperatures $\alpha$ and spring warming $\beta$ as estimated by a generalized additive model. The mean (left) decreases as $\alpha$ and $\beta$ increase. The standard deviation (right) can increase or decrease depending on the relative values of $\alpha$ and $\beta$.

Figure 3: Mean (top) and log standard deviation (bottom) of lilac bloom dates (n = 10,613 from Rosemartin et al. 2015, updated to 2025 using USA NPN records) by winter temperatures $\alpha$ and spring warming $\beta$ as estimated by a generalized additive model. The top two figures show the mean generally decreases as $\alpha$ increases holding $\beta$ constant (left) and $\beta$ increases holding $\alpha$ constant (right). The bottom two figures show the standard deviation can increase or decrease depending on whether $\alpha$ increases with $\beta$ held constant (left) or $\beta$ increases with $\alpha$ held constant (right).

For each record, we estimate winter temperatures $\alpha$ and spring warming $\beta$. We then estimate the conditional mean and log standard deviation functions of the bloom date given $\alpha$ and $\beta$ using a generalized additive model (Hastie and Tibshirani, 1986; Wood, 2017). The two functions are displayed in Figure 2. Select cross sections of these functions are displayed in Figure 3. Note that the bloom date is measured in days since January 1 of the observation year.

These two figures show that increasing either winter temperatures $\alpha$, spring warming $\beta$, or both shift flowering times exactly as predicted by the theory. Increasing $\alpha$ or $\beta$ is associated with earlier mean bloom dates, which generally but not always advance at a decreasing rate.

Increasing $\alpha$ and holding $\beta$ fixed is associated with higher log standard deviation (higher variance), while increasing $\beta$ while holding $\alpha$ fixed is associated with a lower log standard deviation (lower variance). This is the ‘scattered spring’ predicted from the theory: When $\beta$ is fixed, increasing $\alpha$ means more buds bloom in the winter temperatures regime, under which conditions the variance is predicted to be larger. When $\alpha$ is fixed and $\beta$ is increased, the percentage of buds blooming in either regime does not change. Rather the fixed percentage of buds that bloom under spring warming are subject to increased warming, and under these conditions the variance is predicted to be smaller.

We also apply the framework to experimental data from Charrier et al. (2011) in which stems were sampled from 30 walnut trees (Juglans~sp.) trees, chilled at $4^{\circ}$C, and then warmed at various constant temperatures until budburst.

Since temperatures are held constant, this experiment replicates the winter temperatures regime in which the bloom date is approximately normal with mean $\tau/\alpha$ and variance $\sigma^{2}\,\tau/\alpha^{3}$, which we fit using weighted least squares. Figure 4 shows the observed data (points), model fit (solid line), 95% confidence interval (inner grey region, approximately two standard errors), and 95% prediction interval (outer grey region, approximately two standard deviations). The data fit the model well, exhibiting the nonlinear relationship in the mean and variance as predicted by the theory. The relationship cannot be explained by factors like chilling and photoperiod, as these are held fixed.

Figure 4: Mean time until budburst as a function of the average daily temperature in a controlled experiment conducted by Charrier et al. (2011) in which stems were sampled from 30 walnut trees (Juglans~sp.) trees, chilled at $4^{\circ}$C, and then warmed in different environments (x-axis) until budburst. Constant warming reflects the winter temperatures regime. Conditional mean (solid line) estimated using weighted least squares. Inner interval represents the confidence interval (approximately twice the standard error) while the outer interval represents the prediction interval (approximately twice the standard deviation). The figure shows that the mean and variance are nonlinear even after factors like chilling and photoperiod are held fixed. Both decrease with increasing temperatures as predicted by the thermal-sum model.

As a check, we also conduct both analyses using an alternative approach in which we bin the data by quartile of $\alpha$ and $\beta$, and we reach similar conclusions. See Appendix for details. We conclude that the basic thermal-sum model explains the observed shifts in springtime events. In particular, it explains why warming is associated with both decelerating advances and either increasing or decreasing variances (Vitasse et al., 2018; Fu et al., 2015; Stemkovski et al., 2023). We stress that this explanation makes minimal biological assumptions, only that the events are triggered by a thermal-sum, the same assumption made by most growing degree day or spring forcing models (Schwartz, 2024).

The model does not preclude complicated explanations. However, given the fact that the thermal-sum model is parsimonious and nearly two centuries old (Quetelet, 1849), the burden falls on complicated models to justify their complexity. For example, winter chilling and photoperiod are often invoked to explain the observed shifts in springtime events (Fu et al., 2015). These factors are already implicitly accounted for in the thermal-sum model, which allows for idiosyncratic variation. Our results question whether explicitly modeling these factors is necessary and thus justified.

More generally, the theory presented provides opportunities for new experiments and observational insights on which the basic thermal-sum model may be improved. As it stands, the basic thermal-sum model predicts a general desynchronization of spring events as their timing moves closer to the solstice. Absent some compensating mechanism, the timing of biological events across the ecosystem—from crops to forests, and from plants to insects—will in many cases become more variable, with potential cascading consequences for industry and natural ecosystems.

Methods

Our observational data analysis considers all sites in the historical lilac–honeysuckle phenology network (Rosemartin et al., 2015) that recorded the phenophase full bloom and were within 10 miles of a GHCND site with sufficiently complete temperature records. We limit our analysis to common lilac (Syringa vulgaris). We include all observations made from 1955 to 2025 as recorded by the USA National Phenology Network, retrieved using the R package rnpn, resulting in $10,613$ observations.

For each observation, we estimate $\alpha$ using the average daily temperature between January and February. We estimate $\beta$ using the slope of a linear regression model fit by regressing daily temperature on day index over March and April. Temperatures below $0$ are set to $0$, a common base temperature in growing degree day models.

The parameters $\alpha$ and $\beta$ are shown for nine sites in Figures 5-7 in Appendix A.3. Figure 5 shows nine randomly selected site-years and has the same interpretation as Figure 1. Site 14726 in New Mexico (Year 1969) is a good example of a bloom in the spring warming regime since relatively little temperature was accumulated during January and February. The average bloom date was May 8th with a standard deviation of 17 days. Site 14276 in New York (Year 1982) is a good example of a bloom between the winter temperatures and spring warming regimes, which had an average bloom date of March 18 with a standard deviation of 27. Figures 6 and 7 show $\alpha$ and $\beta$ for the nine sites with the longest running records.

We fit a two-parameter Gaussian location–scale model using the gaulss function from the R package mgcv in which we define

$$\mu(\alpha,\beta) =\mathbb{E}[\nu(\tau)\mid\alpha,\beta]$$

$$\sigma(\alpha,\beta) =\mathrm{SD}(\nu(\tau)\mid\alpha,\beta).$$

Both $\mu$ and $\sigma$ are modeled with tensor-product splines with basis dimension k = 10. The fitted functions are shown in Figures 2 and 3.

A. Appendix

A.1 Sketch of Theoretical Results

We follow the notation in Gut (2009). Let $X_{i}=\alpha+\beta i+\epsilon_{i}$ denote the temperature on day $i$, where $\alpha>0$, $\beta\geq 0$, and the $\{\epsilon_{i}\}$ are independent and identically distributed with $\mathbb{E}[\epsilon_{i}]=0$ and $\mathrm{Var}(\epsilon_{i})=\sigma^{2}$. Define the deterministic cumulative sum

$$\xi_{n}=\sum_{i=1}^{n}\mathbb{E}[X_{i}]=\alpha n+\frac{\beta}{2}n(n+1)=\frac{\beta}{2}n^{2}+\beta\gamma n,$$

where $\gamma=\frac{\alpha}{\beta}+\frac{1}{2}$ and $\xi_{n}=\alpha n$ when $\beta=0$.

Denote the sum of the stochastic component $S_{n}=\sum_{i=1}^{n}\epsilon_{i}$. The object of interest is the first-passage time of $Z_{n}=S_{n}+\xi_{n}$,

$$\nu(\tau)=\min\{n\geq 1:\ Z_{n}>\tau\}.$$

In the language of Lai and Siegmund (1977), $\{Z_{n}\}$ is a perturbed random walk. That is, the random walk $\{S_{n}\}$ is perturbed by the deterministic term $\{\xi_{n}\}$. This setup yields two asymptotic regimes, corresponding to the empirical distinction between the winter temperatures regime in which $\beta\approx 0$ and the spring warming regime in which $\beta>0$.

When $\beta=0$, we apply the usual asymptotic argument for the hiting time of a random walk with constant positive drift as $\tau\to\infty$,

$$\nu(\tau)\ \dot{\sim}\ \mathrm{Normal}\!\left(\frac{\tau}{\alpha},\ \frac{\sigma^{2}\,\tau}{\alpha^{3}}\right),$$

see Gut (2009).

When $\beta>0$, the deterministic component $\xi_{n}=\frac{\beta}{2}n^{2}+\beta\gamma n$ grows quadratically so that the threshold $\tau$ is reached near the deterministic crossing time $m(\tau)$ satisfying $\xi_{m(\tau)}\approx\tau$. The solution is

$$m(\tau)=\frac{-\beta\gamma+\sqrt{\beta^{2}\gamma^{2}+2\beta\tau}}{\beta}=\sqrt{\frac{2\tau}{\beta}}-\gamma+o(1).$$

By the functional central limit theorem, the deviation of $\nu(\tau)$ about $m(\tau)$ is approximately normal. Linearization of $Z_{\nu(\tau)}\approx\tau$ around $m(\tau)$ yields

$$\nu(\tau)\ \dot{\sim}\ \mathrm{Normal}\!\left(m(\tau),\ \frac{\sigma^{2}\,m(\tau)}{\bigl(\alpha+\beta\,m(\tau)\bigr)^{2}}\right)\;\approx\;\mathrm{Normal}\!\left(\sqrt{\frac{2t}{\beta}}-\gamma,\ \frac{\sigma^{2}}{\beta^{3/2}\sqrt{2t}}\right),$$

where we have used the fact that $\alpha+\beta m(\tau)\approx\beta m(\tau)\approx\sqrt{2\beta\tau}$ when $\tau$ is large.

It follows that in the spring warming regime, the mean hitting time scales like $\sqrt{\tau}$ (rather than linearly in $\tau$ as in the winter temperatures regime) and the variance shrinks at rate $\tau^{-1/2}$ (rather than growing at $\tau$ as in the winter temperatures regime).

A.2 Additional Tables

In addition to the model-based approach in the main text, we also examine the data using a binning approach. The results are contained in Tables 3-3.

Table 1: Summary data from the Walnut forcing experiment (Charrier et al. 2011).

$\alpha$	$n$	mean	sd
5	30	178	27.28
10	30	88	12.44
15	30	39	11.33
20	30	26	7.67
25	30	20	4.45

Table 2: Mean bloom date (day-of-year) from lilac data (Rosemartin et al. 2015, updated to 2025 using USA NPN records) by average Jan-Feb temperature $\alpha$ (rows) and average Mar-Apr temperature increase $\beta$ (columns).

$\alpha$ $\backslash$ $\beta$	$[-0.28,\,0.07]$	$(0.07,\,0.11]$	$(0.11,\,0.15]$	$(0.15,\,0.68]$
$(0,\,0.7]$	157.70	152.07	147.60	142.20
$(0.7,\,1.9]$	151.12	146.88	141.64	136.15
$(1.9,\,4.7]$	136.30	131.92	127.58	124.94
$(4.7,\,16.4]$	109.56	109.77	107.36	105.71

Table 3: Standard deviation of bloom date (day-of-year) from lilac data (Rosemartin et al. 2015, updated to 2025 using USA NPN records) by winter level $\alpha$ (rows) and spring warming rate $\beta$ (columns).

$\alpha$ $\backslash$ $\beta$	$[-0.28,\,0.07]$	$(0.07,\,0.11]$	$(0.11,\,0.15]$	$(0.15,\,0.68]$
$(0,\,0.7]$	12.80	11.94	10.90	10.95
$(0.7,\,1.9]$	12.76	11.91	12.67	12.99
$(1.9,\,4.7]$	16.68	14.97	14.45	12.55
$(4.7,\,16.4]$	19.39	17.06	15.31	12.35

Table 3 shows data from a controlled experiment on walnut trees (Juglans~sp.) conducted by Charrier et al. (2011). The data examines 30 trees comprising 6 genotypes at 2 locations. In November of the study year, stems were sampled from each tree and cut into 7 centimeter segments containing a single bud. All segments were chilled at $4^{\circ}$C and then transferred to constant ‘forcing’ environments with temperatures set at one of $5,10,15,20,$ or $25^{\circ}$C. For each temperature, the outcome is the response time, the time (in days) until budburst determined from the date at which 50% of buds reached stage 15 of the BBCH scale.

This experiment isolates the winter regime because temperature is held constant over time. To show consistency with this regime, Table 3 reports the mean and standard deviation of the response times across the 30 stems at each forcing temperature. Higher forcing temperatures advance budburst at a decreasing rate (the mean falls with $\alpha$), and the dispersion of budburst timing is smaller at warmer forcing temperatures. Thus, in an experimental setting, the stopped-random-walk approximation captures the basic comparative statics of constant-temperature forcing on both the level and variability of budburst times. This pattern cannot be explained by chinning or photoperiod as these are held constant.

Tables 3-3 show data from the historical lilac–honeysuckle phenology network compiled by Rosemartin et al. (2015), limited to common lilac and updated using data from USA NPN (Switzer et al., 2025) as describd in the main text. We bin sites and years within a grid of baseline temperatures ($\alpha$, rows) and springtime warming temperatures ($\beta$, columns) and calculate the mean and standard deviation of the observed bloom date, as the number of days since January 1. The results are consistent with Figures 2 and 3 in the main text. The mean and variance are nonlinear functions of $\alpha$ and $\beta$. In particular, the variance generally increases when $\alpha$ increases (holding $\beta$ fixed) and decreases when $\beta$ decreases (holding $\alpha$ fixed).

A.3 Additional Figures

Figure 5: Springtime daily midrange temperatures (points) with $\alpha$ (mean from Jan to Feb) and $\beta$ (slope from Mar to Apr) for nine randomly selected sites from lilac-honeysuckle dataset (Rosemartin et al. 2015) as measured by GHCND stations nearby. Temperatures below 0 are set to 0.

Figure 6: Changes in $\alpha$ for the nine sites from lilac-honeysuckle data (Rosemartin et al. 2015) with longest running records as as measured by GHCND stations nearby. Trend line and confidence interval calculated using linear regression.

Figure 7: Changes in $\beta$ for the nine sites from lilac-honeysuckle data (Rosemartin et al. 2015) with longest running records as as measured by GHCND stations nearby. Trend line and confidence interval calculated using linear regression.

A.4 Simulations

We run two simulations to verify the interpretation of our findings. The first checks the asymptotic approximation derived in Section A.1. The second checks the interpretation of the data in Section A.2.

Simulation 1

The first simulation verifies the accuracy of the two asymptotic normal approximations for the hitting time

$$\nu(\tau)=\min\Big\{n\geq 1:\ \sum_{i=1}^{n}X_{i}>\tau\Big\}$$

under the model

$$X_{i}=\alpha+\beta i+\varepsilon_{i},\qquad\varepsilon_{i}\sim\mathrm{Normal}(0,\sigma^{2}).$$

We set $\sigma=20$ and generate daily temperatures until the threshold $\tau$ is reached, recording the corresponding hitting time $\nu(\tau)$. We repeat this procedure $R=10{,}000$ times for thresholds $\tau=\{1000,2000\}$, $\alpha=\{2,4\}$ and $\beta=\{0,0.1\}$. When $\beta=0$, the simulation matches the winter temperatures regime (regime 1), and when $\beta=0.1$, it matches the spring warming regime (regime 2).

We then compare these simulations with the corresponding theoretical approximation. When $\beta=0$, we use the winter temperatures regime approximation

$$\nu(\tau)\ \dot{\sim}\ \mathrm{Normal}\!\left(\frac{\tau}{\alpha},\ \frac{\sigma^{2}\tau}{\alpha^{3}}\right).$$

When $\beta>0$, we use the spring warming regime approximation

$$\nu(\tau)\ \dot{\sim}\ \mathrm{Normal}\!\left(\sqrt{\frac{2\tau}{\beta}}-\frac{\alpha}{\beta}+\frac{1}{2},\ \frac{\sigma^{2}}{\beta^{3/2}\sqrt{2\tau}}\right).$$

For each setting we computed the standardized variable

$$z=\frac{\nu(\tau)-\mu}{\mathrm{sd}}.$$

The simulated distribution of $\nu(\tau)$ is represented in Figure 8 below by the histograms, standardized as described above. Overlaid is the standard normal density. Across both regimes, the standardized histograms align closely with the $\mathrm{Normal}(0,1)$ curve, and the agreement improves as $\tau$ increases, providing a visual confirmation of the derived asymptotic distributions.

Figure 8: Simulations comparing the distribution of bloom dates (histogram) with an asymptotic approximation (solid line) under two regimes: winter temperatures (stationary temperatures, top) and spring warming (increasing temperatures, bottom).

Simulation 2

To verify our interpretation of the data, we simulate daily temperatures

$$X_{i}=\mu_{i}+\varepsilon_{i},\qquad\varepsilon_{i}\sim\mathrm{Normal}(0,\sigma^{2}),$$

with $\sigma=20$, where the deterministic component follows

$$\mu_{i}=\begin{cases}\alpha,&i=1,\dots,90\quad\text{(January--February)},\\[2.0pt] \alpha+\beta(i-90),&i=91,\dots,180\quad\text{(March--April)}.\end{cases}$$

Bloom occurs when cumulative temperature first exceeds a threshold $\tau$, i.e.

$$\nu(\tau)=\min\Big\{n\geq 1:\ \sum_{i=1}^{n}X_{i}>\tau\Big\}.$$

We simulated $R=10{,}000$ independent realizations of $\nu(\tau)$ for each combination of $\alpha\in\{4,8,10\}$, $\beta\in\{0.2,0.4,0.8\}$, and $\tau\in\{1000,2000\}$. We report the mean (top) and standard deviation (bottom) of $\nu(\tau)$ in Tables 4-4, arranged with rows indexing $\alpha$ and columns indexing $\beta$.

The results reproduce the qualitative patterns observed in the lilac data. Mean bloom dates decrease as either $\alpha$ increases or $\beta$ increases (Tables 4 and 4). Within a fixed $\alpha$ row, increasing $\beta$ reduces the standard deviation of bloom dates (Tables 4 and 4), consistent with the prediction that increasing average temperatures diminish the influence of noise on the threshold-crossing time. Holding $\beta$ fixed, the standard deviation is often larger at higher $\alpha$ (particularly for $\beta\in\{0.4,0.8\}$), consistent with the prediction that earlier crossings occur in a flatter portion of the seasonal cycle where microclimate variability has greater influence.

Table 4: Simulated mean and standard deviation of the bloom date $\nu(t)$ across $R=10{,}000$ replicates with thresholds $\tau\in\{1000,2000\}$. Rows index the winter level $\alpha$ and columns index the spring warming rate $\beta$.

(a) mean, $\tau=1000$

$\alpha$ $\backslash$ $\beta$	0.2	0.4	0.8
4	151.66	136.83	124.86
8	114.87	111.27	106.85
10	98.45	97.20	95.72

(b) sd, $\tau=1000$

$\alpha$ $\backslash$ $\beta$	0.2	0.4	0.8
4	15.48	10.42	7.21
8	16.53	13.27	10.50
10	15.94	14.45	12.71

(c) mean, $\tau=2000$

$\alpha$ $\backslash$ $\beta$	0.2	0.4	0.8
4	199.54	171.15	149.16
8	169.81	152.43	137.37
10	156.22	143.40	131.34

(d) sd, $\tau=2000$

$\alpha$ $\backslash$ $\beta$	0.2	0.4	0.8
4	10.96	7.22	4.84
8	10.93	7.50	5.11
10	10.69	7.66	5.36