Skip to main content


Rapid microevolution during recent range expansion to harsh environments



Adaptive evolution is one of the crucial mechanisms for organisms to survive and thrive in new environments. Recent studies suggest that adaptive evolution could rapidly occur in species to respond to novel environments or environmental challenges during range expansion. However, for environmental adaptation, many studies successfully detected phenotypic features associated with local environments, but did not provide ample genetic evidence on microevolutionary dynamics. It is therefore crucial to thoroughly investigate the genetic basis of rapid microevolution in response to environmental changes, in particular on what genes and associated variation are responsible for environmental challenges. Here, we genotyped genome-wide gene-associated microsatellites to detect genetic signatures of rapid microevolution of a marine tunicate invader, Ciona robusta, during recent range expansion to the harsh environment in the Red Sea.


The Red Sea population was significantly differentiated from the other global populations. The genome-wide scan, as well as multiple analytical methods, successfully identified a set of adaptive genes. Interestingly, the allele frequency largely varied at several adaptive loci in the Red Sea population, and we found significant correlations between allele frequency and local environmental factors at these adaptive loci. Furthermore, a set of genes were annotated to get involved in local temperature and salinity adaptation, and the identified adaptive genes may largely contribute to the invasion success to harsh environments.


All the evidence obtained in this study clearly showed that environment-driven selection had left detectable signatures in the genome of Ciona robusta within a few generations. Such a rapid microevolutionary process is largely responsible for the harsh environmental adaptation and therefore contributes to invasion success in different aquatic ecosystems with largely varied environmental factors.


Dissecting evolutionary dynamics of environmental adaptation in organisms is one of the fundamental research themes in evolutionary biology and ecology, particularly in the context of rapid global climate change and distributional shift of organisms in the past five decades [1,2,3]. Species spreading to and subsequently successfully colonizing new regions are, or will be, inevitably suffering from environmental changes/pressures [4, 5]. Adaptive evolution is one of the crucial mechanisms to facilitate species to survive and thrive in novel environments [6,7,8]. Recent studies suggest that adaptive evolution could rapidly occur in species (i.e., rapid microevolution) to respond to novel environments or environmental challenges during range expansion [9,10,11,12].

Rapid microevolution is usually associated with changes of genetic variation in a population, such as changes of allele frequency and generation of new alleles through mutations at functional genes [7, 13, 14]. Multiple studies clearly showed that evolutionary changes could occur within dozens of generations, as those studies detected significant evolutionary signatures using genetic markers in populations living in varied environments [15,16,17]. Available evidence suggests that such a rapid microevolutionary process could be driven by multiple factors and their complex interactions, such as selection (natural or artificial), gene flow, genetic drift and mutations [15, 16, 18, 19]. For example, rapid microevolution occurred in three-spine sticklebacks within less than 50 generations during habitat transitions from marine to freshwater habitats, and many important functional loci and genes associated with adaptive traits were involved in the process of rapid environmental adaptation [20, 21]. However, for environmental adaptation, most studies successfully detected phenotypic features associated with local environments, such as changes of morphological characteristics and life history traits, but did not provide ample genetic evidence on microevolutionary dynamics [15, 22]. Thus, it is crucial to deeply investigate the genetic basis of rapid microevolution in response to environmental adaptation, in particular on what genes and associated variation are responsible for environmental changes/challenges during range expansion [23,24,25].

Biological invasions are promising ‘natural’ experiments that provide opportunities to study rapid microevolutionary processes [10, 26, 27]. Invasive species inevitably face extremely rapid, or even sudden, environmental changes during range expansion, thus providing good materials to deeply investigate how species rapidly adapt to novel environments [26, 28, 29]. The highly invasive ascidian, Ciona robusta, is a notoriously fouling organism with wide salinity tolerance ranging from 12‰ to 40‰ and temperature varying from 3 °C to 30 °C [30, 31]. C. robusta is usually considered to be native to Northwest Pacific and has invaded the coasts of Mediterranean Sea, Atlantic, Pacific and Oceania oceans [31,32,33]. Most invasion events of this species have been recorded in temperate and sub-tropical waters, though a few populations appeared in tropical harbours with transitory states [34, 35]. In different environments, the life cycle characteristics of C. robusta, such as the growth rate, life span and spawning time, vary depending on local environmental factors, particularly the water temperature and salinity [30]. For example, populations in cold and temperate coastal regions could reproduce 2–3 generations per year [30, 36]; however, more than 4 generations/year were observed in warm waters, such as those in warm regions of the Mediterranean Sea [30].

A recent introduction of C. robusta has been reported in the Red Sea, a representative tropical region which is well known for its high temperature and salinity [35]. The first occurrence of C. robusta was detected in the Eilat marina in 2015 [35]. Many environmental factors in the Red Sea, such as water temperature (> 27 °C) and salinity (> 40‰, Table 1), represent the extreme in the distribution ranges of C. robusta reported so far. Genetic analyses based on mitochondrial DNA suggested that the Red Sea population might be introduced from the Mediterranean Sea through the Suez Canal [35]. Since the opening of the Suez Canal in 1869, a large number of non-native species have been introduced from the Red Sea to the Mediterranean Sea [37, 38]. For example, more than 90 fishes have been introduced through the Lessepsian migration, a process of biological invasions from the Red Sea into the Mediterranean Sea [37,38,39]. However, only a handful of opposite-direction introductions from the Mediterranean to the Red Sea (anti-Lessepsian migration) have been recorded so far [40], mainly owing to harsh environmental conditions in the Red Sea, such as high salinity, high temperature and oligotrophic conditions [41, 42]. For the Mediterranean C. robusta populations, the optimal temperature range was 14–23.4 °C [43]; however, the Red Sea population could survive and reproduce in such harsh environmental conditions, indicating a wider environmental tolerance than previously assumed. Thus, the recent successful range expansion provides good materials to study genetic signatures of rapid microevolution, particularly on genetic mechanisms of rapid local adaptation to severe environmental conditions in such a short timescale in the wild.

Table 1 Sampling sites and genetic diversity indices based on genome-wide gene-associated microsatellites of Ciona robusta

In this study, we sampled the recently established population of C. robusta in the Red Sea and genotyped genome-wide gene-associated microsatellites to reveal the genetic signatures of rapid microevolution in the harsh environment. In order to cover various temperature and salinity gradients globally, we included five other populations (two Mediterranean populations, one Atlantic population, two Pacific populations) from Lin et al. [44]. At the genome level, we aimed to (i) investigate the population genetic patterns potentially shaped by local environmental conditions in the Red Sea, and (ii) detect adaptive genetic variation and genes involved in rapid environmental adaptation in the process of range expansion to harsh environments.


Significant difference of environmental factors

Based on the gradients of sea water temperature and salinity (Table 1) and our previous population genetic studies [32, 33, 44], five other global populations from Lin et al. [44] were divided into three groups. Two populations (AM and BL) from the Mediterranean Sea were assigned into one group, the population SA from the Atlantic Ocean formed another group, and the remaining two populations (NMF and GAP) from the Pacific Ocean were the third group. Non-parametric Kruskal-Wallis test of six environmental factors showed statistical difference between the Red Sea population and other three groups (p < 0.01), suggesting significantly different environmental conditions in the Red Sea. The Red Sea population had significantly higher temperature and salinity than the other populations (Additional file 1: Figure S1, S2).

Population genetic diversity and structure

A total of 22 individuals from the Red Sea were confirmed as C. robusta by molecular identification. Based on 146 genome-wide gene-associated microsatellite loci, we identified 975 alleles across all 22 individuals. Microsatellite diversity of the Red Sea population (i.e. AR = 2.524) was significantly lower than that in the other five global populations (Kruskal-Wallis test, p < 0.01), especially two Pacific populations (NMF, AR = 3.937; GAP, AR = 3.789; Table 1). However, we did not detect the signal of recent population bottleneck with the TPM model, and the allele frequency of the Red Sea population exhibited a normal L-shape.

Significant genetic differentiation was detected between the Red Sea population and the others, with the pairwise FST values ranging from 0.053 to 0.215 (Table 2). The pairwise FST values were relatively low between the Red Sea population and two Mediterranean populations (AM, BL). In comparison, the highest genetic differentiation was found between the Red Sea population and one Pacific population collected from New Zealand (NMF, FST = 0.215; Table 2). The Bayesian assignment analyses in STRUCTURE based on all microsatellite loci suggested a two-cluster model (optimal K = 2). Individuals from the Red Sea population were assigned into the same cluster formed by three other populations (two Mediterranean populations: AM, BL; one Atlantic population: SA), whereas two Pacific populations (NMF, GAP; Fig. 1a) were grouped into the other cluster. However, the Red Sea population was clearly separated from the cluster of Mediterranean-Atlantic populations at higher K values (e.g., K = 3, 4, 5; Fig. 1a). Furthermore, the results confirmed that the Red Sea population exhibited genetically different structure when the Red Sea and Mediterranean-Atlantic group was re-analyzed separately (Fig. 1b).

Table 2 Population genetic differentiation (pairwise FST) based on genome-wide gene-associated microsatellites of Ciona robusta
Fig. 1

Bayesian inference population genetic structure of Ciona robusta in STRUCTURE. a K values from 2 to 5 based on the Red Sea population and other five populations; b K values from 2 to 3 based on four populations from the Red Sea, the Mediterranean sea and the Atlantic ocean, respectively. Each genotype is represented by a thin line, with proportional membership in different clusters indicated by different colors. Bold vertical lines separate collection sites, with sites ID as per Table 1

Signatures of rapid local adaptation

Here, we used two genome scan approaches to detect the footprints of selection. Firstly, we identified a total of 32 outlier loci (21.9%) in ARLEQUIN, including 17 loci (11.6%) for directional selection and 15 loci (10.3%) for balancing selection (Additional file 1: Figure S3). Secondly, 18 outlier loci (12.3%) were detected based on the BAYESCAN approach, with 12 loci (8.2%) under directional selection and six loci (4.1%) under balancing selection (Fig. 2a). Interestingly, five out of nine identified outlier loci in BAYESCAN were located on the chromosome 1, suggesting a candidate genomic region under selection (Additional file 1: Figure S4). In total, 14 (9.6%) loci were commonly identified as FST outliers by both approaches (Fig. 2b).

Fig. 2

Microsatellite loci under selection. a, FST-based outliers detected by BAYESCAN, and the solid vertical line represents false discovery rate of 0.05; b, the number of microsatellite loci under selection identified by two approaches (ARLEQUIN and BAYESCAN), as well as the environmental association analysis (matSAM)

In order to further assess environmental influence on outlier loci, spatial analysis method (matSAM) was performed to test the correlation between allelic frequency and environmental factors. Interestingly, a total of 93 alleles (9.5%) at 70 microsatellite loci (47.9%) showed significant association with at least one of the six environmental parameters (P < 0.05), indicating strong effects of local environmental factors on genetic variation. A total of 41 loci were associated with temperature, including 22, 6 and 33 loci for annual average water temperature (AveT), the highest monthly average water temperature (MaxT) and the lowest monthly average water temperature (MinT), respectively. A total of 57 microsatellite loci were detected with the strong correlation with salinity, including 40, 30 and 56 loci for annual average water salinity (AveS), the highest monthly average water salinity (MaxS) and the lowest monthly average water salinity (MinS), respectively. Altogether, 28 loci showed significant associations with both temperature and salinity. Furthermore, a total of 19 loci detected by ARLEQUIN or BAYESCAN approach were significantly associated with environmental factors (Fig. 2b). We used these 19 loci as adaptive outlier loci for further analyses because they were the best candidates associated with adaptive polymorphisms at neighboring genes potentially under environmental selection. Of these 19 adaptive outliers and within their selective sweep windows, a total of 44 genes were identified (Table 3; Additional file 1: Table S1). After GO term enrichment analysis, 29 GO terms were significantly enriched, such as metabolic process, ion transport, and catalytic activity (Fig. 3). In addition, four KEGG pathways were detected, including adherens junction, cell cycle, ubiquitin mediated proteolysis and Wnt signaling pathway. Across all 19 adaptive outlier loci, we detected allele frequency changes at a total of 35 alleles (Fig. 4). The allele frequency of several adaptive loci in the Red Sea population was changed (Fig. 4). For example, the frequency of allele 207 at the locus Cin162 (frequency = 0.81) was much higher in the Red Sea population than that in the other five populations (frequency = 0–0.23). A similar pattern was observed for the allele 270 at the locus Cin153 (frequency = 0.11) in the Red Sea population, which was much lower than that that in the others (frequency = 0.14–0.86). Interestingly, we detected significantly positive or negative correlations between allele frequency and environmental factors at the adaptive loci. For example, at the locus Cin162, the frequency of allele 207 (Cin162–207) was positively correlated with the minimum water temperature (Pearson’s r = 0.87, p = 0.024; Fig. 5). At the locus Cin211, the frequency of allele 168 (Cin211–168) had a positive correlation with the minimum salinity (Pearson’s r = 0.89, p = 0.019; Fig. 5). At the locus Cin153, the frequency of allele 270 (Cin153–270) was negatively correlated with salinity, including AveS (Pearson’s r = − 0.90, p = 0.015), MaxS (Pearson’s r = − 0.84, p = 0.035) and MinS (Pearson’s r = − 0.95, p = 0.004; Fig. 5). These patterns suggest that local environmental factors, such as temperature and salinity, could largely influence the allele frequency of adaptive loci, thus contributing to rapid microevolution.

Table 3 Gene annotation of 19 adaptive outlier loci against the Ciona robusta genome
Fig. 3

Gene ontology (GO) term enrichment analysis. The GO annotation results were based on 44 genes annotated by 19 adaptive outlier loci and within 20 kb selective sweep windows. Gene ontology categories included molecular function, cellular component and biological process. GO categories for each function were sorted by decreasing order of evidence, based on the GO enrichment test P-value

Fig. 4

Heat map of allele frequencies of 35 alleles obtained from 19 adaptive loci. Rows represent specific alleles, and columns represent different populations. Colours represent normalized allele frequencies

Fig. 5

Pearson correlation tests between allele frequencies at adaptive loci and environmental factors. MinT, AveS, MaxS and MinS refer to the lowest monthly average water temperature, annual average water salinity, the highest monthly average salinity, the lowest monthly average salinity, respectively

Most of these adaptive loci were annotated to genes involved in temperature and salinity adaptation (Table 3; Additional file 1: Table S1). For the adaptive locus Cin162, it was annotated as the mitochondrial uncoupling protein 3 gene (UCP3) located on the Chromosome 9 in C. robusta. This gene is involved in metabolic process and especially adaptive thermogenesis in response to cold stress based on GO annotation according to Uniprot or AmiGo2 database (Table 3). The S-phase kinase-associated protein 2 gene (SKP2), which was the best gene hit for the locus Cin95, plays a key role in mitotic cell cycle and protein ubiquitination. Another locus Cin211 was located within the region of the AP-3 complex subunit mu-1gene (AP3M1), and the AP3M1 gene, which regulates the process of protein metabolism and Rab GTPase binding, is involved in salinity stress tolerance. In addition, two genes were detected in the up-stream of the locus Cin74 (Additional file 1: Table S1): one gene was ubiquitin-conjugating enzyme E2C gene (UBE2C), which was reported to get involved in ubiquitin-dependent protein catabolic process under heat stress in Ciona ascidian, and another one was the telomere length regulation protein TEL2 gene (TELO2) with the functions of DNA damage response and Hsp90 protein binding. Indeed, many other genes play crucial roles in various important biological processes related to responses to environmental stresses, such as cell apoptosis, fatty acid metabolic process, and ion transport (Table 3; Fig. 3; Additional file 1: Table S1). Altogether, we detected functional genes associated with these adaptive loci that were involved in the process of rapid adaptation to the harsh environmental conditions in the Red Sea.


Deep insights into the rapid adaptive changes during range expansion will facilitate our understanding on species’ distribution ranges and responses to future environmental changes [3, 15, 18, 45]. Based on the genome-wide gene-associated microsatellites, we investigated genetic signatures of rapid adaptation to harsh environments of the Red Sea in C. robusta. The Red Sea population was significantly differentiated from the other populations (Table 2; Fig. 1), and the genome-wide scan using multiple methods identified a set of adaptive outlier loci (Fig. 2; Additional file 1: Figure S3). Interestingly, the allele frequency at several adaptive loci in the Red Sea population was largely altered (Fig. 4), and we found significant correlations between allele frequency and local environmental factors at these adaptive loci (Fig. 5). Furthermore, many genes were annotated to get involved in local temperature and salinity adaptation, and biological processes and pathways of adaptive genes may underline the effect of rapid adaptive evolution (Table 3 and Additional file 1: Table S1; Fig. 4). All the evidence obtained in this study clearly showed that environment-driven selection had left detectable signatures in the genome of C. robusta within a few generations, thus supporting the hypothesis that local environmental adaptation can rapidly occur, and such a rapid microevolutionary process may largely contribute to invasion success in aquatic ecosystems (i.e., rapid microevolution hypothesis).

Factors for population genetic structure in the Red Sea

For newly established populations, gene flow and some events such as genetic bottleneck, genotype sorting and drift are likely to influence population genetic structure [46,47,48]. Given the increasing shipping activities, especially through the Suez Canal, human-mediated dispersal could facilitate gene flow between populations of the Red Sea and Mediterranean Sea [16, 41, 49]. According to Shenkar et al. [35], mtDNA analyses of C. robusta samples from the same sampling site in 2015 suggested that the Red Sea population might be introduced from the Mediterranean Sea, where C. robusta were detected much earlier in the end of nineteenth century [31, 44]. In our study, the Red Sea population showed significant genetic divergence with the two populations from the Mediterranean Sea (Table 2; Fig. 1). Such a pattern suggests that gene flow might not occur, or local adaptation has largely overridden the effect of gene flow as we detected strong influence of local environmental factors on both loci and associated allele frequency (Figs. 2, 4 and 5). For novel populations with a small number of colonists, they may experience a dramatic decline in genetic diversity due to bottleneck effects, which can contribute to population genetic divergence by retaining a non-random set of genotypes when compared to their original populations [50]. Although we detected reduced population genetic diversity in the Red Sea population when compared with two Mediterranean populations (AM, BL), further analyses showed no signal of recent bottleneck in the Red Sea population, indicating that genetic bottleneck might not be the major factors for the observed genetic divergence of the Red Sea population. Although we could not exclude the potential influence of other processes, such as genotype sorting and genetic drift, multiple lines of evidence in this study clearly supported that natural selection caused by harsh environments had largely influenced the population genetic patterns. Strong selection pressures by environmental challenges, such as global warming, marine pollution and ocean acidification, can result in high population genetic divergence for marine species [4, 51, 52]. Here, we detected significant environmental difference between the Red Sea population and the others. The Red Sea population was firstly detected in 2015 and has experienced ~ 10 generations after introduction. We identified a large number of adaptive loci which exhibited significant correlation with water temperature and salinity, indicating a strong environment-driven selection in the genome of C. robusta in a short timescale. Similar strong and rapid adaptation was also detected in other marine species. For example, significant reproductive isolation was detected within 13 generations among introduced salmon populations in novel environments [53]. Bernardi et al. [16] reported that rapid adaptive evolution promoted the bluespotted cornetfish to rapidly invade over the Mediterranean Sea in less than 7 years, and they further revealed adaptive genomic regions under salinity selection, such as osmoregulatory regions. Results obtained in this study, as well as evidence from related investigations, suggest that rapid local adaptation could rapidly occur within several generations when selection pressures are high.

Sources of genetic variation in rapid environment-driven adaptation

Rapid local adaptation during biological invasions primarily depends on two sources of genetic variation: standing genetic variation and new beneficial mutations [7, 13, 14, 54]. A large number of studies suggested that rapid local adaptation induced by recent natural selection and/or strong selection pressure predominantly relied on standing genetic variation, by which the beneficial alleles were immediately available at higher frequencies compared to new mutations [10, 55]. For the Red Sea population, standing genetic variation could provide preadapted genotypes, such as the genotypes well-suited to high temperature and salinity, to rapidly adapt to the harsh environment in the Red Sea. In addition, Lin et al. [44] revealed that the length of selective sweeps on the genomic regions in C. robusta could be as narrow as 9.8 kb, suggesting the soft selective sweep with standing genetic variation. However, recent studies documented that novel mutations could provide opportunities, especially for small founding populations, to rapidly adapt to the new environments [56]. C. robusta exhibits a high level of mutation rate, accompanied with relatively high fecundity, rapid growth rate and short life history, which may provide new mutations for natural selection [26, 57]. Thus, we could not rule out the effects of new mutations in environment-driven selection, and more than one type of genetic variation might get involved in the rapid local adaptation in the Red Sea population.

Candidate genes for rapid local adaptation

Understanding the function of adaptive genes can provide insights into the genetic basis of rapid microevolution in changing environments [12, 44]. In our study, many candidate genes were identified with important roles in rapid adaptation to temperature and salinity (Table 3 and Additional file 1: Table S1). For example, one strong candidate locus (Cin162) was annotated to mitochondrial uncoupling protein 3 gene (UCP3), and this gene is especially involved in cold stress response. UCP3-mediated effects on adaptive thermogenesis under cold stress have been discovered in marine fish, mainly related to fat metabolism [58, 59]. Interestingly, the allele frequency of this locus (Cin162–207) was detected to be significantly correlated with the minimum temperature (MinT) with much higher allele frequency in Red Sea population than that in the others. This pattern supported that local environments could largely alter allele frequency of functional genes, thus contributing to the rapid microevolution. Another two important genes involved in high temperature adaptation, ubiquitin-conjugating enzyme E2C gene (UBE2C) and telomere length regulation protein TEL2 gene (TELO2), were detected within the selective sweep window of the locus Cin74. For the UBE2C gene, it encodes the protein that mainly plays a role in ubiquitin-dependent protein catabolic process [60, 61]. According to Lopez et al. [61], the ubiquitin-mediated proteolysis pathway, including ubiquitin-conjugating enzyme E2, was enriched by the proteomic analyses of C. intestinalis in high temperature treatment. In addition, this gene can cooperate with Cullin-type ubiquitin ligase to get involved in spermatogenesis and fertilization in C. intestinalis [60, 62]. Another gene TELO2, which can participate in the cellular resistance to DNA damage, especially interactions with the Hsp90 protein, is likely to play a crucial role in coping with thermal stress in the Red Sea population. In addition, we detected a large number of genes involved in salinity adaption. For example, the AP-3 complex subunit mu-1gene (AP3M1) with a much higher allele frequency of the locus Cin211 in the Red Sea population, encodes the AP-3 protein, which is linked to the Golgi region as an adaptor-related protein complex. The ion transport between the trans-Golgi network and the plasma membrane was required for salinity stress adaptation [63, 64]. This gene is also connected to the RABA1 GTPase region, which can get involved in the regulation of ion homeostasis in Arabidopsis [63, 64]. Two genes, histone acetyltransferase p300 gene (EP300) and low-density lipoprotein receptor-related protein 6 gene (LRP6), were categorized into the Wnt signaling pathway. The Wnt signaling pathway in C. intestinalis includes at least three pathways: the canonical pathway, the planar cell polarity (PCP) pathway and the Wnt/Ca2+ pathway based on the KEGG database (the Kyoto Encyclopedia of Genes and Genomes, KEGG). Multiple studies have revealed that the Wnt signaling pathway played a crucial role in many biological processes, such as cell proliferation [65], embryogenesis [66], and biomineralization [67]. A recent proteomic analysis of the Pacific white shrimp suggested that the Wnt pathway could participate in low salinity adaptation [68]. Thus, further transcriptomic and proteomic investigations on C. robusta under high temperature and salinity in laboratory conditions can contribute to our understanding which adaptive genes and related biological processes and pathways play a key role in rapid microevolution.


Although C. robusta is generally considered as a temperate-water species, it has currently successfully established in a Red Sea marina since 2015. Rapid local adaptation is a crucial mechanism for invasive species to cope with novel environments and/or environmental stresses during range expansion. Although it is well-known that phenotypic plasticity is a crucial mechanism for temporal adjustment during range expansion, while microevolution is more important for successful population establishment in a long-term perspective. Interestingly, the results obtained in our study suggest that rapid microevolution could occur within a few generations. We detected a large number of adaptive loci and candidate genes responsible for temperature and salinity adaptation in C. robusta, and these genes are candidate loci for further studying adaptation dynamics based on single genes or gene networks using C. robusta as a model. Our study confirms the rapid microevolution hypothesis that rapid local adaptation is largely responsible for the harsh environmental adaptation and therefore contributes to invasion success in different aquatic ecosystems with largely varied environmental factors.

Materials and methods

Sampling and species identification

C. robusta samples were collected from the Eilat marina (29°33′11.3″ N, 34°57′35.9″ E) on March 2018 (experiencing ~ 10 generations after introduction). All collected individuals were immediately preserved in absolute ethanol at 4 °C for further genetic analyses. Total genomic DNA was extracted using a DNeasy Blood and Tissue Kit (Qiagen). To confirm the species identification based on morphology, one mtDNA fragment cytochrome c oxidase subunit 3–NADH dehydrogenase subunit 1 (COX3-ND1) [33, 69], was used to perform molecular identification. The COX3-ND1 segment was amplified by the primers TX3F and TN1R following the protocol as described in Zhan et al. [33].

Genome-wide gene-associated microsatellites genotyping

Here, we used a total of 146 genome-wide gene-associated microsatellite makers [44, 70] to genotype the Red Sea population. Specially, six microsatellite markers (Cin54, 76, 92, 125, 210, 213) were excluded from the Lin et al. [44] dataset because of poor PCR amplifications for the Red Sea individuals. PCR amplification and microsatellites genotyping were performed based on the protocol as described in Lin et al. [44, 70]. Amplified fragments were separated on an ABI 3730xl automated sequencer with the GeneScan™ 500 LIZ™ internal size standard. Alleles were scored with GeneMapper® software v.4.0 (Applied Biosystems). For the other five populations (two Mediterranean populations: AM, BL; one Atlantic population: SA; two Pacific population: NMF, GAP; Table 1), the datasets were adopted from our previous study [44].

Water temperature and salinity are among the most crucial environmental variables in marine ecosystems, as these factors largely influence survival, development and many physiological processes of marine species. The Red Sea was characterized with a high level of temperature and salinity [35]. Therefore, the average values of six environmental factors (among 1955–2012) associated with sea surface temperature (SST) and salinity were obtained from the National Oceanic and Atmospheric Administration (NOAA;, including three metrics of SST: annual average water temperature (AveT), the highest monthly average water temperature (MaxT), the lowest monthly average water temperature (MinT) and three metrics of salinity: annual average water salinity (AveS), the highest monthly average water salinity (MaxS), the lowest monthly average water salinity (MinS). We used a non-parametric Kruskal-Wallis test to assess differences in temperature and salinity between the Red Sea population and the others in SPSS v.18. In addition, a principal component analysis (PCA) with the ‘prcomp’ package in program R was performed to illustrate the high salinity and temperature in Red Sea population.

Population genetic patterns

Population genetic diversity was measured by allelic richness (AR), observed heterozygosity (HO), expected heterozygosity (HE) and inbreeding coefficient (FIS) using FSTAT v. [71]. The non-parametric Kruskal-Wallis test was used to test the difference in allelic richness (AR) between the Red Sea population and the others. Recent bottleneck effect was tested using program Bottleneck v.1.2.02 [72] based on the two-phase model (TPM) with 90% one-step mutations, and statistical significance was based on 1000 iterations with a one-tailed Wilcoxon test. In addition, we performed the mode shift test, and deviation from L-shaped distributions of allele frequency would suggest recent bottlenecks using Bottleneck program [73].

We assessed population genetic differentiation with pairwise FST values using 10,000 permutations in ARLEQUIN and the significance level was adjusted after sequential Bonferroni correction. To further investigate population genetic structure, we performed the Bayesian clustering in STRUCTURE v.2.3 [74]. For the STRUCTURE analysis, we assessed likelihoods for models with the number of clusters ranging from K = 1 to 6 (the total number of populations), with 100,000 Markov chain Monte Carlo iterations preceded by 100,000 burn-in, and we performed ten independent runs for each K value. In addition, we further conducted Bayesian analyses to assess population assignment between the Red Sea population and Mediterranean-Atlantic population group. The optimal K value was identified with STRUCTURE HARVESTER [75, 76]. The program DISTRUCT [77] was used to visualize the results.

Identification of candidate adaptive loci under selection

In this study, two different population differentiation approaches were used to search for microsatellite loci putatively under selection. Firstly, we screened all microsatellite loci with the fdist2 method under a hierarchal island model in ARLEQUIN [78]. This analysis was simulated based on 1000 demes with 50,000 simulations. Secondly, we identified outlier loci based on a Bayesian approach in BAYESCAN v.2.1 [79] according to default parameter settings. Candidate loci under selection were defined as those with false discovery rate lower than 5% (q-value < 0.05).

In order to identify loci that have strong association with particular environmental variables and further confirm that the outlier loci were the outcome of selection resulting from local environment factors, the environmental association analysis (EAA) was conducted in MatSAM v.2 [80]. MatSAM uses multiple univariate logistic regressions to test the associations between the allelic frequency and environmental variable at particular loci with the Bonferroni correction at a confidence level of 95% based on the cumulated test. Pearson correlations were performed in SPSS v.18 to test the correlation between the allele frequency and environmental factor at the adaptive loci.

Adaptive outlier loci annotation

To annotate gene functions of the putatively adaptive outlier loci, we firstly located these loci on the chromosomes based on the C. robusta KH assembly of Satou et al. [81] at Ensembl ( In addition, we chose 20 kb up- and down-stream of each adaptive outlier loci as a selective sweep window against the C. robusta KH assembly as suggested by Lin et al. [44], and then we obtained the sequences of the selective sweep windows at Ensembl. Furthermore, sequences of candidate microsatellite loci and selective sweep windows were subjected for BLASTN search at the NCBI website to obtain the best gene hits. Finally, UniProt database and AmiGO 2 GO browser [82, 83] were used to conduct gene annotation and then assign gene function based on gene ontology (GO) terms. We performed GO enrichment and KEGG pathways analysis using DAVID 6.7 [84].


  1. 1.

    Salamin N, Wüest RO, Lavergne S, et al. Assessing rapid evolution in a changing environment. Trends Ecol Evol. 2010;25(12):692–8.

  2. 2.

    Sanford E, Kelly MW. Local adaptation in marine invertebrates. Annu Rev Mar Sci. 2011;3:509–35.

  3. 3.

    Briski E, Chan FT, Darling JA, et al. Beyond propagule pressure: importance of selection during the transport stage of biological invasions. Front Ecol Environ. 2018.

  4. 4.

    Hoegh-Guldberg O, Bruno JF. The impact of climate change on the world’s marine ecosystems. Science. 2010;328(5985):1523–8.

  5. 5.

    Moran EV, Alexander JM. Evolutionary responses to global change: lessons from invasive species. Ecol Lett. 2014;17(5):637–49.

  6. 6.

    Lee CE. Evolutionary genetics of invasive species. Trends Ecol Evol. 2002;17(8):386–91.

  7. 7.

    Prentis PJ, Wilson JRU, Dormontt EE, et al. Adaptive evolution in invasive species. Trends Plant Sci. 2008;13(6):288–94.

  8. 8.

    Kingsolver JG, Buckley LB. Evolution of plasticity and adaptive responses to climate change along climate gradients. Proc R Soc B. 2017;284(1860):20170386.

  9. 9.

    Bergland AO, Behrman EL, O'Brien KR, et al. Genomic evidence of rapid and stable adaptive oscillations over seasonal time scales in drosophila. PLoS Genet. 2014;10(11):e1004775.

  10. 10.

    Colautti RI, Lau JA. Contemporary evolution during invasion: evidence for differentiation, natural selection, and local adaptation. Mol Ecol. 2015;24(9):1999–2017.

  11. 11.

    Bernatchez L. On the maintenance of genetic variation and adaptation to environmental change: considerations from population genomics in fishes. J Fish Biol. 2016;89(6):2519–56.

  12. 12.

    Hamilton PB, Rolshausen G, Webster TMU, et al. Adaptive capabilities and fitness consequences associated with pollution exposure in fish. Phil Trans R Soc B. 2017;372(1712):20160042.

  13. 13.

    Bock DG, Caseys C, Cousens RD, et al. What we still don't know about invasion genetics. Mol Ecol. 2015;24(9):2277–97.

  14. 14.

    Peischl S, Excoffier L. Expansion load: recessive mutations and the role of standing genetic variation. Mol Ecol. 2015;24(9):2084–94.

  15. 15.

    Gienapp P, Teplitsky C, Alho JS, et al. Climate change and evolution: disentangling environmental and genetic responses. Mol Ecol. 2008;17(1):167–78.

  16. 16.

    Bernardi G, Azzurro E, Golani D, et al. Genomic signatures of rapid adaptive evolution in the bluespotted cornetfish, a Mediterranean Lessepsian invader. Mol Ecol. 2016;25(14):3384–96.

  17. 17.

    Gillis MK, Walsh MR. Rapid evolution mitigates the ecological consequences of an invasive species (Bythotrephes longimanus) in lakes in Wisconsin. Proc R Soc B. 2017;284(1858):20170814.

  18. 18.

    Hoffmann AA, Willi Y. Detecting genetic responses to environmental change. Nat Rev Genet. 2008;9(6):421.

  19. 19.

    McEvoy PB, Higgs KM, Coombs EM, et al. Evolving while invading: rapid adaptive evolution in juvenile development time for a biological control organism colonizing a high-elevation environment. Evol Appl. 2012;5(5):524–36.

  20. 20.

    Guo B, DeFaveri J, Sotelo G, et al. Population genomic evidence for adaptive differentiation in Baltic Sea three-spined sticklebacks. BMC Biol. 2015;13(1):19.

  21. 21.

    Lescak EA, Bassham SL, Catchen J. Evolution of stickleback in 50 years on earthquake-uplifted islands. Proc Natl Acad Sci U S A. 2015;112(52):E7204–12.

  22. 22.

    Hendry AP, Kinnison MT. An introduction to microevolution: rate, pattern, process. Genetica. 2001;112(1):1–8.

  23. 23.

    Beaumont MA, Balding DJ. Identifying adaptive genetic divergence among populations from genome scans. Mol Ecol. 2004;13(4):969–80.

  24. 24.

    Bentham J, Morris DL, Graham DSC, et al. Genetic association analyses implicate aberrant regulation of innate and adaptive immunity genes in the pathogenesis of systemic lupus erythematosus. Nat Genet. 2015;47(12):1457.

  25. 25.

    Marques DA, Jones FC, Di PF, et al. Experimental evidence for rapid genomic adaptation to a new niche in an adaptive radiation. Nature Ecol Evol. 2018;2(7):1130.

  26. 26.

    Zhan A, Briski E, Bock DG, et al. Ascidians as models for studying invasion success. Mar Biol. 2015;162(12):0.

  27. 27.

    Huang X, Li S, Ni P, et al. Rapid response to changing environments during biological invasions: DNA methylation perspectives. Mol Ecol. 2017;26(23):6621–33.

  28. 28.

    Barshis DJ, Ladner JT, Oliver TA, et al. Genomic basis for coral resilience to climate change. Proc Natl Acad Sci U S A. 2013;110(4):1387–92.

  29. 29.

    Paiva F, Barco A, Chen Y, et al. Is salinity an obstacle for biological invasions? Glob Change Biol. 2018;24:2708–20.

  30. 30.

    Carver CE, Mallet AL, Vercaemer B. Biological synopsis of the solitary tunicate Ciona intestinalis. Dartmouth, Nova Scotia: Bedford institute of Oceanography;2006.

  31. 31.

    Bouchemousse S, Bishop JDD, Viard F. Contrasting global genetic patterns in two biologically similar, widespread and invasive Ciona species (Tunicata, Ascidiacea). Sci Rep. 2016;6:24875.

  32. 32.

    Zhan A, Darling JA, Bock DG, et al. Complex genetic patterns in closely related colonizing invasive species. Ecol Evol. 2012;2(7):1331–46.

  33. 33.

    Zhan A, Macisaac HJ, Cristescu ME. Invasion genetics of the Ciona intestinalis species complex: from regional endemism to global homogeneity. Mol Ecol. 2010;19:4678–94.

  34. 34.

    Monniot C, Monniot F. Additions to the inventory of eastern tropical Atlantic ascidians; arrival of cosmopolitan species. Bull Mar Sci. 1994;54(1):71–93.

  35. 35.

    Shenkar N, Shmuel Y, Huchon D. The invasive ascidian Ciona robusta recorded from a Red Sea marina. Mar Biodivers. 2017:1–4.

  36. 36.

    Millard N. Observations and experiments on fouling organisms in Table Bay harbour, South Africa. Trans Roy Soc S Afr. 1951;33(4):415–46.

  37. 37.

    Coll M, Piroddi C, Steenbeek J, et al. The biodiversity of the Mediterranean Sea: estimates, patterns, and threats. PloS One. 2010;5(8):e11842.

  38. 38.

    Belmaker J, Parravicini V, Kulbicki M. Ecological traits and environmental affinity explain Red Sea fish introduction into the Mediterranean. Glob Change Biol. 2013;19(5):1373–82.

  39. 39.

    Por FD. One hundred years of Suez Canal-a century of Lessepsian migration: retrospect and viewpoints. Syst Zool. 1971;20:138–59.

  40. 40.

    Golani D, Bogorodsky SV. The fishes of the Red Sea–reappraisal and updated checklist. Zootaxa. 2010;2463(1–135).

  41. 41.

    Malaquias MAE, Zamora-Silva A, Vitale D, et al. The Mediterranean Sea as a gateway for invasion of the Red Sea: the case of the indo-West Pacific head-shield slug Chelidonura fulvipunctata Baba, 1938. Aquat Invasions. 2016;11(3):247–55.

  42. 42.

    Boudouresque CF, Perret-Boudouresque M, Verlaque M. Donor and recipient regions for exotic species of marine macrophytes: a case of unidirectional flow, the Mediterranean Sea. Rapports de la Commission Internationale pour la Mer. Méditerranée. 2016;41:426.

  43. 43.

    Caputi L, Crocetta F, Toscano F, et al. Long-term demographic and reproductive trends in Ciona intestinalis sp. a. Mar Ecol. 2015;36(1):118–28.

  44. 44.

    Lin Y, Chen Y, Yi C, et al. Genetic signatures of natural selection in a model invasive ascidian. Sci Rep. 2017;7:44080.

  45. 45.

    Berna L, Alvarez-Valin F. Evolutionary genomics of fast evolving tunicates. Genome Biol Evol. 2014;6(7):1724–38.

  46. 46.

    Ometto L, Li M, Bresadola L, et al. Demographic history, population structure, and local adaptation in alpine populations of Cardamine impatiens and Cardamine resedifolia. PLoS One. 2015;10(5):e0125199.

  47. 47.

    Chen Y, Li S, Lin Y, et al. Population genetic patterns of the solitary tunicate, Molgula manhattensis, in invaded Chinese coasts: large-scale homogeneity but fine-scale heterogeneity. Mar Biodivers. 2017:1–13.

  48. 48.

    Pedersen SH, Ferchaud AL, Bertelsen MS, et al. Low genetic and phenotypic divergence in a contact zone between freshwater and marine sticklebacks: gene flow constrains adaptation. BMC Evol Biol. 2017;17(1):130.

  49. 49.

    Saenz-Agudelo P, Dibattista JD, Piatek MJ, et al. Seascape genetics along environmental gradients in the Arabian peninsula: insights from ddRAD sequencing of anemonefishes. Mol Ecol. 2015;24(24):6241–55.

  50. 50.

    Kolbe JJ, Glor RE, Schettino LR, et al. Genetic variation increases during biological invasion by a Cuban lizard. Nature. 2004;431(7005):177.

  51. 51.

    Kelly MW, Padilla-Gamiño JL, Hofmann GE. Natural variation and the capacity to adapt to ocean acidification in the keystone sea urchin Strongylocentrotus purpuratus. Glob Change Biol. 2013;19(8):2536–46.

  52. 52.

    Nadeau CP, Urban MC, Bridle JR. Climates past, present, and yet-to-come shape climate change vulnerabilities. Trends Ecol Evol. 2017;32:786–800.

  53. 53.

    Hendry AP, Wenburg JK, Bentzen P, et al. Rapid evolution of reproductive isolation in the wild: evidence from introduced salmon. Science. 2000;290(5491):516–8.

  54. 54.

    Barrett RDH, Schluter D. Adaptation from standing genetic variation. Trends Ecol Evol. 2008;23(1):38–44.

  55. 55.

    Hermisson J, Pennings PS. Soft sweeps: molecular population genetics of adaptation from standing genetic variation. Genetics. 2005;169(4):2335–52.

  56. 56.

    Heilbron K, Toll-Riera M, Kojadinovic M, et al. Fitness is strongly influenced by rare mutations of large effect in a microbial mutation accumulation experiment. Genetics. 2014;197:981–90.

  57. 57.

    Tsagkogeorga G, Cahais V, Galtier N. The population genomics of a fast evolver: high levels of diversity, functional constraint, and molecular adaptation in the tunicate Ciona intestinalis. Genome Biol Evol. 2012;4(8):852–61.

  58. 58.

    Mark FC, Lucassen M, Pörtner HO. Thermal sensitivity of uncoupling protein expression in polar and temperate fish. Comp Biochem Physiol Part D: Genomics Proteomics. 2006;1(3):365–74.

  59. 59.

    Bermejo-Nogales A, Calduch-Giner JA, Pérez-Sánchez J. Gene expression survey of mitochondrial uncoupling proteins (UCP1/UCP3) in gilthead sea bream (Sparus aurata L.). J Comp Physiol B. 2010;180(5):685–94.

  60. 60.

    Yokota N, Harada Y, Sawada H. Identification of testis-specific ubiquitin-conjugating enzyme in the ascidian Ciona intestinalis. Mol Reprod Dev. 2010;77(7):640–7.

  61. 61.

    Lopez CE, Sheehan HC, Vierra DA, et al. Proteomic responses to elevated ocean temperature in ovaries of the ascidian Ciona intestinalis. Biol Open. 2017;bio.024786.

  62. 62.

    Rodriguez CI, Stewart CL. Disruption of the ubiquitin ligase HERC4 causes defects in spermatozoon maturation and impaired fertility. Dev Biol. 2007;312:501–8.

  63. 63.

    Asaoka R, Uemura T, Ito J, et al. Arabidopsis RABA1 GTPases are involved in transport between the trans-Golgi network and the plasma membrane, and are required for salinity stress tolerance. Plant J. 2012;73:240–9.

  64. 64.

    Asaoka R, Uemura T, Nishida S, et al. New insights into the role of Arabidopsis RABA1 GTPases in salinity stress tolerance. Plant Signal Behav. 2013;8(9):e25377.

  65. 65.

    Clevers H, Nusse R. Wnt/β-catenin signaling and disease. Cell. 2012;149(6):1192–205.

  66. 66.

    Pellettieri J, Alvarado AS. Cell turnover and adult tissue homeostasis: from humans to planarians. Annu Rev Genet. 2007;41:83–105.

  67. 67.

    Maor-Landaw K, Karako-Lampert S, Ben-Asher HW, et al. Gene expression profiles during short-term heat stress in the red sea coral Stylophora pistillata. Glob Change Biol. 2014;20(10):3026–35.

  68. 68.

    Xu C, Li E, Liu Y, et al. Comparative proteome analysis of the hepatopancreas from the Pacific white shrimp Litopenaeus vannamei under long-term low salinity stress. J Proteome. 2017;162:1–10.

  69. 69.

    Iannelli F, Pesole G, Sordino P, et al. Mitogenomics reveals two cryptic species in Ciona intestinalis. Trends Genet. 2007;23:419–22.

  70. 70.

    Lin Y, Chen Y, Xiong W, et al. Genomewide gene-associated microsatellite markers for the model invasive ascidian, Ciona intestinalis species complex. Mol Ecol Resour. 2016;16(3):784–93.

  71. 71.

    Goudet J. FSTAT, a program to estimate and test gene diversities and fixation indices (version 29.3).2001.

  72. 72.

    Piry S, Luikart G, Cornuet JM. BOTTLENECK: A computer program for detecting recent reductions in the effective population size using allele frequency data. J Hered. 1999;90:502–3.

  73. 73.

    Luikart G, Allendorf FW, Cornuet JM, et al. Distortion of allele frequency distributions provides a test for recent population bottlenecks. J Hered. 1998;89:238–47.

  74. 74.

    Pritchard JK, Stephens M, Donnelly P. Inference of population structure using multilocus genotype data. Genetics. 2000;155:945–59.

  75. 75.

    Evanno G, Regnaut S, Goudet J. Detecting the number of clusters of individuals using the software STRUCTURE: a simulation study. Mol Ecol. 2005;14:2611–20.

  76. 76.

    Earl DA, vonHoldt BM. STRUCTURE HARVESTER: a website and program for visualizing STRUCTURE output and implementing the Evanno method. Conserv Genet Resour. 2012;4(2):359–61.

  77. 77.

    Rosenberg NADISTRUCT. A program for the graphical display of population structure. Mol Ecol Resour. 2004;4:137–8.

  78. 78.

    Excoffier L, Lischer HEL. ARLEQUIN suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and windows. Mol Ecol Resour. 2010;10:564–7.

  79. 79.

    Foll M, Gaggiotti O. A genome-scan method to identify selected loci appropriate for both dominant and codominant markers: a Bayesian perspective. Genetics. 2008;180:977–93.

  80. 80.

    Joost S, Kalbermatten M, Bonin A. Spatial analysis method (SAM): a software tool combining molecular and environmental data to identify candidate loci for selection. Mol Ecol Resour. 2008;8:957–60.

  81. 81.

    Satou Y, Mineta K, Ogasawara M, et al. Improved genome assembly and evidence-based global gene model set for the chordate Ciona intestinalis: new insight into intron and operon populations. Genome Biol. 2008;9:R152.

  82. 82.

    Carbon S, Ireland A, Mungall CJ, et al. AmiGO: online access to ontology and annotation data. Bioinformatics. 2008;25(2):288–9.

  83. 83.

    Mi HY, Muruganujan A, Thomas PD. PANTHER in 2013: modeling the evolution of gene function, and other gene attributes, in the context of phylogenetic trees. Nucleic Acids Res. 2013;41:D377–86.

  84. 84.

    Huang DW, Sherman BT, Lempicki RA. Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nature Protoc. 2009;4(1):44–57.

Download references


Great thanks to Lion Novak and Gal Navon for their assistance with the sample collection in Eilat, the staff of Eilat Marina for their support, and the Israel’s National Monitoring Program of the Gulf of Eilat for providing us temperature and salinity data.

Availability of data and material

Not applicable.

Ethical approval and consent to participate

Not applicable.


This work was supported by the National Natural Science Foundation of China (Nos. 31622011, 31772449) to AZ, the Youth Innovation Promotion Association, Chinese Academy of Sciences (No. 2018054) to SL, the Israel Science Foundation (No. 993/15) to NS.

Author information

AZ and YC conceived this study and wrote the manuscript. NS collected samples from the Red Sea. NS, YC, PN, YL and SL contributed to the lab work and data analysis. All co-authors revised and approved the final manuscript.

Correspondence to Aibin Zhan.

Ethics declarations

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Additional files

Additional file 1:

Figure S1. (a) Temperature and (b) salinity data in six populations. Figure S2. PCA plot on environmental factors for six Ciona robusta populations. Figure S3. Outlier detection in 146 microsatellite loci in ARLEQUIN. Purple line represents 99% confidence intervals; red and green lines represent 95 and 5% confidence intervals, respectively. Figure S4. Manhattan plot showing the distribution of FST-based outliers detected by BAYSECAN across different chromosomes of Ciona robusta KH assembly. The q-value of given locus is the minimum false discovery rate (FDR) at which this locus may become significant. Table S1. Protein MB21D2 (MB21D2 gene). (ZIP 1091 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (, which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Chen, Y., Shenkar, N., Ni, P. et al. Rapid microevolution during recent range expansion to harsh environments. BMC Evol Biol 18, 187 (2018).

Download citation


  • Biological invasion
  • Rapid microevolution
  • Range expansion
  • Invasive species
  • Red Sea
  • Adaptive genes
  • Ciona robusta