This study took place in the Grant Grove Peninsula (GGP), a roughly 6,177-hectare region of Sequoia and Kings Canyon National Parks (SEKI), California, United States ranging from 1,000 to 2,350 m in elevation (Figure 1). In addition to burned areas from the 2015 Rough fire, this region of the southern Sierra Nevada has a rich history of prescribed fire used to manage the forested landscape for habitat structure and resource protection (Swetnam et al., 2009; NPS, 2020).

From 2012 to 2016, California experienced the most severe drought in at least the last century (Griffin and Anchukaitis, 2014), resulting in the death of an estimated 129 million trees in the state (Moore et al., 2018). Large, protected areas, like SEKI, were not spared from these conditions. By 2016, an estimated 20% of all living trees had died in SEKI lower elevation mixed conifer forests (Stephenson et al., 2019), with the magnitude of the effect decreasing with increasing elevation (Paz-Kagan et al., 2017). Much of the mortality was exacerbated by bark beetles, with differential effects depending on the tree and bark beetle species. Large pines suffered particularly high mortality, with 50–70 percent of larger pines dying in some areas (Stephenson et al., 2019). Drought induced tree mortality were effectively moderated by prior burn treatments that reduced tree densities (Restaino et al., 2019; van Mantgem et al., 2016; van Mantgem et al., 2021).

The Rough Fire struck the region at the height of these conditions in 2015. On 31 July 2015, a lightning strike north of Kings Canyon National Park ignited the Rough Fire, ultimately burning 61,324 hectares (613 km²) over 4 months with 3,810 hectares within the park. To assess the severity of the 2015 Rough Fire on the entire study area, multi-spectral satellite imagery from 14 October 2014 (pre-fire) and 9 October 2015 (post-fire) were compared using the satellite image-based summary of RAVG Composite Burn Index (NPS, 2003; USFS, 2015). The Composite Burn Index (CBI) is a standardized fire severity rating based on a composite of effects to the understory vegetation (grass, shrub layers), midstory trees and overstory trees and conveys the degree of burn severity through a comparison with satellite-driven burn severity data (USFS, 2015). Burn severity categories of CBI range from 0 (unchanged/not detected) to 3 (highest severity) (Table 1). In the burned landscape for the Rough Fire high-severity fire accounted for an estimated 31% (Figure 1). The fast growing megafire reached the giant sequoia forests of the GGP on 11 September 2015. While low-severity fire impacts were present in a large portion of the sequoia forest within GGP (73%), several patches of high-severity fire impacts (4%) existed in areas without previous prescribed fires (Figure 1). The most significant impacts to the GGP were on portions of the grove outside the park and not directly part of the study area.

Nine sites were selected within the GGP representing five different fire history categories (Table 1, Supplementary Figure S1): two sites in old growth forest burned by the Rough Fire without a history of prescribed fire (RF1); two sites in old growth forest burned by the Rough Fire with prior prescribed fire (RF1 + Rx); two sites in old growth forest unburned by the Rough Fire but with prior prescribed fire (UB1 + Rx); one site in old growth sequoia forest with no recorded fire history (UB1); and two sites in a second growth sequoia forest with no recorded fire history (UB2).

All nine sites were within giant sequoia-mixed conifer forests, which also included white fir (Abies concolor), red fir (Abies magnifica), sugar pine (Pinus lambertiana), ponderosa pine (Pinus ponderosa), jeffery pine (Pinus jeffreyi), and incense cedar (Calocedrus decurrens) (Haultain et al., 2020). Only the UB2 sites were within a logged area with second-growth giant sequoia trees approximately 100 years old, all other sites were within late succession old growth forests with mature, giant sequoia trees typically more than 1,500 years old.

We placed one SM3 Songmeter (Wildlife Acoustic, Inc.) acoustic recording unit (ARU) at each of the nine study sites (Figure 1, Supplementary Figure S1). Locations for each ARU were determined via random point generation in ArcGIS software (ESRI, Inc.) within each respective fire history categories. Criteria for site selection included within sequoia tree vegetation, 60 m from fire history boundary, at least 25 m from a hiking trail, 40 m from a service road, and 100 m from a meadow or stream and 245 m from another recorder location. The spatial separation of at least 245 m from another site ensured spatial independence for biological sounds.

Each ARU was affixed approximately 2 m above ground to a standing, high-canopied tree (Figure 1C, Supplementary Figure S1). An external A2 acoustic microphone (sensitivity of −4 ± 3dB (0 dB = 1 V/Pa @1 kHz) and dynamic range of 14–94 dB sound pressure level at 0 dB gain) was mounted to a 1 m rod that extended horizontally from the unit (Figure 1C). The ARUs were programmed to record each day for 5.5 h. In the morning, the ARUs recorded for 4 h, starting 1 h before sunrise. In the evening, they recorded 1.5 h before sunset. During Daylight Savings Time (approximately March–November), the ARUs recorded this schedule every day. At other times of the year, we scheduled ARUs to record every other day to extend the recording period in winter months when access to the sites was more difficult. Each recording was saved in an uncompressed WAV audio file format with the following audio settings: gain of 35 dB and sampling rate at 48 kHz.

From the spring through fall, each site was visited weekly to check memory cards, battery life, and equipment damage. In the winter months, site visits were less frequent, ranging from every 2 weeks to once each month depending on snow conditions. Recordings were collected from June 2016 through November 2019 with some variation in total recording days across sites (Table 1).

All audio data were converted to calibrated 1-s one-third octave band sound pressure levels (SPL) with nominal frequencies of 12.5–8,000 Hz using the technical specifications of the ARUs to provide absolute measures of SPL (Merchant et al., 2015) within one-third octave bands. This time-frequency representation has been used effectively in several ecological studies (Buxton et al., 2016; Buxton et al., 2018) and is the standard resolution used in National Park Service standard acoustic monitoring (Lynch et al., 2011). Using previously described methods (Buxton et al., 2018), from the SPLs we calculated the acoustic complexity index (ACI) and a measure of background sound levels in 10-min intervals (R code available: https://irma.nps.gov/DataStore/Reference/Profile/2293158). We used 1.6–8 kHz nominal third-octave frequency bands to represent biological sounds, a frequency range dominated by avian calls. To represent low-frequency background sound (anthropogenic and weather sounds), we used nominal one-third octave bands between 0.315–1.25 kHz. We ultimately used a single index, ACI, rather than a suite of indices because it has shown to be a suitable index to describe ecological conditions, particularly avian diversity (Pieretti et al., 2011; Farina et al., 2015; Buxton et al., 2016; Ross et al., 2018). At each site, all 10-min ACI values were summarized into a daily mean.

The automated data processing routine allowed for rapid processing of audio data when SD cards were retrieved on a weekly basis. The rapid processing allowed us to diagnose data quality issues and provide near real-time assessments of equipment performance and acoustic conditions.

We included weather (temperature and rain) conditions in our analyses to improve interpretations of ACI. At each study site, temperature dataloggers (iButtons; Thermochron iButton, Maxim Integrated, Inc.) were affixed to trees in the forest understory. Each iButton was placed within 2 m of the ARU at approximately the same height above ground. They recorded temperature throughout the day at 30-min intervals. iButtons were retrieved and replaced approximately every 6 weeks, except in 2019 when conditions were too hazardous to visit the sites due to risk of falling branches and trees.

The occurrence of daily rain events was obtained from the Park Ridge weather station located approximately 3 km from the study area (RAWS, 2021). Hourly rain depth was summed for each day and converted to a binary variable indicating presence or absence of rain. We also used this data source for a daily regional air temperature measure.

Daily ACI, transformed by the base 10 logarithm, was the response variable in the models because we were interested in what explanatory variables influenced observed patterns in ACI. This transformation expressed the premise that explanatory variables had multiplicative effects on ACI, so the log transformation resulted in an additive model framework. Fire effects were represented by two explanatory variables: Fire History category (FH, Table 1 for descriptions) and Days Since Rough Fire (dayRF) was contained (06 November 2015). Additional variables were incorporated to explore effects on ACI: seasonal and interannual conditions were represented as Julian day and weather effects were represented by regional air temperature and a binary variable for measurable rainfall on that day. Although we found evidence for site specific temperature improving model fit, we did not include the variable in the final models because iButton (iButtonLink, LLC) data was not collected during the last year of ARU deployments due to unsafe conditions. In some cases, presence of abiotic sounds can influence ACI, so we included a variable for background sound levels below the frequencies of the biotic sounds we sought to measure.

We fitted 192 models we considered to express plausible relationships between these explanatory variables and ACI. Generalized Additive Mixed Models (GAMM) provided the structure for each model, fitted using package mgcv (Core Team, 2017; Wood, 2017). GAMM enabled us to investigate potentially nonlinear relationships between the explanatory variables and ACI, which seemed certain for seasonality and plausible for other factors. Additionally, GAMM framework enabled us to account for temporally autocorrelated data.

The first, simplest model used Site as a random effect. Two more models were generated by adding Days Since Rough Fire as a separate smooth term for each fire history (common to both) plus alternative representations for Julian Day: a simple smooth or distinct smooths for each fire history category. These three options were doubled to six by adding a set of models that included fire history as a factorial variable. These six options were quadrupled to 24 by adding three options to express the effects of air temperature: a shared smooth, distinct smooths by fire history, and distinct smooths by site. These 24 options were quadrupled to 96 by adding three options to express the effects of background noise level: a shared smooth, distinct smooths by fire history, and distinct smooths by site. The final doubling of model formulae was realized by adding a set of models that included the categorical Rain variable. Temporal autocorrelation of daily ACI was represented by the corCAR1 function on date nested within site variables. All models assumed Gaussian errors in fitting the dependent variable (an identity link function was used). Further, we investigated collinearity of explanatory variables a priori and found low correlation among variables.

Multimodel inference in the R statistical environment was used to evaluate which model provided the lowest AICc values: the best mix of prediction accuracy and concision (Burnham et al., 2011; MuMIn v1.43.17; Core Team, 2017). The best-performing model was then used to interpret the patterns in ACI. Partial conditional dependence plots (e.g., Figure 2) affirmed that residual errors were generally unimodal and symmetric about the predicted values.

Using the predicted ACI values from the best performing model, we compared patterns across all fire history categories, between sites burned by Rough Fire, and within the seasonal peak in ACI.

To understand changes in vegetation and habitat structure post-fire, vegetation surveys were conducted at each site in the beginning of the study (June 2016) and repeated at the end of the study (June 2019). Two or three randomly located plots were sampled near each acoustic recorder site to measure features of the vegetation, including ground and canopy cover and fuel load (Table 3). Plot locations were determined through random point generation in ArcGIS software (ESRI, Inc.), under the conditions that each plot was within a 150 m buffer of the study site ARU. Area sampled was based on the nearest 20 large trees [trees ≥10 cm diameter at breast height (DBH)] and nearest 20 small trees (<10 cm DBH) to the plot center. Cover of bare ground, herbaceous plants, and shrubs were determined using three 15 m point-intersect transects run at random directions from the marked plot center. Overstory canopy cover was determined at plot center using a spherical densiometer with four subsamples collected at each cardinal compass direction and then averaged together. Fuel loads, consisting of litter, duff, and coarse woody debris were measured along the three 15 m transects using Brown’s Inventory methods (Brown, 1974).

Our statistical analysis exploring the relationship between ecoacoustic index (ACI) and a suite of plausible explanatory variables revealed strong evidence for variation explained by time-since-fire and fire history (Figure 2) while controlling for effects of season and air temperature. The best performing statistical model using multimodel inference method (Supplementary Table S1) showed a delta AIC of 12 from the next best performing model and therefore this model was used to interpret the relationships between the explanatory variables and ACI (Burnham et al., 2011). The best performing model explained 69% (adjusted r-sq = 0.69) of the observed variation in ACI and included four explanatory variables: day since Rough Fire was conditioned on fire history category, air temperature, season (Julian Day), and site as random effect (Table 2).

Evidence for the influence of fire on ACI was expressed as distinct smooth terms for days since Rough Fire for each fire history category (Table 2; Figure 2). While ACI was similar across all sites in the first spring post fire (∼day 200, Figure 2), the patterns began to deviate in the second spring post fire (∼day 600, Figure 2), with a pronounced decrease at the RF1 category. The RF1 + Rx category showed some variation in ACI after the Rough Fire with lower ACI in the second spring after the fire but then less pronounced fluctuations throughout the remainder of the study (Figure 2). At the sites unburned by Rough Fire, UB1+Rx had low ACI values in the second spring after the fire (∼day 600, Figure 2) and UB1 showed similar, but less pronounced seasonal patterns in ACI throughout the study period (Figure 2) and the only smooth function with low evidence of effect on ACI from days since fire in this fire history category (Table 2). The patterns in UB2 sites showed the strongest seasonal variation expressed by maximum ACI values in early summer and minimum ACI values in mid-winter before a sharp decline beginning around day 1,200 (Figure 2). The UB2 peaks in ACI also occurred earlier in the spring season compared to the other fire history categories (Figure 2).

We found strong evidence that season, modeled as Julian day, was a significant predictor of ACI with higher predicted ACI in late spring months (Table 2; Figure 3A). A comparison of predicted ACI values within the month of June showed differences in magnitude across the fire history categories (Figure 3B). The first June (2016) after the fire had the highest predicted ACI across all sites with a decrease in second season (2017) and evidence for an increase in third season (2018). The predicted ACI in the second spring after the fire (2017) had the greatest differences between fire history categories, especially the RF1 compared to RF1 + Rx (Figure 3). UB1 + Rx also showed a similar decrease in ACI the second spring post fire, similar to RF1. These sites are relatively close to the edge of the Rough Fire (Figure 1). The differences between predicted ACI narrowed in June 2018 and were most similar across all fire history categories. All fire history categories showed lower predicted ACI in June 2019 compared to 2018 (Figure 3), suggesting a regional effect related to drought conditions during the study period, although this was not directly tested. This effect was less pronounced at RF1 + Rx and UB1 fire history categories. Regional air temperature also showed evidence for an effect on AIC (Table 2). AIC was predicted to be lower at high and lower temperatures.

We compared the predicted acoustic complexity at sites burned by the 2015 Rough Fire (RF1 and RF1 + Rx) to evaluate differences in ACI related to prescribed burning prior to wildfire. Across all years, monthly median predicted ACI was higher at the sites with prescribed burning occurring prior to the Rough Fire (RF1 + Rx, Figure 4). The differences in ACI were greatest in the second season after the fire (2017), then the difference between these fire history categories decreased in 2018 and 2019 (Figure 4).

Using the predicted daily ACI in June across all years, we observed evidence for a shift in the seasonal timing peak ACI between RF1 and RF1 + Rx (Table 4). In the first June after the fire (2016) peak day was similar across all sites, but AIC was greater at RF1 + Rx sites. In second and third spring after Rough Fire (2017, 2018) peak day was earlier at RF1 + Rx sites and largest difference in ACI values was in the second spring post fire. In fourth season, peak ACI shifted back to late-June for the RF1 + Rx sites similar to the RF1 sites.

We observed changes in vegetation and habitat structure from the post fire vegetation surveys conducted in spring 2016 and 2019 (Table 4, Supplementary Figure S1) and the changes differed across fire history categories. In general, the vegetation structure on the ground changed to a greater degree at the sites with Rough Fire (RF1 and RF1 + Rx). For example, bare ground was reduced by 50–100% at Rough Fire sites from the growth of shrubs and herbaceous plants, compared to little change or increase at the non-Rough Fire sites (UB1 + Rx, UB1, UB2). Differences in canopy cover varied across sites, yet patterns related to fire history category were less obvious. Canopy decreased by 4% at RF1 and increased at UB2 by 7%. Total fuel load increased at all fire histories plots with the greatest increase at RF1 and the least change at UB2.

Previous studies on avian response trajectories to fire provide important context for interpreting observed changes in ACI related to time-since-fire and fire history. Species abundance is known to increase immediately after low-severity fire (Hutto and Patterson, 2016). In our study, we found strong evidence that the highest ACI occurred in the first spring after the fire and to the greatest degree at the study sites impacted by the Rough Fire with prior prescribed burning (RF1 + Rx, Figure 2), indicating that low-severity burned areas increases ecoacoustic activity. Further, the relatively stable post-fire ACI relationship across the subsequent seasons at RF1 + Rx sites suggested that resident and migratory bird species may have ultimately benefited from low-severity fire; an outcome that has been documented in nearby mixed conifer Sierra Nevada forests (Bagne & Purcell, 2011). The vegetation surveys at the low-severity burned sites (RF1 + Rx) also showed that habitat structure changes likely resulted in increased habitat heterogeneity by reducing forest fuel while maintaining a green, uncharred forest canopy with understory structure (Table 4).

In the absence of prior prescribed fire, the Rough Fire burned at moderate to high-severity at the RF1 sites (Figure 1; Table 1) and we did not see the same response trajectory in ACI compared to the low-severity burned sites (RF1 + Rx). The ACI response pattern at RF1 sites exhibited the most pronounced successive decreases in the first and second years after the fire (Figures 2, 3). By the third year of monitoring (2018) the ACI response pattern at RF1 rebounded and predicted ACI values were most similar to the RF1 + Rx sites (Figure 3), indicating habitat conditions at the moderate to high-severity burned sites improved by the third season. Vegetative surveys in the fourth year (2019) showed a dramatic increase in shrub cover compared to the first season after the fire (Table 4). However, canopy cover was reduced by 4% which may have shifted the bird community compositions to early succession bark foragers and cavity nesting avian species (Stephens et al., 2015). While this study did not investigate shifts in species composition related to post fire conditions, future analysis of the archived acoustic data may further elucidate these observed shifts in ecoacoustic activity in relation to species present near these sites.

The sites unburned by the Rough Fire (UB1 + Rx, UB1, and UB2) offered some context for how fire changed ecoacoustic activity post-wildfire. ACI similarities between RF1 + Rx and UB1 + Rx suggested that overall lower sustained ACI measured at RF1 could be avoided with increased and active prescribed fire management program. However, the unburned sites are all relatively close to where the Rough Fire burned (Figure 1) and our results suggest these sites were impacted by the presence of the nearby Rough Fire, therefore not considered true control sites. Specifically, while the UB1 + Rx had consistent and sustained higher predicted ACI, similar to RF1 + Rx, the response trajectory post-fire was more similar to that of RF1 (Figure 2) with a second spring season post-fire reduction in ecoacoustic activity (Figure 3). These results indicated that wildfire likely changed ecological conditions in nearby unburned forests, even if habitat structure did not change (Table 4).

Prior to the Rough Fire, the fire history at the RF1 sites was like that of the UB1 site (Table 1), both late succession giant sequoia forest. Thus, comparing the RF1 and UB1 response trajectories likely reflected the effect of high-severity fire in late succession forests when fire has been excluded from the system. Compared to RF1, UB1 exhibited higher overall predicted ACI, yet less seasonal ACI variation (Figure 2). The lower ACI variation across season likely related to a lack of early successional bird communities and may even indicate lower species diversity and density, as has been found in other late successional forests in Oregon, United States (Cushman and McGarigal, 2003). Even with reduced habitat structure and ecoacoustic activity at moderate to high-severity burn sites (RF1) when compared to RF1 + Rx, the Rough Fire may have increased species diversity in an otherwise late succession giant sequoia forest when compared to UB1.

Within the Grant Grove peninsula, it was difficult to find late succession sequoia forests that were not previously burned to compare to RF1. The only other options besides one site within UB1 for representing fire suppressed sequoia forests was to include sites in previously logged second generation forests (UB2), which were geographically appropriate but dissimilar in habitat structure (Figure 1; Table 4). Most giant sequoia logging in this area took place in the late 1800 s, and after an early successional herbaceous stage, the logged forests rapidly transitioned to dense thickets of second-generation giant sequoias with very little understory diversity and structure (Biswell et al., 1966). In the absence of disturbance, this condition persisted at UB2, and the resulting ACI suggested strong and consistent seasonal patterns with a pronounced negative response towards the end of the study. The sharp decline in ACI at UB2 towards the end of the study might indicate local movement away from the fire-suppressed, early successional UB2 forests toward the nearby newly recovering areas of the Rough Fire. Further evaluation of species composition across all years is necessary to understand how shifts in ecoacoustic activity may be connected based on responses across the region.

Drought conditions across all monitoring sites likely influenced the observed decrease in ACI in the final year of monitoring (2019, Figure 2, 3) but was not explicitly evaluated in this study. Other studies have documented Sierra Nevada avian population declines from drought and climate change (Roberts et al., 2019). Assessment of multiple and interacting threats to this forest system complicated the observed relationships to fire history presented in this study yet provided some evidence that could be investigated in future studies. To further understand these interacting disturbance regimes on ecological conditions, longer term acoustic monitoring to compare to measures of disturbance is needed.

The daily temporal resolution of the ecoacoustic index showed evidence for phenological shifts in ecoacoustic activity related to fire history category. Daily predicted ACI increased to a local maximum in June in all years at most sites (Figure 3), yet we found that peak activity in June differed between the RF1 and RF1 + Rx sites (Table 3). Later peaks at RF1 sites in 2017 and 2018 suggested that birds may utilize less disturbed habitats before moving to more disturbed habitats across this forest landscape. These changes in phenology through ACI further elucidated post-fire effects across habitats by signaling the onset, magnitude, and timing of seasonal ecoacoustic activity and how these interactions maybe related in a pyrodiverse giant sequoia forest landscape.

Although we documented evidence for acoustic patterns that correspond to ecological changes known to have occurred in the first four years after a fire, longer term studies could address additional suites of questions, such as succession, ecological stability, and fire regime shifts (Fontaine and Kennedy, 2012), as well as interactions with broader regional pattern such as drought. The longer-term ecosystem responses to fire that have been documented (Stephens et al., 2015; Tingley et al., 2016) likely have concomitant ecoacoustic expressions. The modest cost of sustaining a passive acoustic monitoring network and near-real time index of ecoacoustic activity encourages integration into long-term research programs, yet integration into fire management strategies remains challenged by established protocols, budgets, and perceptions of the value added.

Future research that adopts ecoacoustic indices can benefit from the pragmatic lessons presented here. First, co-locating acoustic recorders with other environmental sensors would improve predictions by controlling for micro-climate drivers in ecoacoustic activity. Second, we obtained informative results from a single ecoacoustic index, ACI; however, a multivariate response model employing a suite of ecoacoustic indices may reveal more detailed patterns of the ecological community. The elevated complexity of these analyses will be increasingly justified as the diversity of sampled ecosystems and acoustically active communities increases (Buxton et al., 2018). Specific detectors that identify bird species and call types can merge behavioral information into ecoacoustic trends, enhancing opportunities to identify mechanisms and drivers of observed patterns, like the fire-history related reductions and phenological shifts. Automated techniques for extraction and classification of species-specific acoustic signals are rapidly improving (Balantic and Donovan, 2020; de Oliveira et al., 2020; Knight et al., 2020). Lastly, periodic communication of acoustical summaries and near-real time alerts to anomalous events like biological invasions have the potential to inform resource management decisions and actions as events unfold (Wood et al., 2019). The capacity to archive acoustic data generates future opportunities to extract more detailed information from these data.

Understanding environmental conditions often relies on indices derived from remote sensing to provide simplified trends in conditions. For example, satellite derived drought indices are relied upon for monitoring and evaluating drought severity across landscapes (Lydersen et al., 2016). Climate indices are used to characterize an aspect of ocean circulation pattern (Schwing et al., 2002). The public, government, and researchers routinely utilize these indices as benchmarks to inform their decisions about the status of conditions. The strength of these indices is a simple diagnostic quantity captured over broad spatial and temporal scales; the limitations are the simplification of complex dynamics represented in a single index.

While ecoacoustic indices are known to represent features of the biological community (Bradfer-Lawrence et al., 2019; Metcalf et al., 2020; Stowell and Sueur, 2020), adoption of these methods into larger scale assessments is limited. These acoustic monitoring sites sit within a broader regional landscape, where management agencies routinely use remote sensing indices to better understand conditions. We argue that this predicted ecoacoustic index, when collected at daily resolution and processed nearly in real-time, offers a scalable remote sensing index to monitor ecological condition after wildfire and complements vegetation surveys and burn severity indices. While only evaluated in a giant sequoia forest landscape, other forest types may show similar ecoacoustic trajectories after wildfires.