Long-term variability and density dependence in Hudson River Dreissena populations

This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited. © 2019 The Authors. Freshwater Biology published by John Wiley & Sons Ltd. 1Cary Institute of Ecosystem Studies, Millbrook, NY, USA 2Graham Sustainability Institute, University of Michigan, Ann Arbor, MI, USA 3Kellogg Biological Station, Michigan State University, Hickory Corners, MI, USA 4Department of Environmental Sciences, University of Virginia, Charlottesville, VA, USA


| INTRODUC TI ON
Mussels of the genus Dreissena are among the world's most problematic invasive species, causing large ecological and economic effects in many lakes, reservoirs, and rivers throughout large invaded ranges in North America and western Europe (e.g. Connelly, O'Neill, Knuth, & Brown, 2007;Higgins & Vander Zanden, 2010). The size and nature of these effects depend on the size and other characteristics of Dreissena populations (e.g. Higgins & Vander Zanden, 2010;Strayer, Solomon, Findlay, & Rosi, 2019), so understanding how Dreissena populations vary through time and what controls this variation is important to understanding their impacts on freshwater ecosystems and water infrastructure. Some Dreissena populations fluctuate greatly from year to year or show long-term trends (e.g. Ramcharan, Padilla, & Dodson, 1992;Stanczykowska, 1977;Strayer, Adamovich, et al., 2019), but detailed descriptions of interannual variation still are scarce. A recent comprehensive analysis of the long-term dynamics of Dreissena populations (Strayer, Adamovich, et al., 2019) found just three sites at which populations of adult dreissenids had been sampled annually for at least 20 years (the Hudson River; Oneida Lake, Hetherington et al., 2019; Lake Peipsi, Timm, Kangur, Timm, & Timm, 1996). The paucity of long-term data severely limits our understanding of interannual variation in Dreissena populations and their impacts.
Based on what has been proposed to control Dreissena populations, we offer a general framework for Dreissena demography (Figure 1) to organise research results and identify research gaps. This framework identifies several life stages (e.g. veligers, new settlers, adults), each of which could be subdivided (e.g. larval stages) for more detailed studies, as well as the transitions between these life stages. Each transition may be controlled by several factors, including the abiotic environment (inhospitable temperatures, disturbance), interactions with conspecifics (e.g. competition for food or space, cannibalism), or interactions with other biota (e.g. food quality, predation, disease).
Our paper has two goals. First, we describe the long-term dynamics of Dreissena populations (the zebra mussel Dreissena polymorpha, and the quagga mussel, Dreissena rostriformis, also known as Dreissena bugensis or Dreissena rostriformis bugensis, Stepien et al., 2014) in the freshwater tidal Hudson River. Our 27-year study is one of the longest and most detailed studies on the long-term behaviour of Dreissena populations, and includes information on changes in body size, body condition, and species composition as well as population density. We use these long-term data to test, revise, or extend conclusions made from previous analyses from the Hudson (Strayer & Malcom, 2006Strayer et al., 1996) and elsewhere. Second, in the context of the framework presented in Figure 1, we use long-term field data to try to estimate several key demographic parameters that determine population dynamics, and correlate long-term variation in these parameters with the biotic and abiotic variables that may control Dreissena populations.

| THE S TUDY ARE A
The study area is the freshwater, tidal reach of the Hudson River, which extends from river kilometre (RKM) 248 (measured upriver F I G U R E 1 Conceptual diagram of the demography of Dreissena populations, customised to match the Hudson River sampling programme. Blue italicised text shows possible controls on each of the demographic transitions, which are named in bold from The Battery in Manhattan) at Troy to RKM 99 at Newburgh.
The entire study area is subject to tides with a typical amplitude of 0.8-1.6 m, but sea salt typically is present only downriver of the study area, except during very dry periods (Cooper, Cantelmo, & Newton, 1988). The tides reverse the direction of water flow except during periods of high freshwater runoff and keep the water column from stratifying. The freshwater tidal reach of the Hudson averages 900 m wide and 8.3 m deep, and has a mean annual freshwater discharge of c. 500 m 3 /s. Extensive shallow areas (<3 m deep at low tide) cover c. 15% of the study area, and are colonised by rooted vegetation (chiefly Vallisneria americana). The water is turbid (growing-season Secchi depths usually are 1-2 m), hard (calcium = 25-30 mg/L), and nutrient rich (c. 0.5 mg/L NO 3 -N, c. 30 µg/L PO 4 -P) (Caraco et al., 1997;Strayer, Cole, et al., 2014).
The river bottom is predominately sandy above RKM 150 and muddy below that, but c. 7% of the study area is rocky (Strayer et al., 1996). Summer water temperatures usually reach 25-28°C (Wells & Young, 1992).

| Sample collection and processing
Our earlier publications (Strayer & Malcom, 2006Strayer et al., 1996) described sample collection and processing, as well as the underlying study design. Briefly, in 1993-2017 we used a stratified random design, containing six strata defined by the combination of two substratum types (soft = mud or sand, soft enough to sample with a PONAR grab; hard = rocky and too hard to sample with a PONAR grab) and three geographic regions (upper river = RKM 213-248; middle river = RKM 151-213; lower river = RKM 99-151).
We estimated the area of rocky bottom in the study area by trying to take samples with a petite PONAR grab (15 × 15 cm) at 253 randomly located sites. We classified a site as rocky if we failed to retrieve a sample in five attempts. We thus estimated that 7% of the river bottom in the study area is covered by rocky substrata.
Estimates for 1991-1992 were taken from several sampling programmes, as described by Strayer et al. (1996).
We sampled soft sediments annually from 1993 to 2017 at 44 sites in midsummer (mid-June to early August) using a standard (23 × 23 cm) PONAR grab. Samples were sieved in the field through a 2.8-mm mesh screen, placed on ice, and frozen upon return to the laboratory. We later thawed the samples and picked out all Dreissena, including those attached to shells of living unionid mussels. We hired divers to sample hard sediments from seven sites annually from 1993 to 2017 in early June and early August of each year (August only in 1993-1995 and 2000). Divers haphazardly collected 10 rocks 15-40 cm long from each site, which were then put onto ice and returned to the laboratory and kept refrigerated. Within 1-3 days, we counted all Dreissena >2 mm long. The August collections were intended to sample the population just before the settlement of the new year class, but in some years, early settlement led to large numbers of new settlers. In such cases, we raised the minimum size threshold to 5 mm. This affected our estimates of the density of the total population, but not those of individual year classes (the basis for our demographic analyses-see below), and had little effect on estimates of population biomass or filtration rates. We estimated the area of each rock by tracing its outline and weighing the tracing (this is the projected area of the river bottom covered by the rock, not its surface area-see Bailey, Chase, & Bechtold, 1995).
All mussels from the PONAR samples were preserved in 70% ethanol, and random subsamples (300 mussels/site, if available) from the diver samples were preserved in 70% ethanol. These mussels were later identified to species and measured for shell length using callipers. In addition, we developed shell-length to body mass regressions by taking 30-50 zebra mussels/site, deliberately chosen to represent a wide range of shell lengths, from the diver samples and freezing them. After thawing, we measured shell lengths using callipers, removed the bodies from their shells, dried them overnight at 60°C, and weighed them on an analytical balance (Mettler AE240).
We then ran a log 10 -log 10 regression of dry body mass on shell length for each sampling date and time. If too few mussels were collected from a sampling site, we either used a regression from a nearby site or pooled data with nearby sites to develop a pooled regression.
These regressions were used to estimate the shell-free dry mass of all mussels in the sample, and combined with the mass-filtration rate regression of Kryger and Riisgard (1988) to estimate community filtration rate. We express community filtration rate as m/day (the same as m 3 /m 2 -day) or as the percentage of the water column filtered/day, assuming a mean water depth in the study area of 8.3 m.
Unless stated otherwise, data are expressed as area-weighted river-wide means that apply to the entire study reach (RKM 99-248).
The major exception to this is that demographic analyses are based on samples taken from hard substrata in the midriver . This procedure is consistent with our earlier analyses (Strayer & Malcom, 2006), which focused on this geographically compact and ecologically relatively uniform part of the estuary, which contained well over half of the Dreissena population in the study area.  (Pace, Findlay, & Lints, 1992) has shown that such samples are adequate to represent the wellmixed water column. Samples were preserved in roughly a 1:1 ratio water to preservative for a final concentration of 2% formaldehyde.

| Estimating demographic parameters and their relationships to potential controlling variables
We used our long-term data on abundance of mussels of different size classes (Figures S1 and S2) to estimate fecundity, year-class (= cohort) strength, survivorship, and body condition, which are demographic parameters that underlie population dynamics, and correlate long-term variation in these parameters with various hypothesised intrinsic and extrinsic (= environmental) controls ( Figure 1).

| Realised fecundity
We define realised fecundity (F) as the number of veligers measured in our samples per spawning adult. Conceptually, it is the product of annual egg production per adult (E) and survival of eggs through larval development to the time veligers appear in our samples (S 1 ).
In practice, we estimated F directly from our measured data on densities of adults and veligers. However, we must combine our data with published information to estimate E and S 1 . We made a rough estimate of egg production by Dreissena in the Hudson (E) by developing empirical regressions between egg production and adult shell length, based on the studies of Walz (1978), Sprung (1991), and Stoeckel, Padilla, Schneider, and Rehmann (2004). The latter two studies counted eggs released only in the first episode of spawning, rather than the entire season, so we doubled these estimates following the conclusion of Walz (1978) that about half of the eggs are released in the first spawning episode. We excluded fecundity data from these studies from mussels <12.5 mm long, because such data were few and highly variable. The resulting regression (n = 35, r 2 = 0.66, p < 0.00001) between shell length (SL, in mm) and egg production (E) was: The standard error of estimate of the regression was 0.697, resulting in a correction factor of 1.275 when estimating fecundity in non-ln-transformed units (Sprugel, 1983) We applied this equation to the observed densities and size distributions of Hudson River mussels from June/July samples, assuming a sex ratio of 1:1 and that egg production of mussels <12 mm long was 0. We then calculated S 1 by dividing F (measured in our veliger samples) by this estimate of E.
We used multiple linear regression with model averaging (Burnham & Anderson, 2002) to test whether observed veliger densities (means over 1 May-30 September) were correlated with estimated egg production by adults, filtration rate of large (shell length >20 mm) adults (as a measure of possible cannibalism-see MacIsaac, Lonnee, & Leach, 1995; mean of June and August estimates), water temperature, and freshwater flow (temperature and flow as the mean for 1 May-30 September). We do not have satisfactory data on other factors (quantity and quality of food, predation on eggs and larvae) that might affect F, so they are not included in our models. With respect to food, although we have data on concentrations of phytoplankton chlorophyll a, previous studies (Sprung, 1989;Vanderploeg et al., 1996;Wacker et al., 2002) have shown that algal size and fatty acid composition are more important than bulk chlorophyll in determining larval growth and survival.
We used similar models to predict S 1 from filtration rate of large adults (mean of June and August estimates), water temperature, and freshwater flow (temperature and flow as the mean for 1 May-30 September).

| Year-class strength
We applied finite mixture analysis (Du, 2002) to separate the population into individual year classes and estimate the number of mussels in each year class and their mean shell lengths. Specifically, we used the mixdist package in R (Macdonald, 2018) to analyse data on densities of zebra mussels on hard substrata in midriver. Mixdist fits a set of frequency distributions to the observed size-frequency data; both the type of distribution and the number of distributions are specified as part of the model, and the fit of alternative models can be compared. Mixdist can fit a large number of possible models; we chose the best model from a consistent subset of logical, parsimonious models.
Specifically, for each sampling time, we fit four models (normal and lognormal distributions with both equal standard deviations and equal coefficients of variation) for one to three classes, resulting in a possible total of 12 models (some models did not converge and were therefore excluded). If justified by the data we also fit four-class models, but these four-class models never outperformed one-to three-class models. The goodness of fit of these models is assessed with a χ 2 statistic; we chose the model with the highest p-value (i.e. the least lack of fit) for the χ 2 statistic. Once the best models were chosen, we discarded a few models that suggested biologically improbable dynamics (specifically, large declines in mean shell length of a year class since the previous sampling period) and replaced them with the next-best model.
Numbers of quagga mussels in our samples were too low (see Results) for us to perform similar analyses on that species.
We estimated year-class strength from these models as the number of zebra mussels surviving to June or August in the year after ln E = 2.688 + 3.40 * ln SL (mm) settlement. These two measures of year-class strength (year-class size in June and August) were strongly correlated (r = 0.90), so we constructed a normalised year-class index by combining the June and August data as follows. We calculated a normalised June yearclass strength by dividing the year-class strength in each June by the mean year-class strength over all Junes, did the same for August, then took the average of June and August normalised densities as the year-class index. This index represents the ratio of year-class size in an individual year to the mean year-class size averaged over all years, so a year-class index of four would mean that four times as many zebra mussels as average were present in the summer after settlement.
We used multiple linear regression with model averaging to model the year-class index as a function of Dreissena filtration rate (mean of June and August estimates during the year of settlement), freshwater flow (mean freshwater flow at Green Island between 1 June and 30 September during the summer of settlement), degree-days above 25°C during the summer of settlement, and degree-days below 2°C during the winter after settlement.
We supplemented or replaced linear regression with quantile regression in a few cases in which we expected that the independent variable would bound rather than determine the value of the dependent variable, or observed a triangular scatter of points. We used the quantreg package in R (Kroenker, 2018) to fit a 90% quantile regression, and used the boot option to estimate standard errors and p-values by bootstrapping.

| Survivorship
We used the estimates of year-class size from the finite mixture models to estimate survival of zebra mussels in two ways. First, we plotted log 10 (density) of a single year class over time and fit a line through this plot to estimate the average instantaneous survival rate of that year class (in units of month −1 ; see Figure S3 for an example).
Because it is based on several data points, this estimate of survival is relatively precise (more precise than the following method) and has a readily calculated error (i.e. the error of the slope estimate of the regression). The disadvantage of this method is that the survival estimate is averaged over the long life of the year class (typically 2-5 years), so it provides information only on long-term trends in survival, not on year-to-year variation. Second, we estimated survival of each year class between every pair of sampling times simply as the ratio of population density at time t + 1 to that at time t. This method provides survival estimates that are associated with shorter intervals of time (e.g. August to August, June to August), and so are more readily correlated with variation in environmental variables.
However, individual estimates of population densities have large errors, resulting in imprecise estimates of survivorship (including estimates that are >1). Furthermore, it is not straightforward to calculate errors or confidence limits around these estimates. Because of these problems, we used this second method to estimate survivorship only for year classes with densities >5,000/m 2 12 months after settlement (statistical behaviour of smaller year classes was erratic; Figure S6).
These data limitations meant that we had too few data points (n = 7-9) to support a formal, multivariate model of short-term survival rates. Instead, we simply plotted survival rates against hypothesised controlling factors (density of predatory blue crabs, number of degree-days >25°C, number of degree-days <2°C).

| Shell growth
We took estimates of mean shell length of zebra mussels of each year class at each sampling time from the finite mixture models.
We used these data to analyse shell growth in two ways. First, we used a Walford plot (shell length at age t + 1 against shell length at age t, Walford, 1946) to describe shell growth of each year class from one August until the next. We used multiple linear regression with model averaging to model residuals from the Walford plot as a function of Dreissena filtration rate, phytoplankton biomass (as chlorophyll a), concentration of particulate inorganic matter, and water temperature, using the mean values of these variables from 1 August to 30 September plus 1 May to 31 July (i.e. the presumed period of growth). Second, we estimated daily instantaneous growth rates between the June and August samples, and used multiple linear regression with model averaging to model these growth rates as a function of the same independent variables as just listed (but averaged over the period between the June and August samples) plus June shell length.

| Condition
We calculated a condition index for every zebra mussel that we weighed, as the residual from the log 10 -log 10 regression of shell-free dry mass on shell length (cf. the residual index of Jakob, Marshall, & Uetz, 1996). That is, we ran the body mass-shell length regression for all specimens collected over the course of the study, then calculated the residual from that regression for each individual mussel as a measure of its condition. Positive values of condition indicate that a mussel weighed more than expected, and negative values indicate that it weighed less than expected, as expressed on a log 10 scale. We collected too few quagga mussels to allow us to calculate condition for this species.
As with shell growth, we modelled body condition in three ways.
First, we explored spatial and temporal variability in condition using a mixed effects model including a fixed effect for month (June or August) and random effects for year and for sampling site nested within river section. Second, we used multiple regression with model averaging to model body condition in August as a function of Dreissena filtration rate, phytoplankton biomass (as chlorophyll a), concentrations of particulate inorganic matter, and water temperature, using the mean values of these variables from 1 August to 30 September plus 1 May to 31 July (i.e. the approximate growing season in the year prior to sample collection). Third, we used the change in mussel condition between June and August as the dependent variable, and used the same independent variables as just listed, but averaged over the period between the June and August samples.

| Temporal change in Dreissena populations
Quagga mussels were first detected in the Hudson in 2008, well after the appearance of zebra mussels in 1991, and have constituted a small part (usually 1-10%) of the Dreissena community since then Although soft sediments cover 93% of the study area, the longterm mean percentage of Dreissena living on soft sediments in the Hudson was only 12.5%. Our earlier analysis (Strayer & Malcom, 2006) showed that the proportion of the Hudson's Dreissena population that  Figure S4), even though quagga mussels (which appeared in 2008) can develop large populations on soft sediments in lakes (e.g. Hetherington et al., 2019;Karatayev et al., 2014;Nalepa, Fanslow, & Pothoven, 2010). However, the numbers of Dreissena and the proportion of the overall Dreissena population living on soft sediments did vary greatly from year to year. We do not know the reason for this variation, but it does not seem to be related to the occurrence of high freshwater flows ( Figure S4).
Mean body mass of mussels also varied substantially (5-fold or more) from year to year (Figure 2b). The long-term record of zebra mussel body mass shows both cycling, which is out of phase with the cycles of population density, and a long-term decline (a multiple linear regression of mean zebra mussel body mass in August on Dreissena population density and year has r 2 = 0.40, a negative correlation with p = 0.003 for year, and a negative correlation with p = 0.055 for population density; see Figure S5 for an alternative quantile regression of these data). Quagga mussels have been consistently larger than zebra mussels, and the body masses of the two species have tended to vary synchronously (for 8 years of August samples, r = 0.66, p = 0.08). Quagga mussel body mass in August also tended to be negatively correlated with Dreissena population density (r = −0.61, p = 0.11, n = 8, Figure S5).
A more detailed analysis of the size structures of the Hudson's Dreissena populations shows strong interannual variation in the density of zebra mussels of different shell lengths (Figure 3, see also Figures S1 and S2 for more detail). Specifically, the time courses of densities of smaller mussels (<15 mm long) contain peaks corresponding to the appearance of dominant year classes (see below). The time courses of densities of larger mussels also show these peaks, lagged to account for growth, especially before 2005, but, after this, mussels >20 mm long almost disappeared from the population.
Body condition, defined as the log 10 -transformed residual from the overall shell length-body mass regression (calculated only for F I G U R E 2 Time courses of descriptors of Hudson River Dreissena populations: (a) river-wide mean population density of zebra (black) and quagga (white) mussels; (b) river-wide mean body mass (dry mass [DM], excluding shell) of zebra (black) and quagga (white) mussels; (c) body condition (see text for definition, ±1 SE) of zebra mussels in June (white) and August (black; we did not estimate body condition of quagga mussels); and (d) estimated filtration rate of Dreissena (both species combined)

| Analyses of causes of long-term demographic variation
We estimated long-term variation in several underlying demo-

| Veliger densities and survival from egg to veliger
Veliger densities were positively correlated with freshwater flow, but not with water temperature or estimated egg production (Table 1).

F I G U R E 3
Summary of temporal changes in size structure (shell length) in populations of the two species of Dreissena (zebra mussels in black, quagga mussels in white) in the Hudson River. Densities are weighted river-wide means. See Figures S1 and S2  There was a weak trend for increasing veliger densities over time (r = 0.41, p = 0.10). Although the linear regression analysis (Table 1) provided little evidence for any effect of filtration rates of large Dreissena, the plot of filtration rates versus veliger density suggested a triangular distribution of points (suggesting cannibalism; Figure 5). Estimated survival rates of eggs to the veliger stage were strongly and negatively correlated with filtration rates of large Dreissena (Table 2, Figure 5), but not with water temperature or freshwater flow (Table 2). However, we caution that the strong negative relationship between S 1 and filtration rates of large adults may have arisen, in part, because the dependent variable (veligers/egg) is roughly proportional to the inverse of the adult biomass and the independent variable (filtration rate) is roughly proportional to adult biomass.

| Year-class formation and strength
The year-class index, as measured by the number of mussels recruited to the summer following settlement, varied from 0 (no detectable year class) to 4.5 (i.e. 4.5 times the long-term mean), and showed eight strong year classes over the course of our study (Figure 4a). For reference, a year-class index of 1 corresponds roughly to a density of 10,000/m 2 on midriver rocks, or a river-wide density of 1,000/m

| Survivorship
Year-class survival rates were in the range of 0 to −0.07/month for the 1992-2004 year classes, then decreased to low levels for the [2005][2006][2007][2008][2009] year classes ( Figure 6). Survival of more recent year classes has risen from these low rates, but still tended to be substantially lower than for the 1992-2004 year classes.
Survival estimated over shorter intervals ( Figure 7) followed a similar pattern, with both August-June (overwinter) and June-August (oversummer) rates being highly variable. We could not discern any obvious patterns between this variability and hypothesised controlling factors (density of predatory blue crabs, hot summers, cold winters), with the exception that the lowest oversummer survival rates occurred in the two hottest summers ( Figure 8). However, the number of data points for these analyses was very small (n = 7-9), limiting statistical power.

TA B L E 2
Results of multiple regression models to predict log 10 (S 1 ), the survival rate from egg to veliger. Slope estimates are averaged over all models. The best model had an adjusted R 2 of 0.81

| Shell growth and body condition
We assessed four measures of shell growth or body condition of zebra mussels: the Walford plot of August shell lengths, shell growth of 1-year-olds between our June and August samples, body condition (see Methods for definition) of mussels in August, and change in body condition between June and August. In the Walford plot In general, models for shell growth and body condition had low predictive power, but suggested limitation by quantity and quality of food (Tables 4-7). In cases where variables had summed Akaike weights >0.5, phytoplankton biomass had positive effects on growth or condition, whereas Dreissena filtration rate ( Figure 10) and particulate inorganic matter concentration had negative effects.
Temperature had only weak effects of inconsistent sign.

| D ISCUSS I ON
Our study documented large, persistent interannual variation in mul-  variables: the interannual range in abundance was 79-fold, and the ranges in biomass, community filtration rate, and mean body mass were 129-fold, 119-fold, and 7-fold, respectively. We also found pervasive evidence of negative density dependence as a cause for that variation.

F I G U R E 5 Relationships between filtration rate of
This high interannual variation was a result of two elements: long-term, noncyclic changes in body size, veliger densities, adult survival, and species composition; and quasicyclic variation in abundance and body size. Long-term changes tended to be subtle and difficult to quantify in the face of shorter-term variability.
Nevertheless, we detected a large decline in zebra mussel body size, which was associated with conspicuous changes to the size structure of the population, and drove a long-term decline in community filtration rates. This decline in body size was caused by decreased survival of adults rather than a decline in rates of shell growth.
Survival fell in part as a result of rising mortality from blue crabs (Carlsson et al., 2011), and should have reduced cannibalism rates on veligers, because larger Dreissena are more effective predators on small zooplankton (MacIsaac et al., 1995).
We found only hints of other long-term, noncyclic changes.
There was a trend for decreased body condition over time for August but not June, and veliger densities tended to rise through time. It

F I G U R E 7
Estimates of survival of zebra mussels for large year classes (>5,000/m 2 on hard substrata in midriver at 12 months after settlement). Black lines and symbols apply to mussels in their second year of life (between 1 and 2 years after settlement); the red line and symbols in the upper panel apply to the interval between 1 and 3 years after settlement F I G U R E 8 Survival of 1-year-old zebra mussels between June and August as a function of degree-days >25°C before 1 August. Linear regression: r = −0.63, p = 0.12 F I G U R E 9 Walford plot for zebra mussels based on August samples from hard substrata in midriver. The outlier in the upper left represents the growth of the 1991 year class between August 1991 and August 1992, and was excluded from the statistical analysis reported in Table 4 is notable that we did not find long-term changes in the density of settled mussels of either species or both species combined, the proportion of quagga mussels in the community (after they appeared), or the density or proportion of the community living on soft sediments, contrary to common expectations (e.g. Burlakova, Karatayev, & Padilla, 2006;Heiler et al., 2013;Karatayev, Burlakova, & Padilla, 2015;Strayer & Malcom, 2006).
In addition, the species composition of the community changed over the long term with the appearance of the quagga mussel in 2008.
This species remains a small part (usually 1-10%) of the Hudson's dreissenid community, in contrast to many other sites at which quagga mussels have displaced zebra mussels (Strayer, Adamovich, et al., 2019).
Quasicyclic variation in abundance and body size was conspicuous. These cycles in turn caused cycling in community biomass and filtration rate (Figure 11), although interannual variation in these two variables was muted because the cycles of abundance and body mass were out of phase with one another (i.e. when abundance was high, body mass was low, and vice versa). Our record is TA B L E 4 Results of multiple regression models to predict deviations from the average growth relationship (i.e. residuals from the Walford plot [ Figure 9] of August shell lengths). Slope estimates are averaged over all models. The best model had an adjusted R 2 of 0.02

Summed Akaike weight Slope (SE)
Chlorophyll ( These observations are consistent with our earlier suggestion (Strayer & Malcom, 2006) that strong year classes suppress the recruitment and growth of subsequent year classes, until the size of the initial strong year class decreases to a point at which the negative density-dependence is weak enough to allow the formation of another strong year class. Models show that strong, negative effects of adults on veligers and juveniles can cause population cycling, and that the cycling driven by such mechanisms will persist indefinitely, and even restart if suppressed (Casagrandi, Mari, & Gatto, 2007;Strayer & Malcom, 2006

F I G U R E 11
Comparison of long-term population dynamics of dreissenids (zebra mussels = closed symbols; quagga mussels = open symbols) at three different sites. Biomass units are grams of shell-free dry matter (Hudson River, this study), dry mass including shells (Oneida Lake, New York, Hetherington et al., 2019;Strayer, Adamovich, et al., 2019), or wet mass, including shells (Lake Peipsi, Estonia-Russia, from Timm et al., 1996). Quagga mussels did not occur in Lake Peipsi and year-class strength were low. Nevertheless, low summer survival rates in the warmest summers are consistent with the idea that extended periods of warm temperatures may kill zebra mussels and thus limit populations (White et al., 2015). Young fish eat veligers (Nack, Limburg, & Schmidt, 2015), as presumably do predatory invertebrates, but we do not have any quantitative estimates of predation rates on veligers in the Hudson, and so cannot evaluate this possibility. Consequently, the importance of factors other than negative density-dependence in setting year-class strength in the Hudson is not yet clear.
Both the long-term behaviour and controls of Dreissena populations probably differ across different kinds of ecosystems (Strayer, Adamovich, et al., 2019). We are aware of only three sites where adult Dreissena were sampled annually and quantitatively for at least 20 years ( Figure 11). Detailed information on demographic parameters and controls are not available for these sites, but differences among the sites are obvious. In the Hudson, interannual variability was high, perhaps cyclic, and the arrival of the quagga mussel had only minor effects. Interannual variability in Dreissena biomass in Oneida Lake was much smaller than in the Hudson, and quagga mussels replaced zebra mussels almost one-for-one when they arrived.
Finally, in Lake Peipsi, zebra mussel biomass apparently underwent a step-change between 1964-1979 and 1980-1991. Each of these periods in Lake Peipsi contained a single catastrophic loss of zebra mussel biomass (in 1968 and 1986), but were otherwise relatively stable. We have suggested that Dreissena populations in the Hudson are controlled by cannibalism or competition for food, but not space (see also Strayer et al., 1996). In contrast, Hetherington et al. (2019) suggested that Oneida Lake populations are controlled by competition for space, and Timm et al. (1996) did not discuss what controlled zebra mussel populations in Lake Peipsi. The occasional catastrophic declines in Lake Peipsi would be consistent with either periodically unfavourable environmental conditions (e.g. hypoxia) or epidemics of disease.

Considering the ecological and economic importance of
Dreissena, it is striking how poorly we understand how the conceptual diagram of Figure 1 is translated into population dynamics in different kinds of ecosystems. We do know that survival from eggs to veligers to settled juveniles is low and variable, so that only a small proportion of eggs and juveniles reach the stage of settled juveniles (see also Lewandowski, 1982;Stanczykowska, 1977;Wood, 2013), and that densities of these different life stages are poorly correlated (Hetherington et al., 2019;Jones & Ricciardi, 2014;Nalepa, Wojcik, Fanslow, & Lang, 1995;Strayer, Adamovich, et al., 2019). As in marine invertebrates with free-living larvae (e.g. Fraschetti, Giangrande, Terlizzi, & Boero, 2002), it seems likely that early post-settlement mortality is critical in setting year-class strength (Wood, 2013). However, because the tiny, newly settled larvae are so difficult to study, there are almost no measurements of this mortality under different field conditions, or demonstrations of the importance of hypothesised controls (but see Wood, 2013). Even our understanding of the demography of later stages is unsatisfactory. There are few measurements of demographic transition rates (i.e. F, S 2 -S 4 in Figure 1