Pseudocalanus Copepods in the Bering Sea: Species Identification Intraspecific Diversity and Biogeography

The calanoid copepod genus, Pseudocalanus, is an important secondary producer and widespread throughout the northern hemisphere. While ecologically important, previous work with this genus has been limited due to difficulties in morphometric identification. Here, Pseudocalanus species composition and distribution were examined in the Bering Sea through sequencing of the mitochondrial CO1 gene of 642 individuals from 10 locations. Phylogentic trees of interand intraspecific diversity demonstrated that in the spring of 2010 (March-June), four species of Pseudocalanus co-existed at various locations across the shelf. Species composition in the middle shelf (50-100 m depth) was significantly different from composition on offshore domain (≥200 m depth). Correlation of species composition with salinity indicated that environmental factors influence distribution. Differences between the middle shelf and offshore populations extended to the intraspecific level in two of the four species. Although high haplotypic diversity was observed in all Pseudocalanus spp., a weak structure between the middle shelf domain and offshore stations was only seen in P. acuspes and P. minutus. This indicates that while there may currently or previously have been some barriers, gene flow is possible across the shelf. Evidence of mixed intraspecific and divided interspecific distributions, in addition to correlation with environmental conditions, highlighted the importance of a species-specific approach to future research of this genus and will increase accuracy in the modeling of an ecosystem-wide response to climate change.

The calanoid copepod genus, Pseudocalanus, consists of seven distinct and ecologically important species that are widespread and commonly co-occurring throughout temperate to arctic regions of the Northern Hemisphere. Due to difficulties in morphometric differentiation between species, very few studies have investigated the meso-scale distribution of these copepods, much less the factors influencing those distributions. Here, inter-and intraspecific diversity of Pseudocalanus was examined in the southeastern Bering Sea. The mitochondrial CO1 sequences of 642 individuals from 10 stations across the Bering Sea shelf were determined during the spring of 2010 (March-June). Four distinct species (P. acuspes, P. mimus, P. minutus, and P. newmani) were identified with an average genetic divergence of 14.7% ± 4%. Species were not evenly distributed across the Bering Sea shelf. Stations located in the middle shelf domain (50 -100 m depth) were dominated by P. acuspes and had higher total abundance of Pseudocalanus copepods than stations located in the offshore domain (>200 m depth), which were dominated by P. minutus. A multi-dimensional scaling plot (MDS) of relative species composition data showed two groupings, each of 40% similarity. Based on relative species composition, offshore stations were significantly different from middle shelf stations (ANOSIM p = 0.01). Haplotypic diversity was high in all species and in combination with a negative Tajima's D could indicate population expansion. A weak genetic structure was found between the middle shelf and offshore stations in P. minutus and P. acuspes. This indicates that while there may currently or previously have been some barriers preventing mating and the exchange of genes, gene flow is possible across the shelf. Correlation of species composition with salinity and depth of the chlorophyll maximum suggests that distribution of species across the Bering Sea shelf is a result of environmental selection rather than isolation. Differential distribution of the four species and correlation with environmental conditions highlighted the importance of a speciesspecific approach to future research of this genus and will assist in creating greater accuracy in the modeling of an ecosystem-wide response to climate change.

Introduction:
Copepods belonging to the genus Pseudocalanus are one of the most numerically dominant zooplankton in temperate-boreal waters of the Northern Hemisphere (Corkett and McLaren, 1978). Perhaps because of their high abundance and wide distribution, they serve an important role in the marine food web as secondary producers, directly connecting the primary production of microalgae with higher trophic levels (Coyle and Pinchuk, 2002). Ecologically and economically important species of juvenile fish, such as sprat and herring, rely upon these copepods as a significant prey item (Corkett and McLaren, 1978;Hinrichs et al., 2002;Möllmann and Köster, 2002), and while small in size, they may be more productive in coastal waters than larger copepods (e.g., Calanus spp.) (Corkett and McLaren, 1978).
Despite their importance, difficulties distinguishing between species morphometrically have limited interspecific research, creating a gap in knowledge of factors affecting the abundance and distribution of individual species (Frost, 1989;Bucklin et al., 1998). Within the Pseudocalanus genus, seven species have been recorded throughout the Northern Hemisphere, P. acuspes, P. elongatus, P. major, P. moultoni, P. mimus, P. minutus, and P. newmani (Frost, 1989). As these species are commonly co-occurring (Corkett and McLaren, 1978;Frost, 1989), most studies have focused on the genus as a whole (e.g., Coyle et al., 2011;Coyle et al., 2008;Coyle and Pinchuk, 2002;Vidal, 1980), and physiologies are largely unaccounted for at the species level (Frost 1989). Problems with grouping Pseudocalanus spp. for analysis of an ecosystem arise if individual species respond to unique environmental parameters for growth and survival. For example, there are differences in overwintering strategies among Pseudocalanus spp. (Corkett and McLaren, 1978;McLaren et al., 1989). Species that have a primarily Arctic distribution, like P. minutus and P. acuspes, store large amounts of lipids compared to those with temperate to sub-arctic distributions (McLaren et al., 1989). Thus, different species may have different nutritional values that could affect the diet of juvenile fish.
Although species are difficult to identify morphometrically, there are clear genetic differences at the mitochondrial cytochrome c oxidase I gene (mtCO1) (Bucklin et al., 1998;Bucklin et al., 2003). Interspecific divergence in the mtCOI is estimated to range from 10-20% and intraspecific divergence from 0.5-2% (Bucklin et al., 2003). To date, six of the seven morphologically recognized species of Pseudocalanus have had the mtCOI gene sequenced and found significant divergence appropriate for a species-level classification (Bucklin et al., 2003). These sequences are now available on the international database, GenBank (Holmborn et al., 2010;Bucklin et al., 2003).
With the use of the mtCO1 gene, an increased resolution of distribution patterns of individual species has now been established for Georges Bank, the Baltic Sea, the English Channel, the Black Sea, Svalbard fjords, and Balsfjord, Norway (Bucklin et al., 1998;Bucklin et al., 2001;McGillicuddy and Bucklin, 2002;Unal et al., 2006;Holmborn et al., 2010;Aarbakke et al., 2011). Most of these ecosystems contain only one member of this genus. The exceptions, Georges Bank and the NE Atlantic, have two and three species of Pseudocalanus respectively. By contrast, there are four species of Pseudocalanus in the Bering Sea: the primarily arctic P. acuspes and P. minutus, and the temperate/sub-arctic P. newmani and P. mimus (Frost 1989).
This multi-species system is ideal for field-based studies of ecological differences among Pseudocalanus spp.
The southeastern shelf of the Bering Sea is >500 km wide and accounts for over half the fish and shellfish landings in the United States (National Marine Fisheries Service, http://www.st.nmfs.noaa.gov/commercial-fisheries/fus/fus11/index).
The inter-annual variations in the abundance of Pseudocalanus copepods have been shown to be a key factor determining the success of economically important fisheries such as walleye pollock (Theragra chalcogramma) (Napp et al., 2000;Napp and Hunt, 2001). In years when temperatures above the thermocline are warmer, Pseudocalanus copepods become the primary, if not exclusive prey items for larval fish in the region (Coyle et al., 2008).
Although no fine-scale mapping of Pseudocalanus spp. has been recorded for the Bering Sea, variations in the zooplankton community across the shelf and between warm and cold years indicate that there may be well-defined spatial and interannual distribution patterns for these smaller calanoid copepods (Cooney and Coyle, 1982). Water depth and spring/summer front systems divide the Bering Sea shelf into three domains: the inner (0-50 m), middle (50-100 m) and outer domain (>100 m depth) (Coachman, 1986). A generally weak cross-shelf advection results in distinct communities within these domains (Coachman, 1986). Larger copepods are found in the outer domains and offshore, whereas smaller copepods typically have higher abundances in the inner and middle domains (Coachman, 1986;Cooney and Coyle, 1982). During "cold" years, larger copepods dominate the zooplankton biomass, while smaller copepods are more abundant during "warmer" years (Coyle and Pinchuck, 2002, Napp et al., 2002, Stabeno et al., 2012.
Pseudocalanus spp. have produced a rather mixed interannual signal with varied abundance patterns during "warm" and "cold" years (Coyle et al., 2008. This confusion may be caused by the grouping of multiple species within this genus .
In this study, we used the mtCO1 gene to determine the spatial distribution of

Field Collection
Zooplankton were collected at 10 stations from 2 cruises, the USCGC Polar Sea PSEA1001 and the R/V Thomas G. Thompson TN249, from across the SE Bering Sea shelf between March and June of 2010 (Table 1). Copepods were sampled with a 150 and 333 µm mesh ring net deployed in vertical tows to a depth of 60 m or 1 meter from bottom (if depth <60m) at rates of 15 or 20 m . min -1 . Net tow duration was ≤5 min. Half of the tow contents was filtered through a 150µm sieve and rinsed into 113.4 g jars filled with 95% ethanol (Bucklin et al., 2001). Twelve to twenty-four hours after collection, samples were filtered into fresh ethanol to limit the diluting effects of organism dehydration.
At each station, temperature, salinity, pressure, bottom depth, and fluorescence were collected using a CTD profiler (Seabird 19plus). Extracted chlorophyll a was measured at the surface and chlorophyll maximum using standard fluorometric techniques (Lomas et al. 2012), and ice coverage was recorded and later retrieved from the ship logs.

Morphometrics, DNA extraction, amplification, and sequencing:
Adult female Pseudocalanus were isolated at random from the zooplankton samples (72-100 station -1 ). Each copepod was rinsed in fresh 95% ethanol and photographed. The urosome length, width, depth, and the prosome length were measured from photographs with ImageJ (Abramoff et al., 2004). DNA was extracted using the 96 well plate DNeasy Tissue Kit (Qiagen Inc.) following manufacturer's instruction with an overnight lysis incubation. The 18S rDNA was amplified using eukaryotic universal primers (Medlin et al., 1988) and thermocycling protocol of 94˚C for 30 s, followed by 35 cycles of 94˚C for 30 s, 60˚C for 1 min, and 72˚C for 2 min, followed by 72˚C for 10 min. The mitochondrial COI gene was amplified using the forward primer L1384-CO1 (Machida et al., 2004) and the reverse primer HCO 2198 (Folmer et al., 1994), with thermocycling protocol of 94˚C for 1 min, 35 cycles of 94˚C for 30 s, 49˚C for 30 s, and 72˚C for 50 s, followed by 72˚C for 5 min. Finally, the nITS region (ITS1-5.8S-ITS2) was amplified with the universal primers ITS-4 and ITS-5 (White et al., 1990) under the following thermocycling conditions: 95˚C for 1 min, 35 cycles of 95˚C for 30 s, 50˚C for 30 s, and 72˚C for 50 s, followed by 72˚C for 7 min. PCR amplification was confirmed on a 1% agarose gel with SYBR Green (Applied Biosystems) staining and successful amplicons were purified with ethanol precipitation (APPENDIX B).
Sequencing for the 18S gene and nITS region was performed using the universal eukaryotic primers for single stranded sequencing (Medlin et al., 1988) and universal ITS4 and ITS5 primers (White et al., 1990) respectively. The mtCO1 sequencing was performed in both directions for 62 copepods. The remaining copepods were sequenced only in the forward direction due to full gene coverage with the L1384-CO1 primer (Machida et al. 2004). All sequencing was performed at the Rhode Island Genomics and Sequencing Center with a 3130X1 (ABI Inc.) sequencer.
The mtCO1 gene was sequenced for a total of 642 copepods. Thirty copepods were sequenced to test divergence for the 18S gene, and the nITS region was sequenced for 40 copepods. Copepods were selected for nuclear sequencing to compare with unique sequence groups found with the mtCO1 gene.
DNA polymorphism was examined for the mtCO1 gene by calculating nucleotide diversity (π), haplotype diversity (H), and Tajima's D in DnaSP v5.10.01 (Librado and Rozas, 2009). Statistical significance of Tajima's D was tested by permutations with 1000 replicates. Population structure of individual species was determined across the shelf with a one-way analysis of molecular variance (AMOVA) and 1000 permutations in ARLEQUIN v3.5.1.3 (Excoffier and Lischer, 2010), where populations were defined as sampling stations and groups were defined as middle shelf domain (stations <100 m deep) and offshore domain (stations ≥200 m deep) (Weir and Cockerham 1984, Excoffier et al. 1992, Weir 1996. Thus, for analysis using F-statistics, F ST is defined as the genetic variation within a station relative to the total genetic variation, F SC is the genetic variation within a station relative to the variation between the middle domain and offshore domain, and F CT is the genetic variation between the middle domain and offshore domain relative to the total genetic variation.

Morphometric Analysis
Multivariate analysis of variance (MANOVA) and quadratic discrimination function analyses (QDFA) statistics were used with mtCO1 species assignments to test the effectiveness of size measurements in species identification. Measurements of urosome length, width, depth (EL, EW, ED), and prosome length (PL) were transformed (natural log) and analyzed as prediction variables. Wilks' Lambda was used as a test statistic to evaluate the MANOVA and QDFA analyses. Pillai's trace was additionally included given its reliability with relatively small sample sizes (Olson, 1976). Jacknife values were recorded from the QDFA to estimate variance.

Environmental Correlation
Resemblance matrices for both species abundance data (Bray-Curtis Similarity) and environmental data (Normalized Euclidean Distance for salinity, temperature, depth of the chlorophyll maximum, total chlorophyll, and dissolved oxygen) were created with PRIMER v6 software Ainsworth 1993, Clarke andGorley 2006). Temperature, dissolved oxygen, chlorophyll and salinity were taken from the depth of the chlorophyll maximum where Pseudocalanus spp. may be feeding and the surface layer. Multidimensional Scaling (MDS) plots created with a Bray-Curtis Similarity matrix and overlaid clusters were created for station species composition to visualize station groupings. These groupings were tested using ANOSIM with all possible permutations (210). Abundance data and environmental patterns were compared with the global BEST test (Clarke and Gorley 2006). Within the BEST application, the BIO-ENV procedure was run using the Spearman coefficient (Kendall 1970) as a rank correlation method to search all possible combinations of variables. The BIO-ENV procedure was run for 999 permutations to test statistical significance of resulting correlations.

Interspecific Diversity
Groupings in mtCO1 gene trees generated from UPGMA, Neighbor-Joining, and Maximum Likelihood methods for 493 bp sequences of 632 individuals indicated that four Pseudocalanus spp. were present in the Bering Sea: the temperate P.
newmani and P. mimus, and the arctic P. minutus and P. acuspes (Figure 1). High bootstrap values (>70%) confirmed species groupings. Divergence in the mtCO1 gene among the four Pseudocalanus spp. ranged from 10.7% between P. acuspes and P. mimus to 18.7% between P. minutus and P. newmani. The average overall divergence among all four species was 14.7% ± 4%.
Nuclear divergence was more limited in both the 18S gene and nITS region than in the mtCO1 gene. The 18S rDNA of 62 individuals representing all four mtCO1-defined species were identical along the 1566 bp sequence. Gene trees (UPGMA, NJ, and ML) from the 727 bp length sequences covering ITS1-5.8S-ITS2 of 33 individuals representing four species and divergent haplotypes showed limited divergence with no defined groups or agreement with the mtCO1 sequence groupings ( Figure 2). The average overall divergence in the nITS region was 1.3% ± 1.1%.

Morphometric Analysis
Four hundred and eighty two individuals, with representatives from every station, were photographed and measured. Pillai's trace results showed a significant difference in size (prosome and urosome length, width and depth) between species for both MANOVA and QDFA models with p < 0.001 for both models. The classification success of species based on these measurements was 68% for Pseudocalanus acuspes, 43% for P. mimus, 40% for P. minutus, and 25% for P.
newmani. When the species data was combined, overall success of classification was 43%.

Geographic Distribution
Multiple species (2-4) coexisted at all ten stations in the Bering Sea ( Figure 3). P. minutus dominated the offshore stations and P. acuspes dominated the middle domain stations. Limited numbers of P. mimus were found, but those sampled appeared primarily in the offshore regions. P. newmani was fairly well mixed among all stations. Absolute abundance of adult female Pseudocalanus for TN249 stations agreed with patterns observed in relative abundance data (Figure 4). Stations in the middle domain (50-100 m) had higher total abundance than stations in the offshore domain (≥200 m depth). Due to differences in sampling methods between cruises, abundance data for the PSEA1001 cruise were not used in comparison between stations.
A multidimensional scaling plot (MDS) of relative species composition data showed two groupings of 40% similarity ( Figure 5). One group consisted of stations all located on the middle shelf, in water less than 100 m deep. The other grouping consists of offshore stations (≥200 m depth) and station 156 based on relative species composition. Offshore stations were significantly different from middle shelf stations (ANOSIM, p = 0.01).
One anomaly found in spatial distribution was at station 156, where 2% of all sampled Pseudocalanus spp. were infected with the parasite Hemiarthrus abdominalis. This parasite ranged in size from 300-700 µm in length and was found loosely attached to the dorsal side of live and preserved Pseudocalanus spp.

Intraspecific Diversity
A total of 137 unique mtCO1 haplotypes were identified from 632 individuals.
Among these, 42 haplotypes were classified as haplotypes of P. newmani, 25 as P. minutus, 21 as P. mimus, and 40 as P. acuspes. Nucleotide and haplotype diversity varied among all four species with the two temperate species, P. newmani and P. mimus, having the highest diversity values (Table 2a). All Tajima's D values departed from the standard neutral model (Tajima, 1989). In P. newmani and P. mimus, Tajima's D was negative but not significant (0.10 > p > 0.05). For P. minutus and P.
acuspes, Tajima's D values were negative and significant (p < 0.05), suggesting population expansion in these two species. Results from the AMOVA test showed a greater percentage of variation within populations than between groups, where populations are defined as stations and groups are defined as the middle domain and offshore (Table 2b). The fixation indices, F SC and F ST , were not significant for all four species (p > 0.05), indicating that there is no genetic structure among populations within groups or among populations and among groups. The F CT was small, but significant for P. acuspes and P. minutus. This indicates that a weak genetic structure exists across the shelf in arctic species.

Correlation of Environmental Data with Species Distribution
Environmental data was compiled for the TN249 and PSEA1001 cruises (Table 3). The BIO-ENV procedure showed that optimal correlation (rho) between environmental parameters and species composition resulted from a combination of depth of the chlorophyll maximum and surface salinity. This pairing of environmental parameters had a statistically significant correlation with species distribution (rho = 0.576). The addition of other environmental variables including surface temperature and surface dissolved oxygen to this combination caused an immediate drop in correlation values (Table 4). Total chlorophyll (measured at the chlorophyll maximum) and % ice coverage were not in the top five correlation combinations.
Temperature, salinity, and dissolved oxygen measured at the chlorophyll maximum depth did not register in the top 10 BIO-ENV combinations.

Discussion:
Mitochondrial CO1 sequencing conducted across the shelf demonstrated that four Pseudocalanus spp., the arctic P. acuspes and P. minutus, and the temperate P. mimus and P. newmani, are present in the SE Bering Sea. This array of species matches the morphometric findings of Frost (1989) and is the greatest number of coexisting Pseudocalanus species in any region that has been confirmed with the use of genetic tools. Nuclear divergence of 0% in the 18S gene and 1.3% ± 1.1% in the nITS region suggests that species divergence may be more recent for Pseudocalanus spp. than other copepods. While the 18S gene is generally only useful in identifying between copepod genera (Bucklin et al. 2003), the divergence found in the nITS region of the Pseudocalanus was significantly lower than other copepods like Acartia tonsa, which have a divergence in the nITS region of 39% between clades (Chen and Hare, 2008). An average nITS divergence of 1.3% ± 1.1% for all Pseudocalanus spp.
indicates this is not a useful region for species identification. The mitochondrial CO1 gene was the only marker successful in identifying species.
With the use of the mtCO1 gene, each station was shown to consist of 2-4 coexisting species. The distribution and abundance of these species varied widely across the shelf. There was a greater overall Pseudocalanus abundance in stations in the middle domain than in those beyond the offshore domain, which is consistent with the findings of previous studies that observed spring zooplankton distributions (e.g. Cooney and Coyle, 1982;Smith and Vidal, 1986). At the species level, P. acuspes and P. minutus were dominant but differed in their distribution across the shelf. P.
acuspes was almost exclusively found in the middle domain, and although P. minutus was found at almost all stations, it makes up the largest portion of the species complex in the offshore stations. While P. mimus was also found primarily in offshore stations, the low number of individuals collected does not allow for definitive conclusions about their distribution. This species is thought to have a greater abundance on the shelf in the Gulf of Alaska based on detailed morphometric analyses (Frost, 1989;Napp et al., 2005) and may be present in the Bering Sea as a result of northward current flow across the shelf break. In contrast, P. newmani was fairly evenly distributed throughout the middle domain and offshore stations. Because this study only covered one spring season, it is unclear what seasonal and interannual variations of species presence and distribution patterns exist. The year of 2010 was classified as a " cold year" (Stabeno et al., 2012), which could explain the dominance of the primarily arctic species. These are the only two species that have been shown to store large amounts of lipids during the winter months (McLaren et al., 1989) and thus have been suggested to respond in a similar manner to larger calanoid copepods like Calanus marshallae/glacialis, which increase in abundance during cooler years .
The differences in species composition between the middle domain and offshore stations were statistically supported in similarity groupings. The one exception to the trends in both abundance and distribution was station T156, a middle domain station, which grouped at 40% similarity with offshore stations. Station P3, a station in the same area as T156 but sampled in March rather than June, was grouped with the middle domain stations. One potential explanation for this difference in grouping could be the effect of the parasite Hemiarthrus abdominalis, which was only found at station T156. While this parasite was only loosely attached to the Pseudocalanus urosome and did not appear to be feeding directly on the copepod, the size of the parasite (446 µm ± 72 µm) may have impacted swimming behavior and/or location in the water column, thus impacting feeding and survival. Any preference of Pseudocalanus host species could be a biotic cause of differences in species composition at this location.
Although biotic factors such as parasitic infection, predation and prey availability all have the potential to greatly impact species distribution, the influence of environmental abiotic factors such as temperature and salinity were the focus of this study. A correlation between species composition and the pairing of salinity and the depth of the chlorophyll maximum indicates that the environment may play a key role in determining distribution patterns. Though the relationship to the depth of the chlorophyll maximum is complex and most likely influenced by a combination of environmental conditions, salinity has been previously suggested to influence zooplankton distribution in the Bering Sea. A study conducted by Smith and Vidal (1984) concluded that salinity and lack of advection across the shelf, not temperature, were the critical factors in zooplankton community structure. There are two possible hypotheses that could explain the link between salinity and Pseudocalanus spp.
structure. The first is that the environment is selecting for certain species. Stations that had higher salinity (31-32.5 psu) were composed of greater percentages of P. minutus and P. newmani, while stations with lower salinity (30-31 psu), were dominated by P. acuspes. If the species have different critical salinity ranges, this could explain their differential distribution. In the Baltic Sea, a decreasing success in spring maturation and reproduction in P. elongatus has been linked to decreasing salinity from increased precipitation (Möllmann et al., 2003). A similar relationship may exist between salinity and the four species of the Bering Sea. An alternative hypothesis to explain this link is that salinity is an indicator of water mass, and species composition will be more similar in water masses that are linked by current flow. If this were the case, however, we would expect to see a similar correlation with temperature, which is not found. Although temperature and salinity are conservative at depth, they are subject to evaporation and air heat exchange at the surface where we observed the strongest correlation with species composition. It appears unlikely that the correlation with salinity is entirely related to the water mass. In order to more thoroughly test this hypothesis, data would need to be collected on the source populations of each of these species and how the currents distribute them across the shelf. This would require a species-specific look at developmental stages as well as populations outside of the Bering Sea. With either hypothesis, it does appear that salinity is more closely correlated with species distribution than temperature and other abiotic factors.
The correlation found between salinity and species composition was not observed at the subspecies level. There was, however, a cross-shelf difference in the haplotypic distribution of the Arctic species, P. minutus and P. acuspes. A weak but significant genetic structure between the middle domain stations and the offshore stations demonstrated that there may currently or previously have been some environmental (e.g. salinity) or physical barriers to genetic exchange (e.g., Goetze, 2005;Chen and Hare, 2008;Chen and Hare, 2011;Aarbakke et al., 2011). In combination with negative Tajima's D values, which signify population expansion (Tajima, 1989), it is possible that this structure is a relict of the last glacial period when flow through the Bering Strait was extremely limited. In the last glacial cycle, the Strait was nearly closed at multiple times, with the last closure ~70kyr BP (Hu et al., 2010). One explanation for the weak structure and expanding population of the two Arctic species is that populations of these species were previously isolated, on either side of the Bering Strait, and have only recently (within the last 70kyr BP) been able to exchange genes. This may additionally be the cause of the lower haplotypic and nucleic diversity observed in the arctic relative to the temperate species. The temperate species, P. newmani and P. mimus, which do not have significant structure or Tajima's D, would not have been subject to this isolation given their already more southern distribution. Regardless of whether structure in the Bering Sea is a relict from the last glacial period or represents current conditions, it is clear that there is some genetic exchange within the ecosystem. Existing structure is weak to non-existent in all four species and both haplotypic and nucleic diversity are high relative to Pseudocalanus spp. from other ecosystems. In the Baltic Sea, haplotypic diversity of P. elongatus was 0.024, dramatically lower than even that measured in P. acuspes (H = 0.436). The Bering Sea is clearly a mixing zone for different populations of both arctic and temperate species of Pseudocalanus. The high haplotypic diversity could be a result of the mixing of distinct lineages from different ocean basins (Pacific/Atlantic) in this region and/or the recent exchange of genes following a period of isolation (e.g. Bering Strait closure).
With a better understanding of species distribution and environmental impacts, more accurate modeling may be able to identify the importance of individual species as a prey for larval fish. As the southeastern Bering Sea has previously demonstrated a rapid response to atmospheric perturbations (Napp and Hunt, 2001) and an increase in the frequency of these warm-water years is expected (Coyle and Pinchuck, 2002), studies indicate that Pseudocalanus copepods may become an even more important prey item in the Bering Sea ecosystem (Coyle et al., 2008;Coyle et al., 2011). Future research in the Bering Sea will be able to look at seasonal as well as interannual changes in the distribution patterns established in this study to better predict the impacts of changes in the climate. A better grasp of species-specific abundance and distribution will further assist in the accuracy of fisheries modeling. Pollock have a "higher energy density" in their first winter if they are able to consume mesoplankton with higher lipid levels in a previous cool summer Hunt et al., 2011). If Pseudocalanus spp. with lower lipid content (P. mimus and P. newmani) become the primary Bering Sea species during warm summers, Pollock and other foraging fish species may be negatively impacted and have a lower survival rate the following winter. This study has established the species composition and distribution for a 'cool' spring season as well as possible environmental forcing factors and thus provides a baseline for future studies to investigate seasonal and temporal changes within the Bering Sea.
On a global scale, this research highlights the importance of investigating the role of Pseudocalanus spp. in an ecosystem, from a species-specific perspective. The Bering Sea is only one of several ecosystems that house multiple species within this genus. If salinity is an important determinant in habitat and niche partitioning, changes to the environment might cause a redistribution that would not be seen if the genus were investigated as a whole. Large shifts in the biogeographic distribution of zooplankton have already been observed since the 1960s as a result of climate change (e.g., Beaugrand et al. 2009). It is possible that shifts will occur within a genus allowing some species to move into new territories and others to retreat northward.
Recent advances in molecular protocols have allowed for a more detailed investigation into the present and future of this widespread and ecologically important copepod genus.       Table 4: The top 5 combinations of abiotic factors that most strongly correlate (Spearman's rank correlation coefficient) with the species composition patterns in the SE Bering Sea. An "X" designates which environmental factors were combined to achieve the corresponding correlation value. Abiotic values measured at the chlorophyll maximum depth are not included due to their absence in the BIO-ENV results.