Space‐for‐time is not necessarily a substitution when monitoring the distribution of pelagic fishes in the San Francisco Bay‐Delta

Abstract Occupancy models are often used to analyze long‐term monitoring data to better understand how and why species redistribute across dynamic landscapes while accounting for incomplete capture. However, this approach requires replicate detection/non‐detection data at a sample unit and many long‐term monitoring programs lack temporal replicate surveys. In such cases, it has been suggested that surveying subunits within a larger sample unit may be an efficient substitution (i.e., space‐for‐time substitution). Still, the efficacy of fitting occupancy models using a space‐for‐time substitution has not been fully explored and is likely context dependent. Herein, we fit occupancy models to Delta Smelt (Hypomesus transpacificus) and Longfin Smelt (Spirinchus thaleichthys) catch data collected by two different monitoring programs that use the same sampling gear in the San Francisco Bay‐Delta, USA. We demonstrate how our inferences concerning the distribution of these species changes when using a space‐for‐time substitution. Specifically, we found the probability that a sample unit was occupied was much greater when using a space‐for‐time substitution, presumably due to the change in the spatial scale of our inferences. Furthermore, we observed that as the spatial scale of our inferences increased, our ability to detect environmental effects on system dynamics was obscured, which we suspect is related to the tradeoffs associated with spatial grain and extent. Overall, our findings highlight the importance of considering how the unique characteristics of monitoring programs influences inferences, which has broad implications for how to appropriately leverage existing long‐term monitoring data to understand the distribution of species.

tion/non-detection data at a sample unit and many long-term monitoring programs lack temporal replicate surveys. In such cases, it has been suggested that surveying subunits within a larger sample unit may be an efficient substitution (i.e., space-fortime substitution). Still, the efficacy of fitting occupancy models using a space-fortime substitution has not been fully explored and is likely context dependent. Herein, we fit occupancy models to Delta Smelt (Hypomesus transpacificus) and Longfin Smelt (Spirinchus thaleichthys) catch data collected by two different monitoring programs that use the same sampling gear in the San Francisco Bay-Delta, USA. We demonstrate how our inferences concerning the distribution of these species changes when using a space-for-time substitution. Specifically, we found the probability that a sample unit was occupied was much greater when using a space-for-time substitution, presumably due to the change in the spatial scale of our inferences. Furthermore, we observed that as the spatial scale of our inferences increased, our ability to detect environmental effects on system dynamics was obscured, which we suspect is related to the tradeoffs associated with spatial grain and extent. Overall, our findings highlight the importance of considering how the unique characteristics of monitoring programs influences inferences, which has broad implications for how to appropriately leverage existing long-term monitoring data to understand the distribution of species.

K E Y W O R D S
Delta Smelt, incomplete capture, Longfin Smelt, occupancy model, space-for-time substitution, species distribution model

| INTRODUC TI ON
Natural resource managers invest considerable resources to support long-term monitoring programs to maintain a continuous flow of data on biological systems of interest. Although such data are collected for a variety of purposes, information concerning how and why species redistribute across dynamic landscapes is usually of high interest. However, it is often the case that monitoring data are compromised to some extent by the sampling or observation process. For example, it is widely recognized that our ability to monitor the distribution and abundance of fish and wildlife is complicated by incomplete capture (Bayley & Peterson, 2001;Conroy et al., 2008;Haynes et al., 2013;Pine et al., 2003;Rowe et al., 2019). That is, variation in catch data can be related to both ecological (i.e., true variation in distribution and abundance) and observational (i.e., variation in capture efficiency) processes. Failure to account for incomplete capture complicates our ability to reliably interpret patterns contained in monitoring data because it requires the likely violated assumptions that capture probabilities are constant (or only vary randomly across space and time) and that the magnitude of the biological signal in trends is greater than the noise in the monitoring data (Anderson, 2001;Thompson, 2002).
Occupancy models are probably the most widely used statistical approach to quantify the relationship between species occurrence and environmental conditions while accounting for incomplete capture (MacKenzie et al., 2017). Like all estimators, occupancy models contain a few assumptions: the occupancy state (typically species presence or absence) cannot change at a sample unit across replicate surveys; detection probabilities must be independent across sample units and replicate surveys; heterogeneity in occupancy and detection probability can be explained using covariates; and false-positive detections do not occur. Occupancy models have proven to be an extremely flexible tool as extensions to the original parameterization continue to be rapidly developed to better suite a variety of ecological problems (reviewed in MacKenzie et al., 2017). Still, occupancy models require replicate detection/non-detection data at a sample unit to be fitted (but see Henry et al., 2020). Such data are typically collected by surveying a focal sample unit multiple times, where replicate surveys are far enough apart in time that detection probabilities are independent and close enough in time that the occupancy state at the sample unit is static. However, long-term monitoring programs that were initiated prior to the widespread use of occupancy models typically lack replicate surveys at a sample unit due to the fact that they often come at the cost of not being able to survey a greater number of sample units. In such cases, it has been suggested that managers can potentially treat surveys of multiple subunits within a larger sample unit as replicate surveys when fitting occupancy models (Guillera-Arroita, 2011; Kendall & White, 2009). This space-for-time substitution is attractive because managers can potentially account for incomplete capture without reducing the spatial coverage of monitoring programs, and it potentially allows managers to leverage existing monitoring data that lack temporal replicate surveys at a sample unit in their sampling design. Still, fitting occupancy models to data that implement a space-for-time substitution has not been fully explored (Kéry & Royle, 2016).
The Sacramento-San Joaquin Delta and San Francisco Bay (Bay-Delta) system is the largest estuary on the Pacific Coast.
Management decision making in the Bay-Delta is complex as managers attempt to restore and conserve at-risk fish species while trying to deliver freshwater to meet the needs of the public (Brown et al., 2009;Hanemann & Dyckman, 2009;Moyle et al., 2018). Long-term fish monitoring programs using trawl surveys were initiated as early as the late 1950s to support decision making for various fishes in the Bay-Delta. Over the years, additional monitoring programs that rely on trawl surveys have been initiated by different management agencies to monitor fishes for different objectives. Although fish catch data from these surveys have long been used to inform policy and management decision making, it was not until relatively recent that the effect of the observation process on inferences regarding the distribution and relative abundance of fishes in the Bay-Delta was considered (Goertler et al., 2020;Latour 2016;Mahardja et al., 2017Mahardja et al., , 2021Mitchell et al., 2017Mitchell et al., , 2019Newman 2008;Peterson & Barajas, 2018;Polansky et al., 2019;Thomson et al., 2010). Overall, this collection of work has demonstrated that fish capture efficiency varies markedly across species, time, location, and monitoring program within the Bay-Delta.
Given the above findings, it is clear that the data collected by fish monitoring programs in the Bay-Delta need to be corrected for incomplete capture to properly inform decision making. Like many monitoring programs across the world (Budy et al., 2015;Duarte et al., 2021;Kéry et al., 2010;Sadoti et al., 2013;Van Strien et al., 2013;Whitlock et al., 2020; and many more), occupancy models are increasingly being used to accomplish this for fishes in the Bay-Delta despite sampling designs not explicitly considering occupancy models when many of these monitoring programs were established (Goertler et al., 2020;Mahardja et al., 2017Mahardja et al., , 2021Peterson & Barajas, 2018). However, one of the long-term monitoring programs used to monitor the distribution and relative abundance of Delta Smelt (Hypomesus transpacificus) in the Bay-Delta, the Spring Kodiak Trawl (SKT), lacks temporal replicate surveys .
In 2016, the Enhanced Delta Smelt Monitoring Program (EDSM) was initiated to provide accurate and higher-resolution spatial and temporal data on the distribution and relative abundance of Delta Smelt in the Bay-Delta . It includes temporal replicate surveys at sample units using the same sampling gear as SKT. Our objective was to use these real-world fish monitoring data that deploy the same sampling gear to demonstrate some of the tradeoffs associated with using spatial replicates and occupancy models to monitor variation in the distribution of species. To accomplish this, we fit occupancy models to SKT and EDSM data using both temporal and spatial replicate surveys to estimate factors related to the distribution of the two native pelagic fish species in the Bay-Delta, Delta Smelt and Longfin Smelt (Spirinchus thaleichthys), which are both species of concern and often the focus of the different monitoring programs in the Bay-Delta. We then compared the estimates among the different approaches and datasets.

| Study area
The Bay-Delta system encompasses the area that is bounded by the confluence of the Sacramento and San Joaquin Rivers to the south and east, the San Pablo Bay to the West, and Suisun Bay to the north in California, USA (Figure 1). This region supports more than 500 fish, wildlife and plant species, including several threatened and endangered species. Although the Bay-Delta was once comprised of a series of marsh-wetland complexes in the eastern portion at the confluence of the Sacramento and San Joaquin Rivers, levee construction along stream channels and the construction of islands began in the mid-19th century to, in part, control seasonal flooding and to protect water exports from saltwater intrusion to better meet the increasing demands for freshwater across the State of California (Galloway et al., 1999). Over the years, system-wide changes to the environment have been documented, including a rise in contaminants (Fong et al., 2016), invasion by exotic species (Linville et al., 2002;Underwood et al., 2006), an increase in water clarity, a shift in the distance from the Golden Gate Bridge to the point where the salinity on the bottom is about 2 ppt (locally referred to as X2 and often used as a measure of water quality), and an increase in specific conductance and temperature (Peterson & Barajas, 2018).

| Species background
Delta Smelt are endemic to the Bay-Delta and listed as threatened by the U.S. Fish and Wildlife Service (U.S. Fish and Wildlife Service, 1993) and endangered by the State of California (California Department of Fish and Game, 2008). The species is believed to be pelagic, using primarily open water habitats (Sommer et al., 2007).

| Data collection
The EDSM was initiated in 2016 by the U.S. Fish  of EDSM is to provide accurate and high-resolution data on the distribution and relative abundance of Delta Smelt using two different gears: the Kodiak trawl and 20 mm survey. We did not analyze EDSM data collected by the 20 mm survey for this study because it uses a different sampling gear than the SKT. Surveys using the Kodiak trawl began in late 2016 by dividing the Bay-Delta into strata and using a generalized random-tessellation stratified design (Stevens & Olsen, 2004) with equal probability sampling to select sampling locations. The program began with four strata, but these strata were refined into eight strata in July 2017 and then again into 10 strata in December 2017. After inspection of the data, however, we noticed that nonrandom patterns in location selection seemed to occur prior to July 2017, which may be related to targeted surveys for other purposes (Shawn Acuña, Metropolitan Water District of Southern California, personal communication). Tows for the Kodiak trawl occurred in every month except April and May and used a net that has a 13.9 m 2 mouth opening and a mesh size that ranges from 5.08 cm at the mouth to 0.64 cm at the end. The EDSM implements replicate tows to characterize the presence, absence, and relative abundance of fish at each sampling location. Because of how replicate tows are implemented, a sample unit is the general area around each randomly generated sampling location. Importantly, prior to sampling, the EDSM protocol sets the maximum number of replicate tows that can be taken at a sample location depending on presumed Delta Smelt density and time period of the survey. However, EDSM uses a modified removal design where replicate samples are terminated or tow duration reduced based on the number of Delta Smelt that are captured in the first two tows at a sample location to minimize detrimental effects of sampling on Delta Smelt in the sample unit.
When and how often this occurred is not readily available in the metadata. Each tow was classified based on the amount of debris in the net and the time the net was in the water. We restricted our analyses to data collected from July 2017 through March 2019 and condition 1 tows (i.e., tows that were "good/normal," with no blockage of the net and ≥3 min long; Figure 1). This dataset had 2369 sample units with an average of 4.43 replicate tows (range: 1-9) at each sample unit. Before and after each tow, field crews collected environmental data. Specific conductivity (μs/cm 3 ), water turbidity (NTU), water temperature (°C), and dissolve oxygen (mg/L) at the surface were measured with calibrated meters. Secchi depth (m) and bottom depth (m) at the beginning of the tow also were measured and recorded. At the end of each tow, field crews recorded the tow duration (min), the tow number at the site that day, the volume of water sampled (m 3 ), the time of sampling, the start location for the tow, tow direction (upstream, downstream, or neither), and the number of fish caught at different sizes for each species. Although tide data were recorded, these data were lacking for 16.36% of the tows. Therefore, we downloaded tide height data from the National Oceanic and Atmospheric Administration (NOAA) station at Port Chicago, CA (NOAA. Accessed May 11, 2020, https://tides andcu rrents.noaa.gov/stati onhome.html?id=9415144), and classified tides as flood, ebb, and slack based on changes in tide height. Importantly, we used a threshold of 0.04 m per 15 min change to indicate a change in the tide direction to avoid high turnover in tide classifications due to minor changes in tide height. We used these tide data for all EDSM tows to maintain consistency in the source of the data.
We examined the environmental data for outliers that may have been caused by equipment malfunction and treated these values as missing data per the suggestion of the EDSM data curator (Catherine Johnston, U.S. Fish and Wildlife Service, personal communication).
To compare approaches (i.e., spatial vs. temporal replicates using occupancy models) with these data, our primary analysis of EDSM data used temporal replicate surveys on the same day at a sampling location as the replicate surveys in a sample unit, and our secondary analysis of EDSM data included both spatial and temporal replicate surveys within a subregion on the same day as the replicate surveys in the sample unit (hereafter referred to as EDSM with spatial replicates; Figure 1). Importantly, we opted to include both spatial and temporal replicate surveys within a subregion for the EDSM with spatial replicates analysis to make more direct comparisons with our primary EDSM analysis by having the exact same data included. This Accessed April 15, 2020, http://www.dfg.ca.gov/delta/ proje cts.asp?Proje ctID=SKT). SKT samples fish once a month from January through May at 40 fixed stations (Figure 1). Although supplemental tows have occurred opportunistically over the years, we restricted our analyses to data collected by the formal tows at the fixed sites of the SKT from January 2002 through May 2019. Notably, all sample data were included because tow condition information similar to EDSM are not available for these data. Samples were collected using a 10-min surface tow while deploying almost the identical net used in the EDSM Kodiak trawls. Before and after each tow, field crews collected environmental data. Specific conductivity (μs/cm 3 ), water turbidity (NTU), and water temperature (°C) at the surface were measured with calibrated meters. Secchi depth (m) and bottom depth (m) also were measured and recorded. Tides were recorded as: high (slack), ebb, low (slack), and flood. At the end of each tow, field crews recorded the tow distance (m) and duration (min), the volume of water sampled (m 3 ), the time of sampling, and the number of fish caught for each species. Again, we examined the environmental data for outliers that may have been caused by equipment malfunction and treated these values as missing data. We also treated tow volume measurements as missing data when the SKT metadata noted that there was a malfunction with the flow meter used to calculate tow volume. Similar to the EDSM with spatial replicates analysis, we considered all tows within the same subregion on the same day to be replicate surveys at a sample unit for the SKT data ( Figure 1). This resulted in a dataset with 2031 sample units with an average of 1.71 replicate tows (range: 1-5) at each sample unit. The SKT tows can sometimes be spaced out relatively far in time even if they are conducted in the same subregion on the same day. Therefore, we considered an analysis that implemented a time restriction to only included tow data if the tow was conducted within a similar time span across temporal replicate surveys at a sampling location for EDSM. However, preliminary analyses yielded nearly identical results between the two approaches. Therefore, we did not implement a time restriction to filter the SKT data for our final analyses that we present herein.
We downloaded supplemental environmental data that may influence the distribution and detectability of Delta Smelt and Longfin Smelt in addition to the environmental data collected during trawl surveys. Estimates of X2 were calculated using Dayflow, which is a computer program that was developed to determine daily historical delta boundary hydrology (California Department of Water Resources. Accessed May 11, 2020, https://water.ca.gov/Progr ams/ Envir onmen tal-Servi ces/Compl iance -Monit oring -And-Asses sment/ Dayfl ow-Data). Time since sunrise is a covariate that may capture heterogeneity in detection probabilities across the day. It is similar to the recorded time of each tow, but it may be more biologically relevant since species often adjust behavior in response to daylight and the time of sunrise is variable across the year. We downloaded sunrise times for Rio Vista, CA using the R package StreamMetabolism (Sefick, 2016), which uses maptools based on the NOAA sunrisesunset calculator. We used these sunrise data to calculate time since sunrise for all tows. Finally, how sample units are spatially arranged can have a large influence on the occurrence of fishes in the Bay-Delta in addition to the environmental covariates normally collected and evaluated for these species (Peterson & Barajas, 2018). This can be related to environmental conditions not captured by the other covariates that are correlated with space, such as bottom depth, spatial proximity to other sample units, etc. Thus, we measured distance between the tow location and nearest shoreline to try to account for some of this unmodeled spatial environmental heterogeneity.
For shoreline location, we downloaded the National Hydrography Dataset (NHD) High Resolution spatial data (U.S. Geological Survey. Accessed March 11, 2020, https://www.usgs.gov/core-scien cesyste ms/ngp/natio nal-hydro graph y/natio nal-hydro graph y-datas et?qt-scien ce_suppo rt_page_relat ed_con=0#qt-scien ce_suppo rt_ page_relat ed_con). The NHD represents the water drainage network in the USA with features such as rivers, streams, canals, lakes, ponds, coastline, dams, and stream gages. We downloaded these data for the State of California, clipped the shapefile to the area of interest, dissolved all polygons into a single polygon, and visually inspected the final shapefile to ensure there were no apparent issues. We then calculated the nearest distance between each tow location and the shoreline. Spatial analyses were completed in Esri ArcMap software (version 10.7.1).

| Statistical analyses
We estimated the distribution of Delta Smelt and Longfin Smelt while accounting for incomplete capture using occupancy models (MacKenzie et al., 2002). Thus, our models estimated the probability that a sample unit is occupied by a species (ψ j,k ) and the probability that the species is detected, given that the sample unit is occupied by the species (p i,j,k ), where i, j, and k denote tow, location, and tow set, respectively. Notably, occupancy models do not require that the same number of replicate tows occur in every sample unit.
Furthermore, removal designs similar to the EDSM sampling protocol can be accommodated with this analytical approach, although the precision of the estimates can be less when compared to standard occupancy modeling sampling protocols (MacKenzie & Royle, 2005). Both of the probabilities described above can be related to explanatory variables using logit-linear models .
Note that a central assumption of occupancy models is that the occupancy state across replicate surveys is static. Thus, explanatory variables cannot vary by tow for occupancy probability (ψ), but may vary by tow for detection probability (p). In cases where the explanatory variable was measured for each tow, we used the mean value (ignoring missing tow-specific values) as a covariate for occupancy probability and the tow-specific value as a covariate for detection probability.
We aimed to consider similar explanatory variables across species, datasets, and analyses, but this was not always possible (Table 1). Here, we describe the global model for the primary EDSM analysis first, followed by a description of the adjustments we made for the EDSM with spatial replicates analysis and then the SKT analysis. We evaluated if detection probability (p) was related to tow volume, secchi depth, tow number, time since sunrise, days since EDSM TA B L E 1 Mean, standard deviation (in parentheses), range, and percent missing data (if applicable) for continuous explanatory variables and frequency for categorical explanatory variables included in the occupancy models to estimate the probability that a sample unit is occupied (ψ), and the probability that the species is detected, given that the sample unit is occupied (p)  was initiated, dominant tide direction, and season. We note that days since EDSM was initiated was meant to capture potential learning (i.e., field crews may catch more fish as they gain experience in the field). It is also worth noting that season differed by species based on Merz et al. (2013) and Moyle et al. (2016), and this covariate was meant to represent potential differences in detection probability EDSM with spatial replicates analyses were identical except we did not consider models that had both mean temperature and the quadratic effect of day of year on occupancy probability (ψ) because these variables were highly correlated in these data. We held the tow number to be consistent with the primary EDSM analyses because this covariate was meant to capture trends (either increasing or decreasing) in capture efficiency across replicate tows at a sampling location. The analyses of the SKT data differed slightly than the primary EDSM analyses. For detection probability (p), we no longer incorporated season since all tows occurred in the same season and we replaced day since EDSM was initiated with year since the SKT was initiated. Dissolved oxygen was not measured for SKT surveys, so we did not consider it in the SKT analyses. Furthermore, we never considered day of year and mean temperature in the same model for occupancy probability (ψ) because they were highly correlated in the SKT data. Preliminary analyses indicated the SKT Delta Smelt data were too sparse to incorporate mean temperature on occupancy probability (ψ), so we did not include this variable when analyzing these data. Similarly, preliminary analyses indicated the SKT Longfin Smelt data were too sparse to incorporate mean specific conductivity on occupancy probability (ψ), so we did not include this variable when analyzing these data.
We used a sequential-by-submodel strategy for our model selection procedure due to the large number of explanatory variables we wanted to consider. In particular, we fit every combination of explanatory variables on detection probability (p) while holding occupancy probability (ψ) constant (i.e., an intercept model). We then evaluated submodel importance using Akaike's Information Criterion corrected for small sample size (AIC c ; Burnham & Anderson, 2002), where an important submodel was one that had a ∆AIC c less than five and no uninformative parameters (see Arnold, 2010). We then carried important submodels forward to the next step and repeated the procedure for occupancy probability (ψ). Our final model set was restricted to the submodels that were deemed important at each step of this process. Importantly, we used this approach because it has been shown to recover a substantial portion of the total AIC c model weight, recover the top-ranked model, and reduce the number of models fit by nearly half when compared to fitting every combination of covariates for occupancy analyses (Morin et al., 2020). For all analyses, we standardized continuous covariates to have mean of zero and a standard deviation of one prior to fitting models, and we imputed the mean value (i.e., zero) for all missing values. The respective reference categories for region, season, and tide were Far West, Subadult/Adult (Longfin Smelt) or Subadult (Delta Smelt), and Slack.
We described model parameters by their mean, standard error, and 95% confidence interval. We also calculated odds ratios for model coefficients (Hosmer & Lemeshow 2000) and considered a covariate to be strongly influential if the 95% confidence interval did not overlap zero. We also conducted a bootstrap goodness-of-fit test using the most complex model in the final model set for each analysis.
Specifically, we used the estimates from the most complex model to simulate 10,000 datasets, fit the model to those simulated data, and

| RE SULTS
Although the sequential-by-submodel strategy we used for our model selection procedure resulted in a large number submodel combinations that were evaluated, the number of models in the final model set was relatively low across datasets (Tables 2 and 3). We had competing models for most datasets (i.e., multiple models with ∆AIC c < 2); however, close inspection of the model output verified that the estimated model intercepts and coefficients were nearly identical across all competing models for each analysis. Thus, we based our inferences on the model with the greatest AIC c weight for each analysis (Tables 4 and 5 The primary EDSM analysis estimated Delta Smelt occupancy probability was, on average, very low across regions (Table 4).
Specifically, occupancy probability was, on average, 0.07 in Far West, with the odds of occupancy being 1.14 times higher in North, 10.59 times lower in South, and 1.25 times lower in West. Although, the TA B L E 2 Akaike's Information Criterion corrected for small sample size (AIC c ), change in AIC c (∆AIC c ), model weights (w i ), and number of parameters (k) for the probability that a sample unit is occupied (ψ), and the probability that the species is detected, given that the sample unit is occupied ( (Figure 2). The probability of detecting Delta Smelt given they were present was quite low (Table 4) The EDSM with spatial replicates analysis estimated Delta Smelt occupancy probability was, on average, relatively high (0.62) and, unlike the primary EDSM analysis, did not vary across regions (Table 4).
Occupancy probability still had a strong, positive relationship with mean temperature and a strong, negative relationship with mean specific conductivity; however, now occupancy probability was also found to have a strong, positive relationship with X2. Specifically, the odds of occupancy were 2.53 times higher with every 4.83°C increase in mean temperature, 1.57 times higher with every 8.80 km increase in X2, and 1.86 times lower with every 6997.59 μs/cm 3 increase in mean specific TA B L E 4 Mean, standard error (SE), and 95% confidence interval (CI) for the probability that a sample unit is occupied (ψ), and the probability that the species is detected, given that the sample unit is occupied (p) Note: "na" indicates that the explanatory variable was not considered in the analysis, and "-" indicates that the explanatory variable was included in the analysis but was not in the selected model.
conductivity. Although occupancy now increased with mean secchi depth and decreased with day of year and mean distance to shoreline, these relationships were not strong because the confidence intervals for the coefficient estimates overlapped zero. Although higher for this analysis, the on average occupancy probability decreased late in the year, which was similar to the primary EDSM analysis (Figure 2). Still, contrary to the primary EDSM analysis, no discernable pattern in the on average occupancy probability was found during the season when adults are the dominant life stage present. Overall, patterns in detection probability were similar to before with a few exceptions (  (Figure 2). Here, detection probability was relatively high, with the on average detection probability estimated at 0.22 (Table 4). Detection probability still had a strong, negative relationship with secchi depth and time since sunrise. The odds of detection were 1.28 and 3.25 times lower for every 1.93 h past sunrise and TA B L E 5 Mean, standard error (SE), and 95% confidence interval (CI) for the probability that a sample unit is occupied (ψ), and the probability that the species is detected, given that the sample unit is occupied (p) Note: "na" indicates that the explanatory variable was not considered in the analysis, and "-" indicates that the explanatory variable was included in the analysis but was not in the selected model. 0.41 m increase in secchi depth, respectively. Detection probability now increased with tow number with the odds of detection being 1.12 times higher with every 0.89 tows passed. Detection probability also had a strong, negative relationship with survey year (a covariate not considered in EDSM analyses) with the odds of detection being 2.14 times lower after every 4.99 years. Although tide direction was in the selected model, the effect of this covariate was not strong given the confidence interval for both coefficient estimates overlapped zero.
The primary EDSM analysis estimated the on average occupancy probability for Longfin Smelt was lower than Delta Smelt (Table 5).
Specifically, occupancy probability was 0.01, on average, in Far West.
The odds of occupancy were 9.87 times lower in North, 4.71 times lower in South, and 1.19 times lower in West. It seems regional variation in occupancy was governed by the difference between Far West and North, given the confidence intervals for the other coefficient estimates overlapped zero. Occupancy was positively related to X2 and mean distance to shoreline, but negatively related to mean temperature and mean secchi depth. In particular, the odds of oc- The EDSM with spatial replicates analysis estimated Longfin Smelt occupancy probability was relatively high when compared to the primary EDSM analysis (  SKT in X2, and 7.54 times lower with every 0.50 m increase in mean secchi depth. Although mean dissolved oxygen was now included in the selected model for occupancy probability, the relationship was not strong because the confidence interval for the coefficient estimate overlapped zero. Mean distance to shoreline was no longer in the selected model. Although relatively high compared to the primary EDSM analysis, similar patterns in the on average occupancy probability were found across the year (Figure 2). Detection probability was much lower for this analysis, with the on average detection probability estimated at 0.03 (Table 5) (Table 5).
For this analysis, detection probability was only related with secchi depth. The odds of detection were 17.29 times lower with every 0.41 m increase in secchi depth.

| D ISCUSS I ON
Occupancy models are increasingly used to adjust monitoring data for imperfect detection to more accurately track how and why species redistribute across dynamic landscapes; however, many longterm monitoring programs lack the temporal replicate surveys at a sample unit needed to fit these models using standard approaches.
We analyzed Delta Smelt and Longfin Smelt monitoring data that use temporal and spatial replicate surveys to demonstrate some of the tradeoffs associated with implementing a space-for-time substitution with occupancy models. In doing so, we also evaluated the influence of environmental and sampling conditions that are related to the distribution and detectability of these species of concern.
The probability a sample unit was occupied was much higher when using spatial replicates for both species and both datasets. This finding was expected and is likely related to the change in the spatial scale of our inferences. Subregions represented a sample unit when using spatial replicates, which is considerably larger than the area around a EDSM sampling location. As the spatial extent of a sample unit increases the likelihood a species is present in the sample unit concomitantly increases (i.e., the species area relationship ;Connor & McCoy 1979 (Kendall & White, 2009;MacKenzie et al., 2017).
Importantly, SKT data were sampled without replacement when treating stations as spatial replicate surveys, and we could not retroactively fix this potential issue in the survey data. This may partly explain the differences in the estimates when comparing results from the EDSM and SKT analyses. Thus, it would be worthwhile to modify the existing SKT sampling design to select survey stations within a subregion with replacement if this spatial resolution for occupancy dynamics will be considered further with these data.
As the spatial extent of our sample units increased, our ability to detect environmental effects on the distribution of species became complicated. In particular, the relationship between occupancy and explanatory variables differed by analysis, with some variables no longer having support and others having the opposite effect. Although both monitoring programs use the same sampling gear, it seems plausible that some of these patterns are related to the differences in spatial and temporal coverage between the two monitoring programs. That is, the different surveys may be sampling different segments of the populations. However, differences between the results of EDSM analyses that used the same exact data were also apparent. Thus, we suspect these patterns are also largely related to the tradeoffs associated with spatial grain and extent. As noted earlier, a fundamental assumption of occupancy models is that the occupancy state does not change across replicate surveys. This meant that tow-specific measurements of environmental conditions had to be combined across replicate tows to include as covariates on occupancy probability. However, as the spatial extent of sample units increased the environmental covariates we considered likely displayed more heterogeneity within a sample unit (Wiens, 1989).
Although standard practice, the combined value we used to capture environmental conditions (i.e., the mean) did not capture this heterogeneity and by extension likely were not representative of habitat conditions within a sample unit at the scale of a subregion. Thus, careful consideration needs to be given concerning the heterogeneity of environmental conditions within sample units when determining the appropriate spatial extent of sample units and how to summarize covariate information.
Inferences regarding detection probabilities also differed considerably between spatial and temporal replicate survey designs. When there are temporal replicate surveys within a sample unit, detection probabilities are interpreted as the probability of detecting a species given it occurs in the sample unit. In contrast, when using spatial replicates within a sample unit the detection probabilities estimate the probability of detecting a species at a tow location given they occur within the subregion at that tow location (i.e., availability for capture). The probability of detecting a species in a sample unit depends on three things: the unit is occupied, the probability of capturing at least one individual, and the number of individuals available for capture (Bayley & Peterson, 2001;Royle & Nichols, 2003). Thus, the observed decrease in Delta Smelt detection probabilities through time for SKT could be due to at least three potential mechanisms working individually or in concert. First, the decrease in detection probability could be due to a contraction in the distribution of a species within subregions leading to fewer SKT stations where Delta Smelt are present. Second, the decrease in detection probability could be due to an environmental change over time, such as water becoming clearer, that concomitantly reduces fish capture efficiency with nets.
Finally, the decrease in detection probability could be due to overall decreases in Delta Smelt abundance and thus, fewer fish available for capture. Although we are unable to determine the mechanism(s) responsible, this finding suggests that attempts to adjust historical SKT data with contemporary detection probability models may not be a valid approach.
This leads to the question: are we overly relying on occupancy models for the sake of using them when many long-term monitoring programs were established without considering this analytical approach? The answer to this question depends on the objectives of the monitoring program, the ecology of the species, the characteristics of its existing sampling design, and the willingness of managers to modify the existing sampling design if needed. The SKT and EDSM were both initiated with the aim of monitoring the distribution and relative abundance of Delta Smelt to inform policy and management decision making in the Bay-Delta. It is clear that monitoring data collected in the Bay-Delta are complicated by the observation process and that capture efficiency can vary considerably (Goertler et al., 2020;Latour, 2016;Mahardja et al., 2017Mahardja et al., , 2021Mitchell et al., 2017Mitchell et al., , 2019Newman, 2008;Peterson & Barajas, 2018;Polansky et al., 2019;Thomson et al., 2010;this study). Although trends in occupancy and relative abundance of fishes in the raw catch data for a few species and monitoring programs (with different sampling extents and resolutions) in the Bay-Delta appear unbiased when compared to estimates from multistate occupancy models (Peterson & Barajas, 2018), there are at least two primary reasons why not separating the ecological and observation process is problematic with respect to the goals of SKT, EDSM, and management in the Bay-Delta. First, the point estimates are systematically lower when not accounting for incomplete capture. This means that estimates of distribution and relative abundance of fishes depict a scenario that is much grimmer than reality when not accounting for incomplete capture. Although there seems to be a general sense that erring on the side of being conservative (i.e., estimating the species is worse off than reality) is less problematic when managing at-risk species, management in the Bay-Delta is complex and needs to consider the objectives of multiple stakeholder groups. Being less transparent in the status of species or using less reliable information to inform policy and management decisions can understandably lead to mistrust in the decision-making process, which can have negative ramifications for species conservation over the long term (Conroy & Peterson, 2013;Gregory et al., 2012). Second, the agreement between the trends in raw catch data and occupancy estimates while account- after it is clear that they implement a less than optimal sampling design to reliably track the system states of interest (reviewed in Thompson et al., 1998). The initiation of EDSM represents a potential unforeseen opportunity given the EDSM Kodiak trawl uses the same gear as SKT and EDSM incorporates temporal replicate surveys at a sampling location. In particular, integrated analyses may allow SKT data to be corrected for incomplete capture without modifying the SKT sampling protocol or extent of the sample unit if it can be assumed that detection probabilities do not differ substantially between the monitoring programs (Grabowski et al., 2009). Again, we caution against using this approach for historical SKT data without a proper evaluation of its validity given our results indicate system-wide changes may have occurred such that detection probabilities for historical SKT data are not comparable to detection probabilities for EDSM data. Additionally, occupancy analyses of these historical data would benefit from more integrated analyses that leverage data from the multiple monitoring programs while treating the different sampling gears as different detection methods directly in the model (Nichols et al., 2008). To our knowledge, Mahardja et al. (2021) were the first to try this in the Bay-Delta when analyzing Chinook Salmon detection/nondetection data at the spatial scale of a subregion. Our results suggest this approach requires more careful consideration to identify a smaller spatial scale when the objective is to relate occupancy probability to environmental covariates. Although the appropriate spatial scale is admittedly hard to define, we propose such an approach would likely require some consideration of the following: the spatial and temporal overlap in the monitoring programs to effectively borrow information across datasets; the movement ecology of the focal species to meet model assumptions of closure and independence (or perhaps use estimators that account for dependence in the data); the heterogeneity in environmental conditions to accurately represent habitat within sample units; and the scale at which information is used to inform management decisions to provide information that is directly relevant to managers and policy makers.
A primary objective of many long-term monitoring programs, including the trawl surveys in the Bay-Delta, is to track the distribution and relative abundance of species in relation to changing environmental conditions. N-mixture models have increasingly been used to analyze count data such as these to estimate spatial and temporal patterns in abundance while accounting for imperfect detection; however, these models are extremely sensitive to assumption violations and unmodeled heterogeneity in the monitoring data (Barker et al., 2017;Duarte et al., 2018;Link et al., 2018). This is the primary reason we opted to fit occupancy models in this study. Although occupancy probabilities can be positively correlated with abundance (Clare et al., 2015;Collier et al., 2013;Linden et al., 2017;Mathewson et al., 2012), this relationship is not always supported (Ellis et al., 2013;Stauffer et al., 2021). It is worth highlighting the potential utility of multistate occupancy models to correct monitoring data for imperfect detection while estimating more than two occupancy states (Nichols et al., 2007).
Specifically, occupancy states can be defined in terms of relative abundance: the probability that a sample unit is occupied and the probability that a sample unit is occupied by a large number of individuals (relatively abundant) given that the sample unit is occupied. This approach has proven useful for other systems and species (MacKenzie et al., 2009;Martin et al., 2010;Royle, 2004;Royle & Link, 2005;Whitlock et al., 2020) and for Bay-Delta fishes specifically (Peterson & Barajas, 2018). Still, using multistate occupancy models to estimate occupancy and relative abundance is relatively underutilized in the literature, which is likely related to the lack of an in-depth simulation-based study to evaluate how this estimator performs under various sampling and environmental conditions. We encourage work focused on evaluating how robust relative abundance estimates from multistate occupancy models are to provide some level of credence (or lack thereof) to this approach to estimate the distribution and relative abundance of species from replicate count data.
Information concerning how and why species redistribute across dynamics landscapes is invaluable to help guide policy and management decisions. Thus, it comes as no surprise that managers have invested considerable resources to support long-term monitoring programs. Our results demonstrate how analyzing real-world monitoring data using occupancy models and spatial replicates within a sample unit can influence inferences regarding the distribution and detection probability of species. Overall, we found using spatial replicates within a sample unit was not a substitution at all, but instead changed the scale of our inferences enough that it complicated our ability to relate environmental conditions to the distribution of pelagic fishes in the Bay-Delta. Our work highlights the importance of considering how the unique characteristics of monitoring programs, particularly the spatial extent of sample units, influences our inferences regarding the occurrence and detectability of species. Notably, even in the best of circumstances detecting a declining trend can take multiple years of monitoring, resulting in unnecessarily delayed and reactive management actions. This delay can potentially lead to system states that require a level of investment in resources that is beyond what is available to reverse species declines (reviewed in Grant et al., 2013). Thus, it would be worthwhile to translate results from analyses of monitoring data such as this into predictive models within a decision analytic framework to facilitate proactive management strategies and learning via adaptive management.

ACK N OWLED G M ENTS
We thank J. Kirsch, C. Phillis, and S. Acuña for their helpful com-