Commentary: Demographic response of osprey within the lower Chesapeake Bay to fluctuations in menhaden stock

Atlantic menhaden (Brevooria tyrannus; herein menhaden) support high volume fisheries and are important forage for many fishes, seabirds, and marine mammals. Manager and stakeholder concerns about the impacts of menhaden harvest on ecosystem processes, biodiversity, and productivity of other valuable fisheries have motivated advancements in menhaden assessment modeling and development of an ecosystem management approach (Anstead et al., 2021). Ecological reference points are derived from an ecosystem model including selected fish predators (Chagaris et al., 2020), however, full ecosystem models have quantitatively linked menhaden to seabirds (Buchheister et al., 2017). Watts et al., (2024) summarized valuable long-term demographic and foraging ecology information for osprey (Pandion haliaetus) inhabiting the Mobjack Bay subestuary of lower Chesapeake Bay.From data collected during specific times (May -Jul;1974-1975, 1985, 2006-2007, 2021), the authors reported that osprey reproductive rate (no. surviving young per breeding pair), brood size, fish provisioning, and percentage of diet comprising menhaden declined substantially over the study duration. Particularly concerning was that estimated reproductive rates after 1985 were insufficient to offset adult mortality.These osprey demographic and foraging results warrant consideration of mechanisms that may underlay the documented temporal patterns. One candidate is that menhaden abundance has declined over time thereby reducing availability to osprey, which the authors investigated by regressing mean osprey reproductive rate on coastwide indices of juvenile (young-of-year) abundance (herein coastwide indices). A statistically significant regression was reported (□□ < 0.01, □□ 2 = 0.91) prompting the recommendation "osprey population viability requires that the menhaden stock be restored to 1980s levels". Prior to submitting the paper for publication, the authors made key findings public through a press release (CCB, 2023) that garnered considerable regional (Williams, 2023) and national attention (Hurdle, 2023) resulting in several stakeholder groups demanding an immediate moratorium on menhaden harvest in Chesapeake Bay.The regression analysis assumed the coastwide indices, which are model-based data aggregations from 16 surveys spanning Rhode Island to South Carolina, represent the abundance and availability of menhaden to osprey in Mobjack Bay, Virginia. This was justified by stating the coastwide indices were the only stock-wide empirical abundance metric that fully covered the osprey study period, and through a correlation analysis of those indices with both the Maryland juvenile and the Mid-Atlantic Adult (MAD, age-1 and older) indices. The correlation analysis appears aimed at reconciling discrepancies between the broad spatial scale of the coastwide indices and the localized osprey demographic data, as well as the differing age-classes represented by coastwide and MAD indices, considering that osprey prey on age-1 and older menhaden. The authors apparently reasoned that a significantly positive correlation between the coastwide and Maryland juvenile indices would indicate alignment between broad-scale and regional juvenile abundance patterns. Similarly, a significantly positive correlation between the coastwide and MAD indices would indicate consistency between juvenile and older fish abundance patterns. Therefore, strong positive correlations in these cases would support using the coastwide indices to explain osprey reproductive rate data in Mobjack Bay. However, we argue this reasoning ignores key aspects of menhaden population dynamics and that the data analysis was flawed.(i.)The correlation between the coastwide and Maryland juvenile indices was significant (Spearman Rank, □□ = 63, □□ = 0.61, □□ < 0.05), however, this result is misleading because the Maryland data are included in estimation of the coastwide indices. In fact, the Maryland survey is the only data source contributing to every year of the coastwide time-series, which gives it disproportionate influence and essentially reduces the correlation analysis to a comparison of a variable with itself. Additionally, asynchronous recruitment patterns between menhaden nursery habitats in Southern New England and Chesapeake Bay have been documented (Buchheister et al., 2016), and since the coastwide indices aggregate these divergent patterns, they cannot accurately represent any specific localized area. (ii.)Correlation of the coastwide and MAD indices was significant (Spearman, □□ = 37, □□ = 0.4, □□ < 0.05), however, this analysis is concerning for several reasons. Firstly, the MAD indices began in 1985, yet 33% of the osprey the reproductive rate data are from the 1970s.Assuming the correlation holds prior to 1985 requires extrapolating beyond the temporal scope of the data, which is not recommended. Secondly, a scatter plot of the coastwide and MAD indices showed the data from 1985 and 1986 likely had a strong influence on the estimated correlation coefficient. When these years were excluded, the correlation was no longer statistically significant (Spearman, □□ = 35, □□ = 0.31, □□ = 0.07), indicating an overall weak relationship between the coastwide and MAD indices. Thirdly, because the peak ingress of larval menhaden to Chesapeake Bay occurs in January and February (Lozano and Houde, 2013), these individuals originated from spawning during the previous fall (Latour et al., 2023). Consequently, comparisons between the coastwide and MAD indices should account for this temporal lag. For one-and two-year lags, the correlation coefficients were not statistically significant (Spearman, □□ 1 = 36, □□ 1 = 0.25, □□ 1 = 0.13; □□ 2 = 35, □□ 2 = 0.22, □□ 2 = 0.20). This lack of coherence between metrics of adult and juvenile abundance is not surprising, since attempts to discern a formal stock-recruitment relationship have been unsuccessful (SEDAR, 2020). Collectively, these points show the coastwide indices do not accurately reflect the abundance of age-1 and older fish and therefore are not a reliable proxy for the abundance and availability of menhaden to osprey in Mobjack Bay. (iii.) For the osprey field study, each nest constituted the sampling unit, and several derived statistics were computed from data recorded for each breeding pair. The reproductive rates used in the regression were arithmetic means of counts of young produced per nest. This approach is less than ideal for several reasons. Firstly, analyses should generally be conducted on raw data (e.g., a GLM fitted to counts of surviving young per nest) rather than means because the latter are estimated quantities with error. Regression analyses that treat means as observed data typically ignore this error and thus fail to appropriately characterize uncertainty, which is critical to goodness-of-fit and regression parameter hypothesis tests. Moreover, averaging removes sampling unit replication, which for the osprey study, resulted in a regression model fitted to only □□ = 6 data points. Secondly, simple linear regression assumes the response variable is normally distributed, which is problematic since reproductive rate cannot be negative. Using reproductive rate values (□□ □□ ) inferred from Watts et al., (2024;Figure 4), the 95% data interval (□□̅ ± 1.96 □□□□ □□ ) was computed to be (-0.07, 2.36) which shows the lower limit is negative. While the normal distribution can be used for some nonnegative data types, the osprey reproductive rate data fail the "95% range check" (Limpert and Stahel, 2011) indicating an asymmetric distribution should be considered. Thirdly, and related to the previous point, we contend the regression model relating □□ □□ to the coastwide menhaden indices (□□ □□ ) should be formulated to reflect a multiplicative, lognormal error structure: □□ □□ = (□□ 0 + □□ 1 □□ □□ )□□ □□ □□ such that log(□□ □□ ) = log(□□ 0 + □□ 1 □□ □□ ) + □□ □□ (see Appendix). Application of this form of the model showed a positive relationship (Figure 1A), but the □□ 1 estimate that was on the threshold of statistical significance (□□ = 0.05).Given this result, we conducted a simulation analysis to assess the robustness of the linear regression model. For the □□ □□ values, standard errors (□□□□ □□ ) were also inferred from Watts et al., (2024;Figure 4) and used to parameterize unique lognormal distributions, □□ □□ ~□□□□(□□ □□ , □□ □□ 2 ),where □□ □□ 2 = log ( □□□□ □□ 2 □□ □□ 2 + 1) and □□ □□ = log(□□ □□ ) -1 2 □□ □□ 2 . Using these distributions, 1000 sets of □□ □□ were randomly generated and regressed on the coastwide menhaden indices. Results indicated that only 41.1% of the □□ 1 estimates from the simulated data sets were statistically significant (Figure 1B). All analyses were performed with the R software program (R Core Team, 2024) In summary, while we share concerns about the demographic and foraging trends of osprey in Mobjack Bay, the analyses presented in Watts et al., (2024) do not establish a clear relationship with menhaden abundance and availability. Furthermore, in the absence of a discernable stock-recruitment relationship, managing to "1980s levels" (or any specific level) is quite challenging since menhaden recruitment appears to be shaped by mechanisms associated with interactions among the coastal distribution of spawners, larval transport, climatology, and nursery habitat suitability more than harvest (Buchheister et al., 2016). Moving forward, we strongly encourage concurrent osprey nest and menhaden sampling, followed by the application of appropriate statistical methodologies, to directly assess their ecological linkage, particularly in the broader context of shifting distributions of fish populations in the Northwest Atlantic Ocean.Figure 1. (A) Relationship between osprey annual reproductive rates (□□ □□ ) and coastwide juvenile menhaden abundance indices (□□ □□ ). Point labels are (□□, □□ □□ , □□ □□ , □□□□ □□ ), where □□ □□ and □□□□ □□ were inferred from Watts et al., (2024;Figure 4) and the back-transformed, bias-corrected fitted line with it 95% shaded confidence interval overlaid. Note that Watts et al., (2024) reported □□ 2021 = 0.22 but the correct value is □□ 2021 = 0.36 (ASMFC, 2022), and the equal spacing of coastwide index data on the x-axis of Watts et al., (2024;Figure 4) obscured the true pattern of the data. (B) Histogram of pvalues associated with significance test of □□ 1 from the 1000 simulated regressions with a multiplicative, lognormal error structure where the dashed line denotes □□ = 0.05.