Diversity and distribution of nearshore barnacle cyprids in southern California through the 2015–16 El Niño

Sample collections

All field collections were made under the California Department of Fish and Wildlife scientific collection permit SC-10523. Plankton samples were collected offshore of Bird Rock, La Jolla, California, USA and included five sampling stations at distances of 280, 460, 640, 820, and 1,000 m from shore at water depths of 4, 6, 8, 10, and 12 m, respectively (Fig. 1A). This study site was selected because of the large population of C. fissus barnacles on its shoreline, and because larvae of the predominant barnacle species in the area have been collected in surrounding waters (C. fissus, B. glandula, Pollicipes polymerus; Hoffman, 1989; Pineda, 1991, 1999; Pineda & López, 2002; Tapia & Pineda, 2007; Tapia et al., 2010). Plankton sampling was conducted on a 7.6 m boat, and consisted of 33 collection cruises: nine in spring-summer (May–July) and six in fall-winter (October–December) of 2014, eight in spring-summer (April–July) and four in fall-winter (October–December) of 2015, and six in summer (June–August) of 2016. All cruises were conducted during the day without consideration for lunar phase or tidal cycle. An Ebara 50DWXU6.4S Dominator submersible semivortex pump was used to filter two m³ of seawater in discrete two m depth bins from the surface to the bottom at each of the five sampling stations during each cruise. The pump was oscillated throughout the two m sampling intervals so that the entire water column was sampled at all stations. Each cruise resulted in 20 plankton samples unless sea state conditions or equipment failure interrupted sampling. Seawater was filtered through a 118-μm mesh net to capture barnacle cyprids, and plankton collected on the filter were immediately preserved in 100% ethanol.

Barnacle cyprids were enumerated and identified to species based on published morphological descriptions (Lewis, 1975; Brown & Roughgarden, 1985; Miller et al., 1989; Miller & Roughgarden, 1994; Shanks, 2001) using a dissecting microscope. Additionally, morphologically distinct cyprids were selected for DNA barcoding, preserved individually in ethanol, and refrigerated until molecular analysis was carried out. DNA barcoding was completed to ground truth morphological identifications and classify cyprids with unfamiliar morphology that we were otherwise unable to match to known species. After barcoding, morphological identifications were assessed, and cyprids that were unknown prior to barcoding were compared to reference sequences.

Morphologically distinct adult barnacles were collected from boulders in the intertidal at three locations to serve as a baseline for barnacle species in the area: onshore of the plankton sampling location at Bird Rock, La Jolla, California, USA; north of sampling location at Dike Rock, La Jolla, California, USA; south of sampling location at Alisitos, Baja California, Mexico (Fig. 1B). Collections were conducted on May 7, 2016 at Bird Rock, May 9, 2016 and June 20, 2016 at Dike Rock, and June 22, 2016 at Alisitos. The goal was to find ∼5 representative individuals for each different species, which resulted in 55 barnacles total. All adults were dissected; a small piece of tissue was removed, preserved in ethanol, and refrigerated until molecular analysis was carried out.

Molecular analysis

Molecular analysis methods were slightly modified from the methods outlined by Govindarajan et al. (2015). For cyprids, the whole larva was placed into PCR tubes. For adults, Qiagen DNEasy Blood & Tissue extraction kits were used to extract DNA from the small piece of tissue removed during dissection. The mitochondrial COI gene was amplified using HCO2198 and LCO1490 primers (Folmer et al., 1994) or jgHCO2198 (Geller et al., 2013) and mlCOIintF primers (Leray et al., 2013). The PCR cycle for Folmer primer amplifications was: (a) an initial denaturation at 95 °C for 3 min; (b) 35 cycles of 95 °C 30 s, 48 °C 30 s, 72 °C 60 s; and (c) a final extension at 72 °C for 5 min. A touchdown protocol was used for amplifications with the Geller–Leray primers (Leray et al., 2013). PCR products were visualized with a 1% agarose gel stained with Biotium GelRed, and purified with Qiagen Qiaquick PCR Purification kits. After purification, DNA was qualified and quantified with a NanoDrop One Microvolume UV-Vis spectrophotometer and sequenced in both directions (Eurofins Scientific, Brussels, Belgium). Geneious 9.1.7 (Biomatters Limited, Auckland, New Zealand) was used to assemble and align chromatograms, and sequences were compared to reference sequences in GenBank using Basic Local Alignment Search Tool, accessed on 7 June 2017. All sequences generated in this study were submitted to GenBank (MK496547–MK496616).

Diversity measurements

Species diversity was assessed for each cruise using Hurlbert’s (1971) modification of Sanders (1968) rarefaction index: $$E\left(S_n\right)={\displaystyle \sum _{i=1}^s\left[1-\frac{\left(\begin{array}{c}N-N_i\\ n\end{array}\right)}{\left(\begin{array}{c}N\\ n\end{array}\right)}\right]}$$ where E(S_n) is the expected number of species in a sample of n individuals selected at random from a collection containing N individuals, S species, and N_i individuals in the ith species (Hurlbert, 1971). Rarefaction is used to compare sample collections of different sizes by reducing them to a common size (n), since species richness usually increases with N (Hurlbert, 1971). Mean E(S₂₀) was calculated for each sampling period, and a randomization test (Manly, 1997) was completed to determine if differences in E(S_n) values were significant.

Cyprid distributions

Cyprid concentration (no. m⁻³) at each station was multiplied by station depth to standardize samples to cyprid density (no. m⁻²) to account for vertical dilution of larvae when comparing cyprid concentrations at stations of varying depths. Cyprid density across the five sampling stations for each species was averaged for each cruise. The mean distance from shore (MDS; Hagerty, Reyns & Pineda, 2018; referred to as “average distance offshore” in Shanks & Shearman, 2009) of each cyprid species for each cruise was calculated as: $$\text{MDS}=\frac{\sum \left(\text{no.\hspace{0.17em}}{\text{m}}^{-2}\text{\hspace{0.17em}at\hspace{0.17em}a\hspace{0.17em}given\hspace{0.17em}station\hspace{0.17em}}\times \text{\hspace{0.17em}Distance\hspace{0.17em}from\hspace{0.17em}shore\hspace{0.17em}of\hspace{0.17em}station}\right)}{\sum \left(\text{no.\hspace{0.17em}}{\text{m}}^{-2}\text{\hspace{0.17em}at\hspace{0.17em}all\hspace{0.17em}stations}\right)}$$

A one-way ANOVA was used to test for differences in MDS among cyprid species; ANOVA assumptions were met.

Environmental variables

Temperature and depth measurements were recorded before and after sampling at each station with a SonTek CastAway-CTD. Mean temperature for each cruise was calculated as the average of all CTD temperature measurements, and thermal stratification was calculated for each station. Analysis of 1,251 CTD profiles indicated that thermal stratification is an excellent proxy for density stratification (top-to-bottom temperature difference is correlated with top-to-bottom density difference: r = 0.987). Furthermore, given the relatively widespread deployment of temperature loggers on moorings, our calculations allow comparisons of thermal stratification to be made with other studies that have collected temperature profile data (Hagerty, Reyns & Pineda, 2018; Pineda, Reyns & Lentz, 2018). Salinity was not considered because stratification at our study location is mainly dependent on vertical differences in temperature (Hagerty, 2017).

Thermal stratification was defined as change in temperature m⁻¹ (ΔT °C m⁻¹), as in Hagerty, Reyns & Pineda, 2018, and calculated as: $$\text{Thermal stratification}=\frac{\left(\text{Temperature at surface}-\text{Temperature at bottom}\right)}{\left(\text{Depth of bottom temperature}-\text{Depth of surface temperature}\right)}$$

Zonal stratification (Hagerty, Reyns & Pineda, 2018) was calculated by averaging thermal stratification at the three offshore stations–640, 820, and 1,000 m from shore with bottom depths of 8, 10, and 12 m, respectively—for each cruise. These stations were chosen to represent zonal thermal stratification because they were frequently associated with the highest stratification at our study site.

Cyprid abundance and species diversity

Molecular analysis resulted in 70 successful sequences: 26 cyprids and 44 adults. In all cases, morphological identifications were affirmed by the barcoding results. We found that the Geller–Leray primer set (Geller et al., 2013; Leray et al., 2013) worked better for our samples than the Folmer set (Folmer et al., 1994). During the first round of barcoding, we started with the Folmer set and switched to the Geller–Leray set after achieving little success. A contributing factor may have been contamination with rose Bengal, as some near-bottom samples with heavy sediment content had been stained, and we discovered during our first round of barcoding that DNA from stained cyprids could not be amplified with our procedure (0 of 16 stained cyprids were successfully barcoded). Easton & Thistle (2014) successfully amplified and sequenced DNA from 122 copepods stained with rose Bengal using a Chelex®-based procedure, which is thought to remove any inhibitory compounds from samples. Of the 27 cyprids that did not come from rose Bengal samples during the first round of barcoding, Folmer primers were used for 13 and Geller–Leray primers were used for the remaining 14. DNA was successfully amplified and sequenced for five cyprids using the Geller–Leray method (36%), but only two cyprids were successfully barcoded using the Folmer method (15%). The Folmer primers were not used during the second round of barcoding, when 19 barcodes were obtained from 50 attempted cyprids (38% success rate using only Geller–Leray primer set). Table 1 outlines how many cyprids were barcoded from each sampling period.

Six species of cyprids were identified, belonging to five genera and four families: C. fissus, B. glandula, P. polymerus, Tetraclita rubescens, Megabalanus rosa, and B. trigonus (Figs. 2A–2F). Chthamalus spp., B. glandula, and P. polymerus cyprids were identified morphologically, and were confirmed by molecular analysis (Table 2). C. fissus cyprids were given a tentative morphological identification of Chthamalus spp., because larvae of C. fissus and C. dalli are morphologically identical (Miller et al., 1989). However, C. fissus is the most abundant local barnacle and C. dalli is rare in San Diego, California (Morris, Abbott & Haderlie, 1980). Indeed, all of our barcoded Chthamalus cyprids and adults were identified as C. fissus. Three additional cyprid species, T. rubescens, M. rosa, and B. trigonus, were identified molecularly (Table 2). There were four distinct morphotypes of unknown cyprids, which were not successfully sequenced due to a failure to amplify their DNA, which could be due to low quality or insufficient quantity of DNA (Figs. 2G–2J). Six species of adults were identified, matching the identifications of the known cyprid species other than B. trigonus, and one additional species. Three adults—one collected at Dike Rock and two collected in Alisitos, Baja California—matched an unclassified Sessilia, Balanoidea sp. MM-2014 (GenBank: HG970516.1), collected from a soft coral in the genus Leptogorgia in Florida, USA (Malay & Michonneau, 2014).

Cyprid abundance varied by species, but was highest for most species in spring-summer 2014 and summer 2016 (Figs. 3C–3E). P. polymerus cyprids had relatively similar densities throughout all sampling periods until a large peak in density in summer 2016 (Fig. 3C). In contrast, M. rosa had higher relative abundances during spring-summer periods than during fall-winter periods (Fig. 3E). C. fissus was by far the most abundant cyprid found in samples (maximum density = 346 cyprids m⁻²), followed by P. polymerus (max. = 122 m⁻²) while T. rubescens (max. = 30 m⁻²), B. trigonus (max. = 28 m⁻²), M. rosa (max. = 12 m⁻²), and B. glandula (max. = 12 m⁻²) were found at low densities (Fig. 4). Cyprids were present in relatively large numbers starting in 2014 (average density = 121 cyprids m⁻² in spring-summer 2014), decreased in fall-winter 2014 and remained low throughout 2015 (average density did not exceed 38 cyprids m⁻² during these sampling periods), and increased again in summer 2016 (average density = 152 cyprids m⁻²). The highest density for the dominant species, C. fissus, and the two least abundant species, B. glandula and M. rosa, were collected during spring-summer 2014. The highest densities of P. polymerus, T. rubescens, and B. trigonus occurred in summer 2016 (Fig. 3).

All six identified cyprid species were collected during spring-summer of 2014 and 2015, as well as during the 2016 sampling period, but B. trigonus cyprids were not found in fall-winter of 2014, and T. rubescens cyprids were not collected in fall-winter of 2015 (Figs. 3 and 4). Species assemblage varied by sampling period, but C. fissus cyprids dominated the samples during all five periods (Fig. 4). C. fissus cyprids were collected during all cruises except December 29, 2014 and November 20, 2015, when no cyprids were found in our samples. The relative abundance of C. fissus was noticeably lower compared to other species during spring-summer 2015 and summer 2016, mainly due to an increase of M. rosa and P. polymerus cyprids (Fig. 4). Unknown cyprid species were not collected during fall-winter cruises in 2014 or 2015 (Fig. 4).

Species diversity was highest on May 14, 2014 (E(S₂₀) = 5.92) and lowest on the 6th and 18th of November in 2014, and on the 16th and 20th of July in 2015, when C. fissus was the only species present (Fig. 3B). E(S₂₀) measurements were more variable in spring-summer 2014 and spring-summer 2015 (ranging 1–5.92) than during the other sampling periods (Fig. 3B). Mean E(S₂₀) was comparable for most sampling periods, except in fall-winter 2014, when it was low (Fig. 5).

We conducted a randomization analysis (Manly, 1997) to test for differences in E(S₂₀) means among the five sampling periods: spring-summer 2014, fall-winter 2014, spring-summer 2015, fall-winter 2015, and summer 2016. Random samples of n_ki values were taken without replacement from the vector comprising the entire series of E(S₂₀) (N = 23), where k is sample period, and i is the number of E(S₂₀) values in each sample period. That is, the 23 E(S₂₀) values were randomly redistributed among the five sampling period groups. We then calculated the F one-factor ANOVA statistic, and this was repeated 10,000 times. The approximate p-value is the proportion of F values from the randomized data sets that are larger than the F computed from our data set; p = 0.3932. This indicates there were no significant differences in E(S₂₀) among the sampling periods.

Cyprid horizontal distributions and water column conditions

There were significant differences in MDS among cyprid species (ANOVA results: F = 15.13, p = 0.00016, df = 1), and a Tukey post hoc test was completed for pairwise comparisons. M. rosa and B. glandula had significantly greater MDS values than C. fissus, meaning M. rosa, and B. glandula cyprids were distributed significantly farther from shore than C. fissus (q_s = 225.17, p = 0.000043 and q_s = 172.45, p = 0.0072, respectively). M. rosa was also significantly farther offshore than P. polymerus (q_s = 168.51, p = 0.010). Average MDS for C. fissus cyprids was greater during fall-winter than spring-summer and summer sampling periods (Table 3), which means C. fissus cyprids were distributed closer to shore in spring-summer. In contrast, P. polymerus cyprids were closer to shore on average during fall-winter sampling periods (Table 3). While too few M. rosa cyprids were collected in fall-winter of 2015 to make comparisons, the average MDS value for M. rosa in fall-winter 2014 was lower than that of the spring-summer periods of 2014 and 2015, as well as that of summer 2016. Another notable result was that T. rubescens cyprids were located closer to shore in summer 2016 compared to earlier sampling periods (Table 3).

In general, water temperature was relatively low at the start of spring-summer sampling periods, increased throughout spring-summer, remained high until fall, and then decreased during fall-winter sampling periods (Fig. 3A). The highest temperatures were recorded during fall-winter 2015, and the lowest measurement was during summer 2016. Thermal stratification was high during spring-summer sampling, low during fall-winter sampling, and the highest measurements occurred during the first few cruises of summer 2016 (Fig. 3A). Overall mean temperature for each sampling period ranged from 17.7 to 20.6 °C, and was higher in fall-winter than spring-summer (Fig. 5B). Mean temperature was similar during both sampling periods in 2014 and summer 2016, ranging between 18.7 and 19.0 °C, but decreased to 17.7 °C in spring-summer 2015 and increased to 20.6 °C in fall-winter 2015. Mean zonal stratification ranged from 0.05 to 0.38 °C m⁻¹ and was lowest in fall-winter sampling periods (Fig. 5C). Mean zonal stratification was highest in summer 2016, and was more variable in spring-summer/summer sampling periods than fall-winter sampling periods (Figs. 3A and 5C).