**5. Discussion**

Flow and salinity have been the subject of scientific observation in San Francisco Estuary over more than a century [14,65]. These observations, which have provided a reasonable understanding of contemporary conditions and associated trends in the estuary over the past century, have also supported decision-making related to freshwater flow management in the estuary. However, important knowledge gaps exist that cannot be addressed solely by the available observations. These gaps relate to the fact that large-scale changes in the estuary and its watershed pre-date the observational record by several decades; these gaps also relate to the fact that California is known to have been subject to highly variable climatic conditions over the past millennium. Additional data collection is thus insufficient to provide an understanding of how the system behaved prior to the initiation of large-scale disturbances after the mid-19th century. This work is an attempt to fill these gaps by (i) using an updated set of tree-ring chronologies to represent annual runoff into the Central Valley from the surrounding higher-elevation watersheds over the past millennium and (ii) utilizing a modeling approach to relate runoff to freshwater flow to the estuary and to salinity intrusion in the Delta. This integration of tree-ring based estimates of runoff with models of flow and salinity representing different configurations of the system (pre-development and contemporary) allows for a more nuanced exploration of flow and salinity changes over periods much longer than covered by the instrumented record. This information is of scientific and practical importance because it can help guide decisions related to the restoration of the estuarine ecosystem. These decisions have recently focused on potential changes to flow and salinity management in the estuary [18], decisions with major environmental and economic consequences for California.

The updated tree-ring-based reconstruction shows a mean annual Central Valley runoff (8RI) of approximately 29 BCM, a quantity that is similar to that observed in the contemporary system. The reconstruction also shows large single-year anomalies from the mean, although multi-year anomalies over averaging periods of 5–100 years are minimal. An important observation is that, while high runoff extremes overlap with the instrumented record, their magnitudes are lower than observed. The reconstruction indicates the occurrence of individual years with runoff significantly lower than seen in the gauged record; however, over longer averaging periods, pre-development runoff variations are of similar magnitude to those in the instrumented period and represent broadly similar patterns of wet and dry periods. Our findings are contrary to some prior reconstructions—notably Stine [47] based on the position of tree stumps at Mono Lake and Graham and Hughes [48] focused on Merced River runoff and the Mono Lake Basin inflow—that indicate evidence of more severe droughts over the past millennium than any

seen in the instrumented period. However, our findings are broadly consistent with the work of Meko et al. [20] that focused on a smaller set of tree-ring data from the Sacramento and San Joaquin River basins. Our work confirms that major long-term deviations from the mean runoff in California's Central Valley are often seen in the tree ring record, but the extended drought of the late 1920s and early 1930s compares with the most extreme in the past millennium.

Consistent with the historical pattern of flow variability, California at the end of WY 2020 appears to be in the midst of another severe long-term drought. The Eight River Index, averaged over the preceding 20-years from WY 2020, stands at 25.6 BCM; long-term average runoff of lower magnitude was last seen in the severe drought of the 1920s and 1930s, when 20-year average flows from WY 1931 to 1939 ranged from 22.4 to 25.1 BCM.

Integration of the tree-ring based runoff with downstream flow models shows that, although median freshwater flows into the estuary have not changed significantly, Delta outflow has been more variable in the contemporary period compared to pre-development conditions. Extending to salinity, the contemporary period is associated with greater spring (February–June) salinity intrusion and lesser summer–fall (July–October) salinity intrusion relative to the pre-development period. Both outflow and salinity intrusion are directly affected by contemporary reservoir storage and release patterns, which tends to smoothen the intra-annual extremes that are seen in the pre-development system.

Our work can also be compared with similar research by Stahle et al. [23] who used tree-ring chronologies to relate to salinity intrusion in the Delta. Using a 625-year (1379– 2003) tree ring chronology to reconstruct the February–June average X2 position in San Francisco Estuary, the authors were able to explain 73% of the variance in the observed X2 data over a 1956–2003 calibration period. The authors used their reconstructed salinity record to examine return intervals between single-year X2 extremes and to quantify the frequency of consecutive seasonal maxima and minima over the period of record. [23] recognized that their reconstruction does not mimic pre-development salinity conditions in the estuary. Rather, the authors concede that their X2 time series:

... *provides an estimate for variability in the salinity gradient on interannual-to-decadal time scales, given the present land cover, stream morphology, and regulated flow environment, in response to the range of modern and prehistoric seasonal precipitation totals registered in the tree-ring record over the past 625 years.*

In other words, the authors employed a "level-of-development" methodological approach [22,24,63] and their 625-year X2 reconstruction represents salinity variability under a contemporary level of development. Stahle et al. [23] found the hydroclimatic signal from tree growth to be approximately stationary over the past six centuries; thus, we expect that their estimate would be consistent with the contemporary level X2 time series presented by Gross et al. [22] that utilized a shorter, more recent climate sequence spanning WYs 1922–2003.

In contrast with [23], our work explicitly attempts to reconstruct an extended record of historic salinity conditions in the estuary. Our work relies on modeled relationships between Central Valley runoff, Delta outflow and estuarine salinity under natural and anthropogenically altered hydrologic conditions as they occurred over the period spanning 903–2008 (long record) and 1640–2001 (short record). We expect that our pre-1850 estimates would be consistent with the pre-development level X2 time series presented by Gross et al. [22] similarly, we expect that our post-1944 estimates (following construction of Shasta Dam) would be consistent with the contemporary level X2 time series presented by [22]. Consistent with Gross et al. [22], our work shows that, in spite of a relatively stationary hydroclimatic signal from tree growth (a proxy for Central Valley runoff), the San Francisco Estuary's seasonal salinity pattern has changed due to anthropogenic alterations. Specifically, both studies suggest that the estuary is now more saline in the late winter and spring than it was under pre-development conditions because of contemporary water management. Both studies also show that the contemporary estuary is less saline in the late summer and fall than it was under pre-development conditions because of water management.

Flows in the observed record from WY 1872 onwards display significant change, with a wet period from the late 19th century to the first decade of the 20th century, followed by a severe drought in the late 1920s and early 1930s. As noted earlier, these were drivers of the early water resources engineering activities in California, with a focus on flood control in the late 19th century to be followed by a focus on regulation and storage by the 1920s and beyond. The longer flow reconstruction in this work highlights that this period in the observed record captures the types of extreme shifts that have occurred over the past millennium. These reconstructed data therefore show that widely used level-ofdevelopment modeling approaches, that repeat the instrumental sequence of flows for different estuary and watershed configurations, appear to be an appropriate methodology for representing a range of pre-development conditions.

The agreement between the tree-ring record with the earliest hydrologic observations confirms and draws attention to a period of dramatic hydrologic change during the late 19th and early 20th centuries, a shift driven by natural factors coupled with rapid regional development. Thus, population growth and extensive watershed modification was overlaid on the underlying hydrologic shift from very wet to very dry conditions, which complicates the task of inferring estuarine changes in the early development period (WYs 1850–1911). Putting this hydrologic shift in the context of other anthropogenic drivers is important in understanding how the estuary responded during this early period and in setting salinity targets for estuarine restoration.

Although tree-ring proxies from long-lived species with high sensitivity to drought are a powerful tool for water resources planning in California, some important caveats are noteworthy. As indicated by the reconstruction statistics, tree-ring width is an imperfect proxy for Central Valley precipitation. In our work, uncertainty is greatest in the early part of the tree-ring record because the data network thins. In particular, the tree-ring record suffers from limited blue oak data, a tree which is shorter-lived than several other less moisture-sensitive species. Our Central Valley runoff reconstructions are especially uncertain in extremely wet years, likely due to a weakening in response of tree growth to changes in precipitation in very wet years [19]. The gaged flows themselves may also be more uncertain under high-flow conditions. The phenomenon can lead to problems with non-normality and non-constant variance (as a function of predicted values) or regression residuals. Transformation of flow before regression (employed in this work) can partly help fix such problems with residuals. Despite these efforts, however, error bars are generally wider in wet years than in dry years [53]).

Another aspect of uncertainty in our Central Valley runoff reconstructions that is not necessarily summarized by calibration and validation statistics is the possible lack of detection of runoff variations at very low frequencies (e.g., wavelengths > 200 years); this results from detrending techniques used to standardize ring widths into annual indices of growth. Specifically, climate trends that span periods longer than the time series of ring width from longest-lived individual trees are removed by detrending. Uncertainty associated with reconstructed flows at very low frequencies could be addressed with alternative detrending methods, such as regional curve standardization or age-banding, which rely on estimation of the curve of expected ring width as a function of tree-ring age [102,103]. Such methods require intensive sampling, demand representative age classes over the entire tree-ring records, and may be possible for a small number of tree species in the study area. The tree-ring data set could be extended in the medieval period (and possibly earlier). Such an extension may be possible if remnant preserved wood can be found for species with strong moisture signals—e.g., *Quercus douglasii*, *Psuedotsuga macrocarpa*, and *Pinus balfouriana*.

The flow–salinity models used in this reconstruction generally assume stationary sea level conditions, although sea level and thus salinity intrusion is expected to have changed over the reconstruction periods. An exception is Andrews et al. [26] in which the pre-development model used a (single) lower sea level value than the contemporary model. Additional detailed mechanistic salinity intrusion modeling may help resolve the impact of

continuous sea level changes in the past millennium and allow refinement of the empirical models employed in this work. Although technically feasible, this is generally limited by the cost and computational complexity of running the mechanistic models under a range of observed sea level values.

In conclusion, this work is an important step in integrating tree-ring data with models of flow and salinity in the San Francisco Estuary and its watershed to develop millennialscale ranges of pre-development conditions. Even though refinements are possible, we believe that our findings support regulatory decision making by providing a baseline to inform future flow regulations and restoration activity in the estuary. While it is helpful to have a target reference range, attainment in future years will continue to be a challenge: much in the system remains highly dynamic and will continue to evolve over time, including sea level, precipitation, snowmelt, and runoff patterns, all of which are expected to be affected by climate change, with limited predictability at present.

**Supplementary Materials:** The following are available online at https://www.mdpi.com/article/ 10.3390/w13152139/s1, Supplementary Material A. Metadata of tree-ring chronologies (xlsx file) used in development of Model 1. Two sheets—a table followed by key to columns. Supplementary Material B. Tree-ring Standardization and Flow Reconstruction (pdf file). Descriptions of how treering widths were standardized into site chronologies and how those chronologies were converted to estimated WY flow related to Model 1. Supplementary Material C. Statistics of Single Site Reconstructions (SSRs) (pdf file). Statistics apply to the regression models as calibrated by regression of square-root-transformed annual flow on lagged residual tree-ring chronologies and their squares. Identifies which chronologies are used in Model 1 long record and short record runoff reconstructions. Supplementary Material D. Supplemental figures on outflow and salinity modeling related to Model 3. Supplementary Material E (Referenced in Supplementary Material A). Listing of Annual Reconstructed Flow (tab-separated ASCII text file). Listing of annual observed and reconstructed Central Valley runoff (8RI), along with 50% confidence interval for Model 1 long record and short record runoff reconstructions. Supplementary Material F (Referenced in Supplementary Material A). Compressed tree-ring measurement files. These are compressed ASCII text files of the measured ring widths used to develop the tree-ring chronologies considered for use in the Model 1 annual runoff reconstructions. Individual file names are linked to the site information (metadata in Supplementary Material A) through column "FilePrefix." Supplementary Material G (Referenced in Supplementary Material A). Listing of annual values of residual site tree-ring chronologies. Tab-separated ascii time series matrix of tree-ring chronologies listed in Supplementary Material A. Data columns in same order as chronologies are listed in Supplementary Material A. Row 1 is header that is the column number followed by the SiteCode field from Supplementary Material A. First column is the year. Following rows correspond to tree-ring residual chronologies data for years 899–2016 at each of the 69 sites.

**Author Contributions:** P.H.H. conceptualized the study; D.M.M. and P.H.H. performed the tree ring and salinity modeling, respectively. P.H.H., D.M.M., and S.B.R. synthesized the findings wrote the manuscript. All authors have read and agreed to the published version of the manuscript.

**Funding:** Funding for this work was provided by the Metropolitan Water District of Southern California through Agreement Number 143878 with Tetra Tech Inc. (Pasadena, CA, USA).

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** The data presented in this study are available in Supplementary Materials A–C and E–G.

**Acknowledgments:** The authors would like to acknowledge David Stahle for sharing his historical salinity data and Michelle Goman for her insights on paleosalinity conditions in the estuary.

**Conflicts of Interest:** The authors declare no conflict of interest. The funders reviewed the overall research plan and reviewed a draft of the manuscript; the study design, analyses and interpretation of the results were performed by the authors.
