Plant phenology is the study of the onsets and durations of growth and reproductive cycles such as bud burst, flowering and seed dispersal. When measured over time, phenological data can provide important information regarding plant responses to changes in climate variables, such as temperature and precipitation (Schwartz 2003; Chuine 2010). Phenophase timing shifts are important indicators of climate change, particularly in California's water-limited Mediterranean climate (Lipton et al. 2018). This highlights the importance of monitoring the phenological sensitivity of keystone species, such as oaks, which support high levels of biodiversity, carbon sequestration, and ecosystem resilience (Bargali et al. 2015).

Species have evolved traits to respond to temperature cues and increasing daylength with the early or late onset of bud burst, balancing access to resources, such as light and nutrients, against the risk of frost damage (Chamberlain et al. 2019; Wright et al. 2021). Species with traits for chilling requirements during dormancy are generally better protected from late season frosts by halting bud burst until the chilling requirement has been met (Ford et al. 2016). Additionally, species with strong daylength requirements, which delay bud burst, are less likely to be impacted by frosts, but therefore have limited responses to take advantage of early warming (Chamberlain et al. 2019). Spatial differentiation of phenological responses also occurs within species. For example, populations at higher elevations or latitudes may initiate early onset upon meeting their chilling requirements, while those at lower latitudes, whose chilling requirements may be unmet, exhibit delayed or failed bud-break, or the opposite may occur with lower elevation populations responding to warmer temperatures with earlier onsets than higher elevation populations (Ford et al. 2016). Similarly, the timing of flowering onset and pollen dispersal underlies reproductive success, requiring trade-offs between maximizing phenophase duration and minimizing risk factors to seed and fruit maturation (Koenig et al. 2012).

Oaks are adapted to maximize their use of resources and may exhibit polycyclism, or multiple flushes of leaves or flowers throughout the growing season, responding opportunistically to precipitation and temperature changes. California oaks in particular are adapted to taking advantage of limited and varying amounts of precipitation (Bobinac et al. 2012; Gerst et al. 2017). Onset dates of bud burst and the subsequent length of the growing season are important factors for tree species since they determine the potential amount of height and radial growth available for that year (Askeyev et al. 2005; Ford et al. 2016; Koenig et al. 2021). Multiple flushes may also be responses to herbivory or infectious agents, such as oak moth defoliation or powdery mildew which may also be influenced by climate variables themselves, stimulating trees to produce young leaves to compensate for the loss of large proportions of their leaf canopy (Visser and Holleman 2001).

Previous studies have identified temperature as a driver of oak phenology, with Askeyev et al. (2005) reporting that early spring warming resulted in earlier bud burst onsets in Q. robur L. in Russia, while the role of precipitation has been unclear. Furthermore, Gerst et al. (2017) compared two central and eastern North American oak species (Q. alba L. and Q. rubra L.) with two western species (Q. agrifolia Née and Q. lobata Née) using 2009–2014 USA-NPN citizen science data and identified seasonal temperature variables as being the most significant drivers of leaf and flower flushes for the western oaks, while noting that precipitation was a stronger driver in the eastern species' models while the western species exhibited low sensitivity to precipitation.

There are 20 native Californian oak species, with additional hybrids, occupying overlapping distributions varying with abiotic factors including elevation, aspect and soils (Nixon 2002). High levels of morphological variation occur within and between Quercus species, with the majority of Californian species being from the white oak clade, Quercus section Quercus, or from the red oaks, Quercus section Lobatae (Hipp 2015).

This study aimed to expand our understanding of the phenological responsiveness of five of California's dominant native oaks by using large community-contributed science datasets to assess whether there were significant differences among species' phenophase onsets. Multi-variable models of the drivers of those phenophases were then developed to identify the strongest explanatory climate variables.

Methods

USA National Phenology Network Phenometric Data

Phenology data have been contributed to the California Phenology Project, as part of the USA National Phenology Network (USA-NPN), since 2009. The USA-NPN offers a database whereby professional and community scientists regularly record long-term phenological observations using standardized protocols for over 1200 floral and faunal taxa, including 26 oak species, via Nature's Notebook (Rosemartin et al. 2018). The USA-NPN observation protocols are rigorous and provide data products of phenophase status and abundance of marked plants along with their latitude, longitude and elevation (Rosemartin et al. 2018).

Additionally, the USA-NPN offers climate variables originating from Daymet, as gridded daily surface weather and climate data, to provide spatially explicit data, including measures of daylength and calculations of aggregated values, such as Accumulated Growing Degree Days (AGDD) and seasonal average temperature and precipitation measurements (Thornton et al. 2016; USA-NPN 2020). Table 1 summarizes the geospatial and climate variables used in this study. Additionally, statewide California mean annual temperature and precipitation data for 1902 to 2018 were accessed using NOAA's statewide time series tool (NOAA 2019) to compare recent measures with historical trends (Table 2).

Table 1.

Variables Included in the Datasets from the USA-NPN (2020) and Their Definitions.

Table 2.

Mean Annual Temperature and Precipitation for California 2010–2018 and Comparisons with 1901–2000 Mean Values, Which Were Calculated Using the Statewide Time Series Tool (NOAA 2019).

Firstly, the five oak species with the highest numbers of observations within California (2011–2019) were identified in the USA-NPN database: Q. agrifolia Née., Q. douglasii Hook. & Arn., Q. garryana Douglas ex Hook., Q. lobata Née., and Q. kelloggii Newb. Four of the species were deciduous, with the single evergreen species being Q. agrifolia. The majority of the trees in the study were located in California, although a small number of sites were beyond their normal distribution ranges, including some Q. garryana and Q. kelloggii in Oregon and eight Q. agrifolia in Massachusetts and North Dakota. The onsets of bud burst, flowering (of the male inflorescences, or catkins), and fruiting were selected as the focus of this study, based on the day of the year when a new phenophase was observed for each plant, following clearly defined phenophase descriptions (USA-NPN 2020). Bud burst was observed when a distinct leaf tip emerged from a leaf bud. Flowering was observed with the emergence of flower buds or the presence of inflorescences, and fruiting was noted when juvenile acorns were seen to have developed (USA-NPN 2021).

The five oak species datasets were downloaded from the USA-NPN's database and imported into R-Studio (R version 3.5.1). Observations were selected with an accuracy of seven days, whereby only ‘Yes’ observations preceded by ‘No’ observations within seven days were analyzed; this is the highest level of accuracy offered in identifying phenophase onset dates and reduces the incidence of data entry errors or of onsets missed because of infrequent observations. Partial years of data were discarded and, in order to focus on spring phenophase onsets, each bud burst and flowering dataset was filtered to include only days 1–180 of each year (DOY 1–180) between 2012–2019. For fruiting, days 1–30 of each year were filtered out of datasets to exclude observations before flowering had occurred and which were therefore likely to be non-viable fruit retained from the previous season which had become visible with leaf drop.

One of the species, Q. kelloggii, was problematic regarding fruiting data since these acorns require an 18-month maturation period, rather than the more common one year. Consequently, Q. kelloggii may bear two generations of acorns simultaneously, making it difficult to identify whether observations are for new juvenile acorns or whether onsets of new fruits have been masked by the ongoing presence of the previous season's fruits. Q. kelloggii's fruiting dataset was also very small. Consequently, all fruiting phenophases of ‘fruit’ and ‘ripe fruit’ were combined and included in Q. kelloggii’s analysis.

Results

Fig. 1 shows the locations of each species' observation sites. Most study sites were based throughout California, however some extended northwards since the species distribution ranges of both Q. garryana and Q. kelloggii occur in southern Oregon. Additionally, two sites with eight Q. agrifolia observations in Boston, MA and Jamestown, ND were included in these analyses (USA-NPN 2020). Table 3 shows the geographical ranges of latitude, longitude and elevation for each of the five species with the overall range being 32.9° to 46.9° for latitude, –124.3° to 71.1° for longitude, and 0 to 2053m for elevation.

Fig. 1.

Maps showing locations of study observation sites for each oak species in California, with Q. garryana and Q. kelloggii including sites in Oregon (USA-NPN 2020). Note that eight Q. agrifolia observations were also from two sites in Jamestown, N. Dakota and Boston, Massachusetts not shown on the maps. QUAG: Q. agrifolia, QUDO: Q. douglasii, QUGA: Q. garryana, QUKE: Q. kelloggii, QULO: Q. lobata. Maps created with Natural Earth.

Table 3.

Ranges of Latitude, Longitude, and Elevation for USA-NPN Observation Data.

The earliest mean onsets for all oak species were seen for all phenophases during the peak of the 2012–2016 drought, whereby the earliest bud burst and flowering onsets occurred in the spring of 2015 following a rainy winter and the highest average temperatures during the study period. In 2015, California's average temperature was 1.9°C above the 1901–2000 mean, precipitation was 18.67 cm below the baseline, and the mean bud burst onset of all combined species occurred on day of the year 52, followed by flowering onset on day 62, and fruiting onset on day 138 (Table 2, NOAA 2019).

The mean phenophase onsets (day of the year) between 2012–2019 for the five combined species were 71.302 (SE = 1.344, n = 382) for bud burst, 77.033 (SE = 1.502, n = 322) for flowering, and 156.544 (SE = 3.582, n = 297). Median phenophase onsets for each species were illustrated using box plots (Fig. 2).

Fig. 2.

Box plots of median phenophase onset for five native California oak species, 2012–2019, showing the minimum, first quartile, median, third quartile and maximum for each set of onsets. A. Bud burst. B. Flowering. C. Fruiting. QUAG: Q. agrifolia, QUDO: Q. douglasii, QUGA: Q. garryana, QUKE: Q. kelloggii, QULO: Q. lobata.

No significant differences were identified among bud burst onset means for individual species (F_4,378 = 1.058, P = 0.377,) (Table 4). No significant difference was identified among the mean onsets of the five species for fruiting (F_4,293 = 1.942, P = 0.104). However, significant differences were identified for mean spring flowering onset among the five species (F_4,318 = 4.906, P < 0.001). Quercus. lobata had an earlier onset of flowering than Q. agrifolia (Tukey's HSD Test, P = 0.026) and Q. lobata's onset was earlier than Q. douglasii's (Tukey's HSD test, P = 0.001).

Table 4.

Anova and Tukey HSD Results Testing for Significant Differences (α = 0.05) Among Mean Phenophase Onsets of Five Native California Oak Species Using USA-NPN Datasets Between 2012–2019.

First model outputs for each phenophase identified the primary role of geographic factors, rather than climate variables, in predicting phenophase onset; daylength, latitude, elevation, and longitude had the greatest predictive weights in determining phenophase onsets and producing high R² values. Since this study was focused on climate drivers of phenophases, the geographic variables were removed from subsequent nested models leaving only the climate variables, such as T_max, T_min and seasonal cumulative precipitation values.

Table 5 illustrates the outputs of the best-fitting linear models for the onsets of each phenophase. The strongest explanatory climate variables (Fig. 3) were found to be winter precipitation, mean accumulated precipitation, and winter T_max.

Table 5.

Outputs for the Best-fitting Phenophase Onset Linear Models for Five California Oak Species Using USA-NPN Data, 2012–2019.

Fig. 3.

Coefficient plots for model parameter estimates for climate drivers of California oak phenophases, showing coefficient estimates and their upper (97.5%) and lower (2.5%) confidence intervals. A. Bud burst. B. Flowering. C. Fruiting.

The best-fitting bud burst model with eight climate predictors (winter precipitation, mean accumulated precipitation, winter T_max, summer T_max, summer T_min, spring precipitation, autumn T_min, autumn precipitation) explained 37% of the variance (Table 5). The strongest flowering model included seven predictors (winter precipitation, mean accumulated precipitation, winter T_max, autumn T_max, summer T_max, autumn precipitation, summer T_min) and explained 28% of the variation in phenophase onset. The best fruiting model was less predictive despite including eight explanatory variables (winter precipitation, mean accumulated precipitation, winter T_max, autumn precipitation, autumn T_min, spring T_max, spring precipitation, summer T_min), explaining only 11% of the effect.

Coefficient plots for bud burst and flowering models demonstrate that the temperature and precipitation variables generally have narrow confidence intervals without crossing the zero line (Figs. 3A and 3B). However, in contrast, the fruiting model has wider confidence intervals and three variable coefficients cross the zero line, illustrating the weaker predictive power of the model (Fig. 3C).

Overall, winter precipitation, mean accumulated precipitation and winter T_max were identified as the most significant climate drivers of bud burst, flowering and fruiting phenophases for these five California native oaks.

Discussion

The earliest mean onsets for all five oak species were seen for all phenophases during the peak of the 2012–2016 drought, whereby the earliest bud burst and flowering onsets occurred in the spring of 2015 following a rainy winter and the highest average temperatures during the study period. In 2015, California's average temperature was 1.9°C above the 1901–2000 mean, precipitation was 18.67 cm below the baseline, and the mean bud burst onset of all combined species occurred on day of the year 52, followed by flowering onset on day 62, and fruiting onset on day 138 (Table 2, NOAA 2019). Koenig et al. (2021) similarly reported that Q. lobata exhibited earliest bud burst onset in this same year, 2015, in a 30-year study of this species in central coastal California.

California oaks are adapted to taking advantage of limited and varying amounts of precipitation, as has been demonstrated in this study with the identification of winter precipitation as the most significant climatic driver of spring bud burst, flowering and fruiting. No significant difference was found among the five oak species' bud burst or fruiting onsets (Table 4, Figs. 3A and 3C), suggesting that they are generally adapted to respond to the same climate cues. This effect was seen across both oak clades, the white oaks, Quercus section Quercus, and the red oaks, Quercus section Lobatae, including four deciduous species and the single evergreen species, Q. agrifolia.

Significant differences were, however, found for flowering onset between the species, particularly between Q. lobata and Q. agrifolia and between Q. lobata and Q. douglasii. Quercus kelloggii's median flowering onset was earlier than that of the other oak species, followed by Q. agrifolia, Q. lobata, Q. garryana and then Q. douglasii (Table 4, Fig. 2B). The relative timing of bud burst and flowering varies, with some oak species tending to produce new leaves shortly before staminate flower buds become apparent, while in others leaves and catkins emerge synchronously since terminal buds are frequently mixed, containing both new leaves and male flowers (Oliveira et al. 1994; Keator 1998). The discrepancy between the timing of bud burst and flowering noted here may be attributed to the effect of non-flowering trees on the dataset, since populations may include many non-flowering trees, which may not yet have reached maturity or resulting from other unknown factors, and subsequently will not go on to produce any fruits (Koenig et al. 2012).

Temporal differentiation of flowering and pollen-release periods between species lessens the likelihood of hybridization in mixed species communities, since large amounts of pollen are not being dispersed simultaneously (Cavender-Bares et al. 2017). Hybridization does not occur between Quercus clades, but Q. lobata and Q. douglasii are both in the white oak clade, so that hybridization may occur between these species when their distributions overlap (Hipp et al. 2018). The significant difference between the flowering onsets of these two white oak species may be a genetic adaptation maintaining differentiation between them and reducing the probability of cross pollination and hybridization. Cavender-Bares et al. (2017) attributed temporal segregation in recruitment and regeneration between the lineages as an enhancement towards coexistence between the clades, by maintaining speciation through reproductive isolation, illustrated here by the differing flowering onsets of the white oak Q. lobata and the red oak Q. agrifolia.

However, Q. garryana is also in the white oak clade and has an overlapping flowering phenophase with both white oak species, increasing the probability of hybridization between them. This may explain the occurrence of many hybrid populations of white oaks, including Epling's oaks, Quercus × eplingii, which are hybrids between Q. douglasii and Q. garryana and commonly occur in communities throughout California, such as in blue oak woodland alliances (CNPS 2019).

Model outputs for each phenophase identified the primary role of geographic factors rather than climatic variables in predicting phenology; daylength, latitude, elevation, and longitude had the greatest predictive weights in determining phenophase onsets and producing high R² values. This concurs with a study by Gerst et al. (2017) which compared the phenophases of western (Q. lobata and Q. agrifolia) and eastern US oaks (Q. robur and Q. alba) to assess the effects of species relatedness with geospatial factors and found that species were more responsive to geographic factors than to relatedness between species.

Multiple regression modeling, after removing the geographical variables, identified winter precipitation and mean accumulated precipitation as the most significant climate variables for the onset of all three phenophases, with winter T_max being the most predictive temperature variable. Winters are often the most unpredictable season of the year in California, with swings in precipitation varying from drought conditions to seemingly constant rain or snow. Even without the additional uncertainties of climate change, it is not unexpected that this phenological sensitivity was exhibited by all five species across all phenophases (Swain et al. 2018).

The roles of maximum winter and spring temperatures have been previously identified as strong climatic drivers of the timing of oak phenophases although the role of precipitation has been less clear. For example, Askeyev et al. (2005) reported that winter temperatures in the Volga region of Russia had shifted bud burst 2.8 days earlier for every 1°C of warming for oaks and noted that temperatures in the month before bud burst were the most significant predictors of this phenophase. Recently, Koenig et al. (2021) confirmed that warmer temperatures induced earlier bud burst in Q. lobata, while noting the inverse relationship between higher mean maximum temperatures and rainfall in California. Papper and Ackerly (2021) similarly identified warmer November and March minimum temperatures as drivers of advancing bud burst onsets.

Gerst et al. (2017) demonstrated that traits such as bud, flower, and fruit phenophases are extremely sensitive to temperature, or growing degree days, and that higher temperatures may induce early bud burst and leaf flushing in western tree species such as Q. lobata and Q. agrifolia, while eastern oaks such as Q. robur and Q. alba were more responsive to spatial and temporal precipitation variations. However, a comparison between the 2017 study and this one highlights significant differences between them. For example, Gerst et al. (2017) used a smaller dataset for Q. lobata and Q. agrifolia which included California data between July 2011–2014, with 2012 and 2013 being the coolest years during 2012–2018 and California's 2012 average precipitation was 104.7% of average, while 2013 and 2014 were drought years with approximately 35% and 89% respectively of average precipitation (Table 2, NOAA 2019). Thus, the longer dataset used in this current project (2012–2019) includes the years before, during, and after the 2012–2016 drought with three additional rainy winters. This offers a longer dataset and a more representative sample of California oaks' responses to climate extremes, particularly regarding precipitation, which was identified as the dominant climatic driver of their phenophases in this study.

General limitations of plant phenology studies to predict the relationships between climate and the timing of phenophases, stem from the complex effects of many genetic, environmental and ecological variables, which can confound assumed relationships between climate and phenophases. These include community interactions, herbivory, disease, pests, pollinators, soil, nutrient and water availability, and anthropogenic impacts such as habitat fragmentation and pollution may also be significant (Knapp et al. 2001; Schwartz 2003). High estimates of variance, or R² values, were therefore not expected with these models owing to these complex ecosystem dynamics also influencing phenophase timing.

Strengths of this study include the relatively large datasets (1001 observations) in comparison with many phenological studies, stemming from the observations of community-contributed scientists and the widespread distribution of observation and climate data across Californian sites (in addition to southern Oregon, North Dakota, and Massachusetts) with the inclusion of five native oak species.

Concerns about the USA-NPN datasets and their usage in this study stem from the relatively short duration of the monitoring period. The larger Q. agrifolia and Q. lobata datasets began in 2012, but data collection began in 2014 for the other species, reducing observation numbers, particularly for Q. kelloggii and Q. garryana. Sample size was also an issue in this study, despite the large dataset offered by USA-NPN. Many of the total observations downloaded from the USA-NPN were of phenophases not being assessed here, such as pollen release and leaf senescence, so that some species had only one or two observations of individual phenophases in a single year. Inequities in the weighting of species further complicated analysis of variations in phenophase onset between species, for example, Q. agrifolia comprised 156 of the 297 total fruiting observations, while Q. garryana contributed only 12. Larger datasets with more frequent observations would greatly improve certainty about the differences between phenophases and responsiveness to climate cues. Adapting the initial data analysis parameters by loosening the requirement for only accepting ‘Yes’ observations which had been preceded by a ‘No’ observation within seven days (for example by increasing it to 14 days) would have increased the number of phenophase onset datapoints and thus sample sizes, but may also have increased the inaccuracy of those onset dates.

Additionally, there are some concerns about the accuracy of the data. Many of the observation sites are located at University of California field stations or at scientific or educational organizations, so that many observers are scientists, researchers and educational staff, but many are community-based scientists with varying levels of experience (USA-NPN 2020). Although they are all trained to follow clear protocols, it is difficult to ensure that everyone is following them exactly. Good eyesight combined with careful use of binoculars or hand lenses is essential for the differentiation of whether a leaf bud is truly breaking or whether a catkin flower is open and producing pollen. Male flower buds are also easy to miss during the first days of bud break without experience, particularly when they emerge from mixed terminal buds along with young leaves which may obscure them (Keator 1998). This may lead to the onset of flowering being missed and reported as occurring after leaf bud burst. Lastly, the identification of new fruits is difficult because they are small and are initially difficult to differentiate from axillary buds, especially when high up in a tree. However, the USA-NPN strives to ensure standardization and accuracy by providing training and support; reporting 91% concordance between trained and expert observers and the flagging of conflicting observations (Fucillo et al. 2014).

The frequency of observations, or lack thereof, may also impact whether the onset or duration of a phenophase has been captured, especially for flowering which may be short and easily missed. These limitations may impact the accuracy of the data, particularly for those species, such as Q. kelloggii (n = 36, as compared with n = 541 for Q. agrifolia) with smaller numbers of observations which are indicative of fewer, less experienced observers. Weekend bias has also been postulated in community-based observations, where volunteers may be more likely to contribute to projects on weekends, thus increasing the likelihood of biased onset dates being noted, although Courter et al. (2013) reported that changes in habits seemed to be reducing this tendency.

A final note on possible inaccuracies in observational data concerns hybridization, which is a common feature of oak communities throughout California, particularly where white oak species distributions overlap. Oak hybrids abound while leaf morphology varies highly across individual trees and throughout populations, increasing the difficulty of differentiating between species and hybrids (Hipp 2015). Consequently, it may be assumed that a number of the trees identified by species in the USA-NPN database are actually hybrids, perhaps leading to a greater probability of similar phenophase onsets and durations. This may be a subject for further study.

Further study of the impacts of climate change and extreme weather events on phenophases is critical to understanding the potential impacts on species regeneration and distributions as well as on forest resilience and dynamics, particularly with climate models projecting increasing frequency of droughts and extreme weather events like atmospheric rivers (Chuine 2010). The future of keystone species such as oaks is dependent on their adaptability to changing climate, and to climate extremes in particular. This adaptability rests on their phenological responsiveness, and winter precipitation has been found here to be a main driver of vegetative growth and reproductive potential. Those populations that are able to respond to changing climate variables by bringing forward or delaying phenophases, such as bud burst, flowering or fruit ripening, will be the most resilient and able to expand or maintain their ranges, while those less adaptable species will be unable to survive increasingly rapid rates of change.