Discrete phenotypes are not underpinned by genome-wide genetic differentiation in the squat lobster Munida gregaria (Crustacea: Decapoda: Munididae): a multi-marker study covering the Patagonian shelf
Wang et al. BMC Evolutionary Biology
Discrete phenotypes are not underpinned by genome-wide genetic differentiation in the squat lobster Munida gregaria (Crustacea: Decapoda: Munididae): a multi- marker study covering the Patagonian shelf
Chen Wang 0
Shobhit Agrawal 0
Jürgen Laudien 0
Vreni Häussermann 1 2
Christoph Held 0
0 Alfred Wegener Institute, Helmholtz Center for Polar- and Marine Research , Am Handelshafen 12, 27570 Bremerhaven , Germany
1 Huinay Scientific Field Station , Huinay, Los Lagos , Chile
2 Universidad Católica de Valparaíso, Facultad de Recursos Naturales, Escuela de Ciencias del Mar , Avda. Brasil 2950, Valparaíso , Chile
Background: DNA barcoding has demonstrated that many discrete phenotypes are in fact genetically distinct (pseudo)cryptic species. Genetically identical, isogenic individuals, however, can also express similarly different phenotypes in response to a trigger condition, e.g. in the environment. This alternative explanation to cryptic speciation often remains untested because it requires considerable effort to reject the hypothesis that the observed underlying genetic homogeneity of the different phenotypes may be trivially caused by too slowly evolving molecular markers. The widespread squat lobster Munida gregaria comprises two discrete ecotypes, gregaria s. str. and subrugosa, which were long regarded as different species due to marked differences in morphological, ecological and behavioral traits. We studied the morphometry and genetics of M. gregaria s. l. and tested (1) whether the phenotypic differences remain stable after continental-scale sampling and inclusion of different life stages, (2) and whether each phenotype is underpinned by a specific genotype. Results: A total number of 219 gregaria s. str. and subrugosa individuals from 25 stations encompassing almost entire range in South America were included in morphological and genetic analyses using nine unlinked hypervariable microsatellites and new COI sequences. Results from the PCA and using discriminant functions demonstrated that the morphology of the two forms remains discrete. The mitochondrial data showed a shallow, star-like haplotype network and complete overlap of genetic distances within and among ecotypes. Coalescent-based species delimitation methods, PTP and GMYC, coherently suggested that haplotypes of both ecotypes forms a single species. Although all microsatellite markers possess sufficient genetic variation, AMOVA, PCoA and Bayesian clustering approaches revealed no genetic clusters corresponding to ecotypes or geographic units across the entire South-American distribution. No evidence of isolation-by-distance could be detected for this species in South America. Conclusions: Despite their pronounced bimodal morphologies and different lifestyles, the gregaria s. str. and subrugosa ecotypes form a single, dimorphic species M. gregaria s. l.. Based on adequate geographic coverage and multiple independent polymorphic loci, there is no indication that each phenotype may have a unique genetic basis, leaving phenotypic plasticity or localized genomic islands of speciation as possible explanations.
Phenotypic plasticity; Genetic homogeneity; Squat lobster; Microsatellites; Gene flow
Different species have different morphologies and
lifestyles, which is commonly taken (but not often tested)
to reflect different underlying genotypes. The advent of
affordable DNA sequencing and molecular barcoding
have served to greatly intensify the crosstalk between
molecular and taxonomic disciplines by uncovering a
large number of previously overlooked genotypes [1–3],
many of which could be shown to be associated with
equally overlooked morphotypes that in retrospect were
identified as (pseudo)cryptic species [4–6].
However, the popularity and large number of cryptic
species currently being discovered have led to an
under-appreciation of the notion that sharply distinct
morphotypes are not always the consequence of genetic
differences but can also be invoked from the same
genotype, often called by environmental triggers. The
differences between associated morphotypes and
lifestyles of ecotypes within the same species can be
surprisingly pronounced [7–11].
Proving polyphenism and rejecting cryptic speciation
as an explanation is harder than sequencing a
mitochondrial gene fragment, which may in part explain the
relatively lower number of well-studied cases of
polyphenism. Whilst consistent differences among different
morphotypes in a single mitochondrial marker suffice
to at least flag these clades as candidate cryptic species,
the opposite observation (no consistent differences) is
not a conclusive demonstration of the absence of
genetic differentiation among ecotypes. In order to show
that too slowly evolving markers or other artefacts (e.g.
mito-nuclear discordance ) did not trivially cause
the observed lack of differentiation, considerably more
extensive molecular evidence including multiple
unlinked nuclear loci with sufficiently high substitution
rates is required. Such extensive a posteriori knowledge is
rare (e.g. in the fully sequenced Daphnia pulex [13–15]),
but numerous experimental studies in which the genetic
identity of individuals is known a priori contribute greatly
to our understanding of the importance of polyphenism
and morphological plasticity, e.g. parthenogenetic aphids
[16, 17], marbled crayfish , polyembryonic armadillos
, inbred lines of Drosophila , cloned swine . It
is unclear if the small number of confirmed polyphenism
resulting from similar or identical genetic backgrounds is
a condition truly rare in nature or whether it reflects
mostly a discovery and/or publication bias.
In this paper we investigate the dimorphic squat lobster,
Munida gregaria sensu lato (Fabricius, 1793), which is
currently considered to comprise the ecotypes M. gregaria
sensu stricto Miers 1881 as well as its junior synonym M.
subrugosa Dana, 1852 (see ). For clarity and brevity,
hereafter we refer to them as Munida gregaria s.l.
comprising the two ecotypes gregaria s.str. and subrugosa. In
South America, M. gregaria s. l. occurs in shallow marine
waters off Patagonia, including Tierra del Fuego and the
Falkland Islands/Islas Malvinas, while in the southwestern
Pacific M. gregaria s. l. are reported from off eastern New
Zealand and Tasmania ([22, 23] and references therein).
The taxonomic status of gregaria s. str. and subrugosa
ecotypes has been subject to conflicting interpretations.
Both ecotypes were often regarded as different species
because of morphological differences in adult specimen
(Fig. 1) as well as different behaviors at certain
developmental stages [24–27]. Williams (1973) on the other hand
interpreted gregaria s. str. as a transient, pelagic
ontogenetic stage that would later in life gain the physical features
of subrugosa upon adopting a permanently benthic
lifestyle. Regardless of the taxonomic ramifications, gregaria
s. str. is often found in huge pelagic swarms that
subrugosa lacks [23, 28]. These differences persist even where
both ecotypes co-exist in the same habitat. Nevertheless,
on the basis of a lack of mitochondrial DNA
differentiation , these two ecotypes are currently treated as a
single polymorphic species under the name of M. gregaria
in the most recent taxonomic revision of the family .
But this evidence must be considered insufficient because
the sampled region (Beagle Channel) represents a very
small part of the species’ distribution and the molecular
evidence rest exclusively on two linked mitochondrial
markers (COI and ND1), whereas the results of the only
nuclear marker (ITS-1) had to be excluded from the final
analysis of the only molecular study .
In order to test whether gregaria s. str. and subrugosa
ecotypes correspond to different species of Munida or
represent a single species with variable phenotypes, we
employed multiple independent, fast-evolving nuclear
microsatellite markers  and an expanded set of
mtDNA sequences. The sampled area encompasses nearly
the entire distribution of gregaria s. str. and subrugosa
ecotypes in South America. In addition, we analyzed
morphological differences of both ecotypes and different
ontogenetic stages following the method of  in order to
test if the more complete geographic sampling continues
to support the discrete morphological clusters or if the
boundaries between the two ecotypes vanish under more
complete geographic coverage.
A total number of 219 individuals were used for both
morphological and molecular analyses in this study. These
samples were collected at 25 stations from Patagonia and
off the Falkland Islands/Islas Malvinas ranging from 38° S
to 55° S from the shallow subtidal area to 179 m water
depth by mid-water or bottom trawls (Table 1, Fig. 2).
Sexual maturity was identified as presence of eggs and/or
sexually dimorphic pleopods. In adult males the first two
gregaria s. str.
Fig. 1 General views of gregaria s. str. and subrugosa ecotypes (above) and schematic diagrams showing morphometric measurements (below).
For the dorsal view of both ecotypes, scale bar represents 1 cm. Measurements are made on: ACW, anterior carapace width; RBW, rostrum basis
width; DaW, width of dactylopodus of third maxilliped; PW, width of propodus of the third maxilliped; EL, eyestalk length
pairs of pleopods are modified to form gonopods, the
remaining three pairs are flap-like; in adult females the
first pair is missing and the remaining four pairs are
elongated with long setae for egg-carrying .
Specimens were checked by eye and classified as
gregaria s. str. rather than subrugosa ecotypes based on the
following characteristics: longer eyestalk length (EL),
wider rostrum basis (RBW) and broader and blunter
dactylus of the third maxilliped (DaW) (Fig. 1). These
three morphometric characteristics together with
anterior carapace width (ACW) and width of propodus of the
third maxilliped (PW) were statistically significant in
discriminating ecotypes . We measured our samples
using a Leica MZ-12.5 microscope with intraocular scale
to the nearest 0.1 mm. To determine patterns emerging
from the morphometric measurements of these five
body parts, principal components analysis (PCA) was
plotted using the statistical package PAST 3
(PAlaeontological STatistics, ). Applying discriminant functions
(DF1 and DF2) introduced in , ΔDF (DF1-DF2)
values were calculated based on these measurements
and subsequently plotted using R .
DNA was extracted from ethanol-preserved abdominal
or cheliped muscle tissue using QIAamp DNA Mini Kit
(QIAGEN, Germany). For mtDNA analysis, a region of
the COI gene was amplified using the universal primers
HCO2198 and LCO1490  for 96 individuals
(Table 1). The 10 μl reactions consisted of 0.02U/μl
Hotmaster Taq (5 Prime), 0.2 mM dNTPs, 0.5 μM of
forward and reverse primers, 1 × PCR-buffer and 1 μl
(about 30 ng) of template DNA. PCR was conducted
using an initial denaturation at 94 °C for 2 min,
Table 1 Sampling sites and number of adult and juvenile (in parentheses) gregaria s. str. and subrugosa ecotypes
Station Latitude Longitude gregaria s. str. adults (juveniles) subrugosa adults (juveniles)
Northern Chilean Patagonia (NCP)
Tierra del Fuego archipelago (TdF)
Mar del Plata (MdP)
NmtDNA the number of specimens used with mitochondrial marker, NMSAT the number of specimens used with microsatellites
arefers to one specimen with missing genotype at a certain microsatellite locus
followed by 36 cycles of 94 °C for 20 s, annealing at
47 °C for 20 s, 65 °C for 1 min, and a final extension at
65 °C for 10 min. Size and quality of amplified products
were checked on a 2% agarose gel in TAE buffer, and
then 1 μl of purified PCR product was used for cycle
sequencing with the HCO primer. Sanger sequencing
was conducted on an ABI 3130xl sequencer. Alignment
was done using CODONCODE ALIGNER 4.0
(CodonCode Corp.) and checked for the presence of
ambiguities and stop codons.
DNA polymorphism was examined as haplotype
diversity (HD) and nucleotide diversity (π) for each ecotype
and all samples using DnaSP 5.10 . Genealogical
relationships among haplotypes were inferred using
statistical parsimony implemented in TCS 1.21 .
Pairwise genetic divergences measured as number of
nucleotide differences were calculated within and
between the two ecotypes in MEGA 5.2 . For a better
understanding of the genetic distances and barcoding
gap analysis, we added three congeneric species to our
Isla Gran Malvina
Fig. 2 Sampling sites of gregaria s. str. (solid circle) and subrugosa (open circle) ecotypes. FM, Falklands/Malvinas; TdF, Tierra del Fuego archipelago;
NCP, Northern Chilean Patagonia; MP, Mar del Plata
analysis: M. rutllanti (n = 5; GenBank accession numbers:
JQ306226-JQ306230), M. quadrispina (n = 3; GenBank
accession numbers: DQ882090-DQ882092), and M.
gracilis from our own collection (n = 2; GenBank
accession numbers: KJ544249-KJ544250). Pairwise genetic
distances were calculated within M. gregaria s. l. (pooled
gregaria s. str. and subrugosa) and versus the other three
The COI dataset was analyzed using coalescent based
approaches Poisson tree processes (PTP) model 
and the general mixed Yule coalescent model (GMYC)
[38, 39] for a critical evaluation of species delimitation. As
an outgroup, COI sequences from four congeners, M.
rutllanti (n = 5; GenBank accession numbers:
JQ306226JQ306230), M. quadrispina (n = 2; GenBank accession
numbers: DQ882090 and DQ882092), M. rosula (n = 1;
GenBank accession numbers: AY350994) and M. congesta
(n = 1; GenBank accession numbers: AY350945), were
added. The substitution model that best fits the data was
determined using jModelTest 2.1.5 [40, 41].
PTP does not require an ultrametric tree, as the
transition point between intra- and inter-specific branching
rates is identified using directly the number of nucleotide
substitution . A maximum likelihood (ML) phylogeny
of the COI dataset was reconstructed in RAxML-HPC 8
in CIPRES portal [42, 43], employing a HKY + G model
that was suggested by the corrected Akaike Information
Criterion (AICc) and the Bayesian Information Criterion
(BIC). Nodal support was evaluated using 1000 bootstrap
replicates. The ML phylogenetic tree was used as the
input tree to run PTP species delimitation analysis in
the PTP webserver (http://species.h-its.org/ptp/). We
ran the PTP analysis for 500,000 MCMC generations,
with a thinning value of 100 a burn-in of 10%. Outgroup
taxa were kept since the MCMC chains did not converge
when they were removed.
The GMYC method requires a fully resolved tree with
branch lengths estimates, which was obtained using the
program BEAST 2.4.3 . We used a site-specific HKY
substitution matrix and a gamma distributed model of
among-site rate heterogeneity with four discrete rate
categories. We implemented a strict clock model of 2%
per Myr as suggested for COI sequences in Crustacea
[45, 46], and selected a Yule tree prior. Default values
were used for remaining priors. MCMC analysis was
run for a total of 10 million generations, sampling every
1000 steps. Convergence was assessed by examining the
likelihood plots through time using TRACER 1.6 .
The COI chronogram was then analyzed using the GMYC
package in SPLITS in R (version 3.1.2, www.cran.r-project.
org), using the single threshold approach [38, 39].
In total, 218 individuals were screened for genetic
variation at 11 microsatellite loci that were originally
designed for M. gregaria s. l.  (Table 1). Allele sizes
were binned manually and genotypes were assessed in
GENEMAPPER 4.0 (Applied Biosystems). Null alleles,
stuttering and large allele dropout were tested using
MICROCHECKER . Because of too many missing
data and possible null alleles, locus Mgr63 and Mgr105
were excluded from subsequent analyses. Genetic
diversity within each ecotype was summarized as allelic
richness (Ar) in FSTAT 22.214.171.124  using the rarefaction
approach, which was also used to determine the
number of private alleles using standardized sample sizes in
ADZE 1.0 . Detection of linkage disequilibrium
between loci and deviations from Hardy-Weinberg
equilibrium (HWE) per ecotype were performed using
GENEPOP 4.2 . All loci were tested for
positive/diversifying or balancing selection using LOSITAN ,
which simulates an expected distribution of FST as a
function of expected heterozygosity under an island
model of migration. The statistical power of this set of
microsatellite loci to detect significant genetic
differentiation between populations/ecotypes was tested with
POWSIM 4.1  using both Chi-square (χ2) and Fisher’s
exact test analysis. Various levels of differentiation
(measured as FST in the range from 0.001 to 0.01) were
determined by combining different effective population
size (Ne) and times since divergence (t). In addition,
POWSIM allows calculating type I error probability,
which is the probability of rejecting the null hypothesis
of genetic homogeneity although it was true by drawing
the alleles directly from the base population (t = 0).
The genetic differentiation among the three major
sampling areas, i.e., FM, NCP and TdF (see Table 1), was
assessed for each ecotype separately with AMOVA in
ARLEQUIN 3.5 . Both FST and RST estimators were
calculated over all nine loci with 1000 permutations. To
provide a visual representation of species separation and
potential subdivision, Principal Coordinate Analysis
(PCoA) was performed in GENALEX 6.5 .
Bayesian assignment tests were used to evaluate the
level of genetic clustering. We used STRUCTURE 2.3.4
 first without giving any prior population
information, letting K range from 1 to 5. We also checked
whether individuals could be assigned correctly to
clusters if the number of ecotypes was given a priori (K = 2).
Both conditions were run with the correlated allele
frequencies option under the non-admixture model, i.e.
under the assumption that there is no gene flow between
ecotypes, as well as under an admixture model, i.e.
allowing limited introgression between clusters. Twenty
runs with 200,000 Markov chain Monte Carlo (MCMC)
iterations after a burn-in period of 25,000 steps were
carried out for each K. The results were uploaded onto
STRUCTURE HARVESTER  and K was determined
using the ad hoc statistic ΔK , as well as mean
estimates of posterior probability L(K) . Results from
the 20 replicates of the most likely value for K were
averaged using the software CLUMPP 1.1.2  and the
output was visualized using DISTRUCT 1.1 .
Since L(K) does not always provide the correct number
of clusters and the ΔK statistic cannot evaluate K = 1 or
the largest value chosen for K , we also applied
STRUCTURAMA 2.0 , which can directly estimate
the number of clusters in which a sample can be
subdivided. We allowed the number of populations to be a
random variable following a Dirichlet process prior, ran
the MCMC analysis for 1,000,000 cycles, sampled every
100th cycle, and discarded the first 400 samples as
The impact of isolation by distance (IBD) across
southern South America on genetic differentiation was
estimated by Mantel tests as implemented in IBD Web
Service 3.23 . For this purpose, only the geographic
position but not the ecotype of the samples
(mitochondrial and microsatellite data) were used (Table 1; Fig. 2).
Pairwise FST values for mitochondrial data and pairwise
(δμ)2 genetic distance  for microsatellites were
obtained in ARLEQUIN 3.5, where spatial distances were
calculated using the Geographic Distance Matrix
Generator 1.2.3 (http://biodiversityinformatics.amnh.org/open_
source/gdmg/index.php). Geographical distances were
log-transformed to account for two-dimensional habitat
distribution , and the significance of the slope of
the reduced major axis (RMA) regression was assessed
by 30,000 randomizations.
PCA comparison of five key morphometric
characteristics  revealed clearly distinct groups corresponding
to ecotype and age, with the first principal component
explaining 92.48% of the variation. Samples of adult and
juvenile gregaria s. str. formed two isolated groups, both
of which were clearly distinct from the subrugosa
samples. The subrugosa individuals comprise the adult and
juvenile sub-groups that overlap in part (Fig. 3). The
discriminant functions with these five morphometric
characteristics yielded result coherent with the PCA.
Samples of adult gregaria s. str. and juvenile gregaria s.
str. clustered separately, while due to the limited number
of subrugosa juveniles (n = 12) it is hard to ascertain
whether juvenile and adult subrugosa show statistically
significant differences (Fig. 4). Results of morphological
analyses ascertain that even among samples from the
entire South American distribution gregaria s. str. and
subrugosa are distinct ecotypes at different ontogenetic
stages with discrete morphological traits rather than
forming the extremes of a continuous distribution.
Mitochondrial COI sequence variation
A total number of 96 COI sequences from 61 gregaria s.
str. and 35 subrugosa individuals were obtained as an
alignment of 618 bp (GenBank accession numbers:
KJ544251 - KJ544346). These sequences were collapsed
into 30 different haplotypes, possessing 29 variable
(segregating) sites, of which 11 were parsimony-informative.
The subrugosa ecotype showed slightly higher genetic
diversity than gregaria s. str. (Additional file 1: Table S1).
In the 206 codons of the alignment, the 29 variable sites
were all synonymous substitutions and no stop codons
Genealogical relationships among haplotypes showed a
very shallow, star-like structure. The most common
haplotype (n = 59) was shared by both ecotypes as well
gregaria s. str._adult
gregaria s. str._juvenile
Component 1 (92.48% of Variation)
Fig. 3 Principal component analysis (PCA) biplot for the morphometric measurements of all the samples. Dashed lines show vectors from the five
representative body characteristics that are statistically significant in discriminating ecotypes. The 95% concentration ellipses are given for juvenile
and adult gregaria s. str. as well as juvenile and adult subrugosa. A dashed 95% concentration ellipse represents all subrugosa individuals. Solid
arrows indicate suggested interpretation as ontogenetic transition within ecotype (horizontal) and morphological discreteness between ecotypes
(vertical; see discussion)
Morphological differentiation within M. gregaria s. l.
as by all the sampling regions, which differed from the
other haplotypes in 1 to 3 mutational steps (Additional
file 2: Figure S1).
Extent of intraspecific and interspecific COI divergence
The mean number of differences among sequences
within each ecotype was 0.809 for gregaria s. str. and
1.408 for subrugosa, between ecotypes it was 1.107.
The plotted pairwise genetic distances show complete
overlap of distributions within and between ecotypes
(Fig. 5). The maximal number of differences was six
base pairs and no barcoding gap between ecotypes
could be identified. By contrast, pronounced barcoding
gaps exist between M. gregaria s. l. and each of the
three other Munida species. These interspecific
distances were at least ten times larger than distances
between gregaria s. str. and subrugosa ecotypes (Fig. 5).
The results of COI data including samples from almost
the entire South American distribution provide no
evidence of genetic separation between gregaria s. str. and
PTP and GMYC species delimitation
The ultrametric tree obtained in BEAST was used to
illustrate the delimitation of putative species recognized
by the different approaches conducted with the COI data
(Fig. 6). Both PTP and GMYC analyses detected 5
candidate species corresponding to the current taxonomic
units, that is, the four outgroup Munida species and a
single species M. gregaria including all the sequences of
gregaria s. str. and subrugosa.
Overall, the number of alleles per locus ranged between
six (Mgr46) and 38 (Mgr60) with an average of 14.9. For
each locus the number of alleles, allelic size range, allelic
richness as well as observed and expected heterozygosities
per population and ecotype are reported (Additional file 3:
Table S2). The mean number of private alleles per locus
when sample size was standardized was slightly higher for
subrugosa (0.198 ± 0.120) compared to gregaria s. str.
(0.111 ± 0.065), both of which were very low, reflecting a
high degree of allelic sharing between the two ecotypes.
Compared to the limited variation among COI sequences,
the nine microsatellite loci exhibited broader allelic
ranges and orders of magnitude higher allelic variation
(Additional file 4: Figure S2). All loci showed no linkage
disequilibrium. Locus Mgr90 showed significantly higher
heterozygosity than expected, but since excluding Mgr90
had only minor effect on the results, it was kept in this
study. The power test suggested that our microsatellite
dataset was sensitive enough to detect very weak genetic
differentiation (FST = 0.005) in probabilities of close to
100% using both chi-square (χ2) and Fisher’s exact test
(Additional file 5: Figure S3). The FST outlier analysis
showed that none of the loci was under potential selection
at 95% confidence level thus they were regarded as neutral
in our interpretation of the results.
Genetic differentiation and individual assignment inferred
Hierarchical AMOVA showed almost all the genetic
variance distributed among individuals within sampling
within (gregaria s. s.tr + subrugosa)
(gregaria s. str. + subrugosa) vs M. rutllanti
(gregaria s. str. + subrugosa) vs M. quadrispina
(gregaria s. str. + subrugosa) vs M. gracilis
Pairwise genetic distance (number of base differences)
Fig. 5 Estimates of divergence between sequences among the two ecotypes and three other Munida species. Pairwise comparisons were
performed using number of base differences
areas (Table 2), thus corroborating the weak geographic
structure in the distribution of mitochondrial haplotypes.
PCoA showed approximately equal distributions along
the first three axes, which accounted for 20.04, 19.56
and 18.76% of the total genetic variance, respectively
(Fig. 7). This result indicates that there is no single
factor (ecotype or other) that would dominate the
distribution of total genetic variance for the high-resolution
Bayesian cluster analyses with STRUCTURE suggested
the best K was 1 according to average log probability
(L(K)), but K = 2 was indicated by the highest statistic ΔK.
This is because the change in log probability does not
account for the smallest and largest K. Even under K = 2,
each individual possessed a roughly equal probability of
being assigned to the first versus the second cluster, which
indicates that all individuals belong to one single group
(Fig. 8). Congruent distribution of posterior probabilities
for each individual was obtained given the a priori
assumed number of putative populations (=ecotypes).
The STRUCTURAMA analysis corroborated the
inferences of the STRUCTURE analysis. The sampled
individuals all belonged to one group with a posterior
probability of 1. Eventually, both Bayesian analyses
showed no correlation of our microsatellite data with
either ecotypes or geographical units.
An absence of genetic structure over the entire
distribution range in South America was also found in the IBD
tests. Based on both mitochondrial and microsatellite data
sets, Mantel tests indicated no significant correlation (for
COI, r = 0.0608, P = 0.229; for microsatellites, r = −0.005,
P = 0.508) between genetic and log-transformed
geographic distances (Additional file 6: Figure S4).
Stability of morphological dimorphism in M. gregaria s. l.
In theory, populations belonging to a so-called ‘ring species’
might appear sharply distinct in an area of secondary
overlap, but appear more gradually changing in morphology or
genetics through areas of their distribution that have been
Fig. 6 Phylogenetic relationships for M. gregaria s. l. and outgroup taxa. Ultrametric phylogenetic tree inferred from COI in BEAST species. Scale
axis showed ages in millions of years (Ma). Species delimitation scenarios obtained from different methods are indicated in columns to the right.
Numbers at nodes/tips are Bayesian support values in PTP model and AIC-based support values in GMYC model for the delimited species
Table 2 Hierarchical analysis of molecular variance based on microsatellites for both ecotypes
Source of variation
Among sampling areas 3
Within sampling areas
FCT = −0.00103
FSC = 0.00437
FST = 0.0334
RCT = −0.00448
RSC = 0.00451
RST = 0.00005
PCoA Axis 1 (20.04% of Variation)
Fig. 7 Principal coordinate analysis (PCoA) of 218 individuals of gregaria s. str. and subrugosa ecotypes based on microsatellites
more continuously inhabited (see [65, 66] and references
therein). Inadvertently sampling only in the zone of
secondary overlap might therefore create the incorrect impression
of discrete morphotypes or genotypes when populations
with intermediary morphotypes remain unsampled.
Increasing the sampling area from a single location in
the Beagle Channel  to a continental scale, our data
suggest that the boundary between two morphotypes
(gregaria s. str. and subrugosa) is nonetheless not blurred
across the South American shelf (Figs. 3 and 4).
The expanded morphometric analysis further suggests an
ontogenetic dimension in the morphometry. It may be
expected that the gap between adults and Northern Chilean
Patagonia (NCP) juveniles in gregaria s. str. might be closed
by inclusion of juveniles from other populations and reveal
a continuous ontogenetic transition as can already be found
in subrugosa (Fig. 3). The discreteness of the subrugosa and
gregaria s. str. morphotypes, however, is not a sampling
artefact and stable with respect to a more representative
sampling scheme as well as inclusion of different life stages.
gregaria s. str.
Fig. 8 Individual probabilities of cluster assignments from the software STRUCTURE. The most likely number of clusters K = 2 is shown for
gregaria s. str. (n = 161) and subrugosa (n = 57) using nine microsatellite loci. Each vertical line represents the probabilities for a single
individual to be assigned to one of the clusters (K)
Since the proposal that phenotype and genotype form
two fundamental different levels of biological
abstractions , untangling the relationship between
phenotypes and the underlying genotypes has long been
challenging and intriguing. The advent of molecular
techniques has greatly fostered studies of
phenotypegenotype interaction, especially in the wake of helped
discovery of (pseudo)cryptic genetic divergence whereas
corresponding phenotypes appeared identical. Such
unexpected genetic diversity, which was later often
corroborated by other independent evidences from
morphology , breeding behavior  or multiple,
independent and informative nuclear markers , has
become an important supplement for the phenotypic
identification of an organism to species or sub-species
level in taxonomic practice [71–73].
In other cases, however, molecular marker-based
examination found no genetic differentiation matching
discrete phenotypes, which is exemplified by the present
Munida gregaria case. Nonetheless, the lack of
differentiation at a single marker is insufficient to extrapolate to
the entire genome, especially in view of the different
inheritance in the mitochondrial and nuclear genomes [74,
75]. A previous molecular study used only mitochondrial
evidence and found no consistent genetic differentiation
associated with each ecotype  but failed to
demonstrate genetic homogeneity in the nuclear genome. The
only nuclear locus (ITS 1) was excluded from the final
analysis in  due to conflicting information and
possible paralogy of sequences. The inference of genetic
homogeneity in  thus rested exclusively on two fully
linked mitochondrial markers, COI and ND1 (the third
mitochondrial marker 16S yielded identical sequences
among all individuals). In the absence of
recombination, mitochondrial genes are vulnerable to
introgressive hybridization, sex-biased dispersal, incomplete
lineage sorting and heteroplasmy [12, 76–79]. The
determination of a ‘barcoding gap’ (i.e., significant
difference between inter- and intraspecific variation) may
fail in case of close phylogenetic relationship or recent
However, the shortcomings of previous analyses 
were addressed by our more expansive sampling and
the inclusion of multiple unlinked microsatellites, thus
suggesting that the distinct phenotypes in M. gregaria
s. l. are not caused by different genotypes.
A case of phenotypic plasticity
A common caveat to marker-based population genetic
studies in case of no differentiation detected among
populations (i.e., different phenotypes in this case) is
that there may be still unsampled isolated regions of
differentiation within genome. Such ‘genomic islands of
differentiation’ [83, 84] are usually associated with genes
under divergent selection, whilst selectively neutral
markers are not involved [85–87]. This alternative is
hard to falsify and might be true for any marker-based
study in organisms with incompletely known genomes.
Adaptive divergence associated with certain selected
genes has been demonstrated in the presence of gene
flow [88–90]. The availability of genome-wide
sequencing may help identify such individual genes, if they exist
indeed, contributing to the phenotypic differentiation
between the two ecotypes.
Except for the possibility of ‘genomic islands’
underpinning different phenotypes, the different ecotypes
within M. gregaria s. l. are then strongly suggestive of
phenotypic plasticity. The exact nature of a trigger that
determines which of the morphotypes will be expressed
is unknown at present. In similar examples from
parthenogenetic Daphnia and aphids, sharply distinct
morphotypes arise from the same genetic background
[16, 17, 91] and in some examples the environmental
triggers controlling which phenotype is preferentially
expressed are known. The sex of offspring from the
same clutch was found to be determined by temperature
among various gonochoristic organisms (those having
separate sexes), e.g. in invertebrates [92, 93], fishes [7, 94, 95],
turtles  and crocodilians . Dramatically different
morphologies can be expressed in presence or absence of
predators in Daphnia water fleas [98–100], barnacle
Chthamalus fissus , whereas little genetic
correspondence is involved in the predator-induced morphological
changes [13, 14, 101].
Although genomic islands of speciation cannot be
completely ruled out, some anecdotal evidence suggests
that one or several as yet unknown environmental
factors may be involved in the determination of Munida
ecotypes. In its South American distribution, Munida
gregaria s. l. occupies extensive latitudinal distribution
along both coasts of Patagonia and wide bathymetry
(from water surface down to 1137 m recorded for
subrugosa [102, 103]), which involves a strong gradient of
environmental conditions (temperature, salinity, oxygen
concentration and food resources). In some species the
feeding performance and diet composition during larval
phases can induce development into different
morphotypes or sex reversal [104, 105]. Since gregaria s. str. and
subrugosa differ in feeding habit as deposit feeders and
actively swimming planktonic feeder, respectively ,
changes in environmental food composition may affect
the metamorphosis of M. gregaria s. l. in an adaptive
fashion, favouring its development into one ecotype
rather than the other. Long-term observations of the
proportion of both gregaria s. str. and subrugosa ecotypes
in the Beagle Channel and San Jorge Gulf demonstrate
the existence of ecotypes is patchy and not stable over
time (see  and references therein). Recent
hydroacoustical evidence postulates that major pelagic swarms
of gregaria s. str. on the Argentine continental shelf are
associated with productive areas such as frontal zones
that vary considerably in spatial and temporal scales
, implying the availability of phytoplankton in
frontal zones might favor the expression of gregaria s.
Heterochrony, which is generally defined as a
developmental change in relation to size and shape in the
timing or rate of ontogenic events (see review in ),
might be a possible mechanism involved in the
observed plasticity in M. gregaria s. l.. Heterochronic
process such as paedomorphic plasticity was postulated
in a widespread squat lobster in the Pacific of South
America, Pleuroncodes monodon . A clear
boundary exists in its distribution where to the north it is a
smaller, pelagic form and to the south it is a larger,
benthic form. Like gregaria s. str. and subrugosa, these
two forms showed no mitochondrial DNA
differentiation either. A similar developmental variation might
be involved in M. gregaria s. l., since the population
from San Jorge Gulf was shown to have faster growth
rate and earlier reproductive investment in its early life
history than the southern populations from Beagle
Channel and Strait of Magellan .
Whether or not the two ecotypes spring from a
genetically entirely homogenous background or whether small
localized genomic islands associated with each exist, our
data have made it abundantly clear that the simple,
perhaps too simple, model of a genome-wide 1:1
relationship between the genotype and an associated phenotype
(however ill equipped we may be to recognize the latter)
does not apply to the Munida gregaria case.
A next-generation sequencing approach  with
higher number of loci and vastly improved coverage of
the genome is a promising way to determine if islands of
genetic differentiation associated with the ecotypes are
involved or if the trigger determining the expression of
one or the other morphotype from identical genotypes
may be independent of genetic differentiation and under
the control of an extrinsic factor.
Based on extensive sampling of the species’ distribution
in South America and using nine independent
polymorphic nuclear microsatellite loci in addition to new
mitochondrial COI sequences, we were able to show
that the lack of genetic differentiation between distinct
gregaria s. str. and subrugosa ecotypes is not an artefact
due to insufficient genomic and/or geographic
sampling or slowly evolving markers. Instead they are likely
expressed from a single underlying genotype although
two largely identical genotypes with interspersed
localized genomic islands of differentiation cannot be
fully ruled out without a more complete coverage of
the genome. Morphological tests affirmed the
boundaries between the two ecotypes were not blurred with
continental-scale geographic sampling, and remain
stable despite an ontogenetic dimension in the data.
These findings corroborate the current taxonomic view
of M. gregaria s. l. (Fabricius, 1793) as a single,
dimorphic species, thus demonstrating a pattern very
unlike cryptic speciation commonly found in DNA
taxonomy and DNA barcoding studies. Our study also
emphasizes the necessity of incorporating
complementary nuclear multi-locus markers in studies aiming at
taxonomy and genotype-phenotype relationship, in
view of the increasing numbers of reported discordance
between mtDNA and nuclear DNA. M. gregaria is
developing into a model affording deeper insights into
the phenotype-genotype relationship, environmental
control of ontogeny and ultimately into the process of
Additional file 1: Table S1. Genetic diversity of COI sequences per
ecotype. (DOCX 15 kb)
Additional file 2: Figure S1. Haplotype genealogies for gregaria s. str.
and subrugosa ecotypes. Each branch represents one substitution; small
filled circles represent hypothetical, unsampled haplotypes. Radii reflect
number of individuals that share a particular haplotype. FM, Falklands/
Malvinas; NCP, Northern Chilean Patagonia; TdF, Tierra del Fuego
archipelago; MdP, Mar del Plata. (EPS 450 kb)
Additional file 3: Table S2. Diversity indices of nine microsatellite loci
for the two ecotypes. Reported are number of alleles nA, fragment size
range, observed heterozygosity HO, expected heterozygosity HE and
allelic richness Ar. Significant deviation from Hardy-Weinberg equilibruim
(P < 0.05, based on 10,000 permutations) after Bonferroni correction
were labeled in bold. (DOCX 19 kb)
Additional file 4: Figure S2. Allelic frequencies (in logarithm) of the nine
microsatellite loci polymorphism. (EPS 1291 kb)
Additional file 5: Figure S3. Tests of statistical power for microsatellite
data set as inferred with POWSIM 4.1. (EPS 1184 kb)
Additional file 6: Figure S4. Correlation between genetic distances
(Pairwise FST values for COI and pairwise (δμ)2 values for microsatellites)
and log-transformed geographical distances for mitochondrial and
microsatellite data for specimens from 25 sampling sites listed in
Table 1. (EPS 1500 kb)
Additional file 7: Data S1. Cytochrome c oxidase subunit I (COI) sequence
alignment file. (NEXUS 59 kb)
Additional file 8: Data S2. Microsatellite data file for each ecotype under
three major sampling areas. (TXT 17 kb)
We are very grateful to Erika Mutschke (Universidad de Magallanes, Punta
Arenas), Felipe González (Reserva Añihue, Chile), Kareen Schnabel (NIWA,
New Zealand) and Fernando L. Mantelatto (University of São Paulo, Brazil)
for providing samples. We thank the Reserva Añihue for equipment support
and their hospitality. We also thank Andrea Eschbach for technical assistance
in the lab and Florian Leese for his help in the field.
C.W. was supported by the Chinese Scholarship Council (CSC grant Nr.
2009633009). This is publication no. 148 of Huinay Scientific Field Station.
Munida around the Falklands/Malvinas were collected by C.H. during the
ICEFISH expedition supported by NSF grant OPP 01–32032 to H. William
Detrich (Northeastern University, Boston, USA).
Availability of data and materials
All DNA sequences from this study are available on GenBank, accession
numbers KJ544249 - KJ544346. Alignment of COI sequences from individuals
of both ecotypes is provided in Additional file 7: Data S1. Fragment lengths
of nine microsatellite loci are provided in Additional file 8: Data S2.
CH and CW conceived and designed the study. The field work was carried
out by CH, CW, JL and VH. CW carried out the experiments and the analyses.
SA and CH provided guidance in data analysis. CW and CH wrote the
manuscript, SA contributed helpful comments. All authors read and
approved the final manuscript.
The authors declare that they have no competing interests.
Consent for publication
1. Knowlton N. Molecular genetic analyses of species boundaries in the sea . Hydrobiologia . 2000 ; 420 ( 1 ): 73 - 90 .
2. Hebert PDN , Penton EH , Burns JM , Janzen DH , Hallwachs W. Ten species in one: DNA barcoding reveals cryptic species in the neotropical skipper butterfly Astraptes fulgerator . Proc Natl Acad Sci U S A . 2004 ; 101 ( 41 ): 14812 - 7 .
3. Bickford D , Lohman DJ , Sodhi NS , Ng PK , Meier R , Winker K , Ingram KK , Das I. Cryptic species as a window on diversity and conservation . Trends Ecol Evol . 2007 ; 22 ( 3 ): 148 - 55 .
4. Saez AG , Lozano E. Body doubles . Nature . 2005 ; 433 ( 7022 ): 111 .
5. Krabbe K , Leese F , Mayer C , Tollrian R , Held C. Cryptic mitochondrial lineages in the widespread pycnogonid Colossendeis megalonyx Hoek, 1881 from Antarctic and Subantarctic waters . Polar Biol . 2009 ; 33 ( 3 ): 281 - 92 .
6. Vanelslander B , Créach V , Vanormelingen P , Ernst A , Chepurnov VA , Sahan E , Muyzer G , Stal LJ , Vyverman W , Sabbe K. Ecological differentiation between sympatric pseudocryptic species in the estuarine benthic diatom Navicula phyllepta (bacillariophyceae) . J Phycol . 2009 ; 45 ( 6 ): 1278 - 89 .
7. Römer U , Beisenherz W. Environmental determination of sex in Apistogrammai (Cichlidae) and two other freshwater fishes (Teleostei) . J Fish Biol . 1996 ; 48 ( 4 ): 714 - 25 .
8. Padilla DK . Inducible phenotypic plasticity of the radula in Lacuna (Gastropoda: Littorinidae) . Veliger . 1998 ; 41 ( 2 ): 201 - 4 .
9. Jarrett JN . Predator-induced defense in the barnacle . J Crustacean Biol . 2009 ; 29 ( 3 ): 329 - 33 .
10. Zakas C , Wares JP . Consequences of a poecilogonous life history for genetic structure in coastal populations of the polychaete Streblospio benedicti . Mol Ecol . 2012 ; 21 ( 22 ): 5447 - 60 .
11. Harris KDM , Bartlett NJ , Lloyd VK . Daphnia as an emerging epigenetic model organism . Genet Res Int . 2012 ; 2012 :8.
12. Toews DPL , Brelsford A. The biogeography of mitochondrial and nuclear discordance in animals . Mol Ecol . 2012 ; 21 ( 16 ): 3907 - 30 .
13. Spanier KI , Leese F , Mayer C , Colbourne JK , Gilbert D , Pfrender ME , Tollrian R. Predator-induced defences in Daphnia pulex: Selection and evaluation of internal reference genes for gene expression studies with real-time PCR . BMC Mol Biol . 2010 ; 11 ( 1 ): 1 - 11 .
14. Tollrian R , Leese F. Ecological genomics: steps towards unraveling the genetic basis of inducible defenses in Daphnia . BMC Biol . 2010 ; 8 ( 1 ): 1 - 4 .
15. Colbourne JK , Pfrender ME , Gilbert D , Thomas WK , Tucker A , Oakley TH , Tokishita S , Aerts A , Arnold GJ , Basu MK , et al. The ecoresponsive genome of Daphnia pulex. Science . 2011 ; 331 ( 6017 ): 555 - 61 .
16. Müller CB , Williams IS , Hardie J. The role of nutrition, crowding and interspecific interactions in the development of winged aphids . Ecol Entomol . 2001 ; 26 ( 3 ): 330 - 40 .
17. Dombrovsky A , Arthaud L , Ledger TN , Tares S , Robichon A. Profiling the repertoire of phenotypes influenced by environmental cues that occur during asexual reproduction . Genome Res . 2009 ; 19 ( 11 ): 2052 - 63 .
18. Vogt G , Huber M , Thiemann M , van den Boogaart G , Schmitz OJ , Schubart CD . Production of different phenotypes from the same genotype in the same environment by developmental variation . J Exp Biol . 2008 ; 211 ( 4 ): 510 - 23 .
19. Loughry WJ , McDonough CM . Phenotypic variablity within and between litters of nine-banded armadillos . Southeastern Naturalist . 2002 ; 1 ( 3 ): 287 - 98 .
20. Astauroff BL . Analyse der erblichen Störungsfälle der bilateralen Symmetrie . Z Indukt Abstamm Vererbungsl . 1930 ; 55 ( 1 ): 183 - 262 .
21. Archer GS , Dindot S , Friend TH , Walker S , Zaunbrecher G , Lawhorn B , Piedrahita JA . Hierarchical phenotypic and epigenetic variation in cloned swine . Biol Reprod . 2003 ; 69 ( 2 ): 430 - 6 .
22. Baba K , Macpherson E , Poore GC , Ahyong ST , Bermudez A , Cabezas P , Lin C-W , Nizinski M , Rodrigues C , Schnabel KE . Catalogue of squat lobsters of the world (Crustacea: Decapoda: Anomura-families Chirostylidae, Galatheidae and Kiwaidae) . Zootaxa . 1905 ; 2008 : 1 - 220 .
23. Zeldis JR . Ecology of Munida gregaria (Decapoda, Anomura): distribution and abundance, population dynamics and fisheries . Mar Ecol Prog Ser Oldendorf . 1985 ; 22 ( 1 ): 77 - 99 .
24. Retamal MA . Catálogo ilustrado de los crustáceos decápodos de Chile . Gayana Zoología . 1981 ; 44 : 1 - 110 .
25. Boschi E , Fischbach C , Iorio M. Catálogo ilustrado de los crustáceos estomatópodos y decápodos marinos de Argentina [Illustrated catalogue of stomatopod and decapod crustraceans of the marine littoral of Argentina] . Frente Marítimo . 1992 ; 10 , sec. A: 7 - 94 .
26. Hendrickx ME . The temperate species of the genus Munida Leach (Crustacea, Decapoda, Galatheidae) in the east Pacific, with the description of a new species and additional records for tropical-subtropical species . Bulletin de l'Institut royal des Sciences naturelles de Belgique . 2003 ; 73 : 115 - 36 .
27. Tapella F , Lovrich GA . Morphological differences between 'subrugosa' and 'gregaria' morphs of adult Munida (Decapoda: Anomura: Galatheidae) from the Beagle Channel, southern South America . J Mar Biol Assoc UK . 2006 ; 86 ( 5 ): 1149 - 55 .
28. Tapella F , Lovrich GA , Romero MC , Thatje S. Reproductive biology of the crab Munida subrugosa (Decapoda: Anomura: Galatheidae) in the Beagle Channel , Argentina . J Mar Biol Assoc UK . 2002 ; 82 ( 4 ): 589 - 95 .
29. Pérez-Barros P , D'Amato ME , Guzmán NV , Lovrich GA . Taxonomic status of two South American sympatric squat lobsters, Munida gregaria and Munida subrugosa (Crustacea: Decapoda: Galatheidae), challenged by DNA sequence information . Biol J Linn Soc . 2008 ; 94 ( 2 ): 421 - 34 .
30. Wang C , Held C. Isolation and characterization of 11 microsatellite markers from the squat lobster Munida gregaria (Decapoda: Galatheidae) around the Falkland Islands/Islas Malvinas . Conserv Genet Resour . 2015 ; 7 ( 1 ): 147 - 9 .
31. Hammer Ø , Harper DAT , Ryan PD . PAST: Paleontological Statistics Software Package for Education and Data Analysis . Palaeontol Electronica . 2001 ; 4 ( 1 ):9. http://folk.uio.no/ohammer/past.
32. R Core Team . R: A language and environment for statistical computing . Vienna: R Foundation for Statistical Computing; 2015 .
33. Folmer O , Black M , Hoeh W , Lutz R , Vrijenhoek R. DNA primers for amplification of mitochondrial cytochrome c oxidase subunit I from diverse metazoan invertebrates . Mol Mar Biol Biotechnol . 1994 ; 3 ( 5 ): 294 - 9 .
34. Librado P , Rozas J. DnaSP v5: a software for comprehensive analysis of DNA polymorphism data . Bioinformatics . 2009 ; 25 ( 11 ): 1451 - 2 .
35. Clement M , Posada D , Crandall KA . TCS: a computer program to estimate gene genealogies . Mol Ecol . 2000 ; 9 ( 10 ): 1657 - 9 .
36. Tamura K , Peterson D , Peterson N , Stecher G , Nei M , Kumar S. MEGA5: molecular evolutionary genetics analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods . Mol Biol Evol . 2011 ; 28 ( 10 ): 2731 - 9 .
37. Zhang J , Kapli P , Pavlidis P , Stamatakis A. A general species delimitation method with applications to phylogenetic placements . Bioinformatics . 2013 ; 29 ( 22 ): 2869 - 76 .
38. Pons J , Barraclough TG , Gomez-Zurita J , Cardoso A , Duran DP , Hazell S , Kamoun S , Sumlin WD , Vogler AP . Sequence-based species delimitation for the DNA taxonomy of undescribed insects . Systematic Biol . 2006 ; 55 ( 4 ): 595 - 609 .
39. Monaghan MT , Wild R , Elliot M , Fujisawa T , Balke M , Inward DJG , Lees DC , Ranaivosolo R , Eggleton P , Barraclough TG , et al. Accelerated species inventory on madagascar using coalescent-based models of species delineation . Systematic Biol . 2009 ; 58 ( 3 ): 298 - 311 .
40. Darriba D , Taboada GL , Doallo R , Posada D. jModelTest 2: more models, new heuristics and parallel computing . Nat Methods . 2012 ; 9 ( 8 ): 772 .
41. Guindon S , Gascuel O. A simple, fast, and accurate algorithm to estimate large phylogenies by maximum likelihood . Systematic Biol . 2003 ; 52 ( 5 ): 696 - 704 .
42. Miller MA , Pfeiffer W , Schwartz T. Creating the CIPRES Science Gateway for inference of large phylogenetic trees . In: Gateway Computing Environments Workshop (GCE). New Orleans: 2010 . p. 1 - 8 .
43. Stamatakis A. RAxML version 8: a tool for phylogenetic analysis and postanalysis of large phylogenies . Bioinformatics . 2014 ; 30 ( 9 ): 1312 - 3 .
44. Bouckaert R , Heled J , Kühnert D , Vaughan T , Wu C-H , Xie D , Suchard MA , Rambaut A , Drummond AJ . BEAST 2: A software platform for bayesian evolutionary analysis . Plos Comput Biol . 2014 ; 10 ( 4 ): e1003537 .
45. Schubart CD , Diesel R , Hedges SB . Rapid evolution to terrestrial life in Jamaican crabs . Nature . 1998 ; 393 ( 6683 ): 363 - 5 .
46. Wares JP , Cunningham CW . Phylogeography and historical ecology of the north Atlantic intertidal . Evolution . 2001 ; 55 ( 12 ): 2455 - 69 .
47. Drummond A , Rambaut A. BEAST : Bayesian evolutionary analysis by sampling trees . BMC Evol Biol . 2007 ; 7 : 214 .
48. Van Oosterhout C , Hutchinson WF , Wills DPM , Shipley P. Micro- Checker: software for identifying and correcting genotyping errors in microsatellite data . Mol Ecol Notes . 2004 ; 4 ( 3 ): 535 - 8 .
49. Goudet J. FSTAT, a program to estimate and test gene diversities and fixation indices (version 2.9.3) . Lausanne: Lausanne University; 2001 .
50. Szpiech ZA , Jakobsson M , Rosenberg NA . ADZE: a rarefaction approach for counting alleles private to combinations of populations . Bioinformatics . 2008 ; 24 ( 21 ): 2498 - 504 .
51. Raymond M , Rousset F. An exact test for population differentiation . Evolution . 1995 ; 49 ( 6 ): 1280 - 3 .
52. Antao T , Lopes A , Lopes R , Beja-Pereira A , Luikart G . LOSITAN: a workbench to detect molecular adaptation based on a Fst-outlier method . Bmc Bioinformatics . 2008 ; 9 ( 1 ): 323 .
53. Ryman N , Palm S. POWSIM: a computer program for assessing statistical power when testing for genetic differentiation . Mol Ecol Notes . 2006 ; 6 ( 3 ): 600 - 2 .
54. Excoffier L , Lischer HE . Arlequin suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and Windows . Mol Ecol Resour . 2010 ; 10 ( 3 ): 564 - 7 .
55. Peakall R , Smouse PE . GENALEX 6: genetic analysis in Excel. Population genetic software for teaching and research . Mol Ecol Notes . 2006 ; 6 ( 1 ): 288 - 95 .
56. Pritchard JK , Stephens M , Donnelly P. Inference of population structure using multilocus genotype data . Genetics . 2000 ; 155 ( 2 ): 945 - 59 .
57. Earl D , von Holdt B. STRUCTURE HARVESTER: a website and program for visualizing STRUCTURE output and implementing the Evanno method . Conserv Genet Resour . 2012 ; 4 ( 2 ): 359 - 61 .
58. Evanno G , Regnaut S , Goudet J. Detecting the number of clusters of individuals using the software STRUCTURE: a simulation study . Mol Ecol . 2005 ; 14 ( 8 ): 2611 - 20 .
59. Jakobsson M , Rosenberg NA . CLUMPP: a cluster matching and permutation program for dealing with label switching and multimodality in analysis of population structure . Bioinformatics . 2007 ; 23 ( 14 ): 1801 - 6 .
60. Rosenberg NA . DISTRUCT: a program for the graphical display of population structure . Mol Ecol Notes . 2004 ; 4 ( 1 ): 137 - 8 .
61. Huelsenbeck JP , Andolfatto P. Inference of population structure under a Dirichlet process model . Genetics . 2007 ; 175 ( 4 ): 1787 - 802 .
62. Jensen J , Bohonak A , Kelley S. Isolation by distance, web service . BMC Genet . 2005 ; 6 ( 1 ): 13 .
63. Goldstein DB , Ruiz Linares A , Cavalli-Sforza LL , Feldman MW . Genetic absolute dating based on microsatellites and the origin of modern humans . Proc Natl Acad Sci . 1995 ; 92 ( 15 ): 6723 - 7 .
64. Rousset F. Genetic differentiation and estimation of gene flow from Fstatistics under isolation by distance . Genetics . 1997 ; 145 ( 4 ): 1219 - 28 .
65. Liebers D , De Knijff P , Helbig AJ . The herring gull complex is not a ring species . Proc R Soc Lond B Biol Sci . 2004 ; 271 ( 1542 ): 893 - 901 .
66. Irwin DE , Bensch S , Irwin JH , Price TD . Speciation by distance in a ring species . Science . 2005 ; 307 ( 5708 ): 414 - 6 .
67. Johannsen W. The genotype conception of heredity . Am Nat . 1911 ; 45 ( 531 ): 129 - 59 .
68. Mobley KB , Lussetti D , Johansson F , Englund G , Bokma F. Morphological and genetic divergence in Swedish postglacial stickleback (Pungitius pungitius) populations . BMC Evol Biol . 2011 ; 11 ( 1 ): 1 - 18 .
69. Gómez A , Wright PJ , Lunt DH , Cancino JM , Carvalho GR , Hughes RN . Mating trials validate the use of DNA barcoding to reveal cryptic speciation of a marine bryozoan taxon . Proc R Soc B Biol Sci . 2007 ; 274 ( 1607 ): 199 - 207 .
70. Leese F , Held C. Identification and characterization of microsatellites from the Antarctic isopod Ceratoserolis trilobitoides: nuclear evidence for cryptic species . Conserv Genet . 2008 ; 9 ( 5 ): 1369 - 72 .
71. Hebert PDN , Gregory TR . The Promise of DNA Barcoding for Taxonomy . Systematic Biol . 2005 ; 54 ( 5 ): 852 - 9 .
72. Hajibabaei M , Singer GA , Hebert PD , Hickey DA . DNA barcoding: how it complements taxonomy, molecular phylogenetics and population genetics . Trends in genetics : TIG . 2007 ; 23 ( 4 ): 167 - 72 .
73. Valentini A , Pompanon F , Taberlet P. DNA barcoding for ecologists . Trends Ecol Evol . 2009 ; 24 ( 2 ): 110 - 7 .
74. Galtier N , Nabholz B , Glémin S , Hurst GDD. Mitochondrial DNA as a marker of molecular diversity: a reappraisal . Mol Ecol . 2009 ; 18 ( 22 ): 4541 - 50 .
75. Fisher-Reid MC , Wiens JJ . What are the consequences of combining nuclear and mitochondrial data for phylogenetic analysis? Lessons from Plethodon salamanders and 13 other vertebrate clades . BMC Evol Biol . 2011 ; 11 ( 1 ): 300 .
76. Ballard JWO , Whitlock MC . The incomplete natural history of mitochondria . Mol Ecol . 2004 ; 13 ( 4 ): 729 - 44 .
77. Ebach MC , Holdrege C. DNA barcoding is no substitute for taxonomy . Nature . 2005 ; 434 ( 7034 ): 697 .
78. Rubinoff D , Cameron S , Will K. A genomic perspective on the shortcomings of mitochondrial DNA for “barcoding” identification . Journal of Heredity . 2006 ; 97 ( 6 ): 581 - 94 .
79. Wares JP . Natural distributions of mitochondrial sequence diversity support new null hypotheses . Evolution . 2010 ; 64 ( 4 ): 1136 - 42 .
80. Wiemers M , Fiedler K. Does the DNA barcoding gap exist?-a case study in blue butterflies (Lepidoptera: Lycaenidae) . Front Zool . 2007 ; 4 ( 8 ): 1 - 16 .
81. Moritz C , Cicero C. DNA barcoding: Promise and pitfalls . Plos Biol . 2004 ; 2 ( 10 ): e354 .
82. Meyer CP , Paulay G . DNA barcoding: error rates based on comprehensive sampling . Plos Biol . 2005 ; 3 ( 12 ): e422 .
83. Turner TL , Hahn MW , Nuzhdin SV . Genomic islands of speciation in Anopheles gambiae . Plos Biol . 2005 ; 3 ( 9 ): e285 .
84. Harr B. Genomic islands of differentiation between house mouse subspecies . Genome Res . 2006 ; 16 ( 6 ): 730 - 7 .
85. Wu C-I , Ting C-T. Genes and speciation. Nat Rev Genet . 2004 ; 5 ( 2 ): 114 - 22 .
86. Gavrilets S , Vose A. Dynamic patterns of adaptive radiation . Proc Natl Acad Sci USA . 2005 ; 102 ( 50 ): 18040 - 5 .
87. Nosil P , Funk DJ , Ortiz-Barrientos D. Divergent selection and heterogeneous genomic divergence . Mol Ecol . 2009 ; 18 ( 3 ): 375 - 402 .
88. Hemmer-Hansen J , Nielsen EE , Frydenberg J , Loeschcke V. Adaptive divergence in a high gene flow environment: Hsc70 variation in the European flounder ( Platichthys flesus L.). Heredity . 2007 ; 99 ( 6 ): 592 - 600 .
89. Larmuseau MHD , Huyse T , Vancampenhout K , Van Houdt JKJ , Volckaert FAM . High molecular diversity in the rhodopsin gene in closely related goby fishes: A role for visual pigments in adaptive speciation? Mol Phylogenet Evol . 2010 ; 55 ( 2 ): 689 - 98 .
90. Limborg MT , Helyar SJ , De Bruyn M , Taylor MI , Nielsen EE , Ogden ROB , Carvalho GR , Consortium FPT , Bekkevold D. Environmental selection on transcriptome-derived SNPs in a high gene flow marine fish, the Atlantic herring (Clupea harengus) . Mol Ecol . 2012 ; 21 ( 15 ): 3686 - 703 .
91. Stibor H , Lüning J. Predator-induced phenotypic variation in the pattern of growth and reproduction in Daphnia hyalina (Crustacea: Cladocera) . Funct Ecol . 1994 ; 8 ( 1 ): 97 - 101 .
92. Cline TW . A sex-specific, temperature-sensitive maternal effect of the daughterless mutation of Drosophila melanogaster . Genetics . 1976 ; 84 ( 4 ): 723 - 42 .
93. Klass M , Wolf N , Hirsh D. Development of the male reproductive system and sexual transformation in the nematode Caenorhabditis elegans . Dev Biol . 1976 ; 52 ( 1 ): 1 - 18 .
94. Baras E , Prignon C , Gohoungo G , Méalard C. Phenotypic sex differentiation of blue tilapia under constant and fluctuating thermal regimes and its adaptive and evolutionary implications . J Fish Biol . 2000 ; 57 ( 1 ): 210 - 23 .
95. Koumoundouros G , Pavlidis M , Anezaki L , Kokkari C , Sterioti A , Divanach P , Kentouri M. Temperature sex determination in the European sea bass , Dicentrarchus labrax (L. , 1758 ) (Teleostei, Perciformes, Moronidae): Critical sensitive ontogenetic phase . J Exp Zool . 2002 ; 292 ( 6 ): 573 - 9 .
96. Bull J , Vogt R. Temperature-dependent sex determination in turtles . Science . 1979 ; 206 ( 4423 ): 1186 - 8 .
97. Lang JW , Andrews HV . Temperature-dependent sex determination in crocodilians . J Exp Zool . 1994 ; 270 ( 1 ): 28 - 44 .
98. Green J. The distribution and variation of Daphnia lumholtzi (Crustacea: Cladocera) in relation to fish predation in Lake Albert , East Africa . J Zool . 1967 ; 151 ( 1 ): 181 - 97 .
99. Parejko K , Dodson SI . The evolutionary ecology of an antipredator reaction norm: Daphnia pulex and Chaoborus americanus . Evolution . 1991 ; 45 ( 7 ): 1665 - 74 .
100. Tollrian R , Harvell CD . The ecology and evolution of inducible defenses . New Jersey : Princeton University Press ; 1999 .
101. Jarrett JN . Inter-population variation in shell morphology of the barnacle . J Crustacean Biol . 2008 ; 28 ( 1 ): 16 - 20 .
102. Retamal MA. Los Decápodos de Chile. 1994 .
103. Arntz WE , Gorny M , Soto R , Lardies MA , Retamal M , Wehrtmann IS . Species composition and distribution of decapod crustaceans in the waters off Patagonia and Tierra del Fuego, South America . Sci Mar . 1999 ; 63 Suppl 1 : 303 - 14 .
104. Parsons KJ , Robinson BW . Foraging performance of diet-induced morphotypes in pumpkinseed sunfish (Lepomis gibbosus) favours resource polymorphism . J Evolution Biol . 2007 ; 20 ( 2 ): 673 - 84 .
105. Zupo V. Effect of microalgal food on the sex reversal of Hippolyte inermis (Crustacea: Decapoda) . Mar Ecol Prog Ser . 2000 ; 201 : 251 - 9 .
106. Romero MC , Lovrich GA , Tapella F , Thatje S. Feeding ecology of the crab Munida subrugosa (Decapoda: Anomura: Galatheidae) in the Beagle Channel , Argentina . J Mar Biol Assoc UK . 2004 ; 84 ( 2 ): 359 - 65 .
107. Diez MJ , Cabreira AG , Madirolas A , Lovrich GA . Hydroacoustical evidence of the expansion of pelagic swarms of Munida gregaria (Decapoda, Munididae) in the Beagle Channel and the Argentine Patagonian Shelf, and its relationship with habitat features . J Sea Res . 2016 ; 114 : 1 - 12 .
108. Smith KK . Heterochrony revisited: the evolution of developmental sequences . Biol J Linn Soc . 2001 ; 73 ( 2 ): 169 - 86 .
109. Haye P , Salinas P , Acuna E , Poulin E. Heterochronic phenotypic plasticity with lack of genetic differentiation in the southeastern pacific squat lobster pleuroncodes monodon . Evol Dev . 2010 ; 12 : 627 - 33 .
110. Varisco M , Vinuesa JH . Growth and reproduction investment of the young of the year of the squat lobster Munida gregaria (Crustacea: Anomura) in the Patagonian coast . Sci Mar . 2015 ; 79 ( 3 ): 345 - 53 .
111. Martin SH , Dasmahapatra KK , Nadeau NJ , Salazar C , Walters JR , Simpson F , Blaxter M , Manica A , Mallet J , Jiggins CD . Genome-wide evidence for speciation with gene flow in Heliconius butterflies . Genome Res . 2013 ; 23 ( 11 ): 1817 - 28 .