Cryptic phenotypic plasticity in populations of the North American freshwater gastropod, Pleurocera semicarinata

Abstract

Background

Phenotypic plasticity is termed ‘cryptic’ when it becomes so extreme as to prompt an (erroneous) hypothesis of speciation. Populations of pleurocerid snails nominally identified as Pleurocera (‘Goniobasis’ or ‘Elimia’) livescens are common elements of the macrobenthos in rivers and on lakeshores from Vermont through northern Ohio and northern Indiana to Wisconsin, USA. Small streams in southern Ohio, southern Indiana, and Kentucky are inhabited by Pleurocera semicarinata, and larger rivers by Lithasia obovata. The three nominal species differ only in qualitative aspects of shell morphology - P. semicarinata demonstrating a slender shell with a small body whorl, L. obovata a robust shell with a large body whorl, and P. livescens intermediate.

Results

Here I use allozyme electrophoresis to estimate genetic divergence over 11 enzyme loci among three populations of P. livescens, two populations of P. semicarinata, and one population of L. obovata sampled across 650 km of their combined range. With each of these six populations I sample a control population of Pleurocera canaliculata, demonstrating that the genetic divergence among the six known conspecific populations is comparable to that observed among the six study populations.

Conclusions

The specific nomina L. obovata and P. livescens would appear to be junior synonyms of P. semicarinata. The shell morphological differences by which these taxa have heretofore been distinguished would appear to result from ecophenotypic plasticity, driven perhaps by predation, substrate, or current.

Background

In recent years, a great deal of research interest has been directed toward freshwater gastropods as model organisms for the study of phenotypic plasticity (e.g., Auld and Relyea [2011]; Brönmark et al. [2012]). Most research has focused on the basommatophorans, which by virtue of their rapid generation times and ease of culture adapt well to experimentation. Reared in the presence of crushing predators, basommatophoran snails tend to develop shells with broad body whorls (Lakowitz et al. [2008]; Brönmark et al. [2011]), while in the presence of predators extracting snails from their apertures, body whorls tend to develop more narrowly (DeWitt [1998]; Langerhans and DeWitt [2002]; Hoverman et al. [2005]). The phenotypic effects of temperature, current, and environmental calcium on shell morphology have also been well documented in laboratory populations of basommatophoran snails (Britton and McMahon [2004]; Rundle et al. [2004]; Dillon and Herman [2009]).

Research on phenotypic plasticity in the longer-lived, more slowly growing caenogastropod snails has tended to focus on systematic observations in the field, typically calibrated with independent genetic controls. Dillon ([2011]) documented a striking tendency for the shells of Pleurocera clavaeformis to develop a broader and more ‘robust’ morphology as a function of river size in East Tennessee, suggesting strongly that the shell characters upon which the distinction between the genera Pleurocera, Goniobasis, and Elimia have been based are attributable to ecophenotypic plasticity. Dillon et al. ([2013]) used geometric morphometrics to document an identical phenomenon in populations of Pleurocera canaliculata sampled in streams of varying size from New York to Alabama, USA. Both of these studies involved populations previously assigned to multiple specific nomina by 19th-century taxonomists, and both were controlled using genetic similarities calculated over approximately ten polymorphic enzyme loci. Dillon and colleagues coined the term ‘cryptic phenotypic plasticity’ to describe the phenomenon demonstrated by these diverse pleurocerid populations - interpopulation morphological variance so extreme as to prompt an (erroneous) hypothesis of speciation.

The purpose of the present study is to extend the findings of Dillon ([2011]) and Dillon et al. ([2013]) to include freshwater gastropod populations currently identified as Pleurocera semicarinata, Pleurocera livescens, and Lithasia obovata. ‘Melania’ (now Pleurocera) semicarinata was described from ‘Kentucky’ by Thomas Say in 1829, its range now generally understood to include small streams through most of Kentucky, southern Ohio and southern Indiana (Goodrich [1940]). ‘Melania’ (now Lithasia) obovata was also described by Say ([1829]), just a few paragraphs subsequent to his description of M. semicarinata, from ‘Kentucky River,’ and ‘also at the Falls of the Ohio.’ Its range is similar to that of P. semicarinata, although restricted to large rivers (Goodrich [1940]). ‘Melania’ (now Pleurocera) livescens was described from ‘Lake Erie, New York’ by Menke ([1830]). Populations of P. livescens inhabit streams, rivers, and lake shores from Wisconsin through northern Indiana and northern Ohio to Ontario and Vermont (Goodrich [1939], [1940]; Clarke [1981]).

Pleurocera semicarinata and P. livescens are distinguishable only by qualitative characters of the shell (Figure 1), populations of the former tending to demonstrate higher spires, and those of the latter typically displaying broader body whorls. L. obovata is more distinctive, bearing a strikingly robust shell with a much-expanded body whorl, although juvenile L. obovata may be indistinguishable from juvenile Pleurocera inhabiting small streams nearby (Figure 2).

Figure 1

Example shells from the three nominal species. Left to right: Pleurocera semicarinata from Kentucky population SK, Pleurocera livescens from New York population LN, and Lithasia obovata from Indiana/Kentucky population OF.

Figure 2

The shells of four individual L. obovata selected from population OF to show a growth series.

Snails bearing shells of a morphology intermediate between that of typical P. semicarinata and P. livescens are common in natural populations. Thus, workers in the field have conventionally followed Goodrich ([1940]), assigning the nomen ‘semicarinata’ to populations in the southern half of the range and the nomen ‘livescens’ to populations in the north. In Indiana, for example, Goodrich and van der Schalie ([1944]) reported the range of P. livescens to be ‘the Wabash River, the Maumee River, and the St. Joseph River of Lake Michigan; probably also in that part of Lake Michigan which borders Indiana.’ They gave the range of P. semicarinata as ‘the two forks of the White River, the Whitewater River, the upper parts of the Big Blue River and probably other streams in southern Indiana.’ This division of the state corresponds roughly to the glacial maximum; P. livescens inhabiting the more lentic environments of the north and P. semicarinata the more lotic south.

Dunithan et al. ([2012]) sampled 39 populations of Pleurocera from Indiana and selected 191 shells varying imperceptibly from those bearing high spires and small body whorls to those with low spires and broad body whorls. The relative warp axis describing most of this variation was significantly correlated with ten environmental variables (of 17 tested), most strongly latitude, longitude, temperature, flow, conductivity, and substrate. Shells sampled from populations in the more lotic lower latitudes tended to be more slender. Although Dunithan and colleagues identified all 39 of their study populations as P. livescens, the 16 most southern populations they sampled would have been identified as P. semicarinata by Goodrich and van der Schalie ([1944]).

In the present study, I use allozyme electrophoresis to test the hypothesis that populations historically referred to Pleurocera semicarinata, P. livescens, and Lithasia obovata are conspecific, the nomen P. semicarinata (Say [1829]) having priority for the single species represented. Over a span of 30 years, the electrophoretic resolution of allozyme variation has demonstrated great utility as a tool for the measurement of genetic divergence among pleurocerid populations throughout eastern North America, with calibration against morphological (Dillon [1984a]), chromosomal (Dillon [1991]), and mtDNA sequence divergence (Dillon and Frankis [2004]). The levels of divergence expected both within and among conspecific pleurocerid populations have been established, and the identification of distinct biological species is typically unambiguous (Chambers [1980]; Dillon and Davis [1980]; Dillon and Reed [2002]; Dillon and Robinson [2011]). But given the extended distances across which populations of P. semicarinata, P. livescens, and L. obovata are understood to range, from Vermont west to Illinois and south to Kentucky, calibration of the observed levels of genetic divergence against some external control would seem desirable. Thus, I also sample a control population of P. canaliculata along with each study population, comparing the genetic similarities observed among the set of known conspecifics to a set of unknowns.

Methods

Because there are no criteria by which nominal P. semicarinata and nominal P. livescens can be unambiguously distinguished, six study populations were selected by prior references in the published literature: three P. livescens, two P. semicarinata, and one L. obovata. Each study population was sympatric with a population of P. canaliculata, five of which were previously analyzed by Dillon et al. ([2013]): AN, AM, AI, AK, and CF. Also included in the present analysis was population SV of P. semicarinata from Virginia, which has served as a reference standard for allozyme studies since the work of Dillon and Davis ([1980]). Sample sites are mapped in Figure 3.

Figure 3

Eastern and Midwestern states of the USA, showing selected rivers and sample sites. Filled circles indicate populations of Pleurocera semicarinata, squares (nominal) P. livescens, and the open circle is (nominal) Lithasia obovata.

Figure 4

Nei's ([1978]) unbiased genetic identities. Based on 11 allozyme loci resolved across the 13 populations of pleurocerids sampled for the present study, with a neighbor-joining tree based on Cavalli-Sforza and Edwards ([1967]) chord distances.

Topotypic population LN of P. livescens and control population AN of P. canaliculata were sampled together near the mouth of Silver Creek, approximately 800 m upstream from Lake Erie, the type locality of ‘Melania’ livescens as given by Menke ([1830]). Study population LM of P. livescens and control population AM of P. canaliculata were the well-studied ‘Station 2’ populations of Dazo ([1965]) and Dillon ([1991]). Study population LI of P. livescens and control population AI of P. canaliculata were sampled together from one of the Wabash drainage study sites of Pyron et al. ([2008]), specifically listed as P. livescens habitat by Goodrich and van der Schalie ([1944]). Study population SO of P. semicarinata and control population AO of P. canaliculata were sampled together from the Scioto drainage, specifically listed a P. semicarinata habitat by Goodrich ([1940]). Topotypic population SK of P. semicarinata and control population AK of P. canaliculata were sampled together from ‘Kentucky,’ which was the type locality of Thomas Say's ([1829]) ‘Melania’ semicarinata, no more specific data given. Topotypic population OF of Lithasia obovata and control population CF of P. canaliculata were sampled together approximately 60 km upstream from the falls of the Ohio, which was given as a type locality for ‘Melania’ obovata by Say ([1829]). Reference population SV of P. semicarinata was population ‘PINE’ in the allozyme study of Dillon and Davis ([1980]) and population ‘Gs’ in the mtDNA sequence study of Dillon and Frankis ([2004]). Detailed locality data for all sites are given in Table 1, and voucher specimens for all populations deposited in the Academy of Natural Sciences of Philadelphia.

Table 1 Locality data

The allozyme analysis reported here was conducted simultaneously with the analysis of Dillon et al. ([2013]), using identical techniques. At least 30 individuals from each pleurocerid population were returned alive to the laboratory, where they were cracked and frozen in Tris tissue buffer for electrophoretic analysis. Techniques and apparatus for horizontal starch gel electrophoretic resolution of allozyme variation in homogenates of molluscan tissues are detailed in Dillon ([1992]), along with recipes for all the buffers and stains employed here. Allozyme bands interpretable as the products of codominant genes segregating in Mendelian fashion were resolved at 11 loci using 10 enzyme stains as in Dillon ([2011]):

The Tris Cit 6 buffer (buffer XIII of Shaw and Prasad [1970]) was used to resolve 6-phosphogluconate dehydrogenase (6PGD), octopine dehydrogenase (OPDH), and isocitrate dehydrogenase (2 loci, the cathodal IDHF and the anodal IDHS). A Poulik ([1957]) discontinuous buffer system was employed for glucose-phosphate isomerase (GPI), sorbitol dehydrogenase (SDH), and octopine dehydrogenase (a second time). The TEB8 buffer system (buffer III of Shaw and Prasad [1970]) was used to analyze phosphoglucomutase (PGM - the strong, fast locus only), xanthine dehydrogenase (XDH), and mannose phosphate isomerase (MPI). A TEB9.1 buffer (Dillon and Davis [1980]) was used for octanol dehydrogenase (OLDH), esterases (EST1 - the strong, slow locus only), and xanthine dehydrogenase (a second time).

Mendelian inheritance of allozyme phenotype has been confirmed for GPI, OPDH, and EST1 by Dillon ([1986]) and for 6PGD by Chambers ([1980]). Putative allelic designations for each zone of allozyme activity were assigned by setting the most common band in reference population SV to ‘100,’ and naming all other alleles by the mobility of their allozymes (in millimeters) relative to this standard.

Gene frequencies and mean direct-count heterozygosities (the unbiased estimate of Nei ([1978])) were calculated using Biosys version 1.7 (Swofford and Selander [1981]). Because large numbers of alleles were resolved at some loci, sample sizes dictated that genotypes be pooled into three classes before testing for Hardy-Weinberg equilibrium: homozygotes for the most common allele, common/rare heterozygotes, and rare homozygotes together with other heterozygotes. Yates-corrected chi-square statistics were then employed for this purpose. I calculated matrices of Nei's ([1978]) unbiased genetic identity and distance, as well as Cavalli-Sforza and Edwards ([1967]) chord distance. As chord distances are Pythagorean in Euclidean space (Wright [1978]), they were used to calculate a neighbor-joining tree with Phylip v3.65 program NEIGHBOR (Felsenstein [2004]).

Results

No variation was detected across all 13 populations examined at two loci: ISDHS and SDH. Over the remaining data set of 13 populations × 11 variable loci resolved, 50 loci were polymorphic by the 95% criterion. Chi square goodness-of-fit tests to Hardy-Weinberg expectation returned only two statistics significant at the nominal 0.05 level, a result clearly attributable to type 1 error. Allele frequencies and direct-count heterozygosities at all variable loci across the seven populations newly surveyed here are available as an Additional file 1. See Dillon et al. ([2013]) for allele frequency data from the five P. canaliculata controls previously surveyed (AN, AM, AI, AK, and CF) and reference population SV.

Values of Nei's unbiased genetic identity among the six study populations of nominal P. semicarinata, P. livescens, and L. obovata ranged from 0.633 to 0.967, very similar to the range of 0.661 to 1.00 observed within the six control populations of P. canaliculata (Figure 4). Genetic identities between the six study populations and the six controls ranged from 0.672 to 0.348.

The neighbor-joining tree based on Cavalli-Sforza and Edwards chord distances showed two discrete clusters, the six study populations clustered in a fashion quite similar to the six controls. The maximum interpopulation chord distance within the P. canaliculata controls was the 0.529 observed between Kentucky population AK and New York population AN, greater than the 0.516 maximum chord distance observed between study populations SK from Kentucky and LM from Michigan.

Discussion

Genetic divergence at allozyme-encoding loci between conspecific pleurocerid populations is a function of both environmental difference and geographic distance (Dillon [1984b]). Sampled from upper Tennessee River drainages at a maximum separation of approximately 250 km, a set of nine P. clavaeformis populations returned a median Nei's Genetic Identity I = 0.90 with a range of 0.77 to 1.00 (Dillon [2011]). Chambers ([1980]) reported a median I = 0.84 (range 0.73 to 0.96) among eight populations of Pleurocera floridensis sampled from Florida, where the distance between the most distant pair of populations was approximately 350 km. A pair of similar studies each involving eight Pleurocera catenaria populations from the Atlantic Piedmont returned identical medians of I = 0.88 (range 0.78 to 0.99, over 350 km) in Georgia (Dillon and Robinson [2011]) and I = 0.88 (range 0.72 to 1.00, over 400 km) in North and South Carolina (Dillon and Reed [2002]).

The values of Nei's genetic identity shaded in Figure 4, median I = 0.78 (range 0.66 to 1.00) within the six control populations of P. canaliculata and median I = 0.81 (range 0.63 to 0.97) within the six study populations of nominal P. semicarinata, P. livescens, and L. obovata, are both somewhat lower than those published for pleurocerid populations further east, but certainly within expectation for populations sampled across a range extending approximately 650 km (LI to LN). And the striking similarity of the divergence levels seen in the study and control groups is consistent with the hypothesis that all six study populations are conspecific, the nomen P. semicarinata (Say [1829]) having priority over L. obovata (Say [1829]) and P. livescens (Menke [1830]) for the combined group.

The situation with P. semicarinata, P. livescens, and L. obovata seems identical to that demonstrated by Dillon ([2011]) for nominal P. clavaeformis, Elimia acutocarinata, and Pleurocera unciale, and that demonstrated by Dillon et al. ([2013]) for nominal P. canaliculata, Pleurocera acuta, and Pleurocera pyrenellum. In all three of these sets of populations, the interpopulation phenotypic variance in shell morphology is so great that 19th-century taxonomists erroneously divided a single species into three specific nomina, effectively rendering phenotypic plasticity among them cryptic. The present demonstration that freshwater snail populations previously assigned to the nomen ‘L. obovata’ are conspecific with P. semicarinata may be the most striking case of cryptic phenotypic plasticity yet documented in freshwater snails.

Yet in hindsight, taxonomists of earlier generations were not entirely deceived. Wiebe ([1926]) suggested that the relatively broad body whorl demonstrated by P. livescens in Lake Erie might be a simple consequence of a larger foot, which might in turn represent a local adaptation to heavy wave action. And in 1934 to 1941, Calvin Goodrich authored an eight-paper series entitled ‘Studies of the Gastropod Family Pleuroceridae,’ documenting many striking examples of intraspecific variation in pleurocerid shell morphology and calling into question its heritability. Goodrich ([1934]) observed, for example, ‘that low-spired forms of L. obovata occur near the mouth of Green River, that high-spired forms occur in the upper reaches and tributaries of the river, and that, in a general way, there is a progressive alteration of form downstream to upstream.’ He attributed the greater shell ‘obesity’ of pleurocerid populations inhabiting larger rivers to physical disturbance in the ‘rapid and sometimes tumultuous’ downstream environments, suggesting that ‘variations in relative obesity appear to represent environmental responses rather than genetic differences.’ Goodrich ([1937]) extended his observations on shell ‘obesity’ to populations of Pleurocera sampled throughout the South and Midwest, proposing as a ‘natural law’ that ‘Pleurocera increases in relative diameter in a downstream direction’ and again adding that ‘findings point to environmental influences or controls.’

Sixty years later, Urabe ([1998]) confirmed a significant correlation between shell shape and microhabitat choice in populations of the Japanese pleurocerid Semisulcospira reiniana, snails active in the most rapid currents indeed demonstrating a significantly larger body whorl and lower spire. Urabe's controlled breeding studies also supported Goodrich's speculation that the heritable component of such shell shape variation may be negligible.

Most of the relevant research published in recent years has focused, however, on the response of shell phenotype to the effects of predation (Langerhans and DeWitt [2004]; Covich [2010]; Hoverman and Relyea [2011]). Crayfish can be important predators of snails in smaller rivers and streams, for example, inserting their claws into the apertures of tougher shells and cracking the lip in the style of a can opener. Krist ([2002]) reared P. livescens from Indiana in the effluent from crayfish feeding on conspecific snails, demonstrating a significant narrowing of the body whorl and lengthening of the spire compared to controls not exposed to such predation. Thus, it seems possible that the P. semicarinata populations inhabiting smaller streams may be responding to greater predation pressure exerted by crayfish, while those of the larger rivers may be responding to whole-shell crushing predators such as drum.

The complete range of P. semicarinata, as properly understood, seems to extend across twelve US states and two Canadian provinces, from northern Vermont and southern Quebec west to Wisconsin and Iowa, south through Kentucky. The older field observations and more recent experimental results can be combined with the morphometric survey of Dunithan et al. ([2012]) to yield an overall model for shell form in P. semicarinata across this broad range. Waters in the south, below the glacial maximum, are more lotic in their character, with a well-defined gradient of current, substrate, and predator community from headwaters to the big rivers. Crayfish predation in small streams (possibly in combination with other factors) induces the typical semicarinata form in small streams in the south, and crushing predators in large rivers induce the obovata form. Waters in the northern half of the range, above the glacial maximum, are more lentic in their character and the upstream-downstream gradient less well-defined. Thus, the intermediate P. livescens morphology prevails.

Across all 13 populations surveyed, the greatest genetic divergence was demonstrated by reference population SV. Although identified as P. semicarinata on the basis of shell morphology, the reference standard from Virginia proved slightly more genetically similar to control P. canaliculata than to any nominal conspecific. The maximum value of Nei's genetic identity between reference SV and any other population of P. semicarinata was 0.389, clearly in the range one might expect for a distinct species (Chambers [1980]; Dillon and Reed [2002]; Dillon and Robinson [2011]).

The shell morphology displayed by reference population SV is typical of Pleurocera common and widespread throughout the Great Valley of Virginia in the upper New, Roanoke, James, and Shenandoah drainages (Dillon and Davis [1980]). This range is apparently disjunct from the primary range of P. semicarinata as it is now understood, however, as no nominal conspecifics are known to inhabit the waters of central West Virginia or Maryland. Whether the Virginia populations indeed represent an isolated pocket of P. semicarinata, or may rather constitute a previously unrecognized species cryptic under P. semicarinata, would be a fertile area for future research.

Conclusions

The evidence presented here suggests that populations of North American freshwater prosobranch snails previously identified as Pleurocera semicarinata, Pleurocera livescens, and Lithasia obovata are conspecific, the shell characters by which these taxa have heretofore been distinguished likely resulting from ecophenotypic plasticity. Possible drivers for such phenotypic variation in shell morphology include predation, substrate, and flow.