Fewer and farther between: changes in the timing of Longfin Smelt (Spirinchus thaleichthys) movements in the San Francisco Estuary

Abundance of estuarine fish species has declined globally. In the San Francisco Estuary (SFE), long-term monitoring documented declines of many species including the anadromous species Longfin Smelt (Spirinchus thaleichthys). To improve management and recovery planning, we identified patterns in the timing, seasonal occupancy, and distribution of Longfin Smelt in a monitoring study (San Francisco Bay Study) for five regions of the SFE using a generalized additive model. We then investigated the year-toyear variability in the shape of the seasonal relationships using functional data analysis (FDA). FDA separated the variability due to population size from variability due to differences in occupancy timing. We found that Longfin Smelt have a consistent seasonal distribution pattern, that two trawl types were needed to accurately describe the pattern, and that the pattern is largely consistent with the hypothesized conceptual model. After accounting for variability in occupancy due to year-class strength, the timing of occupancy has shifted in three regions. The most variable period for the upstream regions Suisun Bay and Confluence was age-0 summer and for the downstream region Central Bay, was age-0 late fall. This manifested as a recent delay in the typical fall re-occupation of upstream regions, reducing Longfin Smelt abundance as calculated by another monitoring study (Fall Midwater Trawl); thus, a portion of recent reductions in Fall Midwater Trawl abundance of Longfin Smelt result from changes in behavior rather than a decline in abundance. The presence of multiple monitoring surveys allowed analysis of distribution from one data set to interpret patterns in abundance of another. Future investigations will examine environmental conditions as covariates during these periods and could improve our understanding of what conditions contribute to the shifting occupancy timing of Longfin Smelt, and possibly provide insight into the long-term quality of the San Francisco Estuary as habitat.


Introduction
Declines in the abundances of estuarine fish species have been documented globally and these declines have been linked to a variety of causes. For example, in estuaries in southern Australia hydrologic extremes caused by drought resulted in reduced catch of several estuarine species (Gillson et al. 2009). The loss of eelgrass bed habitat in the estuaries of Massachusetts was associated with reduced fish diversity and biomass of the fish community (Hughes et al. 2002). In southeast Louisiana, a long-term decline in apex predators, particularly Bull Sharks (Carcharhinus leucas) and Alligator Gar (Atractosteus spatula), was associated with declining water quality, habitat destruction, and overfishing (O'Connell et al. 2007). In Brazil, increases in precipitation into an estuary produced declines in euryhaline species and replacement with freshwater species (Garcia et al. 2001). A review of declines in estuarine fish species concluded that overfishing and environmental degradation were two of the main drivers of change (Whitfield 1996).
Like other estuaries, the San Francisco Estuary (SFE) has experienced declines in many species of fish, including both endemic (e.g., Delta Smelt, Hypomesus transpacificus, and Splittail, Pogonichthys macrolepidotus) and wide-spread (e.g., Steelhead Trout, Oncorhynchus mykiss, and Longfin Smelt, Spirinchus thaleichthys) species of management concern (Moyle et al. 2015. The SFE is the largest estuary on the west coast of North America and is the terminus of a watershed that drains approximately 40% of the area of California (Jassby and Cloern 2000). The SFE serves as a primary source of water for much of the state, providing drinking water for over 25 million people and irrigation water for 750,000 acres of farmland (CDWR 2015). Water export pumping stations operated by State and Federal water operators, which started in 1950s, remove up to 65% of the freshwater that flows into the estuary (Arthur et al. 1996;Sommer et al. 2007). Water operations have also reduced the seasonable variability in salinity patterns and freshwater flows through the system (Jassby et al. 1995) and have directly entrained millions of fishes and other organisms (Arthur et al. 1996;Grimaldo et al. 2009). Early in the development of the water export facilities, managers recognized the potential for detrimental effects and initiated studies to document and assess resources at risk (Erkkila et al. 1950;Kelley 1966;Turner and Kelley 1966). Over time, long-term monitoring was implemented to track changes and the effects of water management (Chadwick 1964;Stevens 1977;Kimmerer and Orsi 1996;Lehman 1996, Schemel et al. 1996Honey et al. 2004).
The suite of long-term monitoring studies has documented a substantial loss of pelagic productivity in the SFE in the late 1980s (Kimmerer and Orsi 1996;Kimmerer 2002), but in the early 2000s, the abundance of several fish species collapsed, in a phenomenon referred to as the Pelagic Organism Decline (POD; Sommer et al. 2007). Among the POD species whose abundance collapsed is the Longfin Smelt (Spirinchus thaleichthys), once one of the most abundant pelagic fishes in the estuary (Orsi 1999). Freshwater outflow (or a correlated variable) has been and continues to be considered the driving factor in Longfin Smelt abundance in the SFE (Stevens and Miller 1983;Jassby et al. 1995;Kimmerer 2002;Nobriga and Rosenfield 2016;Tamburello et al. 2018), though in the years after the introduction of the overbite clam, P. amurensis, declines in copepod Eurytemora affinis and the virtual loss of the mysid Neomysis mercedis as principal food sources had a negative effect on Longfin Smelt abundance (Kimmerer and Oris 1996;Orsi and Mecum 1996;Kimmerer 2002;Mac Nally et al. 2010; Thomson et al. 2010); adult stock size also exerts a strong influence on abundance of age-0 fish (Nobriga and Rosenfield 2016).
To document these declines, diagnose their causes, and measure the effectiveness of interventions, researchers use data from monitoring programs to calculate abundance metrics. In the SFE, Longfin Smelt abundance is monitored primarily by two long-term monitoring surveys: the Fall Midwater Trawl Survey (FMWT) and the San Francisco Bay Study (Bay Study; Honey et al. 2004;Rosenfield and Baxter 2007). The FMWT samples only upstream regions, and has been frequently used to identify drivers of the POD, however it has proven ineffective for explaining Longfin Smelt abundance. The Bay Study samples throughout the entire range of the Longfin Smelt year-round, so may provide a more comprehensive view of longfin smelt use of the SFE (Figure 1).
Effective ecological management requires methods for data analysis that can connect the time series into a coherent story about population trends as well as the factors that contribute to those trends. We developed a conceptual model of Longfin Smelt regional distribution within the SFE throughout its life cycle ( Figure 2). The conceptual model was developed to describe potential changes in distribution prior to developing an analytical approach to analyzing the data for evidence of such changes; thus, the conceptual model serves as a hypothesis. The conceptual model indicates that age-0 Longfin Smelt recruit to the gear in spring in all regions except South Bay ( Figure 2); South Bay recruits emigrate prior to size of counting in Bay Study gear. Presence in upstream regions (West Delta and Suisun Bay) declines in summer through early fall, before recovering in late fall and winter in all regions. Presence is high in all regions until the following spring when the pattern of decline and recovery repeats itself such that fish are again present throughout the estuary by late fall. Age-2 fish are not present in the estuary after the spring spawning season. This conceptual model does not show the periodic presence of soon-to-be age-3 fish in December of their age-2 year ( Figure 2).
The goal of this paper is to develop analytical means to test the conceptual model and the hypothesis that Longfin Smelt distributions and timings of region use have changed over time. We used a datacentered approach, meaning that we developed a statistical model to identify spatial and temporal patterns in the presence of Longfin Smelt and then compared the geographic shifts indicated by that statistical model to our conceptual model. Specifically, we asked: 1. Do patterns in presence of Longfin Smelt support the patterns hypothesized by the conceptual model? 2. Does the recent pattern of occupancy differ from the historical pattern spatially? 3. Does the recent pattern of occupancy differ from the historical pattern temporally? 4. If so, what are the implications for interpretation of recent abundance trends?
The outcome of our approach is a set of functions that describe seasonal and long-term fluctuations in distribution by region, which together describe how Longfin Smelt use the SFE. The functional data approach allowed us to separate the effects of abundance from the effects of shifts in timing on the probability of observing Longfin Smelt in different regions of the SFE.

Study species
Longfin Smelt (Spirinchus thaleichthys) is a small (≤150 mm FL), facultatively-anadromous, pelagic fish belonging to the family Osmeridae (true smelts). It is native to coastal waters of the eastern Pacific Ocean and bays and estuaries, ranging from Monterey Bay to Hinchinbrook Island, Prince William Sound, Alaska (Moyle 2002) and possibly to Bristol Bay, Alaska (Gilbert 1895, Dryfoos 1965, Garwood 2017. Longfin Smelt tolerate salinities ranging from fresh to marine (euryhaline). Both anadromous and resident populations exist at various locations, and the San Francisco Estuary population is anadromous (Dryfoos 1965;Moulton 1974;Rosenfield and Baxter 2007). The Longfin Smelt of the San Francisco Estuary represent the southern-most reproductive population and appear genetically distinct from other populations, although there is evidence of gene flow from the SFE population northward to the Humboldt Bay and Washington populations (Garwood 2017;Saglam et al. in press). Based in part on its precipitous abundance decline, Longfin Smelt was listed as threatened under the California Endangered Species act in 2009 (OAL 2010) and the USFSW determined that protection of Longfin Smelt under the U. S. ESA is warranted (U.S. Office of the Federal Register 2012).

Study area
The SFE is comprised of four, large, open-water embayments (South Bay, Central Bay, San Pablo Bay, and Suisun Bay); the confluence of the Sacramento and San Joaquin rivers (Confluence) also described as the West Delta; and the Sacramento-San Joaquin Delta, which is a network of tidal channels to the east (upstream) of the Confluence; and remnants of the emergent wetlands that surround these aquatic habitats ( Figure 1). We refer to the embayments and confluence as "regions" in our descriptions and analyses. Water in the SFE generally flows from east to west and it becomes more saline as it flows towards the ocean (Kimmerer 2004). The outlet for the SFE is the Golden Gate, which connects Central Bay to the Pacific Ocean. The SFE experiences a Mediterranean climate: cold wet winters and springs with rain in the valleys and snow in the Sierra Nevada Mountains followed by hot and dry summers and falls. Historically, river discharge resulted from low elevation rain runoff in winter and snow melt in spring and early summer followed by declining flows through the late fall.

Data
We used data collected by the California Department of Fish and Wildlife's San Francisco Bay Study (Bay Study) monitoring program (San Francisco Bay Study 2017). The Bay Study is one element of the suite of long-term monitoring studies conducted by the Interagency Ecological Program for the San Francisco Estuary. The Bay Study samples fish and collects environmental data at fixed sampling locations throughout the estuary once per month, year-round. We used temperature from Bay Study sampling, described below, and flow data from California Department of Water Resources Dayflow dataset (https://data.ca.gov/dataset/dayflow) to investigate monthly trends in those measures within the estuary. Water measurements (±0.1 °C) were collected once at each sampling station each month. Daily mean Delta outflow from Dayflow was grouped into months for the period January 1988 through December 2015 and portrayed in a box plot. We summarized seasonal patterns in these data using boxplots for pooled monthly data from 1988-2015.
At each sampling location, the Bay Study collects fish using two different gears: an otter trawl catches fish at the bottom of the water column and a midwater trawl collects an integrated sample throughout the water column. The Bay Study added stations over time, introducing three subsequent series of stations to the "series 1" or original 35 stations: Series 2 stations (n=7) were added in 1986 to provide a more balanced number of channel stations (> 7 m depth) and shoal stations in each embayment; series 3 and 4 added stations within the Delta in 1991 and 1994 to provide better coverage of the distributions of estuarine organisms believe to be extending upstream into the Delta in a drought. For this analysis, we used series 1 and series 2, which are the two longest-running groups of sampling sites. Although stations in series 1 represent those originally sampled starting in 1980, we used data from 1986 -2015 to include data from both series and because the period captured our period of interest in the early 2000s. Beginning in 2016, there Bay Study experienced a protracted period of missed samples by one or both gear types, so we did not use more recent data. Data for our analyses included counts by age class of Longfin Smelt by gear type at each station.
Presence or absence of Longfin Smelt was recorded for each gear type at each individual sampling station (3 -11 stations per region) for each age-class ( Figure 1). We recognize that gear retention likely negatively affected spring and summer detection of age-0 Longfin Smelt (c.f., Mitchell et al. 2017Mitchell et al. , 2019 and that mesh size and additional debris and fish collected by the Bay Study OT likely improved summer detection relative to the MWT (Mitchell et al. 2017); however, no Longfin Smelt retention information is available for the MWT nor is any retention information available for Bay Study OT, so this analysis ignores selectivity effects in these analyses. Age class records were labeled by cohort year (the year when fish hatched and were detected as age-0 fish) to facilitate consistent predictions across calendar years as fish change ages (i.e., to prevent jumps in predicted probability of catch from December to January). This allowed the model to account for different patterns in movement for different life stages. We assign January as month 1 of life for all cohorts. Months were also labeled consecutively through the life cycle to reflect potentially changing use patterns of the age classes through time (age-0: months 1-12; age-1: 13-24, age-2+: 25-36). Data (presence/absence in the MWT and OT) from individual stations were used in the models, rather than aggregated values for regions. This provides some variability in presence/absence within each month so that a seasonal trend can be evaluated for each region.
Although we did not use FMWT Longfin Smelt data directly in our modeling, we describe it briefly because of its importance in both describing long-term trends of various fishes, and also because of its use in recent analyses to better describe and attribute causes to those trends. The FMWT samples stations ranging from western San Pablo Bay upstream through Suisun Bay and part of Suisun Marsh and throughout the Sacramento-San Joaquin Delta monthly from September through December (Rosenfield and Baxter 2007). Each month the crew conducts a single MWT tow (using approximately the same gear and methods as Bay Study MWT sampling) at each of the 100 stations used to calculate abundance indices ( Figure 1); additional stations have been added to the sampling panel over time (n=22), but are not used in current abundance index calculation. The 100 index stations were grouped by proximity into 14 regions used in index calculation. Each region has an associated regional volume weighting factor. The annual abundance index is calculated the sum of September through December monthly indices; monthly indices are calculated as the sum of the products of mean regional catch per tow and the regional weighting factor. Because they are consistently obtained, the FMWT index values are assumed to provide a relative approximation of population trends over time.

Modeling distribution
The probability of presence of Longfin Smelt in each gear type and region was modeled using a Generalized Additive Model (GAM) with a logit link. A GAM for probability of presence is similar to a logistic regression model, but with additional flexibility. The probability of presence was fitted as a function of a seasonal trend (month), long-term trend (year) and the change in seasonal trends over time (interaction of month and year) and each of these three terms was allowed to vary by region. Additional parametric intercept terms were added to the model to account for differences between the gear types in each region (interaction of gear and region and the associated main effects).
All analyses were conducted using R (version 3.3.2, "Sincere Pumpkin Patch;" R Core Team 2016). GAMs were fitted using the gam function in package mgcv (version 1.8-16; Wood 2011). Smooths were fitted using thin-plate regression splines (Wood 2003), and tensor products (Wood 2006a). Smooths for the long-term trend and seasonal trends used thin-plate regression splines. Tensor product interactions (ti) were used for the change in seasonal trends over time, rather than a full tensor product smooth (te), because we included individual smooths for year and month separately in the model. The documentation for "smooth.terms" in package mgcv has a thorough description of the types of smoothing terms that can be chosen and rationale for using each.
The number of knots (k) for the smooth terms in these models needed to be carefully limited to avoid over-fitting. Fewster et al. (2000) contains an excellent discussion of considerations for choosing the number of knots for smooth terms in GAMs. The models we fit intentionally ignore environmental covariates that might affect Longfin Smelt distributions within a season to examine general patterns in timing and distributions. It could be tempting to increase the number of knots to improve model fit, but the improved fit would artificially attribute the patterns caused by environmental factors to the temporal trends. Future work will include the effects of environmental covariates on abundance of Longfin Smelt. Thus, these models were not intended to explain all the variability in timing. Following advice in Fewster et al. (2000) for how to choose an appropriate level of smoothing, we plotted the smooth terms for month and year class over time, shared across gear types, at various values of k. We found that k values approximately equal to 0.3 times the number of x values were sufficient to produce a time series that described the major features of the data without excess noise.

Functional Data Analysis
Functional Data Analysis (FDA) comprises a suite of techniques for extracting information about values that vary smoothly over a continuum. These functions could take the form of a curve that varies along one explanatory variable or a surface that varies with the combined effects of two explanatory variables. The predicted probabilities that are generated by the GAM are a type of functional data. The predicted probabilities for each region vary on two temporal dimensions and we can generate smooth curves to show how the probability of occupancy varies on a seasonal scale or the long-term scale across years. We used the function fdata (package fda.usc; Febrero-Bande and de la Fuente 2012) to convert the predicted probabilities into functional data objects that could be analyzed under the FDA framework. Throughout this paper we will use the terms "regions" and "bays" interchangeably to refer to specific embayments within the SFE.
The characteristics of the shape of the seasonal curve provide information about when Longfin Smelt are most likely to occupy a particular region. Areas with higher probability indicate times during the life cycle when Longfin Smelt use a particular region. The width of the peak indicates how long Longfin Smelt reside in that region. The general pattern of Longfin Smelt occupancy is of interest for quantitatively depicting the conceptual model and substantiating how the species uses SFE. Another interest is when and where variation occurs in the typical pattern that could indicate changes in occupancy patterns, such as a shift from a historical pattern to a modern pattern.
The probability of detecting Longfin Smelt in a region at a particular month varies from year to year. This variability consists of two components: amplitude variation and phase variation. Amplitude variation is the variability in the height of the curves, or the magnitude of the probability. Phase variation is the lateral displacements or deformation of curves that are not explained by variability in amplitude. In the case of predicted probabilities of Longfin Smelt presence, phase variation indicates a shift in occupancy to earlier or later in the life cycle. Such changes can influence measures of abundance if species are no longer occupying historical habitat within the monitoring survey's sampling range for the same number of months each year.
To determine a typical pattern for reference, the curves for predicted probabilities of individual yearclasses need to be aligned and a functional average calculated. Two curves are perfectly aligned if they only vary in amplitude. To align two curves, a "warping function," is used that transforms the curves to make them as similar as possible. Aligning curves decouples phase and amplitude variability without losing information because the warping function contains the phase variability, and the amplitude variability is the variability that remains (Sangalli et al. 2010). The aligned functions are then averaged to produce a curve that describes the expected pattern of occupancy for each region. The variability provides information about how consistent the expected seasonal pattern is from year to year and values of the warping functions at each month provide information about whether there has been a directional shift over time. We used the function time_warping (package fdasrvf; Tucker 2017) to align the seasonal curves for each region and produce a functional average. The time_warping function aligns functions using the elastic square-root slope framework. Methods for comparing individual years to the mean functions are available in the code appendix. We provide the relative values of amplitude and phase variance by region and gear type to show the major sources of variation.
To identify timing shifts in regional occupancy, we first plotted phase variance by month through the life cycle for each gear type and identified the month of peak variance, and noted variance of proximal months. We then plotted values from the warping functions for these peak months separately by gear type to look for consistent, directional shifts over time. We chose to only plot functions from months with the peak variability because although other months have similar patterns, the directional changes are most pronounced where variability is highest.

Environmental variable plots
Monthly patterns in temperature and outflow reflect the Mediterranean climate of the SFE. Water temperatures are coolest in January and reach a peak during July, August, and September ( Figure 5). Summer temperatures commonly exceed 20°C in the upper estuary. Outflow is lowest from July to November and higher but more variable during December through May ( Figure 5).

Smoothing (GAMs)
We found significant parametric effects of bay and gear type as well as their interaction ( Table 2). The model with the optimal levels of smoothing to meet the goals of the study had 10 degrees of freedom for the seasonal trend and 8 degrees of freedom for the long-term trend (k = 11, 9, respectively; Table 1, Figure 3 and Figure 4). These values allow the relationship between probability of catching Longfin Smelt and the two trends in time to capture the major features of the relationships without overfitting.
Smooths for the long-term trend, seasonal trend, and their interaction were significant for each region ( Table 2). The equation for the optimal model is where each terms represent a distinct set of smooth functions, with k knots, and where terms represent parametric intercepts for two gears (i) and five regions (j).
Smooths for the long-term trends indicate that the probability of catching Longfin Smelt fluctuates over the long-term, but the trend is overwhelmingly decreasing for all regions since the late 1990s ( Figure 4; graphs for other bays are in appendix B). The significant smooth for the interaction of the long-term and seasonal trend indicates that the seasonal pattern of occupancy changes over time. The most pronounced changes in the seasonal pattern of occupancy over time were in Central Bay, Suisun Bay, and the Confluence.
The optimal model explained slightly less than one third of the variability in the presence/absence data (deviance explained = 26.8%). The model achieved full convergence at 11 iterations. There is some evidence of under-prediction: there are low probability values where Longfin Smelt were caught. This indicates that the model makes conservative predictions.
The effect of gear type on the probability of catching Longfin Smelt varied by bay. The otter trawl was more likely to catch Longfin Smelt in South Bay and Central Bay, whereas the midwater trawl was more likely to catch Longfin Smelt in Suisun Bay and the Confluence, and catching Longfin Smelt was equally probable in San Pablo Bay ( Figure 6). Information from both gear types is necessary to describe the distribution and timing of Longfin Smelt in the SFE because of regional differences in catch probabilities (i.e., detection effectiveness of the gear types). The predicted distributions for each gear type are somewhat complementary.

Functional Data Analysis
General features of the mean probability curves -Mean functions for the predicted probability of Longfin Smelt occupancy depict when Longfin Smelt are most likely to be observed in each region over the course of the life cycle ( Figure 6). We plotted separate regional functions for each of the gear types because the model results indicated significant differences between the gears ( Figure 6). The patterns of presence for each of the regions fit the life history narrative represented by the conceptual model. Although there are some notable differences, there are also some patterns that occur in both the MWT and OT (Figure 6). For each gear type, the shapes of the curves are similar within regions, but the difference appears to be that the magnitude of the peaks differs by gear type. Thus, the general timing of occupancy is confirmed by both gear types, particularly broad regional presence during winters and relative absence during late summer (i.e., August, month 20), but magnitude (abundance) near bottom (OT) or within the water column (MWT) varies by region and season ( Figure 6).
For both gears, probability of presence peaked in winter at the end of the first year of life (Dec or Jan; months 12-13) in all bays except Central Bay, which peaked in summer (July, month 7) in both gears ( Figure 6). An earlier, lesser peak also occurred in San Pablo Bay in summer (July, month 7) for both gears; though of different magnitude, this peak closely matched that of Central Bay timing-wise. These early lesser peaks, particularly in the OT, indicate that age-0 Longfin Smelt are moving downstream in substantial numbers prior to reaching the 40 mm size of counting (age-0 Longfin Smelt begin reaching 40 mm in Apr and have mostly achieved 40 mm by August; Baxter 1999). Other, lesser peaks occur the following winter (December-February, months 24-26) when mature fish move up estuary to spawn ( Figure 6). For every bay, these later peaks are smaller (lower amplitude) than the first peak, reflecting a decline in year-class abundance between the end of their first and second years of life, and possibly an incomplete return of newly age-2 fish (i.e., the small fraction of fish destined to spawn at age 3).
There are three major nadirs in probability of observing Longfin Smelt. The probability of observing Longfin Smelt is essentially zero for months 1-3. This is an artifact of the dataset: the fish are too small to be caught reliably in the nets and no fish < 40 mm were included in the dataset. This is consistent with other studies because typically no fish < 40 mm are counted from these gear types. There is a second nadir in summer of their second year, around months 16-17 to 21-23 (Apr-May to Sep-Nov of age-1). This occurs when the conceptual model predicts that these age-1 fish leave the upper estuary (West Delta and Suisun Bay, and South Bay) to move downstream and at least some go out to the coastal ocean (Rosenfield and Baxter 2007;CDFG 2009). The third nadir is in summer of their third year, around months 27 -34 ( Figure 6). This corresponds to the timing when spent adults typically die postspawning or leave the SFE. San Pablo Bay and the Confluence exhibit a small increase again at months 35 & 36. These likely reflect delayed-maturity of soon-to-be age-3 fish returning to spawn. Earlier in the life cycle, San Pablo Bay exhibits a minor decline in probability centered on month 9, while probability in Central Bay begins a prolonged decline to the following summer ( Figure 6). This coincides with the conceptual model's narrative that Longfin Smelt move out of the Delta and downstream as water temperature increases in summer and early fall, prior to re-occupying the entire SFE beginning in late fall.
Agreement with the conceptual model -The conceptual model predicts that some Longfin Smelt are present initially from the West Delta downstream into Central Bay beginning in month 4, and that presence generally declines through summer in these upstream regions. Our models detect Longfin Smelt presence beginning in April and across all regions other than South Bay. The seasonal decline in upstream regions described in the conceptual model is mostly obscured in our predictive models by effects of initially poor but improving gear retention. Most larvae and small juveniles are too small in winter and spring to be retained or counted; some age-0 Longfin Smelt reach 40 mm in Apr, and most surpass 40 mm by August (month 8; Baxter 1999). Later, increased presence in all regions in late fall and winter obscures any summer-fall declines, with the exception of a shallow dip in San Pablo Bay. On average, the peaks for Longfin Smelt occupancy are around months 11 and 25 for the Confluence ( Figure  6), similar to the conceptual model ( Figure 2). For Suisun Bay, the mean functions have strong peaks around month 13, but only a small increase around month 5 when the conceptual model indicates they should be present. For age-1 (13-23 months) and age-2 fish (≥month 24), the mean function is consistent with the conceptual model. Finally, the mean functions for the Confluence indicate that age-0 Longfin Smelt are present months before and age-1 Longfin Smelt are present slightly before the December abundance peak the conceptual model predicts. For San Pablo Bay, the mean functions for predicted probability indicate that Longfin Smelt are most likely to be caught around months 6, 13 and 25. These peaks are consistent with the conceptual model. The shallow dip centered at month 9 for San Pablo Bay is also consistent with the conceptual model prediction that Longfin Smelt are relatively absent from San Pablo Bay during months 8 -10. The narrowness of this window of absences is a likely cause of the shallowness of the dip. The mean function for San Pablo Bay also indicates that Longfin Smelt are not present during the age-1 summer (months 17 -22), as predicted by the conceptual model. The conceptual model predicts that some Longfin Smelt occur in Central SF Bay from month 4 -25, except for a period in months 18 -20 ( Figure 2). This pattern is somewhat reflected in the mean function for the MWT and OT, although the nadir is brief for most regions and centered at month 20 or 21 for each ( Figure 6).
Variability -Variation from the mean predicted probability function can provide insights into changes in the distribution of Longfin Smelt in the SFE. The considerable variability in timing of occupancy for the age-0 Longfin Smelt in summer and fall in both Suisun Bay and the Confluence (Figure 7) suggest that the observed pattern in any year may not follow the average predicted pattern. The amplitude variance is at least an order of magnitude greater than the phase shift variance (Table 3). The fact that phase variability is small suggests that the seasonal patterns of occupancy are fairly consistent through time; nonetheless, phase variance is distinctly higher in Suisun Bay and the West Delta relative to other regions, particularly for the MWT (Table 3). The relatively higher amplitude variability suggests that annual abundance likely controls how strongly that signal appears in the Bay Study data from year to year.
The pattern and height of regional peaks across years of study roughly corresponds to changes in population indices from monitoring programs: we are more likely to detect fish when more fish are present ( Figure 4 and Figure 8). This is a well-known issue with modeling species occupancy. There was also some variability in the magnitude of peak occupancy from year to year. The variance in the values of the warping functions shows the months when variability in the predicted probability functions cannot be accounted for by changes in amplitude alone (Figure 7). These represent phase shifts presence and indicate that fish are changing their use of a particular region in a way that is different from the overall pattern observed from year to year. For both gear types, the predicted probability curves for South SF Bay were perfectly aligned for all years because there was no significant cohort*month interaction in the GAM. Figure 7 illustrates this lack of variability as a variance of 0 for all months for South SF Bay. San Pablo Bay also had little phase variation, indicating that the pattern of occupancy there is also mostly driven by fluctuations in population size. Central SF Bay, Suisun Bay, and the West Delta all had clear peaks in phase variation relatively early in the life cycle. For the West Delta and Suisun Bay, the peak in phase variation was in the summer (Jun-Jul, months 6-7) and fall (Aug, month 8) of the age-0 year, respectively. For Central SF Bay, the peak in phase variation was in winter around months 12 and 13, and extended a couple months later.
The timing and direction of phase shifts was investigated for peak variance months (Figure 7), and the temporal shift patterns were plotted relative to those peak variance months for each gear and region ( Figure 9). The West Delta shows the highest variance and shifts in timing of age-0 summer (month 7 for MWT, month 7 for OT) occupancy (Figure 7). The shift in timing gradually changes from two months earlier than average to two and a half months later than average over the course of the time series. The shift from around average timing to later happened rapidly from about 2005 to 2015 (MWT) and 2008 to 2014 (OT; Figure 9). The highest variance in occupancy and thus shifts in timing of occupancy in Suisun Bay occurred during the age-0 summer (month 8 for MWT, month 8 for OT) varies from two months earlier to two months later than the mean function ( Figure 9). Major shifts in timing occurred in the midto late-1990s, when the direction changed from later than average to earlier than average, and in the mid-2000s, when the pattern reversed. The late-early-late pattern coincides with the general pattern of low-high-low population indices (Figure 8).
For the Central Bay, variance peaks in the winter (Figure 7; month 13 for MWT, month 12 for OT) provide evidence of shifts in timing for age-0 Longfin Smelt. The amount of phase variation for Central Bay was much higher in the OT than for the MWT. Through the mid-1990s, Longfin Smelt MWT-based occupancy of Central Bay shifted from later than average to earlier by 1996 and this shift to earlier timing was maintained through the end of the study period (Figure 7). Similarly, OT based occupancy of Central Bay in December shifted from 2.5 months later in 1987 to one month earlier by 2000, where it remained (Figure 7). These three examples represent the maximum variances identified among several months of increased variances for each region and gear; patterns of shifts for those increased-variance months not shown generally matched those that were shown.
Patterns in the monthly index values for the FMWT reflect shifts in timing of fall movements, as indicated by the functional data analysis of the SFBS data (cf. Figure 10 and Figure 9). Specifically, the September index values are low when the August pattern for Suisun Bay is later than average (Figure 9). During the period when the August pattern for Suisun Bay is earlier than average, the monthly FMWT index values appear to be higher in general, and the September index values tend not to be the lowest.

Discussion
The conceptual model of Longfin Smelt life history depicting the consistent features of occupancy ( Figure 2) generally agrees with the historical pattern of Longfin Smelt distribution; however, there is evidence that the timing of occupancy, which we interpret as movements, has changed in recent years. There was also evidence of a long-term decline in occupancy, as indicated by a decline in amplitude and in predicted peak probabilities of detection, which corresponds to declines in abundance reported by other researchers (Rosenfield and Baxter 2007;Sommer et al. 2007;Thomson et al. 2010). Some periods of the life history are more variable than others in terms of distribution and timing. These periods may be important for understanding when Longfin Smelt are reacting to a changing environment and for deciding when management actions may have the most impact on increasing Longfin Smelt abundance. Finally, such changes in distribution timing may also influence our interpretation of current abundance indices.

Conceptual model
Regional patterns of occupancy of Longfin Smelt in the SFE are predictable based on year and season alone: increases during the cool, often wet winters and springs, and decreases during warm, dry summers and fall. We interpret the seasonal increases and decreases in regional occupancy early in the life cycle as reflecting recruitment to the gear (i.e., fish growing to sufficient size to be retained) and then, later in the life cycle, as movements into and out of the respective regions as well as into and out of the estuary. The estuary-wide declines in probability of capture that occur during the summer (months 15-21, 27-31) cannot be attributed to mortality alone because probabilities increase again in the fall. The patterns we derived from Bay Study data continue to support the conclusion that Longfin Smelt in the SFE are anadromous and that beginning at age 0 but more commonly at age 1, they most likely spend the summer months outside the SFE in the coastal ocean and in cool marine salinities in Central Bay (Rosenfield and Baxter 2007).
Timing of movements within and into/out of the estuary may reflect large-scale changes in environmental factors (e.g., water temperature) as well as changes in habitat suitability with ontogeny (e.g., salinity tolerance increases) that also follow predictable seasonal patterns. Our model was formulated to address broad patterns in Longfin Smelt distributions and timing of movements with the knowledge that future analyses, including addition of environmental variables, will explore sources of this variability. Environmental factors, particularly variability in temperature and outflow, clearly influence patterns of occupancy, particularly early in life (Dege and Brown 2004;Jeffries et al. 2016;Grimaldo et al. 2017). The general pattern in the occupancy curves for age-0 Longfin Smelt is explained to a large extent by recruitment to the gear, which likely occurs fully in their first fall based on retention of Delta Smelt a similarly sized and shaped fish to Longfin Smelt (Mitchell et al. 2017) and by general temperature patterns, which decline fall into winter throughout the estuary (Orsi 1999b). Distribution varies most for larvae and small juveniles whose distribution location is dependent upon outflow in two ways. First, outflow determines the location of the boundary between fresh and brackish water within the upper estuary (Kimmerer and Monismith 1993;Kimmerer 2004), which is reflected in the changing distributions of mature fish (CDFG 2009) and perhaps spawn locations. Second, outflow determines influences upper estuary circulation patterns (Kimmerer 2004) and thus the distribution of larvae immediately after hatching (Baxter 1999;Dege and Brown 2004;Rosenfield & Baxter 2007). In high outflow years substantial numbers of larvae and small juveniles are distributed in San Pablo Bay; in low outflow years, most are found in Suisun Bay and upstream (Baxter 1999;Grimaldo et al. 2017). Also the earliest hatchlings tend to recruit the farthest downstream as reflected by the earliest increases in predicted probabilities in San Pablo Bay for both gears (and Central Bay for the OT) increasing prior to those of upstream regions ( Figure 6).

Temporal variation: long-term declines
The probability of catching Longfin Smelt has varied considerably through the time frame we analyzed. Most of the variability in individual years can be accounted for in amplitude variation resulting from changes in year-class abundance. Decreasing probability of capture in every region over the long-term reflects the declining abundance of Longfin Smelt in the SFE (Rosenfield and Baxter 2007;Sommer et al. 2007;Thomson et al. 2010;. Declines in Longfin Smelt have been attributed to declines and changes in outflow (Stevens and Miller 1983;Jassby et al. 1995;Kimmerer 2002;Thomson et al. 2010;Tamburello et al. 20018) and subsequent effects on juvenile recruitment and adult abundance (Nobriga and Rosenfield 2016), the invasion of the Amur River clam (also known as the Overbite Clam; Potamocorbula amurensis) in 1986 (Kimmerer 2002;Thomson et al. 2010) and the Pelagic Organism Decline, which began in the early 2000s (Sommer et al. 2007; Thomson et al. 2010). The recent prolonged drought (2012)(2013)(2014)(2015) also contributed to the decline of Longfin Smelt abundance.
Similarly, natural declines in abundance with age resulting in a lower seasonal presence of age-1 Longfin Smelt in the estuary in their second fall and winter ( Figure 6, months 22-28) contributed to a lack of substantial phase variability during that period (Figure 7) even though a greater proportion of age-1 Longfin Smelt (mature fish) appear to be present in the narrow time frame from late December through mid-January than in the past (Weeks 24-25; Baxter unpublished data).

Temporal variation: seasonal shifts
In addition to declines in abundance and probability of capture, there have also been shifts in the timing (phase variation) of regional occupancy for Longfin Smelt. The early part of the life cycle is the time when most of the variability in timing occurs. Variance due to phase shifts is low compared to variance due to amplitude shifts, which indicates that the seasonal pattern is strong and that shifts in timing are relatively lesser adjustments to predicted probabilities, when compared with the larger effect of abundance. Nevertheless, these phase shifts reveal important changes in Longfin Smelt behavior patterns. Two such shifts in timing are: one, that age-0 Longfin Smelt now appear to leave the upper estuary sooner in the summer (inferred from decreasing presence in other sampling gear: https://wildlife.ca.gov/Conservation/Delta/20mm-Survey) and two, that they return later in the fall (observed; Figure 9). The net effect of these shifts in timing is that Longfin Smelt spend more time in the lower estuary and possibly the coastal ocean than they did historically. The summer-fall pattern of occupancy for the age-0 in the West Delta shifted steadily later from the mid-1980s to 2015. Seasonally increasing age-0 Longfin Smelt recruitment through summer and fall, as well as increasing occupancy suggest that Longfin Smelt return from the lower estuary and ocean to re-occupy regions of the SFE as they exhibit more favorable conditions.
For Suisun Bay, phase shifts indicate that summer and fall age-0 occupancy patterns oscillated from later than average (briefly) to earlier and back to later (Figure 9). The earlier timings of summer and fall patterns in the mid-1990s to mid-2000s ( Figure 9) match up with a period of relatively high predicted probabilities ( Figure 6) and abundance values (Figure 8) associated with a period of numerous relatively high outflow years. Although strictly correlative, this suggests a consistent pattern: high outflow winters and springs are associated with higher abundance and enhanced summer and fall occupancy of Suisun Bay. Wet winters and springs may lead to weather that creates habitat conditions in summer and fall that result in increased occupancy for Suisun Bay earlier in the year than occurs on average. During the summer and fall, Delta outflows are typically low for the SFE and it is common for temperatures to approach and exceed thermal maxima for Longfin Smelt (Jeffries et al. 2016), especially in the upper estuary and South San Francisco Bay. This suggests that high temperature, rather than flow, may be the main factor that limits distribution of Longfin Smelt during the summer and fall. By summer Longfin Smelt are no longer limited to brackish water; they can rear throughout the estuary, but high temperatures appear to limit their distribution.
Late-year age-0 Longfin Smelt (and new year age-1 Longfin Smelt) occupancy in Central Bay consistently shifted earlier, indicating decreased use of this region in the winter period beginning at the start of our analysis period and culminating in 1996 in the MWT, but not until about 2000 in the OT (Figure 9). During these later years (>2000), occupancy shifted from relatively high values to a declining trajectory in those months-in the early 1990s to those same values and trajectory in November since the late 1990s or 2000, depending on gear (Figure 9). The timing in peak variation for age-0 in Central Bay was in December (OT) and age-1 in January (MWT), but the period of relatively high variability for both gears encompasses the age-0 fall and age-1 winter and spring, which is normally the end of the period when Longfin Smelt have returned to the estuary from the ocean and are well-distributed throughout the estuary having moved inland. The patterns normally observed during the winter in Central Bay have shifted three to four months earlier than they were in the late 1980s. Since 2000, the OT occupancy pattern beginning in November should depict a steeper decline, suggesting that Longfin Smelt moving upstream are not using Central Bay as much during the winter as they had in the past. It may be that the high degree of variability in the timing of Longfin Smelt use of Central Bay reflects an overall variability in use between the marine habitats of Central Bay and the coastal ocean.

Implications for interpretation of recent abundance trends
The decline of several fish species in the SFE, including Longfin Smelt, has been well-documented (Sommer et al. 2007;Thomson et al. 2010), leading to increasing concern for the Longfin Smelt population as a whole . Our model, based on Bay Study data, revealed shifts in occupancy timing of Suisun Bay and the West Delta to later in fall months that could result in a reduction in apparent FMWT Longfin Smelt abundance beyond what would be expected by an already declining population. The spatial extent of the FMWT covers the San Pablo Bay, Suisun Bay and West Delta regions of the Bay Study. The annual FMWT index is calculated as the sum of four monthly indices (September-December), each calculated as the sum of regional mean catch-per-tow weighted by a regional volume (Rosenfield and Baxter 2007). Delayed occupancy of these regions by one month or more in fall would be reflected as a decline in apparent abundance in the FMWT because of the way the index is calculated. However, if this change in timing of fall distribution strongly affects the annual index, it would be evident as different relative magnitudes in the monthly indices (i.e., relative magnitude of Sep and Oct abundance should be lower). This appears to be the case (Figure 8). Thus, we conclude that in recent years, age-0 Longfin Smelt partially resided downstream of the FMWT sampling regions until later in the sampling period. Returning later in the season reduces opportunities for capture, reducing apparent abundance. It is not clear whether presence later in the sampling period compensated for the absence early in the sampling period, but patterns of presence indicate that age-0 Longfin Smelt returning in fall and early winter remain in upstream regions throughout the winter and early spring, and similarly, maturing age-1 fish remain until spawning is complete providing repeated opportunities for capture each age class.
A previous study using data from the Bay Study and the Suisun Marsh Study found coincident declines in abundance that correlated with those of the FMWT index and concluded the FMWT decline was not caused by a shift in temporal or spatial distribution patterns (Rosenfield and Baxter 2007). However, these authors did not examine regional distribution in detail and the shifts in timing and distribution may only become apparent with the longer (more recent) time series that we used in our study. The previous study used data through 2004, and our results agree that it would be difficult to determine a temporal shift if the dataset were truncated to that year, though declines in Suisun Bay and the West Delta were well underway by 2004 (Figure 9). It is also likely that our regional analysis methods made the shift in timing more apparent. We used a functional data analysis framework to separate the variability introduced by abundance from the variability introduced by shifts in timing. Since abundance introduces more variability in this dataset, it could be that this source of variation masked the effect of timing shifts in the previous study.
The potential for shifts in occupancy timing to affect monitoring survey results, especially those that monitor a specific period of months each year, highlights the value of using more than one study and gear types to monitor each species life stage throughout the Estuary. In this case, different surveys act as both a check on and a complement to each other when we interpret trends in fish abundance and distribution. The FMWT does not cover the entire range for older life stages of Longfin Smelt, but as long as the species uses the estuary consistently from year to year, the relative abundance index calculated should still track consistent relative trends in abundance change, because the relative abundance index will also maintain its relative relationship to the true population. Conversely, if a species shifts its distribution relative to survey sampling range over time (or changes timing of regional occupation), a greater or lesser proportion of the population will be available for capture and will change the relationship between the sampled population (and the abundance index calculated) and the true population size. Data from the Bay Study presented here indicate such a distribution shift and provide additional context for interpreting the FMWT index for Longfin Smelt and recent analyses of abundance trends (i.e., Thomson et al 2010); however, this analysis cannot estimate the magnitude of the change in abundance. Regional occupancy patterns of the different gears used by the Bay Study also illustrate this point. The patterns in Longfin Smelt occupancy that are evident in the Bay Study OT and MWT are somewhat complementary in that both gear types are necessary to illustrate the full picture of Longfin Smelt life history patterns and use of the SFE.

Management Implications
Management practices have long accounted for the life history and generalized distribution of species, but information about the main effects and flexibilities in timing can also inform effective management practices. Increased winter-spring freshwater outflow has been shown to increase production in Longfin Smelt (Stevens and Miller 1983;Jassby et al. 1995;Rosenfield and Baxter 2007;Thomson et al. 2010;Nobriga and Rosenfield 2016). Conditions, particularly outflow, early in the first year of life are the most important for setting the trajectory of the cohort (e.g., Nobriga and Rosenfield 2016). In addition to effect on abundance, the distribution of larval Longfin Smelt is also strongly affected by flow (Baxter 1999;Dege and Brown 2004). Current fish monitoring surveys that target larval and small juvenile Longfin Smelt do not sample the full geographic range of Longfin Smelt larvae (https://wildlife.ca.gov/Conservation/Delta/Smelt-Larva-Survey; https://wildlife.ca.gov/Conservation/Delta/20mm-Survey), particularly during years of persistent high outflow when larvae and small juveniles are known to strongly inhabit all of San Pablo Bay, occasionally Central and South San Francisco Bays (Baxter 1999;Grimaldo et al. 2020;Lewis et al. 2020). Larval and young juvenile Longfin Smelt are too small to be reliably sampled by the Bay Study OT and MWT gears until the summer, but their initial summer distribution likely reflects in large part where they were conveyed as larvae during the peak flow season (e.g., Baxter 1999; Dege and Brown 2004;Grimaldo et al. 2020). Our results also indicate that even for a pelagic species such as Longfin Smelt, using data from two different gear types provides a fuller picture of spatial distribution. Using the Bay Study midwater trawl gear alone, our model would have underestimated the use of the Central Bay by young Longfin Smelt.
After Longfin Smelt reach the second half of their second year, when most fish approach spawning age (months 20-23), their seasonal movement patterns vary little from year to year ( Figure 7) and abundance in the SFE largely determines the ability of monitoring programs to detect them within their regional habitats. If Longfin Smelt were present in the estuary during this period, managers could use this predictability to target management actions to minimize mortality and improve growth; however, like Rosenfield and Baxter (2007), we interpret low detections in late-summer and fall as evidence that maturing Longfin Smelt inhabit the coastal ocean. For example, managing conditions in Suisun Bay to maximize food availability or reduce temperature in the fall may attract adults to potentially stressful temperatures or benefit adult survival through spawning (CNRA 2016); some of this work is being attempted and investigated for another declining species, Delta Smelt (https://resources.ca.gov/Initiatives/Delta-Smelt-Resiliency-Strategy). Also, discontinuing dredging and sand mining when and where Longfin Smelt may stage and spawn could improve spawning success; this has been recognized as important (USACE, USEPA, BCDC, SFBRWQCB 2001). However, the reduced abundance of older fish relative to younger life stages and resulting reduced occupancy levels makes detecting shifts in timing difficult.

Future Research
Our study aimed to describe general patterns in occupancy and as such, it does not include environmental variables such as salinity, temperature, and food availability that would likely better inform our modeling (i.e., explain more variance). Including this additional information could help inform predictive models of Longfin Smelt distribution, given potential alternative management actions. As described previously, the larval and early juvenile distributions of Longfin Smelt are influenced by freshwater outflow and water temperature (e.g., Dege and Brown 2004;Jeffries et al. 2016), so there is no reason to believe that distributions of older Longfin Smelt would not be influenced by the same and additional environmental factors. Future research could build on our results by using a functional data analysis framework to investigate the effects of various environmental factors on the geographic and depth distribution and timing of Longfin Smelt occupancy.