Opening the black box of winter: Full‐year dynamics of crustacean zooplankton along a nearshore depth gradient in a large lake

Past studies of zooplankton seasonality in large temperate lakes have often neglected the winter period. Winter conditions are rapidly changing (e.g., reduced ice cover extent and duration, altered thermal and mixing regimes) in northern lakes, making it important to fill the existing winter knowledge gap. In this study, we sampled five stations in Lake Superior across a nearshore depth gradient through the full year to assess the phenology of crustacean zooplankton communities and the effect of environmental drivers on them. Across stations, zooplankton densities were the lowest in winter (0.9 ± 0.6 Ind. L−1) and highest in summer (14.2 ± 15.1 Ind. L−1). Zooplankton abundances and community composition were less seasonally variable at deeper stations compared to shallower and more terrestrially affected regions. Cladocerans were the dominant taxonomic group during the summer across all stations, while cyclopoid and calanoid copepods were more important during the fall, winter, and spring. Among feeding groups, herbivores were most abundant in summer while omnivores and carnivores dominated in winter. We found that water temperature and food availability were the main drivers of total zooplankton densities through the year and during the cold seasons, but the effect of these factors varied among the main taxonomic groups. Our study demonstrates seasonal and spatial variation in crustacean zooplankton and environmental parameters, with the highest fluctuation at shallower stations. This study offers new information on seasonal crustacean zooplankton dynamics and contributes to understanding the effects of climate change on large lake ecosystems.

ecological processes in winter can be more dynamic than traditionally believed and that conditions during the cold periods of the year are important for lake ecosystem function throughout the year (Salonen et al. 2009;Hampton et al. 2017;Hébert et al. 2021).
To date, most winter research has focused on relatively small lakes. Large lakes have received less attention during winter, mainly because of the logistical, instrumental, and safety challenges of cold weather research on these systems . However, there are several reasons why the winter dynamics of large lakes require dedicated study. Large northern lakes contain much of the world's accessible surface freshwater and provide ecosystem services to millions of people (Sterner et al. 2020). The size of large lakes results in many important differences in their physical, chemical, and biological processes compared to smaller systems, making it difficult to extrapolate information from smaller systems to the world's largest lakes (Tilzer 1990;Ozersky et al. 2021). Finally, there is evidence that large, ice-covered lakes are experience more rapid warming than their smaller counterparts (O'Reilly et al. 2015). Thus, more comprehensive full-year studies of large lakes are needed in order to understand and predict the impacts of ongoing climate change on these important ecosystems.
Crustacean zooplankton plays important roles in the trophic ecology of lakes (Richardson 2008) and can also be used as indicators of environmental change (Jeppesen et al. 2011). It is well known that zooplankton abundance and community composition vary seasonally within lakes and also among lakes during the ice-free period (Sommer et al. 2012;Leech et al. 2018). The drivers of seasonal and among-lake variation in zsooplankton abundance and community composition during the open water season are many and include water temperature, food availability, lake morphometry, and top-down forces (Dodson et al. 2005;Yurista et al. 2009;Straile 2015;Paquette et al. 2022). Less is known about under-ice zooplankton communities and the environmental factors that shape them. Existing winter studies show that, while less abundant than during the open water period, under-ice zooplankton communities remain active during winter, evolve through the ice cover period, and vary among lakes (Grosbois et al. 2017;Hampton et al. 2017;Shchapov et al. 2021). Surveys across lakes show large variation in the abundance of winter zooplankton, often with higher densities in more productive lakes (Hampton et al. 2017;Kalinowska et al. 2019;Shchapov et al. 2021), suggesting food availability plays a role in maintaining winter zooplankton populations. Community composition analysis shows consistent shifts away from cladocerans and toward calanoid and cyclopoid copepods during the ice cover period. The increased relative abundance of copepods is likely explained by their ability to accumulate lipids (Grosbois et al. 2017) and higher flexibility in their feeding habits compared to many obligate filter-feeding cladocerans (Santer 1993;Sommer and Sommer 2006). Given the large winter-to-summer variation in the forces that are known to structure zooplankton communities and the ongoing changes in winter climate, integrating winter into our understanding of zooplankton ecology is an important priority.
To understand how zooplankton communities change seasonally in regions experiencing different winter and summer conditions, we studied the zooplankton communities of Lake Superior through the year at five nearshore stations across a depth gradient. Our objectives were to assess seasonal variations in biotic and abiotic conditions important to zooplankton across our study region and determine the patterns and drivers of full-year seasonal succession of crustacean zooplankton communities. We hypothesized that: (1) biotic and abiotic conditions will fluctuate more at shallow, nearshore stations compared to deeper stations, due to former's more dynamic thermal regimes and greater effect of seasonally variable terrestrial inputs; (2) as observed elsewhere (e.g., Hampton et al. 2017), zooplankton abundances will be lower during winter compared to other seasons across all stations; however, (3) zooplankton abundance and community composition at deeper stations will be less seasonally variable than at shallow stations because of more stable environmental conditions in deeper locations; (4) food availability and temperature will be the main drivers for densities of both total zooplankton and for densities of the main taxonomic groups (cladocerans, calanoids, and cyclopoids) across all stations and seasons.

Sampling locations
Five nearshore stations in the western arm of Lake Superior (Fig. 1, Table 1) were sampled from February 2017 through May 2018. Winter (December 21-March 20) and spring (March 21-June 21) were prioritized and sampled at least 3-5 times to capture the seasonal development of biotic and abiotic conditions across stations. Summer (June 22-September 21) sampling was carried out once or twice at each site shortly after complete water stratification. Fall (September 22-December 20) sampling was done during the fall overturn. The study stations varied in physical, biological, and chemical parameters and represented a gradient of depth, terrestrial influence, and winter conditions (Table 1). Duluth Harbor is a shallow station located in the channel between St. Louis River and Lake Superior, making it the most terrestrially and anthropogenically affected station in our study. Chequamegon Bay was the second shallowest sampling station, followed by the Duluth Entry station, which is often affected by runoff from the steep hills of the City of Duluth and outflow of the St. Louis River. The other two stations, Madeline Island and McQuade Harbor, were the deepest and least terrestrially affected stations in our study. Madeline Island, Duluth Harbor, and Chequamegon Bay were ice-covered during both study winters, while McQuade Harbor and Duluth Entry were ice free. While only 750 m offshore, the McQuade Harbor station is located along an exposed and steeply sloping shoreline

Shchapov and Ozersky
Winter zooplankton of Lake Superior and is the most similar of our stations to open Lake Superior in thermal regime and water chemistry.

Sample collection
Samples were collected from the ice or from a small boat during ice-free periods. Ice thickness and snow depth were measured during winter at stations with ice cover. We used a YSI EXO2 multiparameter sonde (YSI, Yellow Spring, OH, USA) to record depth profiles of water temperature, conductivity, dissolved oxygen, total dissolved solids, and fluorescent dissolved organic matter. Water column transparency and light attenuation were estimated (Table S1) using a LI-COR meter equipped with a cosine LI-192 sensor (LI-COR Biosciences., NE, Lincoln, USA). Light attenuation (k d ) by ice and snow cover was measured by recording photosynthetically available radiation (PAR) values at the ice/snow surface and by submerging and placing the light sensor against the undersurface of the ice cover. Light measurements through the water column were made at 0.5 m increments from the water surface to 2-15 m depth and were used to determine water k d coefficients. Light measurements were used to calculate the euphotic depth (Z eu ) for each station by accounting for water column, ice (combined with snow) k d values (Table S1). Z eu was calculated as the depth at which mean daily PAR equals 1 μmol m À2 s À1 according to Silsbe et al. (2016) and Shchapov et al. (2021) by using water column, ice, and snow k d values in combination with MODIS-Aqua satellite data on daily PAR at the earth surface for our study locations (https://earthdata. nasa.gov).
We used a 3.7-liter Van Dorn sampler to collect water from the surface, various depths through the water column, and 0.5 m above the lake bottom. Water was stored in acid-washed 2-liter bottles until return to the laboratory and later used for chlorophyll a (Chl a), dissolved organic carbon (DOC), particulate organic carbon (POC), particulate organic nitrogen (PON), and total phosphorous (TP) analyses. We used a 64 μm mesh net with 0.5 m mouth diameter to collect crustacean zooplankton samples from ca. 1.5 m above the lake bottom to the lake surface. Zooplankton were fixed in 90% ethyl alcohol, then transferred to 70% ethyl alcohol and stored until counting and identification.

Lab analyses
We estimated phytoplankton biomass as Chl a; water from each depth was filtered through 0.2 μm nitrocellulose filters, which were soaked in 90% acetone for an 18 h extraction in the dark. Extracts were analyzed on a Turner Designs 10-AU fluorometer (Turner Design, Sunnyvale, CA, USA) using 436 nm for excitation and 680 nm wavelength for emission. DOC samples were filtered through precombusted Whatman GF/F filters into precombusted 40 mL amber glass vials and analyzed on a Shimadzu TOC-V autoanalyzer (Shimadzu Co., Kyoto, Japan). POC and PON samples were filtered onto precombusted Whatman GF/F filters and kept at À20 C. Then filters were dried at 60 C and analyzed using elemental analysis-isotope ratio mass spectrometry (Finigan Delta Plus XP, Thermo Electron Corporation, Bremen, Germany; EA-IRMS) at the Large Lakes Observatory. NIST (RM8548ammonium sulfate) and in-house calibration standards such as acetanilide, B-2153 soil, B-2159 flour, and caffeine were used repeatedly after every 10 samples.
Subsamples of preserved zooplankton were counted and identified in a Bogorov chamber under an Olympus SZH10 stereoscopic microscope. We used Balcer et al. (1984) for taxonomic identification and feeding group classification. At least 300 individuals were counted and identified in each sample. Our identification was limited to the species level for adult copepods and cladocerans, while copepodites were distinguished only as cyclopoids or calanoids. Nauplii were counted without taxonomic identification. Adult individuals were separated into three main taxonomic (calanoids, cyclopoids, and cladocerans) and feeding groups (herbivores, omnivores, and predators). We noted the number of copepods carrying egg sacks and of cladocerans with eggs or embryos in their brood pouches. Zooplankton abundances are reported as individuals per liter (Ind. L À1 ).

Data analysis
Depth-integrated values for environmental parameters (Chl a, DOC, POC, PON, and TP) were calculated by using trapezoidal integration (Shchapov et al. 2021). Seasonal changes in zooplankton communities were visualized by nonmetric multidimensional scaling (NMDS) using the R vegan package (Oksanen et al. 2018), on forth-root transformed zooplankton abundance. Nauplii were not included in the total zooplankton, main, and taxonomic groups abundances due to relatively high densities. Seasonal differences in zooplankton community composition were assessed using permutational multivariate ANOVA (PERMANOVA) using the vegan package.
The vegan package was also used for similarity percentage analysis (SIMPER) to identify species that contributed most to dissimilarity between seasons.
To find potential environmental parameters that explain total zooplankton, calanoid, cyclopoid, and cladoceran abundances (the first three categories included adults and copepodites) across stations and seasons, we used multivariate linear regression (MLR) analysis. For the full models we used surface temperature, PON concentration, and Z eu as predictors. Surface temperature is an important physical factor for zooplankton community structure and seasonal dynamics, while PON can indicate food availability (bottom-up effect) for zooplankton. We used Z eu as a metric of the light environment, and hence potential for predation on zooplankton by visual predators. MLR analysis was conducted in the R car package (Fox et al. 2012). To meet the assumption of normality, total zooplankton and cyclopoid abundance data were log10 transformed, while calanoid and cladoceran abundance data were log10 + 0.001 transformed. Surface temperature and PON were log10 transformed, and Z eu was square-root transformed. Surface temperature was used as a polynomial term (squared) due to curvature in the relationship between variables for total zooplankton, cyclopoid, and cladoceran abundance models. Collinearity among predictor variables was explored before running MLR analysis by using a correlation matrix plot ( Fig. S1; Schloerke et al. 2018) together with variance inflation factor values. The relationship between the dependent variable and each predictor parameter was evaluated using added variable plots (Figs. S2-S7). To find the most parsimonious models for zooplankton abundance, full models went through backward manual selection, until all predictors remained significant. Further, the best fitted models were identified by adjusted R 2 , F-test, and Akaike's information criterion (AIC) values. Along with models representing full seasonal data, we ran MLR analysis restricted to the colder water period defined by water surface temperature < 4 C. MLR and all other statistical analysis were performed in the R software environment (version 3.6.2; R Core Team 2017).

Environmental factors
Surface water temperature, light availability, winter ice cover, and snow conditions differed throughout the sampling period across sampled locations (Fig. 2, Table 1 and Table S1). Across all stations, average water surface temperatures were 1 AE 1.3 (SD), 4 AE 3.5, 16.2 AE 3.6, and 3.9 AE 2.4 C in winter, spring, summer, and fall, respectively. The Duluth Harbor station had the highest surface temperature in summer (20.2 C) and, together with Chequamegon Bay, warmed up faster than other locations (to ≥ 10 C by May; Fig. 2). Other stations' maximum summer surface temperatures ranged between 16 C and 18 C.
Three of our sampling stations were ice-covered, with measured ice thickness 12-43 cm in 2017 and 39-62 cm in 2018. Snow depth over the ice ranged 0-7.8 cm (3 cm average) in winter 2017, and 2.8-12.8 cm (5 cm average) in winter 2018. The maximum snow cover was recorded for the Madeline Island station (12.8 cm) in February 2018.
Light conditions across seasons and stations were assessed as euphotic depth (Z eu ) ( Fig. 2 We observed seasonal changes in depth-integrated water chemistry parameters across all five stations (Fig. 3). Overall, winter depth-integrated Chl a, DOC, POC, PON, and TP values were lower than in other seasons. Duluth Harbor had higher, and more seasonally variable, concentrations of Chl a, DOC, POC, PON, and TP than other stations. Chequamegon Bay was the second most productive station based on Chl a, POC and PON concentrations. Duluth Entry, Madeline Island and McQuade Harbor were similarly unproductive and showed muted seasonality in most water chemistry parameters.

Zooplankton abundance
Total zooplankton (adults and copepodites) and nauplii densities varied across seasons and stations (Fig. 4). Zooplankton densities were highest in summer at all stations, but peak abundances varied dramatically among stations. The highest  densities were also found in summer for all stations, except for Duluth Entry with a peak in May (Fig. 4). Like total zooplankton values, nauplii densities were the higher for shallow Duluth Harbor (13.6 Ind. L À1 ) and Chequamegon Bay (10.9 Ind. L À1 ) compared to deeper, less terrestrially influenced stations (Fig. 4). Minimum nauplii values were similar (≤ 0.1 Ind. L À1 ) among all stations and were recorded in December or March for different stations. Calanoid and cyclopoid copepodites were most abundant in summer with more than two-fold higher values for Duluth Harbor and Chequamegon  Bay compared to other stations. The minimum abundances of calanoid and cyclopoid copepodites were similar and recorded during winter 2017 and 2018.

Main taxonomic and feeding groups
The densities and relative abundance of the main zooplankton taxonomic groups changed across seasons and stations (Figs. 4 and 5). Absolute densities of adult calanoids were relatively similar through the sampling period across all stations, with slightly higher values at the end of summer or fall for Madeline Island and Duluth Entry stations. Cyclopoids had similar seasonal changes in density across stations, with peaks in summer and fall (Fig. 5). Cladoceran densities were the highest at the shallow Duluth Harbor and Chequamegon Bay stations during summer and were low at the deeper stations throughout the year.
Calanoids had the highest relative abundance in fall and winter 2018, with the lowest values recorded in summer across all stations (Fig. 5). However, Duluth Entry and McQuade Harbor calanoid relative abundances were still high in the beginning of summer (> 50%). The highest relative abundances of calanoids for deeper stations Duluth Entry, Madeline Island, and McQuade Harbor were recorded during winter and spring and ranged 40%-70%. The more terrestrially affected and shallower stations Duluth Harbor and Chequamegon Bay showed high relative abundance of cyclopoids during winter and spring. The lowest relative abundances of cyclopoids were recorded during summer across all stations (Fig. 5). Summer zooplankton across all stations was dominated mostly by cladocerans, especially in Duluth Harbor and Chequamegon Bay. abundance of cladocerans was high during winter only at Duluth Harbor (57%).
Relative abundance of main feeding groups showed that herbivores were more abundant in summer (51.5 AE 38.9%) than winter (7.2 AE 18.3%) across all stations, except Duluth Entry (Fig. 5). At Duluth Entry, the main feeding groups were dominated by omnivores (59.4 AE 14.6%) and predators (37.4 AE 15.6%) throughout the year, with a slight increase of relative herbivore abundance during summer (10.3 AE 10.5%). Herbivore percentage abundance was low during fall, winter, and spring at deeper Madeline Island and McQuade Harbor stations and ranged between 0% and 3%. At the shallower Duluth Harbor and Chequamegon Bay stations, herbivore relative abundance was high in winter 2017 (68%) and spring 2018 (84%), respectively. Omnivore and predator relative abundances were high for all stations during winters and transition periods. Those two groups comprised about 100% of the total adult zooplankton population during winter and spring in Duluth Entry, Madeline Island, and McQuade Harbor stations.

Zooplankton species composition
We identified 15 species in summer, 14 species in fall, 11 in spring, and 9 in winter across all stations. The highest number of species for Duluth Harbor (9 spp.) was found in fall, while the highest numbers for Chequamegon Bay (9 spp.), McQuade Harbor (8 spp.), Madeline Island (7 spp.), and Duluth Entry (6 spp.) were found in summer. The lowest number of species was similar across stations (2 to 3 spp.) and was recorded in winter and spring.
We performed NMDS to visually assess the differences in zooplankton communities among seasons across sampled stations and their correlates (Fig. 6). Surface temperature, Chl a, and PON were significantly correlated with seasonal changes in zooplankton community structure. Community structure differed significantly across all seasons (PERMANOVA p = 0.0001) and stations (PERMANOVA p = 0.0002). We found that communities winter vs. summer (p = 0.0001), spring vs. summer (p = 0.0002), and summer vs. fall (p = 0.003) were significantly different across all stations. Winter vs. fall (p = 0.08), winter vs. spring (p = 0.6), and spring vs. fall (p = 0.1) were not significantly different. Among stations, deeper Duluth Entry, Madeline Island, and McQuade Harbor were more different from Table 2. Summary of the best models based on surface temperature, PON (particulate organic nitrogen), and euphotic depth for total zooplankton, calanoid, cyclopoid, and cladoceran abundances including the best models only for the winter. shallower Duluth Harbor and Chequamegon Bay stations according to PERMANOVA (p < 0.05). SIMPER analysis of zooplankton communities across all stations identified several cladoceran and cyclopoid species as the main drivers of winter-tosummer and spring-to-summer community dissimilarity (Table S2). Three cladoceran species were responsible for most of the summer-to-fall dissimilarly. A combination of calanoids, cladocerans, and cyclopoids drove winter-to-spring and fall-towinter dissimilarities.

Best models
Zooplankton and environmental parameters MLR analysis showed that surface temperature explained more than 53% of the variation in total zooplankton abundance across all stations and seasons ( Table 2). The best model of calanoid abundance across the full year included surface temperature and PON, which together explained 39% of variation. Cyclopoid abundance was significantly predicted only by surface temperature, which explained 39% of abundance variation. The best model of cladoceran abundance included surface temperature and PON. Together, those factors explained 77% of the variation, with surface temperature having the stronger effect among the two variables (Table 2). We found only two significant (p < 0.05) models for cold-water season zooplankton abundances (water surface temperature ≤ 4 C; Table 2). The model for calanoids included only PON (with a negative slope) and explained 23% of the abundance variation. The best model for cladocerans also included only PON (but with a positive slope) and explained more than 50% variation in abundance.

Discussion
We investigated changes in zooplankton communities across a depth gradient of Lake Superior throughout the year, focusing on the understudied winter and transition seasons. Zooplankton abundance was higher in summer compared to winter across all stations. However, the magnitude of differences between summer and winter was much higher at shallower than deeper stations. Across stations, zooplankton abundances and community composition were more similar in winter than other seasons, especially summer. Calanoid and cyclopoid copepods dominated the zooplankton community in winter, spring, and fall, while cladocerans were more abundant and dominant during summer. Omnivorous and predatory strategies were more common among zooplankton during winter and transition seasons, while herbivory was more common across all stations during the summer. Temperature and food availability were important drivers of total zooplankton abundances, with some important differences among main taxonomic groups. These results are relevant to understanding how changes in climate impact the lower trophic levels and food web structure of large, spatially heterogeneous lakes.

Seasonal changes of environmental conditions
Zooplankton densities and community composition are affected by many factors (Pinel-Alloul et al. 1995;Dodson et al. 2005). Our study showed large seasonal variations in environmental conditions across all stations. In partial agreement with our first hypothesis, shallow stations displayed larger variation in most abiotic and biotic conditions between summer and winter. Light can indirectly affect zooplankton as a key driver of phytoplankton production, but can also enhance predation on zooplankton by planktivorous fish (Williamson et al. 2011). Light conditions (measured here as euphotic depth) showed large seasonal variations across most stations. Euphotic depth seasonality was relatively muted at Duluth Harbor, our shallowest station, where turbidity from the St. Louis River contributed to high light attenuation throughout the year. As expected, ice and snow attenuated light penetration during winter, and ice-free stations had greater euphotic depths during the winter-spring period compared to stations with ice and snow cover. Nonetheless, at most times during our winter observations, the euphotic depth (which we defined as the depth of minimum net primary production requirements; Silsbe et al. 2016) at icecovered stations extended 0.4-14 m into the water column, suggesting that some phytoplankton growth is possible at our sites through most of the year.
Along with light, water temperature is an important physical driver for phytoplankton and zooplankton. Overall, surface temperature had similar seasonal patterns across all stations, with peaks in summer and low values in winter. However, shallow stations warmed faster and reached higher maximum summer temperatures than deeper stations. Average surface water temperatures are increasing in many northern lakes due to climate change (O'Reilly et al. 2015;Sharma et al. 2021), and observations suggest that nearshore regions are warming especially fast in the Great Lakes (Mason et al. 2016). Faster nearshore warming in large lakes in the future may exacerbate the summer-period spatial differences in the abundance and community composition of zooplankton that we observed along our depth gradient (see below).
We measured seston quantity as the depth-integrated Chl a (a proxy for phytoplankton biomass), POC, and PON concentrations. Phytoplankton is an important part of the diets of herbivorous and omnivorous zooplankton and is widely used to indicate marine and lacustrine trophic status (Bell and Kalff 2001;Kasprzak et al. 2008). Many previous studies showed that phytoplankton biomass and production diminish during winter (e.g., Hampton et al. 2017;Shchapov et al. 2021). However, high under ice phytoplankton biomass and production has been documented in many systems including Lakes Baikal, Michigan, and Erie and smaller lakes in Poland, Japan, and the USA (Maeda and Ichimura 1973;Vanderploeg et al. 1992;Mackay et al. 2006;Twiss et al. 2012;Kalinowska et al. 2019;Bramburger et al. 2022). While all of our stations showed relatively low seston concentrations during winter, there were important differences among them. The shallower Duluth Harbor and Chequamegon Bay stations had higher winter period Chl a, PON, and POC concentrations than deeper sites. Deeper sites (especially Duluth Entry and Madeline Island) also had less seasonal variation in seston concentrations than shallow stations. This suggests that more productive and shallow regions of large lakes can maintain relatively high phytoplankton concentration during winter, an observation supported from studies of smaller lakes across productivity gradients (Hampton et al. 2017;Shchapov et al. 2021).
Although phytoplankton is a key component in aquatic food webs, its contribution to the seston pool can sometimes be small since seston is a mixture of autochthonous and allochthonous particles (Hessen et al. 2003). Chl a, POC, and PON all followed similar seasonal patterns at our study sites with generally low values in winter. This suggests that seston POC and PON had phytoplankton origin throughout most of the year, which is also supported by the C : N ratios of seston (Table S1). Generally, low C : N ratios (4-8 M) indicate high algal contribution to the seston pool, while high ratios (> 10) suggest a higher proportion of terrestrial matter (Thornton and McManus 1994). In our study, average seston C : N ratios did not display large seasonal variation and corresponded to algal origins (6.2-8.2) across all stations and throughout the year with two exceptions. The Duluth and McQuade Harbor sites had relatively high seston C : N ratios (11.5) in fall 2017, following a large storm that likely contributed large quantities of terrestrial material.

Seasonal succession of crustacean zooplankton
Most of the spring, summer, and fall total zooplankton densities we observed are similar to values previously reported for Lake Superior (Watson and Wilson 1978;Barbiero et al. 2001). However, summer zooplankton densities in shallow Duluth Harbor and Chequamegon Bay were more similar to values reported by Barbiero et al. (2001) for the warmer and more productive Lakes Ontario and Erie. Total zooplankton densities were lower in winter compared to other seasons (especially summer) across our study stations, in agreement with our second hypothesis. Results also supported our third hypothesis, showing that total zooplankton densities were much more seasonally variable at shallow than deep stations. Along the depth gradient (Duluth Harbor, Chequamegon Bay, Duluth Entry, Madeline Island, McQuade Harbor), summerwinter densities varied 329-, 53-, 8-, 10-, and 33-fold, respectively. The lower zooplankton densities we observed in winter agree with other studies (Hampton et al. 2017;Kalinowska et al. 2019;Karpowicz 2020, Shchapov et al. 2021) where winter zooplankton were reported to be 1.2to 38-fold lower compared to summer values. The reason for this large variability of summer-to-winter differences in zooplankton densities across systems is an open question, but our results suggest that shallow and productive systems that warm rapidly after ice-off may display larger summer-to-winter variations than deeper and more oligotrophic ones.
The seasonality of zooplankton abundance and community structure depends in part on the seasonality of reproduction (Vargas et al. 2006). In our study, we used egg-and embryocarrying adults, nauplii, and copepodite abundances as indicators of zooplankton reproduction. Egg-carrying copepods in Chequamegon Bay, Duluth Entry, and McQuade Harbor stations reached their maximum densities at the end of winterbeginning of spring. We found that copepod nauplii densities increased from spring to their peak during summer and were generally followed by an increase in copepodites densities for all stations. This is likely driven by the rise in water temperature and food availability toward summer, creating favorable conditions for zooplankton development. For example, our shallower, faster-warming, and more productive stations had three-to four-fold higher densities of nauplii and gravid cladocerans in summer compared to deeper stations. Egg-and embryo-carrying cladocerans were absent in all deeper stations during winter but were present in low densities in shallower Duluth Harbor and Chequamegon Bay stations. Generally, cladocerans (especially gravid cladocerans) are absent during winter due to cold temperatures but have been reported in systems with high phytoplankton biomass (Slusarczyk 2009;Shchapov et al. 2021).
Overall, our results showed the dominance of calanoid and cyclopoid copepods during winter and transition periods, while cladocerans were more abundant during summer. Spring and fall codominance of calanoid and cyclopoid groups was previously reported for Lake Superior and can be explained by the bivoltine life cycle of some cyclopoid and calanoid species (Link et al. 2004;Pawlowski et al. 2018), along with their ability to store lipids and use them during winter months (Grosbois et al. 2017;Hébert et al. 2021). Moreover, cyclopoid and calanoid copepods can change their feeding strategies in cold environments with low primary production (Santer 1993). For instance, cyclopoids can switch from omnivory to raptorial feeding and consume other zooplankton when phytoplankton biomass is low (Balcer et al. 1984). The omnivorous calanoid Limnocalanus macrurus was a common species in our winter samples, and are believed to be a key predator of copepod nauplii in the Great Lakes during winter (Warren 1985). Indeed, our feeding group results showed that predators and omnivores were the dominant groups during winter, spring, and fall seasons across all stations except Duluth Harbor in winter 2017. This suggests that predatory and omnivorous zooplankton may be better adapted to cold conditions due to lower reliance on phytoplankton compared to herbivorous zooplankton. In contrast, we found shallow stations strictly dominated by herbivores in summer, while deeper stations showed the codominance of omnivores and predators, with herbivore abundance peaking by the end of summer. This delay in herbivore density increase may be explained by lower water temperatures, weak water column stratification, and reduced phytoplankton biomass at our deeper stations.
In our study, herbivores were mostly represented by cladocerans. Many cladocerans require water temperature above 12 C for optimal growth and reproduction (Gillooly and Dodson 2000). Densities of Daphnia spp. in our study increased at the end of spring when water temperature exceeded 10 C. In total, across all stations, we found three Daphnia species most of which were present during summer and the early fall. In summer, cladoceran relative abundance varied between 50 and 90% of total zooplankton across all our stations, which is higher than reported for Lake Superior previously (10%-20%; Watson and Wilson 1978;Sprules and Jin 1990;Barbiero et al. 2001). The observed higher values are likely due to our focus on relatively shallow, nearshore stations. Another possible explanation for the relatively high cladoceran densities we observed as compared to historic values may be rising water temperatures (Link et al. 2004). For instance, long-term observations in Lake Baikal showed a three-fold increase in summer cladoceran population since 1946 in association with warming (Hampton et al. 2008). Considering the predicted rise in water temperature of Lake Superior and other large systems, cladocerans may become a more common and dominant member of future zooplankton communities in large and cold lakes.
Seasonally, the most evident differences in diversity and community composition were observed between summer and fall and summer and winter zooplankton communities. In comparison, zooplankton community composition and structure were relatively similar across all sites during winter. The sharpest decline in the number of species coincided with the decrease of surface water temperature to below 3 C during the fall-winter transition. We found seven species in summer that were completely absent or below detection during fall and winter across all stations. Together, those seven zooplankton species represented between 30 and 97% of summer total abundances for shallower Duluth Harbor and Chequamegon Bay while for deeper Duluth Entry, Madeline Island, and McQuade Harbor stations they contributed between 3 and 35% of summer abundances. Of all recorded species, only L. macrurus, Diacyclops thomasi, Leptodiaptomus sicilis, Acanthocyclops vernalis, Bosmina longirostris, and Tropocyclops prasinus mexicanus were present year-around. Among stations, Duluth Harbor showed the highest winter diversity, with eight species observed in winter 2017. Out of those eight species, five were present only in Duluth Harbor and not found in other stations. Interestingly, up to 50% of 2017 winter zooplankton in Duluth Harbor was composed of the small cladoceran B. longirostris, whereas in 2018 only the copepods D. thomasi and L. sicilis were present. This may be due to better conditions for zooplankton during winter 2017 compared to winter 2018 (less snow and thus higher light levels and food concentrations). Chequamegon Bay, Duluth Entry, Madeline Island, and McQuade Harbor all had similar winter species composition. There, the copepods L. macrurus, D. thomasi, and L. sicilis were the most abundant in both winters, contributing between 88 and 100% of the total adult zooplankton population. Our results generally agree with those of other studies, which show typical zooplankton dominance by calanoid and cyclopoid species during the winter (Hampton et al. 2008;Grosbois et al. 2017;Jensen 2019;Shchapov et al. 2021). For example, Twiss et al. (2012) reported that copepods contributed more than 90% of total zooplankton density under the ice in Lake Erie in winter, but no more than 50% in summer (Barbiero et al. 2001). In Lake Baikal, winter zooplankton communities are usually dominated entirely by the calanoid Epischura baicalensis (Afanasyeva 1998). Similarly, studies of small eutrophic lakes in Finland and Poland by Ventelä et al. (1998) and Kalinowska and Karpowicz (2020) showed winter dominance of calanoid and cyclopoid species.

Drivers of zooplankton abundance and community structure
The effects of environmental drivers on zooplankton abundance and community structure during winter are not well studied. Open-water season studies show that factors such as temperature, light, mixing, nutrient limitation, and top-down controls combine to determine among-lake variation in zooplankton abundance, community composition, and seasonal succession (McCauley and Kalff 1981;Shuter and Ing 1997;Leech et al. 2018;Paquette et al. 2022). While we were unable to evaluate the role of all these factors, multiple linear regression and NMDS analysis identified water temperature and food availability as the most important drivers of zooplankton abundance and community composition across the full year in our study. While this finding is in apparent agreement with our fourth hypothesis, we saw some unexpected differences among the main taxonomic groups. Total zooplankton and cladoceran (which comprised a large part of total abundances during the summer) across the full year were strongly related to water temperature, with peak densities coinciding with peaks in summer water temperatures. As shown in previous studies, water temperature can be a cue for zooplankton reproduction (Gillooly and Dodson 2000;Vandekerkhove et al. 2005) and is important for growth and the seasonal dynamics of zooplankton (Shuter and Ing 1997;Straile 2015). Interestingly, temperature was much less predictive of cyclopoid and calanoid copepod densities. The relatively weak effect of temperature on cyclopoid and calanoids in our study can potentially be explained by their adaptation to cold conditions (Albers et al. 1996), complex life cycle (Reid and Williamson 2010), and reliance on lipid-rich food availability prior winter (Grosbois et al. 2017).
We used PON concentrations as a proxy for zooplankton food availability. Our results showed that PON concentrations (along with temperature) were a significant predictor of cladoceran abundance across the full year, with a positive relationship between PON and abundance. PON was also an important predictor of cladoceran abundances (explaining 50% of variation) during the cold part of the year, again with a positive relationship. PON did not emerge as a significant predictor of total zooplankton abundance or the abundance of cyclopoids. Interestingly, the relationship between seston PON and calanoid abundance was negative throughout the year and also during winter. PON is a coarse metric of food availability, and is probably a much better measure of resources for filter-feeding herbivorous cladocerans than for omnivorous and predatory copepods, which are more selective feeders and can feed on wider range of food types, including other zooplankton (Balcer et al. 1984;Becker et al. 2004;Sommer and Sommer 2006). It is possible that the negative relationships with PON we observed are the result of confounding or indirect interactions we were unable to account for. A subject of particular interest may be effects of fish predation on the relatively sparse but lipid-rich (e.g., Grosbois et al. 2017) zooplankton of ice-covered lakes. Link et al. (1995) showed active winter feeding by Lake Herring on zooplankton in Lake Superior, with strong selection for calanoid and against cyclopoid copepods. We attempted to include the top-down effect of visual predators on zooplankton abundance by including euphotic depth in our models. Euphotic depth did not emerge as a significant predictor, and at best can only be considered as an indirect metric of potential fish predation. Further studies evaluating the impact of broad environmental drivers on zooplankton through the full annual cycle will help to better predict what impacts loss of ice, warming waters, and other anthropogenic changes will have on lake zooplankton.

Data availability statement
Original data for this work are available via the Data Repository for University of Minnesota (https://doi.org/10.13020/ 852y-7p95).