Nesting phenology of migratory songbirds in an eastern Canadian boreal forest, 1996–2020

. The migration phenology of many bird species has changed over the past few decades, but whether such changes lead to changes in the nesting phenology remains little known. Studying bird nesting in the boreal forest comes with challenges because of the large size of this biome. We evaluated songbird nesting phenology for the past 25 yrs in a boreal forest in eastern Canada, Forêt Montmorency. We used the observation of food transport in adults as an index of parental status, considering the imperfect detection of this status through hierarchical models of site occupation. We estimated annual phenology as the Julian date of the inflection point of the logistic fit of proportion of sites with parental activity as a function of Julian date. Contrary to expectations related to the advance of spring migration in North America, models did not show an advancement in the nesting season. Models showed that passerines can move their nesting date back or forward by 1 to 9 d. Models suggested that short-distance migrants delayed their nesting date by 2 wks against 1 mo for long-distance migrants. These results show the capacity of songbirds to adjust their nesting time and remind us of the value of regional studies when we are interested in reproductive phenology.


INTRODUCTION
Climate change observed in recent decades is thought to affect the phenology of bird migration and breeding, with consequences on population viability (Murphy-Klassen et al. 2005, Jones and Cresswell 2010, Mayor et al. 2017).For example, nesting birds are expected to synchronize peak demand for food by nestlings with peak food abundance, often composed of insect larvae (Baker 1939) subjected to short-term changes in weather.Phenological changes in migration can disrupt the synchrony between peaks in food availability and physiological demand for egg production (Perrins 1970) or nestling growth (Thomas et al. 2001) in the case of late-arriving birds.Conversely, early arrival enhanced by weather along migratory routes can result in birds experiencing cold conditions and scarce food, leading to mortality, especially for insectivorous birds (Ahola et al. 2004, Newton 2007).This temporal mismatch (sensu Bourret et al. 2015) is likely to have negative consequences on reproduction, as exemplified in the Great Tit (Parus major; Visser et al. 2006).And as discussed by Bourret et al. (2015), the mismatch that can be observed between the demand for food by offspring and the peak of resources does not necessarily mean a "mistiming," the result of a non-optimal decision by the birds.This is all the truer, for example, if the food is not restrictive and the net benefits in fitness peak at an abundance of prey lower than the prey's peak abundance.
Climate change has a strong geographic component (Intergovernmental Panel on Climate Change (IPCC) 2007) and thus may lead to regional variation in migration and temporal mismatch between food supply and demand (Visser et al. 2011).Such regional variation has been shown in North America (Marra et al. 2005, Murphy-Klassen et al. 2005, Horton et al. 2020).Marra et al. (2005) showed an advancement of spring arrival dates by 1 d per 1°C increase in temperature during April and May in North American long-distance migratory birds.Studies have pointed to the impact of temperature (Both et al. 2006), precipitation (McDermott and DeGroote 2016), and to a lesser extent, growing degree-days (GDD) on the nesting phenology (Courter 2012, Sockman andCourter 2018).Torti and Dunn (2005) linked higher spring temperatures to earlier nesting dates in six North American bird species.McDermott and DeGroote (2016), with a long-term data set, showed early nesting in 13 out of 21 passerine species in Pennsylvania in years of higher mean daily spring temperatures.Jones and Cresswell (2010), using the difference in temperature increase between the non-breeding area and the breeding area as an index of asynchrony, found that North American species experiencing a mismatch were more likely to be in decline.
Migrants that arrive earlier in breeding areas and experience colder conditions and less abundant food resources have higher rates of mortality, especially insectivorous birds (Newton 2007).This challenge has been exemplified by work on the Pied Flycatcher, (Ficedula hypoleuca (Ahola et al. 2004).According to Ahola et al. (2004), favorable weather during the spring migration is associated with early arrival at the flycatcher's nesting grounds in Finland.But even though the climate remained relatively stable in its nesting ground, this flycatcher population did not change its breeding date, leading to an extension of the prebreeding period, i.e., the period spent at the breeding ground before reproduction onset.Changes in reproductive phenology are likely to have negative consequences on reproduction, as shown on brood size and body mass of Great Tit fledglings (Visser et al. 2006).Zaifman et al. (2017) reported changes in migration timing in three regions of North America (Alaska, Maine, and South Carolina) due to changing temperatures in more long-distance migrants compared with short-distance migrants.Halupka and Halupka (2017), in a meta-analysis of 65 long-term studies concerning 54 species from the Northern Hemisphere, showed that within the last 45 yrs, sedentary species and short-distance migrants prolonged their breeding seasons more than longdistance migrants.
Measuring the productivity of bird populations is notoriously difficult, especially over large landscapes and species assemblages (Gunn et al. 2000).One way to do this is by using nest records collected by volunteers (Rousseu and Drolet 2017).To counteract the low precision of volunteer nest records, Rousseu and Drolet (2017) used back-calculations, which enabled an estimation of the start and end of the nesting period.Another way to address the high cost of direct nest monitoring is to use the proportion of adults raising young as a proxy for true reproductive success, while considering imperfect detection of parental behavior (Corbani et al. 2014).This method offers the possibility of estimating productivity as a proportion of birds with broods over large areas and thus working at the bird community level (Corbani et al. 2014).Given that the proportion of birds with broods from the onset to the peak of a nesting season increases in a sigmoid fashion, it becomes possible to quantify the timing of breeding as the date corresponding to the inflection point (Hurlbert and Liang 2012).
Following numerous studies, the IPCC reports that, in general, climate change is pushing forward spring migration in the Northern Hemisphere (Settele et al. 2014) but does not make a claim about possible changes in the nesting phenology.As changes in migration timing appear to vary between long-distance vs. short-distance migrants (Halupka and Halupka 2017;A. Desrochers, A. Florea, and P.-A. Dumas, unpublished manuscript), one would also expect that changes in nesting phenology would depend on these groups.
In this study, we document the reproductive phenology of migratory songbirds between 1996 and 2020 in a 412 km² boreal forest in eastern Canada.We tested whether (1) annual mean of GDD changed monotonically, (2) birds advanced their nesting season, (3) but to a different extent depending on migration distance, and (4) annual changes in phenology were associated with GDD.

STUDY AREA AND METHODS
The observers collected data at Forêt Montmorency, located in the Laurentians, 70 km north of Quebec City, Quebec, Canada over 20 yrs between 1996 and 2020.Forêt Montmorency is characterized by a latitudinal gradient from south to north over a length of 44 km, with a succession of balsam fir (Abies balsamea)/yellow birch (Betula alleghaniensis), balsam fir/white birch (Betula papyrifera) and black spruce (Picea mariana) forests in the north.Between 1981 and 2010, the average annual temperature was 0.5°C, and the precipitation was 1583.1 mm, of which nearly 40% fell in the form of snow (Environment Canada 2021).The area is characterized by rugged terrain, with an altitude varying between 600 and 1000 m (Darveau et al. 1997).Balsam fir, the dominant tree species, is accompanied by black and white (P.glauca) spruce and white birch (Darveau et al. 1997).

Data collection in the field
Data were collected during 7,106 bird surveys and at 3,466 stations (Table 1), with an annual mean of 173 stations, most of them being concentrated in the former extent of the Forêt Montmorency (Fig. 1).Originally limited to 66 km², the Forêt Montmorency was extended to 412 km² in 2014, which has led to more frequent surveys in the extended Forêt Montmorency post-2013 than in previous years.The sampling took place during the nesting period, from the end of May to the end of July, between 5:00 and 10:00 EDT and in favorable weather, i.e., low or no rain and winds.Most (55%) stations were either distributed along forest roads, paths, and trails, others were located within stands (for more details, see also Drolet et al. 1999).
Each bird survey consisted of observing for birds and recording species, sex, and distance to each visually detected individual, as well as whether they were carrying food, fecal sacs, or nothing.Sampling events lasted 5-30 min, the variation due to projects with different sampling protocols conducted over the 25-yr period, but the observers were stationary for all of them (i.e., at a fixed point in space).In addition, two methods were used to inventory birds at Forêt Montmorency over the 25 yrs: passive or with playbacks of Blackcapped Chickadee (Poecille atricapillus) mobbing calls with 5-watt sound amplifiers (Gunn et al. 2000).From 2007, only the latter method, i.e., with playback, intended to increase visual contact with birds, was used.Playback had an impact on the number of birds seen (slope = 726.3,p < 0.05; we   . use the term "slope" to refer to regression estimates-that is why we used them-but playbacks had no significant effect on the proportion of birds observed with food (slope = 3.604, p = 0.83).
We visited each site two to three times within a maximum of three consecutive days to allow the estimation of the detection probability of parental status, assuming a fixed parental status (closure assumption, Rota et al. 2009).We excluded 1995We excluded , 1999We excluded , 2008We excluded , and 2019 from the analysis because of missing repeated visits to the same survey sites.
We used the atlas of breeding birds of Quebec (Robert et al. 2019) to classify short-and long-distance migrants.The meteorological data, i.e., GDD, came from Environment Canada, "Forêt Montmorency" station (1996)(1997)(1998)(1999)(2000)(2001)(2002) and "Forêt Montmorency RCS" station (2003-).Growing degree-days accumulate whenever the daily mean temperature is above a specified threshold temperature.As we wanted to address the effect of temperatures both on plant and insect growth, a base temperature of 5°C was used to calculate GDD (Environment Canada 2021).

Data preparation
We excluded from the analysis waterbirds, species only heard or in flight, including raptors.We also excluded Fringillids, as they do not carry food for offspring conspicuously in their bills.An individual was considered as breeding only if it was carrying food or a fecal sac, as well as if it was in distress or if its nest was found.Distress cases included only birds displaying behaviors such as the "broken wing," alarm calling birds were not included.

Statistical modeling
As adults do not spend all their time carrying food or fecal sacs, the detection of their parental status is imperfect.To account for imperfect detection of parental behavior, we used the MacKenzie two-state site occupancy model (MacKenzie et al. 2002).It is based on the likelihood of the presence and detection probabilities given the detection histories, i.e., vectors of 1s and 0s corresponding to the detection and non-detection of parental behavior on every survey, respectively.The MacKenzie occupation model was initially developed to estimate the proportion of sites occupied by a species when the probability of species detection is less than 1 (MacKenzie et al. 2002).This model normally involves two random processes: occupation and detection.In our case, we substituted the observation of the species with the observation of parental behavior.Thus, the probability of nestlings occurring at a given site was assessed, as well as the probability of detecting at least one instance of parental behavior during a visit, conditional on the site having nestlings.
Given the very short (3 d) maximum span across visits at the same site, we assume that the parental state of birds remained constant across visits, i.e., that the parental state was "closed" (sensu Rota et al. 2009).To avoid confusing the non-detection of the feeding of the young with the non-detection of the species itself, when no visual observation of the species was made, missing data were assigned.Given that birds were not individually marked, we were unable to infer parental status at the individual level.Such an inference would be unreasonable because individuals from adjacent territories often found themselves together when attracted by recordings (A.Desrochers, S. Boukherroub, personal observation).Instead, we collapsed species-specific data into two species groups: long-distance migrants (LDM) and short-distance migrants (SDM), and scored as "with nestling(s)" any survey that yielded at least one bird of the species group showing parental behavior.
We estimated the probabilities of the presence of adults raising young as a function of Julian date (number of days since 1 January), year, GDD, and latitude.We modeled the year as a quantitative effect because we were interested in whether nesting dates advanced or receded over the 25-yr period.We opted for GDD from among other climatic variables (mean daily temperature and precipitation) because of the role of accumulated heat available for plant emergence (Häkkinen et al. 1998), the growth and maturation of insects (Herms 2004, Philips et al. 2012), the largest feeding category among nesting birds in the region considered (Morse 1995).The annual cumulative GDD up to 20 May was used instead of the daily GDD to avoid the correlation between the latter and the Julian date.Latitude was included in some models (Table 2) because of its possible influence on phenology, given the large extent of the study area and the shifting of sampling sites between years.All models had an intercept of zero, assuming that on 20 May, the number of birds with nestlings was equal to zero.We modeled the probability of detection as a function of the number of birds seen during each survey.As the objective was to assess parental status, the number of covariates was kept to a minimum.We compared eight candidate models, all with the number of birds seen as a detection covariate, but each with a different combination of year, Julian date, GDD, latitude, and interactions (Table 2).We selected the best models based on the Akaike information criterion (AICc; Burnham and Anderson 2002).Growing degree-days accumulate whenever the daily mean temperature is above a specified threshold temperature.As our interest was on the effect of temperatures both on plant and insect growth, we used a GDD calculated with a base temperature of 5°C.
We estimated annual phenology as the Julian date of the inflection point (hereafter, "inflection date") of the logistic fit of proportion of sites with parental activity as a function of Julian date.The inflection point of a logistic curve corresponds to the X value, from which the slope stops increasing and starts decreasing, i.e., with a second derivative of zero.The said X value corresponds to half of the maximum of the curve.This is the same metric as that used by Hurlbert and Liang (2012) to measure migration phenology.We used this value rather than the date corresponding to the peak proportion of birds with parental activity because, in some years, the peak was not reached before the end of the field season.
The fit of the global model, i.e., including all the explanatory variables for each parameter, was tested using the MacKenzie and Bailey (2004) randomization test (n = 1,000 resamples), comparing the frequency of the detection histories observed with those predicted by the model.All models were fitted with the occu function of the "unmarked" package (Fiske and Chandler 2011) in R software version 3.6.2(R Core Team 2019).The alpha level error rate was set to 5%, with correction for multiple tests, i.e., the model selection.The model averaging including first models with Delta AICc under 4.0 and fit test were performed using "mb.gof.test" and "modavg" functions of the AICcmodavg package, respectively.

Trends in growing degree-days
Since 1996, GDD, the annual mean or values at Julian day 170 varied substantially among years (Fig. 2), but with no clear trend (slope = .08,p = 0.21 and slope = .08,p = 0.89, respectively).

Phenology trends by species groups
An annual mean of 1,387 birds and 94 parental behaviors were reported from 36 migratory species, of which 31 were longdistance migrants (Table 3).During the 16 sampling years included in the analysis, most of the earliest dates of observed parental provisioning for each species occurred in June (Table 3).
For long-distance migratory species, the model that performed the best among the competing models was the most complex (Table 2) and provided a good fit to the data (MacKenzie and Bailey test, p = 0.15).For short-distance migrants, six models had a Delta AICc under 4.0 (Table 2); therefore, we used model averaging to calculate parameter estimates for the latter group of species.The fit of the most complex model (highest number of parameters) provided a good fit to the data (MacKenzie and Bailey 2004 test, p = 0.01).
For long-distance migrants, adults with parental activity were observed at 371 out of 1,228 sites, corresponding to 30% of sites.
For short-distance migrants, adults with parental activity were observed at 201 out of 924 sites, corresponding to 22% of sites.
After accounting for imperfect detection, those naïve estimates led to estimated proportions of sites with parental activity approaching 100% late in the nesting season (Fig. 3).
As one would expect, for the two categories of species group studied, the proportion of adults with parental activity increased from the start of the sampling period until its end (Fig. 3).Year was negatively correlated with the proportion of sites with parental activity at a given Julian date, for both groups (Table 4), which reflected a delay rather than a decrease in nesting success (Fig. 3).For long-distance migrants, the positive interaction year/ Julian date indicated that the proportion of sites with parental activity increases progressively more quickly over the years.The inflection date was successfully determined for all 16 yrs included in analysis except for 1996 for LDM, and 1996 and 2000 for SDM.
The inflection date was later in more recent years for both longdistance migrants (slope = 1.66 ±0.10, p < 0.001 ; Fig. 4) and short-distance migrants (slope = 2.38 ±0.25, p < 0.001).Birds were more likely to have active nests at high GDD and latitude at a given Julian date in the case of long-distance migrants (Table 4).Unsurprisingly, the probability of detecting parental behavior Fig. 3. Shifts in the timing of nesting, expressed as the increasing proportion of sites with nestlings, at Forêt Montmorency, between 1996 and 2020.Predicted values generated from the best performing models (see Table 2), with GDD and latitudes set to their respective means.Dots represent inflection dates.
increased according to the number of birds seen (Table 4).For all species combined, it varied between 0.02 and 0.99 for longdistance migrants and between 0.08 and 0.99 for short-distance migrants.How to explain the delay in nesting given the recent literature on the advancement of spring migration?It seems that the advance in migration is not as ubiquitous as it appears at first glance (Knudsen et al. 2011;A. Desrochers, A. Florea, and P.-A. Dumas, unpublished manuscript).A recent paper based on weather radar data claims that spring migration has advanced by several days, if not weeks, in North America since the 1970s (Horton et al. 2020).However, an examination of ca.4,000,000 birder records spanning 50 yrs and 152 species in the province of Quebec, led to a different conclusion (A.Desrochers, A. Florea, and P.-A.

DISCUSSION
Dumas, unpublished manuscript).In fact, the latter study based on first and median migration dates, and accounting for temporal trends in sampling effort and species abundances, showed that long-distance migrants tend to arrive later in recent decades than they did in the 1970s, whereas short-distance migrants did not tend to arrive earlier or later than they did in the 1970s.These observations are consistent with our conclusions on the delay of the breeding date in certain of the long-distance migrants studied, such as the Swainson's Thrush (Catharus ustulatus), the Blueheaded Vireo (Vireo olivaceus) and the Nashville Warbler (Leiotlhypis ruficapilla), by a delay of the spring arrival date.
The delay in nesting phenology may be related to local conditions.The presence of the GDD in the model, usually used as an index of environmental phenology (Laaksonen et al. 2006), such as the bursting of buds (Herms 2004) and emergence of insects (Philips et al. 2012) and plants (Haekkinen et al. 1998), or to compare its effect on the initiation of reproduction between species (Sockman and Courter 2018), showed that this variable plays a role in nesting date in our study area.Nevertheless, the effect of GDD as an explanatory variable was minimal compared with the year.In the long-distance migrants nesting model, the year coefficient was 19.5 times the GDD coefficient.
The delay in nesting phenology may be due in part to changes in species composition in the data over the period 1996-2020, with late-nesting species possibly more prominent in the data in recent years.We are unable to provide a species-specific assessment of phenological change because of insufficient data in most species studied.However, we found no correlation between median arrival date (details in A. Desrochers, A. Florea, and P.-A.Dumas, unpublished manuscript) and species trend as determined from a Poisson regression of species counts from our data (r = 0.11, p = 0.53).Furthermore, none of the correlations between year and the beginning, duration, and end of the survey were significant.
Our results, which showed the presence of nestlings at 50% of the sites in Forêt Montmorency, between the beginning of June and the end of July, are consistent with those of Rousseu and Drolet (2017).In fact, in the massifs of Lac Jacques-Cartier, the nesting period modeled (from laying to the departure of the young) by Rousseu and Drolet (2017) for short-distance migrants can be as early as 6 May and end as late as 9 August.Assuming a mean egglaying/incubation period of 15 d, the earliest date for hatching would be 19 May.At this date, more than 60% of the nests of short-distance migrants were considered with nestlings.In addition, it is possible to estimate for this group of species that the period of feeding young would extend from 19 May to 9 August, with a median value on 29 June, a date that is in the period of nesting where more than 80% of the parents were raising nestlings (Rousseu and Drolet 2017).Thus, the probability of detecting food transport is expected to be highest from 23 May to 20 July.In short, the models established by Rousseu and Drolet (2017) confirm that the dates of 50% of nests with young are not too early for short-distance migrants.For long-distance migrants, the median value of the period of feeding young falls on 2 July, which also validates the values we obtained (feeding period in this case from 26 May to 7 August).
On average, the probability of detecting parental provisioning was low, which emphasizes the need to consider the imperfect detection when estimating parental states based on behavior.We Table 4.Estimated parameters on the logit scale and their standard error of parental status and probability of detection.Raw parameters are shown for long-distance migrants (Delta AICc = 0) and model-averaged estimates and their 95% confidence limits are shown for short-distance migrants (Delta AICc < 4).have also seen that this same detection was dependent on variables linked to the conditions prevailing during the visit, hence the importance of taking them into account.On the other hand, although some authors point out the existence of interspecific heterogeneity in response to mobbing calls (Rae et al. 2015) as well as the correlation between response and the distance between the count survey and the bird nest (Doran et al. 2005), the use of playback of chickadee mobbing calls is justified when the increase of visual contact is desirable (Rae et al. 2015) and when the focus is on the bird community rather than the individual, as was the case in the present study.
Our study could be criticized for being conservative as it obscures differences among single species, not to mention individual pairs.We were constrained to merge multiple species because of the low numbers of food provisioning cases observed, which led to frequent model convergence problems when run with single species.Furthermore, site-based models are presumably less sensitive to the condition of site closure to changes in parental status between the first and the last visit compared with models at the species scale, as a predator is less likely to prey upon all nests at a given site than it is to change the parental status of a particular species.Another gap in our framework is that we did not study the end of the breeding season, which deprives us of information on the duration of nesting and its possible fluctuations and so making our definition of nesting phenology at Forêt Montmorency partial.
The magnitude of our estimates of delay in nesting phenology is uncertain, but nevertheless contradicts the expectation that the nesting season advanced over the last 25 yrs.This apparent contradiction may stem from the assumption that spring and early summer temperatures have recently risen in the study area, which may not be the case.Future studies should try to get a more complete picture of the nesting phenology, particularly the ending of the nesting season.Furthermore, our results raise queries about the possible impact of the delay in nesting season on its duration, on post-nesting stages (i.e.molting), or even on the phenology of autumn migration.A species-specific approach could help to determine whether later reproduction could be a conservation concern.Ultimately, investigating the relationship between local climate conditions, specifically GDD, and the phenological patterns of plants, invertebrates, and birds could provide a valuable approach to identifying any potential mismatch with periods of maximum food availability.
This study is, to our knowledge, the first long-term, in-depth, and community-wide analysis of bird nesting phenology in the boreal forest.A major conservation implication of this study is the demonstration of the ability of boreal birds to substantially shift their nesting phenology within a matter of decades, if not years.This finding is especially noteworthy given that the breeding season is relatively short in the boreal latitudes (Baker 1939, Wyndham 1986).Another implication would be on the influence on survey protocol design to make sure it is flexible enough to adapt to changes in phenology.This study reminds us of the value of regional studies when we are interested in reproductive phenology, given that climate and other environmental change are highly heterogeneous from one region to another.

Fig. 1 .
Fig. 1.Locations of the 3466 bird survey stations at Forêt Montmorency (Quebec, Canada) done during the breeding season (June-July) over 20 yrs.

Fig. 2 .
Fig. 2. Annual mean and cumulated GDD on 20 May at ForêtMontmorency, between 1996 and 2020.Growing degree-days accumulate whenever the daily mean temperature is above a specified threshold temperature.As our interest was on the effect of temperatures both on plant and insect growth, we used a GDD calculated with a base temperature of 5°C.
Contrary to expectations based on migration phenology and nesting dates in many parts of the world(McDermott and DeGroote 2016, Hallfors et al. 2020, Horton et al. 2020), songbirds have been nesting progressively later but in a more synchronous manner at Forêt Montmorency during the last quarter century.This trend was not accompanied by an observable change in the proportion of sites with parental activity, which approached 100% in most years.The magnitude of the delay in the nesting season, 2-4 wks in 25 yrs, appears unsustainable and may be part of a long-term cycle.

Fig. 4 .
Fig. 4. Shift in the timing of breeding, as determined by the date corresponding to the inflection point in the logistic growth of the proportion of sites with nestlings, against Julian date.Data from Forêt Montmorency, 1996-2020.LDM: longdistance migrants ; SDM: short-distance migrants.The line stands for the regression line between the inflection points and years, dots stand for the inflexion point, and the bars for the confidence levels.

Table 1 .
Number of stations and bird surveys realized each year at Forêt Montmorency between 1996 and 2020 with first (earliest) and last date of sampling (latest), total number of birds seen, and of those, the number seen transporting food.

Table 2 .
Model selection for long-distance migrant and shortdistance migrant species.All models included a detection component, with the number of birds seen as a covariate.Jul stands for Julian day, GDD for growing degree-days, lat for latitude.

Table 3 .
Total number of observations with and without parental provisioning for species for which parental provisioning was observed at sites at Forêt Montmorency.
Intercept parameters were excluded for Phi and included for p in all models. †