Evidence of Range Shifts in Riparian Plant Assemblages in Response to Multidecadal Streamflow Declines

Riparian corridors are thought to form hydrological refugia that may buffer species and communities against regional climate changes. In regions facing a warming and drying climate, however, the hydrological regime driving riparian communities is also under threat. We examined recruitment in response to streamflow declines for species inhabiting the riparian zone in southwest Western Australia, testing the extent to which the riparian system has buffered riparian communities from the drying climate. We stratified 49 vegetation transects across the >600 mm per annum regional rainfall gradient encompassed by the Warren River Catchment. Local hydrological conditions were estimated over two 10-year periods; 1980–1989, and 2001–2010, to quantify changes in the flood regime. Mixed effects models tested the relationship between rainfall and flooding on the relative frequency of immature to mature individuals of 17 species of trees and shrubs common to the riparian zones. At the low-rainfall extent of their geographic range, the relative frequency of immature riparian species decreased with declining flow, whereas at the high-rainfall extent of their geographic range the relative frequency of immature individuals increased with declining flow. These results suggest that the geographic ranges of riparian species may be contracting at the low-rainfall margin of their range, while at the high-rainfall margin of their geographic range, reduced flooding regimes appear to be opening up new habitat suitable for recruitment and narrowing the river corridor. No such patterns were observed in upland species, suggesting the river may be buffering upland species. We discuss these findings and their implications for ongoing management and species conservation in a region projected to face further, significant rainfall declines.


INTRODUCTION
For many species, survival over the coming decades will depend on their ability to adapt to the new climatic conditions in situ, or shift geographic range to maintain their climatic optimum (Parmesan, 2006;Aitken et al., 2008). In stark contrast to mobile organisms where analyses of distributional shifts have been shown to match climatic shifts (Chen et al., 2011;Lenoir et al., 2020), sessile organisms, particularly those with longer generation times, can be much more constrained in their responses to climate change (Lenoir and Svenning, 2015; but see Ström et al., 2011).
In a strict sense, determination of range shifts in response to climate change requires temporally replicated data (Bertrand et al., 2011;Feeley et al., 2011;Telwala et al., 2013). In lieu of such datasets, range shifts in plant species have been inferred by examining the skew in abundance distributions (Breshears et al., 2008;Murphy et al., 2010), or by exploiting the long generation times and comparing the distribution of seedlings relative to the adult population (Lenoir et al., 2009;Zhu et al., 2012Zhu et al., , 2014Fei et al., 2017). The assumption is made that the range inhabited by new recruits is representative of the optimal climatic envelope within current climate space, while the distribution of adults represents a climate envelope characteristic of former conditions (Lenoir et al., 2009).
Species range expansion is typically observed as the establishment of seedlings beyond the former adult range (Galiano et al., 2010;Vitasse et al., 2012). While range contraction can manifest as recruitment failure at range margins (Zhu et al., 2012;Bell et al., 2014), adult mortality events are typically taken as more 'conclusive' evidence. These mortality events can be evident for long-lived species following catastrophic disturbance events (Allen et al., 2010;Brouwers et al., 2013;Matusick et al., 2013;Harris et al., 2018), or as more gradual declines in growth and crown dieback (Stella et al., 2013). Juvenile recruitment within an existing range, however, can be much more sensitive to incremental changes in environmental conditions and provide an early indication of a site becoming unsuitable (Lloret et al., 2009;Bell et al., 2014;Garssen et al., 2014).
There is a growing body of evidence that topographic and hydrological features in the landscape may reduce exposure to regional climate changes and buffer organisms from broader environmental changes (Dobrowski, 2011;Lenoir et al., 2017;McLaughlin et al., 2017). Riparian zones are predicted to buffer the suboptimal climatic conditions in regions experiencing a warming and drying climate, affording species more time to adapt to the new environmental conditions. The buffering effect of riparian systems in the face of climate warming hinges on the resilience of the hydrological processes under changed climatic conditions (Davis et al., 2019). In regions with warming and drying climates, the hydrological regime driving water availability in riparian systems is also under threat (Barron et al., 2012), placing specialist riparian species under stress from both hydrological conditions and increasing competition from upland species (Garssen et al., 2014).
Here, we test the ecological impacts of multidecadal streamflow declines on the riparian plant communities of southwest Western Australia (SWWA). We examine the relative frequency of immature versus mature individuals of riparian and upland species inhabiting the riparian zones in response to recent changes in local hydrological conditions. We hypothesize, first, that riparian species that are restricted to near-channel habitats, due to their higher sensitivity to surface water availability, will show a contracted geographic range of immature individuals relative to the observed geographic range of the adult population. Second, we hypothesize that less moisture-demanding upland species will be buffered from regional declines if they inhabit the riparian zone, and therefore the distribution of immature individuals is likely to match the distribution of mature individuals for upland species.

Study System
The Warren River, and its major tributaries the Tone River and Murrin Brook of the SWWA are cumulatively about 275 km in length ( Figure 1A) and transect a shallow topographical gradient to a maximum elevation of 385 m asl (Geoscience Australia, 2016). Native vegetation, in either reserves or indigenous forestry, is the dominant land cover across approximately two-thirds of the catchment ( Figure 1A) (DPIRD, 2020). A few subcatchments in the lower Warren are important agricultural regions and withhold water in dammed gullies for irrigation. In these sub-catchments, water abstraction limits are set to ensure a sustainable ecological flow (Department of Water, 2012).
Mean annual temperatures vary little across the catchment (between 14.3 and 15.7 • C). By contrast, mean annual rainfall ranges from over 1200 mm per annum (mm pa) at the mouth of the river, to less than 550 mm pa in the headwaters (between 1901and 1960, Bureau of Meteorology, 2010. A significant decline in the frequency and magnitude of wet weather systems has been observed in SWWA since the 1970's (Hope et al., 2006), resulting in a 16% decline in rainfall, and culminating in reductions of over 50% in surface runoff in some rivers and water storage dams (Petrone et al., 2010). This trend in declining rainfall is predicted to continue under modeled future climate scenarios (Barron et al., 2012;Silberstein et al., 2012;Hope et al., 2015). As the major climatic driver of vegetation types across the region, rainfall declines are predicted to shift optimal climatic envelopes for vegetation communities (Hamer et al., 2015;Ramalho et al., 2017).

Topography
Accurate spatial quantification of topography and vegetation structure running the length of the Warren River was obtained using an aerial LiDAR (light detecting and ranging) survey. A 500 m wide band, spanning the length of the Warren and its tributaries was carried out from the 13th to the 16th of January 2015 by AAM Geospatial Pty Ltd from a fixed wing aircraft using a Q780 laser system with a pulse rate frequency of 180 kHz. The laser returns had a horizontal accuracy of 0.55 m and vertical accuracy of 0.30 m, and were supplied as point clouds, comprised of x: longitude, y: latitude, and z: elevation meters above sea level. The point clouds were classified algorithmically by AAM into ground, infrastructure and vegetation points, and a 1 × 1 m resolution digital ground model (DGM) was interpolated from the ground points. The LiDAR survey was conducted during summer when the river had ceased to flow and much of the upper half of the catchment was dry. This enabled us to capture elevation across the dry river bed or water level of permanent water bodies.

Sampling Design
Vegetation sampling sites were stratified by rainfall isohyet, defining five strata, ≤600, 600-800, 800-1000, 1000-1200 and >1200 mm ( Figure 1A). Within each stratum, 20 potential survey locations, spaced at least 1 km apart and randomly assigned to the true left or true right bank, were randomly generated  (Department of Water and Environmental Regulation, 2016), and the coverage of native vegetation (DPIRD, 2020), on basemaps sourced from Geoscience Australia (2016). (B) The layout of vegetation transects, demonstrated with T19, where transects were rectified using aerial imagery, field photographs and digital elevation models. The elevation of each tree was calculated relative to transect origin. Forest structure was quantified within buffer zones of 2.5 m for individual plants, and 100 m for transects. Note that the river in the canopy surface model is white, where no vegetation points were recorded.
in ArcGIS 10.3.1 (ESRI Inc.). The feasibility of sampling at each of these locations was assessed during site visits, aiming to survey 10 sites per zone. Sites were rejected if: (1) the area was disturbed by human infrastructure such as roads or bridges or the site was visibly impacted by (2) herbicide use; or (3) recent fire. Vegetation survey sites were classified into flood plains or steep banks. Once five sites of either landform had been selected within a zone, all further sites of that landform were rejected to ensure surveys covered a representative range of geomorphic habitats. Difficult access (e.g., steep granite cliffs) within the 1000-1200 mm pa zone meant that only eight sites were sampled in this zone. An additional site was sampled in the <600 mm pa zone to increase the total area sampled across the narrow riparian zones of the upper catchment. A total of n = 49 sites (≤ 600 mm, n = 11; 600-800 mm, n = 10; 800-1000 mm, n = 10; 1000-1200 mm, n = 8; > 1200 mm, n = 10) were surveyed once each, during one of two consecutive summers, December 2013 to April 2014, and November 2014 to May 2015.

Vegetation Sampling
A 10 m wide transect was run at each site from the water's edge (the edge of the dry river bed in the upper catchment) to the end of the riparian zone. Transect length was determined by site geomorphology, and a change in the dominant vegetation type, varying in length from 5 to 90 m ( Figure 1B). The coordinates at the ends of each transect were marked using a GPS unit (GPSMAP R 62s, Garmin) and the location of each plant rooted within the transect recorded. The position of each plant relative to the transect line was recorded to the nearest 0.5 m. All trees and shrubs were identified to species level following the nomenclature of the Western Australian Herbarium (Western Australian Herbarium FloraBase, 1998).
As the majority of the trees and woody shrubs in the region (predominantly Fabaceae, Myrtaceae, and Proteaceae; Table 1) retain a woody capsule/fruit after seed set, each plant was searched for the presence of fruit or flowers to assess whether an individual was reproductively immature or mature, as a binary response. Agonis flexuosa was the only species observed to sucker and reproduce vegetatively, but without a means to identify an immature individual decisively as the product of sexual or vegetative reproduction both were included and the results treated with added caution. Astartea leptophylla and, less often, Eucalyptus rudis were observed to layer where a branch had been damaged, or coppice from the base of mature individuals, and these were considered part of the parent plant. In these low-power systems (Pettit and Froend, 2001b), layering and coppicing seemed to serve as an adaptation to extend the life of a single diseased or damaged individual, rather than as a means to replicate (H. White personal observation). A binary response of reproductive status was used as opposed to size-class structure based on heights or stem diameters for two reasons: first, to estimate broad 'recruitment' trends and reduce the temporal bias of the single time point sample; and second, a size-based age classification could not be consistently applied across the length of the catchment due to the variation in growth rates imposed by the rainfall gradient.
The transect coordinates were spatially rectified to the DGM, LiDAR vegetation point clouds and field photographs in ArcGIS 10.3.1 to account for field GPS error. The position of each plant within a transect was spatially adjusted to the rectified transect position to obtain corrected geographic coordinates. The absolute elevation (m asl) and elevation relative to the transect origin (i.e., above base flow) for each plant was extracted from the DGM and used to calculate hydrological parameters.

Streamflow
To investigate species distribution patterns in relation to surface water and flood regime, parameters describing ecologically important aspects of flow were estimated using the DGM and maximum daily stage height data obtained from the West Australian Department of Water and Environmental Regulation (DWER) for four gauging stations situated in the main channel along the Warren and Tone Rivers ( Figure 1A, Hiller Road, Bullilup, Wheatley Farm and Barker Road Crossing, Department of Water and Environmental Regulation, 2016). Past and recent flow conditions were estimated by selecting two 10-year periods. 1980 to 1989 represents a past condition and was selected as the first decade where three of the gauges on the main channel were in place (Barker Road was established April 1966, Wheatley Farm in 1970and Bullilup in 1978. The second period, 2001 to 2010, was selected to represent recent conditions based on the availability of data from the Hiller Road gauge station, which was established in June 2000, but was then offline from 2013 to 2015 (Figure 2). Unfortunately, continuous flow records are not available for periods predating the 1970s 'step decline' in precipitation, so the selected periods are likely to underestimate overall flow reduction. However, there is likely to have been a lag period between rainfall change and subsequent ecological impacts of flow reduction, therefore we assume that the 1980s period reflects 'low' impacts of flow reduction, while the 2000s period reflects 'high' impacts of flow reduction. Importantly, further significant shifts have been observed in streamflow since the 1980s (Figure 2; Petrone et al., 2010).
To measure the change in the flooding regime between these two periods, linear models were constructed modeling the maximum daily stage height above baseflow at the four gauge stations, as a function of elevation (m asl) for each day of the two 10-year periods (see Supplementary Material 1 for detailed methods). Then, taking the lowest absolute elevation of each vegetation sampling transect (i.e., the dry river channel or summer water level as baseflow) the estimated maximum stage height was interpolated (or extrapolated) at each transect elevation to generate a time-series for each transect site (Supplementary Figure S1.2).
These time-series were used to estimate the frequency and duration of inundation events experienced by each plant based on its elevation to the nearest 0.1 m from 0.5 m above baseflow. Frequency (F) describes the estimated frequency of inundation presented as a value between 0 and 1: where 0 indicates that the individual plant was not inundated during the 10 years period; and 1 indicates that the water level equaled or exceeded the elevation at least once a year for each year of the 10 years period. Duration (D) indicates the mean number of days a plant is estimated to have been inundated per year over the 10-year period, calculated as the mean number of days per calendar year that the water level was estimated to equal or exceed the elevation of an individual plant. The differences in frequency and duration between the historical and recent rainfall periods were calculated, and the resultant differences, F (change in frequency) and D (change in duration), were retained alongside recent frequency and recent duration, respectively, as predictors in statistical models.

Forest Structure
To control for variation in surrounding land use and microclimate on seedling establishment (Davis et al., 2019), we used the LiDAR vegetation point clouds to quantify habitat structure using LAStools (Isenburg, 2017) and ArcGIS 10.3.1. The forest structure was described at two scales, the landscape scale was described by an area encompassing a transect and a 100 m buffer, and the microclimate of individual plants was approximated using a 2.5 m buffer around the point location of each tree or shrub ( Figure 1B). Metrics describing forest structure were calculated from ground normalized point clouds and a 1 × 1 m resolution canopy surface model (CSM) describing the maximum canopy height in each pixel (see Supplementary Material 1 for detailed methods). At the transect scale, we obtained the following metrics for each transect and its surrounding 100 m buffer: the maximum point height, and the laser penetration rates through six vertical height strata: the penetration rate to 24 m; penetration through 24 to 16 m; 16 to 8 m; 8 to 3 m;  Australia, 2021), and Flora of the South West (Wheeler et al., 2002). *Rainfall range estimated from species records obtained from the Atlas of Living Australia (2021). Rainfall range presents the 5th and 95th percentiles of mean annual precipitation (Precipitation -Annual bio12). † Traits observed by the authors in the field. ? Trait uncertain.
3 to 0.5 m and penetration to ground level (≤0.5 m) from the point clouds. A further four metrics were obtained from the CSM: the range, mean, coefficient of variation (CV) and variance (var) in maximum canopy height across each transect and buffer zone. At the scale of each individual, we obtained the laser penetration rates through the six vertical height strata, and the mean and the maximum height from the CSM around each tree.
To manage collinearities and reduce the number of predictors, a principal component analyses (PCA) was run for each of the transect and individual variable sets in the 'vegan' package (Version 2.4-2; Oksanen et al., 2017). The resultant PC1 and PC2 axes, T_PC1 and T_PC2 (accounting for 73 and 16% of the variation, respectively) at the transect level, and I_PC1 and I_PC2 (accounting for 52 and 18% of the variation, respectively) at the individual tree level, were incorporated as covariates in the following analyses (see Supplementary Material 1 for detailed methods).

Statistical Analysis
We selected 17 species of trees and woody shrubs that were sufficiently abundant (>50 individuals) to investigate the effects of rainfall and shifts in hydrological regime on recruitment. Of the 17 species, six were classed as obligate riparian species, four as facultative riparian species and seven were classed as upland species (Table 1 and Figure 3). As information on the flood ecology of many of these species is limited to just the original taxonomic descriptions in several cases, classification of life history classes was necessarily based on limited published information supplemented by field observations of the authors, adapting definitions by Rood et al. (2010).
A generalized linear mixed model (GLMM) was fitted for each species using a binomial distribution and a logit link (Bolker et al., 2009) in the 'lme4' package (Version, 1.1-23; Bates et al., 2015) in R version 4.0.2 (R Core Team, 2020). The response variable indicates the probability (the 'log odds') that a given plant is immature (coded as '1') rather than mature (coded as '0'), as predicted by the set of modeled environmental covariates described below. For clarity of graphical display, model predictions calculated using the predict.lme function in lme4 were back-transformed on to a 'proportion' scale, rather than log-odds ratios. Values approaching 1 indicate a largely immature population, with few adult individuals and many recruits. Conversely, values approaching 0 indicate a largely mature population, with few recruits, and a low probability of encountering an immature individual ( Figure 4A).
The variables describing transect level (T_PC1 and T_PC2) and individual level (I_PC1 and I_PC2) variation in forest structure were included as fixed covariates to control for variation in land use, and microclimatic conditions independent of the hydrological parameters. To account for non-independence of individuals sampled within a transect, a random intercept was included for transect identity.
Hydrological conditions were defined using combinations of five fixed predictors F, F, D, D, and historic rainfall, and their interactions (Supplementary Material 2). Rainfall was included as historic rainfall only, rather than a measure of change in rainfall, because the two variables were highly correlated, and historic rainfall is more easily interpreted (Supplementary Material 1). Moreover, F and D were also highly correlated (Pearson's r = 0.78, p < 0.05), indicating individuals that were regularly inundated were also inundated for longer durations annually. Rather than discarding one set of the collinear parameters, a two-phase model selection approach was used to identify the parameter set which best fitted each species.
Initially, global models were constructed to test the effects of estimated flood duration (D, D, rainfall inclusive of all twoand three-way interactions) and estimated flood frequency (F, F, rainfall inclusive of all two-and three-way interactions) separately. The global models were simplified using model selection procedures comparing Akaike Information Criterion for small sample sizes (AICc) in the 'MuMIn' package (Version 1.15.6; Barton, 2016). For each parameter set, the most parsimonious model within 2 AIC units of the top model was selected as the 'best fit' model (Arnold, 2010). Subsequently, the AIC of the best flood frequency model was compared to the best model fitted for flood duration, and the final model was taken as the model with the lowest AIC out of either model set. A complete list of candidate models and resultant AIC values are presented in Supplementary Material 2. Prior to analysis, all of the continuous predictors and covariates were centered and scaled by 2 standard deviations (Gelman, 2008). Models were assessed for over-dispersion, however no adjustment was necessary. Model fit was assessed using the Nakagawa and Schielzeth (2013) R 2 approach.

Riparian Species
Reliable models were obtained for three of the six obligate riparian species examined, Astartea leptophylla, Eucalyptus rudis, and Melaleuca rhaphiophylla (Table 2A). Of these, the relative frequency of immature to mature individuals only differed significantly along the examined hydrological gradients in A. leptophylla (Table 2A and Figure 4B). The relative frequency of immature A. leptophylla was greater within regions experiencing a greater decline in annual flood frequency ( F , Table 2A and Figure 4B), particularly in regions with medium to high rainfall. Populations at the low rainfall extent were less likely to be immature individuals ( Figure 4B); no immature A. leptophylla were surveyed in regions drier than 850 mm pa, despite mature individuals ranging out to 640 mm pa (Figure 3). The relative frequency of immature M. rhaphiophylla differed with estimated duration of flooding, change in duration ( D) and their interaction, but with a total of 13 immature individuals of the 144 plants surveyed, the model was not statistically significant ( Table 2A). The modeled co-variates did not explain variation in the relative frequency of immature to mature individuals in the common and wide-ranging E. rudis. Models failed to converge for the remaining three obligate riparian species, Melaleuca cuticularis, Melaleuca viminea, and Taxandria juniperina due to the narrow range of observed variation in responses to predictors (Supplementary Material 2).
The four facultative riparian species, Banksia seminuda, Agonis flexuosa, Hakea oleifolia and to a lesser extent Callistachys lanceolata, revealed relationships with differing aspects of the estimated flow regime and the rainfall gradient ( Table 2A). The relative frequency of immature B. seminuda was higher under high rainfall conditions, and immature individuals were more likely to be found in sites with low flood frequency across its entire range (Table 2A and Figure 4F). While the relative frequency of immature Agonis flexuosa, including both new recruits and immature individuals suckering off nearby plants, was higher under high regional rainfall, there were also significant interaction effects between rainfall, recent flood frequency and its change, indicating that high rainfall alone is insufficient to predict the occurrence of immature individuals (Table 2A and Figures 4C,D). As observed for A. leptophylla, TABLE 2 | Generalized linear mixed effects models testing the relative frequency of immature to mature individuals of (A) riparian and (B) upland species along the riparian zones of the Warren River as a function of mean annual rainfall (Rn), and either inundation duration (D) and change in duration ( D) or flood frequency (F) and change in flood frequency ( F). Variation in forest structure is described at transect and individual level as the covariates, T_PC1 and T_PC2, and I_PC1 and I_PC2, respectively. The proportional change in variance (PCV) for the random effect (transect identity) is calculated between the null and final models. The Akaike Information Criterion (AICc) is a relative measure of goodness of fit scaled to the number of parameters in the model. R 2 GLMM(m) is the marginal variance explained by all fixed factors and R 2 GLMM(c) is the conditional variance explained by both fixed and random factors (Nakagawa and Schielzeth, 2013). NA indicates a term was not tested due to collinearities within the fixed predictor set. In species where the model fit was not better than the null model, results are shown for the null model only. Model coefficients highlighted in bold indicate significant predictors. Note, some models failed to converge due to insufficient variation within the tested environmental variables or age classes, and are therefore not presented for Melaleuca cuticularis, M. viminea, and Taxandria juniperina. †denotes an obligate riparian species. A complete list of tested models is provided in Supplementary Material 2  the relative frequency of immature A. flexuosa individuals under high rainfall increased in floodplain zones that were undergoing declines in flood frequency, but only where sites still experienced a high probability of flooding (Table 2A and Figure 4C). By contrast, sites experiencing declining flood frequency at the low rainfall extent of the geographic range, particularly where flood frequency had declined to less than 0.1 (i.e., flooded once in the 10 years period), had low relative frequency of immature A. flexuosa (Table 2A and Figures 4C,D). For Hakea oleifolia, the riparian zones that underwent the greatest declines in flood duration had the lowest relative frequency of immature H. oleifolia, regardless of rainfall zone or recent duration (Table 2A and Figure 4E). For Callistachys lanceolata, although the best fit models detected differences in the relative frequency of immature individuals with rainfall, flood frequency and change in flood frequency, there were very low frequencies of immature individuals (just 9 of 64 individuals) recorded in total (Table 2A), which severely limited the power of the analysis.

Upland Species
Recruitment varied significantly in relation to the examined hydrological or rainfall gradients in just two of the seven upland species, Trymalium odoratissimum subsp. trifidum and Melaleuca incana (Table 2B, Figure 5, and Supplementary Material 2). The relative frequency of immature individuals was higher in frequently flooded regions for both species (Figures 3, 5). In T. odoratissimum subsp. trifidum, this result was largely driven by the presence of 82 (of a total of 87) immature individuals recorded within a single, low rainfall transect, all of which were present at sites more frequently inundated than observed in the adult population ( Figure 5A).

DISCUSSION
Recruitment failure at a species range margin can be indicative of a disconnect between the geographic and climatic ranges of a species and an advanced warning of an impending range shift. Utilizing one of the world's most striking, geographically stratified rainfall gradients, that has undergone one of the greatest observed declines in rainfall, we tested the effect of streamflow decline on the riparian plant species in SWWA. We show that the relative frequencies of immature versus mature individuals of a number of riparian species differed significantly with the magnitude of divergence from the historical hydrological regime. At the drier, low rainfall, margins of species geographic ranges, declines in streamflow were a key driver of reduction in the frequency of immature individuals, indicative of recruitment failure and impending range contraction at the geographic range margins. At the higher rainfall margins of species geographic ranges, however, immature abundance increased in response to streamflow declines suggesting that riparian communities are expanding into habitats where historically they may have been excluded by high inundation regimes. In contrast to riparian species, the majority of the upland species examined showed little recruitment response to changing hydrological gradients or regional rainfall gradients. This consistency in recruitment could indicate that the river may be stabilizing recruitment processes across the current distributions of upland species despite regional rainfall declines. Here, we discuss these findings and their implications for ongoing management FIGURE 4 | Variation in the relative frequency of immature to mature individuals of riparian tree and shrub species modeled as a function of their hydrological conditions. A positive slope indicates an increasing likelihood of immature individuals being found at sites experiencing a greater decline in flooding; conversely, a negative slope indicates a decreasing likelihood of immature individuals being found at sites experiencing a greater decline in flooding. Hydrological parameters examined include mean annual rainfall (mm pa, percentiles), recent flood frequency (for the period 2001 to 2010), and the change in flood frequency (between two the ten-year periods 1980 to 1989 and 2001 to 2010). Flood frequency is presented as the probability that individuals are likely to be inundated at least once in any one calendar year, where 1 indicates annual flooding and 0 indicates individuals were never flooded. The fitted lines (±95% confidence intervals) represent the 10th, 50th, and the 90th percentiles of model predictions from binomial generalized linear mixed effects models (see methods for details). Note that percentiles for each predictor vary between species because each species is distributed over a distinct rainfall or hydrological range.  in a region projected to face further, significant rainfall declines.

Geographic Shifts in Average Environmental Conditions for Riparian Species
The declines in streamflow observed over the past 30 years have resulted in a marked change in recruitment for a number of riparian species in response to declining frequency and duration of inundation. Declines in flood frequency interacted with rainfall for many species. The probability of immature occurrence was lower at the low rainfall extent of a species geographic range in the riparian species Agonis flexuosa, Astartea leptophylla, and B. seminuda. Although statistically insignificant, it is worth noting that in the sampled populations of a further species, riparian obligate M. rhaphiophylla and facultative riparian species, C. lanceolata, immature plants were observed at estimated relative frequencies of just 9 and 14% respectively. Cumulatively, the results presented here demonstrate lower occurrences of immature individuals at the drier extent of these species geographic ranges, which could precede a contraction in range or selection for a shift in their climatic optima.
The major assumption in examining distributions at a single time point to deduce range mismatch between immature and mature populations, is that differences are indicative of a shift in environmental conditions away from climatic optima rather than the natural divergence between the recruitment niche and the mature niche (Grubb, 1977;Máliš et al., 2016). This phenomena is particularly apparent in riparian systems, where early establishment is highly dependent on surface and shallow soil water until root systems gain access to permanent groundwater sources (Mahoney and Rood, 1998;Capon et al., 2016). Moreover, mature vegetation has the potential to significantly alter its own flow regime over its lifetime by redirecting currents and altering depositional processes (Corenblit et al., 2007;Osterkamp and Hupp, 2010). Here, by including estimates of recent frequency and duration of inundation as independent parameters from the observed changes over time, our results strongly suggest that it is the change in streamflow rather than (or in addition to) the absolute streamflow driving the range mismatch.
The relationship between occurrence of immatures and streamflow decline suggests species are failing to recruit under current conditions, be it through sexual reproduction or through vegetative reproduction as in the case of Agonis flexuosa. While in many riparian systems, a high variability in disturbance and flood regime are known to drive recruitment events, the high predictability and low interannual variance in the flood regime in these river systems has been suggested to drive a stable and continuous recruitment regime (Pettit and Froend, 2001b;Pettit et al., 2001). An examination of the reproductive phenology of the two obligate canopy trees, M. rhaphiophylla and E. rudis, indicates that phenological coupling has evolved to time seed fall with drying weather and the retreat of the river channel each summer (Pettit and Froend, 2001b). For most of the riparian species in this region, seed is stored in the canopy ( Table 1) and release is triggered by environmental cues, namely, the desiccation of the seed capsule. In M. rhaphiophylla, this results in a continuous low level of seed fall throughout the year, which is similar in E. rudis except that E. rudis has a peak during the early summer coinciding with the summer decline in stream flow, exposing suitable germination sites (Pettit and Froend, 2001b). Astartea leptophylla is the only obligate riparian species which is not reported to store seed for any length of time, it flowers over summer and produces seed quickly (Rye, 2013). The reduction in winter rainfall, and subsequent stream flow in the SWWA has not been even; instead, the greatest change observed has been a reduction in rainfall in autumn and early winter which acts to extend the length of the summer period by delaying the onset of wet winter conditions Hope et al., 2015). The observed decline in recruitment in the drier extent of species ranges could be due to drying conditions restricting seedling establishment, attributable to the lower rainfall itself, or to the greater intermittency of surface waters at the lower rainfall sites (Stromberg et al., 2005;. In this case, the changing hydrological conditions could be acting differently for different species. For spring seeding species, a longer summer season could reduce survival in the first few summers before the roots have access to the water table, whereas for summer and autumn seeding species there may be a lack of suitable recruitment sites with an extension to the summer dry season. As coarse estimates of hydrology are used here, further investigation is required to fully understand the demography and recruitment requirements of these species, particularly for A. leptophylla where seed fall occurs during summer and autumn and does not seem to be triggered by environmental cues. The assumption here is that the change in environment is acting on seedling establishment, but equally a reduction in seed production or viability could also be driving the decline. Longer periods without surface waters could reduce seed production from stressed mature plants (Jensen et al., 2008), or lower reproductive success could also be attributed to poor pollination success from fragmentation, particularly in the upper catchment (Hopley and Byrne, 2018), changes to biotic interactions, or a number of these factors acting synergistically.

Evidence for Narrowing of the Riparian Corridor
In contrast to the reductions observed in immature abundance under lower rainfall conditions, declining flood frequencies under higher rainfall conditions increased the relative frequency of immature to mature individuals in the examined riparian species. Previous examination of size-class distributions on the neighboring Blackwood River found no differences in the distributions of E. rudis or M. rhaphiophylla seedlings and mature trees in relation to elevation across the floodplain, and determined that these low-power river systems do not experience the destructive forces which can restrict succession on more powerful rivers (Pettit and Froend, 2001b) and can lead to a naturally high abundance of new recruits and saplings on the high flood frequency sites. Increases in riparian seedling abundances, as well as vegetation density and cover, are known to result from flow reduction due to damming or water extraction (Shafroth et al., 2002;Gordon and Meentemeyer, 2006;Kingsford, 2016), particularly within facultative riparian species (Rood et al., 2010). The initial increase in vegetation cover post-damming, is principally attributed to increases in germination sites with declining flood waters, i.e., moist, damp sediments, as well as a reduction in the erosive flows seasonally clearing establishing seedlings (Mahoney and Rood, 1998;Taylor et al., 1999;Johnson, 2000;Polzin and Rood, 2006). The increases in the frequency of immature individuals observed in the areas of greatest deficit, is a strong indicator that the riparian corridor may be beginning to narrow; a repeat survey of selected sites in the future would confirm this.

Stability in the Upland Populations
In five of the seven upland species examined, the distribution of immature and mature individuals did not differ with regard to metrics describing aspects of streamflow or regional rainfall. Notwithstanding the hydrological parameters, rainfall is considered one of, if not the most important, abiotic determinants of species distribution in the region (Hopper and Gioia, 2004). Consistency in the relative frequency of immature to mature individuals across the rainfall gradient indicates stable range margins within the riparian zones. In species where the surveyed area included the easternmost limits of their distribution, such as the upland shrubs L. propinquus, L. obovatus subsp. revolutus, and H. elliptica (Western Australian Herbarium), these results may be indicative of the river buffering species from regional rainfall declines observed to date (Reside et al., 2014;McLaughlin et al., 2017). The ability of the riparian corridors to extend the range of a species outside normal climatic conditions has been demonstrated in SWWA. Trymalium odoratissimum subsp. trifidum (formerly known as T. floribundum; Hancock et al., 1996) is considered an upland species in the Warren Region, in lower rainfall forests to the north of the study region, however the northern subspecies T. odoratissimum subsp. odoratissimum is restricted to gullies (Hancock et al., 1996). Further examination of the relative frequencies of immature to mature across the non-riparian extent of their ranges is critical to understanding range change in the upland species.

Implications for Climate Change
Over the past 30 years, the riparian vegetation of the Warren Catchment has been subjected to reductions in mean estimated duration of inundation of up to 27 days per year and sites are becoming inundated over fewer winters. The results presented here demonstrate that these declines are affecting recruitment in the common riparian species. At the drier extent of species ranges, declines have resulted in a lower relative frequency of immature to mature individuals, potentially presenting early warnings of a longitudinal range contraction. Meanwhile at the higher rainfall extent of the catchment, increasing frequencies of immature plants on riparian floodplains experiencing declines in inundation frequency indicate the expansion of the riparian vegetation into areas that were previously uninhabitable, and potentially narrowing of the river channel. Downscaled climate models for the SWWA project declines of between 5 and 75 fewer flow-days per year by 2030, on top of the deficits already observed (Barron et al., 2012). Given the apparent shifts in climatic optima already observed here, further flow reductions are likely to significantly impact the riparian vegetation. How these impacts will manifest remains to be seen, but the results presented indicate that a geographical contraction of range of the riparian species is likely. The majority of the riparian species (both facultative and obligate) did not show an upper rainfall limit to their distribution, i.e., all species were observed within the lower reaches of the river, so there is almost no potential for compensatory range expansion.
Despite these observations, we do not anticipate a complete collapse of the riparian flora. First, the projections for summer rainfall are highly uncertain (Hope et al., 2015), but have the potential to drive selection away from the stable continuous recruitment pattern that currently exists to one which may be more episodic and dependent on environmental extremes. Sporadic summer storms could ease drought conditions for the first, critical summer of seeding establishment. As a number of these species exhibit a serotiny with varying degrees of environmental plasticity in seed release, a boom-bust style recruitment pattern could become more important in the future. For example, wide spread flooding resulting from a cyclonic depression in 1978 has been attributed with a mass recruitment event of E. rudis and M. rhaphiophylla (Pettit et al., 2001;Bureau of Meteorology, 2021). Second, even with significant reductions in river flows, the river is unlikely to cease to flow completely (Barron et al., 2012) thus habitat will continue to be available for riparian species, where soils are too saturated for upland species to survive, albeit over a smaller geographic range. Instead, we might expect to see a compositional shift to a greater proportion of mesic, facultative and upland species as further reductions in the inhibitory high flow events are observed (Froend and Sommer, 2010;Garssen et al., 2014). The results presented here suggest that in this river system, and likely others where rainfall is low in the upper catchment, riparian zones may have a limited capacity to buffer species from climate change induced range shifts. Rather, those species already dependent on the riparian zones as refugia in arid landscapes are actually at greater risk. A complementary study in a system with a contrasting rainfall gradient in which the majority of rain falls in the headwaters would provide a valuable insight into the ability of riparian zones to buffer species range movements under differing precipitation regimes. Understanding which features of the landscape have the potential to provide refugia in a changing climate is critical to advance our ability to predict range shifts and identify important sites for protection.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
HW, JS, and RD conceived and designed the study, contributed to the preparation of the manuscript. HW collected the data. HW and RD analyzed the data. All the authors contributed to the article and approved the submitted version.

FUNDING
Funding for the LiDAR survey and project costs was provided by the Warren Catchments Council through the Australian Government's Clean Energy Future -Biodiversity Fund. HW was supported by CSIRO, The Warren Catchments Council, and the Australian Government Research Training Program and the University of Western Australia.