Skip to main content

Reproductive protein evolution in two cryptic species of marine chordate



Reproductive character displacement (RCD) is a common and taxonomically widespread pattern. In marine broadcast spawning organisms, behavioral and mechanical isolation are absent and prezygotic barriers between species often operate only during the fertilization process. Such barriers are usually a consequence of differences in the way in which sperm and egg proteins interact, so RCD can be manifest as faster evolution of these proteins between species in sympatry than allopatry. Rapid evolution of these proteins often appears to be a consequence of positive (directional) selection. Here, we identify a set of candidate gamete recognition proteins (GRPs) in the ascidian Ciona intestinalis and showed that these GRPs evolve more rapidly than control proteins (those not involved in gamete recognition). Choosing a subset of these gamete recognition proteins that show evidence of positive selection (CIPRO37.40.1, CIPRO60.5.1, CIPRO100.7.1), we then directly test the RCD hypothesis by comparing divergence (omega) and polymorphism (McDonald-Kreitman, Tajima's D, Fu and Li's D and F, Fay and Wu's H) statistics in sympatric and allopatric populations of two distinct forms of C. intestinalis (Types A and B) between which there are strong post-zygotic barriers.


Candidate gamete recognition proteins from two lineages of C. intestinalis (Type A and B) are evolving more rapidly than control proteins, consistent with patterns seen in insects and mammals. However, ω (dN/dS) is not significantly different between the sympatric and allopatric populations, and none of the polymorphism statistics show significant differences between sympatric and allopatric populations.


Enhanced prezygotic isolation in sympatry has become a well-known feature of gamete recognition proteins in marine broadcast spawners. But in most cases the evolutionary process or processes responsible for this pattern have not been identified. Although gamete recognition proteins in C. intestinalis do appear to evolve more rapidly, on average, than proteins with other functions, rates of evolution are not different in allopatric and sympatric populations of the two reproductively isolated forms. That sympatry is probably human-mediated, and therefore recent, may explain the absence of RCD.


Reproductive isolation between incipient species is of particular relevance to the process of speciation. Reproductive character displacement - 'the pattern of greater divergence of a (prezygotic) isolating trait in areas of sympatry between closely related taxa than in areas of allopatry' [1, 2] is a common and taxonomically widespread pattern which is of great interest when studying reproductive isolation [3]. Evidence for RCD comes from groups as diverse as fungi [4, 5], plants [6], insects [7, 8] mollusks [9], fish [10, 11] and amphibians [12]. However, the majority of RCD examples come from the Drosophila literature [1317].

The study of RCD has historically been tied to the process of reinforcement, the evolution of prezygotic isolation resulting from selection against hybrid individuals [1315]. More recently, however, workers have emphasized that RCD can be caused by other factors, including ecological variables [18, 19]. But even where selection has been shown to play a role in RCD [20, 21], the specific action of this selection remains unknown.

RCD has been documented in several taxa of marine broadcast spawning organisms [18]. The study of RCD is often more tractable in broadcast spawners than in internal fertilizers. Broadcast spawning individuals simultaneously release sperm and eggs into the water column, where fertilization occurs. Therefore, behavioral and mechanical isolation are absent in broadcast spawning organisms and if two sympatric broadcast spawning species are not temporally isolated [22], prezygotic barriers between these two species occur only during the fertilization process. Such barriers are usually a consequence of differences in the way in which sperm and egg proteins interact, so RCD is often manifest as faster evolution of these proteins between species in sympatry than allopatry [23, 24]. Rapid evolution of these proteins appears to proceed through the action of positive (directional) selection, increasing the proportion of nonsynonymous substitutions (those that cause amino acid changes) relative to synonymous substitutions (the dn/ds ratio) [25].

New information has allowed us to use the broadcast spawning ascidian Ciona intestinalis, a longtime model in developmental biology, to study the early stages of RCD. C. intestinalis comprises two distinct and divergent lineages, now termed Type A and Type B [2628], Type A is thought to be native to the Northwestern Pacific Ocean and to have invaded the Eastern Pacific Ocean, the Mediterranean Sea, the Atlantic coast of South Africa, and the Black Sea [29, 30]. Type B, native to the Northern Atlantic Ocean [31, 32], has invaded the Western Atlantic Ocean [27]. The ranges of Type A and B overlap in the English Channel and the Atlantic coast of France, where limited introgression occurs between Type A and B individuals from these sympatric populations [26]. Hybrids from crosses between allopatric individuals are sterile or inviable in the laboratory [28], but the evidence for introgression in nature shows that these two lineages are not completely reproductively isolated. Nothing is yet known about pre or postzygotic barriers in sympatry.

Here, we first identify a set of candidate gamete recognition proteins (GRPs) from C. intestinalis using proteomic and bioinformatic techniques and test whether these proteins evolve more rapidly than control proteins. Then, choosing a subset of these proteins that show evidence of positive selection, we directly test the RCD hypothesis by comparing divergence and polymorphism statistics in sympatric and allopatric populations of Type A and B C. intestinalis. We ask whether signatures of positive selection are stronger in sympatric than allopatric populations.


Identification of candidate GRPs from sperm: proteomics

161 sperm proteins were identified by LC-MS/MS (Liquid Chromatography-Mass Spectrometry/Mass Spectrometry). Each of these proteins was subsequently identified using one or more of the following sections of the CIPRO database ( C iona intestinalis Protein): the descriptive summary available for many proteins, the Pfam Domain search, and the BlastP search, or the GO (Gene Ontology) program. 144 of the proteins are unlikely to be GRPs; they are likely involved in the movement or metabolism of the sperm (e.g. actin, dynein, myosin, tektin, α-tubulin, ATP-synthase, creatine kinase, enolase, malate dehydrogenase). The identities of these 144 proteins are available from the authors. Of the remaining 17 proteins, seven were likely GRPs, and 10 could not be identified as similar to any known proteins. Of these 17 proteins, four could not be analyzed because primers based on the Type A genome sequence would not amplify the homolog from Type B individuals. An additional gene could not be analyzed because the gene encoding the protein did not have significant tblastn hits to the Type A genome (and therefore no primers could be designed). We selected the remaining 12 proteins for further analysis, 3 of which were likely GRPs.

Identification of candidate GRPs proteins from sperm: bioinformatics

We also used a bioinformatic approach to identify potential GRPs. First, we accessed the functional classifications of Type A testis ESTs sequenced by [33] and selected 25 ESTs that might code for GRPs. We then located the genes corresponding to these ESTs with the KOG (EuKaryotic Orthologous Groups) tool provided on JGI's (Joint Genome Institute) C. intestinalis genome browser These genes were then searched against the CIPRO database using blastx to identify resulting protein matches. 19 of these proteins were determined to be GRP candidates, but 10 failed to amplify and/or sequence in Type B individuals; the remaining nine were selected for further study.

Second, we located every protein in the CIPRO database that was identified as being expressed only in testes. We chose a subset of 10 of these proteins that were likely GRPs, of which nine could be amplified from Type B cDNA. In total, we selected 12 candidate GRPs identified proteomically and 18 identified bioinformatically. We also selected 9 control proteins (not involved in the fertilization process) from the proteomic analyses to compare with the putative GRPs. These control proteins, CIPRO129.24.1, CIPRO183.42.1, CIPRO2.134.1, CIPRO32.59.1, CIPRO53.35.1, CIPRO57.34.1, CIPRO68.37.1, CIPRO81.38.1 and CIPRO963.2.1, were selected because Pfam Domain and BlastP Searches in CIPRO identified them as proteins involved in basic cellular processes rather than fertilization (e.g. tubulin, NADH dehydrogenase, aconitate hydratase).

Comparison of dN/dSvalues between candidate GRPs and control proteins

The dN/dS values for candidate GRPs are significantly higher than the dN/dS values for control proteins (Figure 1, p = 0.004891 using a one-tailed Mann-Whitney U Test).

Figure 1

Pairwise omega (d N /d S ) values. This graph shows pairwise dN/dS values (Type A vs. Type B C. intestinalis) for 30 candidate gamete recognition proteins (GRPs) and 9 control proteins. The line is dN/dS = 0.5.

However, PAML analyses assume that dS values are constant across a sequence. If some sites across the sequence have unusually low dS values, a high dN/dS value could be inferred in the absence of positive selection [34]. Thus, significantly higher dN/dS values for candidate gamete recognition genes than for control genes could be the result of either higher dN or lower dS values in the candidate recognition genes, and only higher dN values provide evidence of positive selection.

To address this issue, we performed two-tailed Mann-Whitney U test in R (the Shapiro-Wilk test found evidence for non-normality), comparing dN values in candidate gamete recognition vs. control genes, and dS values in candidate gamete recognition vs. control genes. dN values were significantly different between candidate gamete recognition and control genes (p = 0.002), whereas dS values were not (p = 0.269). These tests are consistent with the assumption that dN, rather than dS, is driving the observed pattern.

Figure 1 also shows that two proteins have a dN/dS ratio greater than 0.5, which is the value above which we consider positive selection likely to be occurring when conservative pairwise dN/dS comparisons are used. A study by [35] across many different taxa showed that if a pairwise comparison yielded a dN/dS ratio of greater than 0.5, evidence for selection was subsequently observed when more sensitive site-specific tests were used. The rationale here is that dN/dS < 1 can still indicate positive selection at some subset of the amino acid residues within a protein, although much of the protein may be subject to constraint. Of the two proteins with dN/dS > 0.5, CIPRO37.40.1 was identified from the CIPRO database as having testis-only expression, and contains domains similar to ricin-type beta-trefoil lectin domains (dN/dS = 0.618, dN = 0.054, dS = 0.087).

CIPRO100.7.1 was identified proteomically, and its function is unknown. CIPRO100.7.1 is a large protein (1,225 amino acids), and was therefore divided into three sections for sequencing. Sections 2, 3, and all 3 sections analyzed together had dN/dS values less than 0.5, so only Section 1 was analyzed in subsequent sympatric vs. allopatric comparisons. Section 1 of CIPRO100.7.1 had a dN/dS equal to 0.531, a dN equal to 0.050, and a dS equal to 0.095. One protein, CIPRO60.5.1, has a dN/dS ratio lower than 0.5, but was identified in the CIPRO database as a metalloproteinase and has a GO biological function of "sperm binding to zona pellucida" (dN/dS = 0.366, dN = 0.052, dS = 0.142). We chose these three proteins for the sympatric vs. allopatric comparisons, either because they showed evidence of positive selection (in the case of CIPRO37.40.1 and CIPRO100.7.1), or because their putative function was so clearly related to gamete recognition (in the case of CIPRO60.5.1).

Sympatric vs. allopatric divergence analyses

The results of the divergence analyses are shown in Table 1. For all of the candidate GRPs and control proteins, the 95% credible interval of the distribution of the differences between sympatric and allopatric ω (dN/dS) values included zero. This indicates that omega values were not significantly different between sympatry and allopatry for any of the proteins examined. Not enough variation was present in the sympatric Type A alleles for CIPRO60.5.1 or mtCOI to obtain reliable omega values, so the sympatric Type A vs. allopatric Type A comparison could not be performed for these proteins.

Table 1 Sympatric vs. allopatric comparisons of omega values

Sympatric vs. allopatric polymorphism analyses

The summary statistics are shown in Table 2. No consistent differences between sympatric and allopatric Type A or between sympatric and allopatric Type B can be discerned for any of these statistics for any of these candidate genes. Table 3 presents the results of the McDonald-Kreitman tests. Fixed nonsynonymous substitutions are not more common in sympatric than allopatric comparisons; these tests provide no evidence for positive selection on the genes encoding these three candidate GRPs in sympatry. None of the statistics for which significance was determined by coalescent simulation (D, D*, F*, H) showed significant differences between sympatry and allopatry for any of the genes examined, as the confidence intervals for sympatric populations always contained the mean values for allopatric populations, and vice versa (Table 4 for Tajima's D, Table 5 for D* and F*, Table 6 for H).

Table 2 Summary Statistics.
Table 3 Results of the McDonald Kreitman Tests for all genes.
Table 4 Tajima's D statistics for all genes
Table 5 Fu and Li's D* and F* statistics for all genes
Table 6 Fay and Wu's H test for all genes.


Comparison of dN/dSvalues between candidate GRPs and control proteins

Candidate GRPs in C. intestinalis are evolving more rapidly than control proteins, and this pattern is likely driven by substitutions at nonsynonymous sites. Rapid evolution has been documented at specific GRPs (i.e. lysin, VERL and bindin) in marine broadcast spawners, and dN/dS values are lower for mtCOI than the GRPs lysin and VERL for green and pink abalone [36]. The pattern we see in C. intestinalis has also been documented in insects and mammals (e.g. butterflies: [37], field crickets: [38], mouse and human: [39], primates: [40]). This study suggests that a pattern of faster evolution in reproductive proteins may apply to external as well as internal fertilizers.

One explanation for faster evolution in candidate GRPs than control proteins is sperm competition, which occurs in Ciona as it does in internal fertilizers. Selection could be acting on any proteins that determine how quickly sperm fertilize eggs: proteins involved in metabolism, motility, binding, penetration, etc. However, as Figure 1 shows that candidate GRPs are evolving more rapidly than sperm proteins that are not candidate GRPs (control proteins), proteins directly involved in sperm-egg interactions are more likely to be experiencing directional selection than those involved in facilitating sperm access to the egg.

Rapid evolution of sperm GRPs might also result from sexual conflict, which occurs when the optimal outcomes of fertilization are different for sperm and eggs. For sperm, the optimal outcome is fertilization of an egg as quickly as possible. But in many taxa, fertilization of eggs by multiple sperm (polyspermy) results in developmental defects. Therefore, the optimal outcome of fertilization for an egg may be slower fertilization, to avoid polyspermy [25]. Ascidians like C. intestinalis often live in close proximity to many conspecific individuals [41], and an individual usually sends sperm into the water column before eggs (ascidians are hermaphrodites). So eggs are usually released into a vast amount of sperm spawned from many neighbors, making the risk of polyspermy very high. Perhaps in response to this risk, ascidians have evolved two separate blocks to polyspermy, whereas many other marine broadcast spawners have a single block [41]. Given these effective polyspermy blocks, sexual conflict resulting from polyspermy is not likely to be major driver of GRP evolution in ascidians.

Reinforcement could be driving the enhanced prezygotic isolation in sympatry. We know that hybrids between allopatric populations Type A and B C. intestinalis are sterile or inviable in the laboratory, so selection could favor GRPs that allow sperm to preferentially bind to conspecific eggs. However, reinforcing selection could be counteracted by gene flow from allopatric populations [42], especially when secondary contact is recent. Pairwise FST calculations between sympatric and allopatric Type A, and sympatric and allopatric Type B populations show that these populations are not significantly differentiated at any of the three genes encoding the GRPs (data not shown); migration may therefore be occurring between sympatric and allopatric populations.

Lastly, egg surface proteins could be changing rapidly to prevent pathogens from entering the egg. If the same proteins involved in preventing microbial attack are involved in sperm/egg recognition, this could lead to the rapid evolution of sperm proteins to keep up with the ever-changing egg proteins.

Evolution of candidate GRPs in Ciona intestinalis- no evidence for RCD

Our data provide no evidence that positive selection is enhanced in sympatry, and if these candidate GRPs are involved in prezygotic isolation, we have no evidence for enhanced prezygotic isolation. The polymorphism statistics likewise give no indication that RCD is occurring in these three proteins.

We cannot necessarily conclude from lack of evidence for RCD in CIPRO37.40.1, CIPRO60.5.1 and CIPRO100.7.1 that RCD is not occurring in this system. If enhanced prezygotic isolation between Type A and B does exist, there are several reasons why we might not have detected it in this study. First, primers for candidate GRPs were developed from the Type A genomic sequence and were used to amplify and sequence both Type A and B individuals (the Type B genome has not been sequenced). But Type A and B are substantially divergent (p-distances: 0.124 at mtCOI, 0.035 to 0.116 for six nuclear loci; [43]), which could explain why 15 genes encoding GRP candidates could not be successfully amplified and/or sequenced in Type B individuals. It is possible that the genes that could not be amplified and/or sequenced (and were therefore excluded from the analyses) encode proteins that are evolving more rapidly between Type A and B than those that were included in the analyses. If this is the case, we may have missed proteins for which dN/dS values are >1, proteins that would have been included in the sympatric vs. allopatric tests of RCD.

It is, of course, also possible that the three proteins that have high dN/dS values are themselves not involved in gamete recognition. While some aspects of the fertilization process in solitary ascidians such as C. intestinalis are well-characterized [44, 45], the genes and corresponding proteins responsible for species-specificity have not been identified, as they have in other marine broadcast spawners [23, 24].

Two sperm proteins that interact directly with the egg have been identified in C. intestinalis, but it is not known whether these proteins are involved in gamete recognition. The first protein is α L-fucosidase, which binds to the vitelline coat of the egg [46]. Five of our candidate GRPs (CIPRO187.4.1, CIPRO19.75.1, CIPRO33.15.1, CIPRO552.7.1, and CIPRO58.12.1) had domains also found in α L-fucosidases, but none of these proteins were expressed in testis tissue, based on expression data in CIPRO. CIPRO187.4.1 could not be amplified from Type B, and the other four genes did not have dN/dS ratios > 0.5. A second protein involved in sperm-egg interactions is a chymotrypsin-like enzyme that may dissolve the vitelline coat of the egg [47]. However, the amino acid sequence for this protein is not available and we identified dozens of chymotrypsin-like proteins in the genome.

Temporal isolation of Types A and B could also contribute to prezygotic barriers. We do not know whether Type A and B release gametes at the same time of day or in the same season in the English Channel. Since some gene flow has occurred [26], spawning must at least partially overlap.

Gene flow between sympatric and allopatric populations could be obscuring RCD if it is occurring in this system. Pairwise FST calculations show that sympatric and allopatric populations of each type are not significantly differentiated, so gene flow between sympatric and allopatric populations is a possibility.

Finally, because secondary contact between Type A and B C. intestinalis in the English Channel may well be recent, RCD may not have yet taken place. Type B is native to Northern Europe and presumably a long-time resident of the English Channel. We do not know when Type A invaded the English Channel. The first published record of Type A in this area was in 2007 [26], but as Type A and B were only recognized in 2005 [48], Type A living in this area prior to 2005 would not have been distinguished from the native Type B. However, the introduction of Type A was likely human-mediated [26], and therefore a recent invasion on an evolutionary timescale. Evidence for RCD has been found in several systems where secondary contact is relatively recent. For example, RCD in mate choice is observed when limnetic and benthic sticklebacks co-occur in glacial lakes in British Columbia [20]. Similarly, RCD has been documented in the Mus musculus and Mus domesticus hybrid zone, which is thought to represent secondary contact during the Neolithic (since 9500 BCE) [49].

Positive selection on GRPs in other marine broadcast spawners: RCD

Some of the most rapidly evolving proteins yet discovered are GRPs in marine broadcast spawners (e.g. bindin in sea urchins, lysin in abalone and mussels). In sea urchins, the bindin protein facilitates sperm attachment to the egg and fusion of sperm and egg [50]. In three genera of sea urchin that include sympatric species (Echinometra, Heliocidaris, and Strongylocentrotus), regions of bindin show evidence of positive selection ([51] and references therein). In Arbacia, Lytechinus and Tripneustes, genera that do not contain sympatric species, bindin shows no evidence of positive selection ([51] and references therein). This pattern is consistent with RCD on bindin in sea urchins.

Stronger evidence for RCD comes from a study of Echinometra oblonga, which has populations that are sympatric and allopatric with Echinometra species C [23]. Substantial divergence in bindin alleles between E. oblonga and E. sp. C. occurs where the two species are sympatric, but not where they are allopatric [23].

In abalone and mussels, sperm proteins known as lysins are involved in dissolution of the egg vitelline envelope, enabling the sperm to enter the egg. The best-characterized lysins are in the abalone genus Haliotis. An early study of 19 sympatric Haliotis species (19 sympatric and 1 allopatric species) found many pairwise comparisons with dN/dS values > 1 [52]. A later study of 25 species corroborated the pairwise results of [52] and also used maximum likelihood models of codon substitution to identify lineage and site-specific evidence of positive selection [53]. Lineages containing sympatric or closely related species usually had dN/dS values > 1, whereas lineages with distantly related allopatric species always had dN/dS values < 1, a pattern consistent with RCD [53]. But the authors also note a dN/dS value > 1 for the two branches separating a group of Japanese species from two groups of Californian species; this speciation event was likely allopatric.

In the mussel Mytilus galloprovincialis, two divergent clades of Lysin-M7 have been found: G and GD [24]. Evidence of positive selection is seen between G and GD, and within GD [24]. The divergence between the two clades is the result of rapid evolution in the GD clade, and GD alleles are found at higher frequencies in sympatric populations of M. galloprovincialis (where it hybridizes with Mytilus edulis) than in allopatric populations [24].


Enhanced prezygotic isolation in sympatry has become a well-known feature of GRPs in marine broadcast spawners. But in most cases the evolutionary process or processes responsible for this pattern have not been identified. Differentiating between the processes that can lead to patterns of RCD will provide important insights into the process of speciation in marine broadcast spawners.


Identification of candidate GRPs from sperm: proteomics experiment

Sperm was collected from Type A individuals living in Santa Barbara, CA, filtered through a 70 μm nylon cell strainer (BD Biosciences) by centrifuging for 3 minutes at 3,000 rpm, and stored dry at -80°C. Sperm samples from several different individuals were later pooled and diluted 5-fold in phosphate buffer; the concentration of this dilution was determined to be 915 μg/ml. 500 μl of this diluted sperm was shipped to the University of Victoria Genome BC Proteomics Centre for the experiments described below.

9.5 M urea, 50 mM NH4HCO3 and 0.2% SDS were added to the sample, which was then sonicated. The proteins then underwent disulphide reduction and sulphydryl alkylation (200 mM DTT and 200 mM iodoacetamide) and were digested overnight at 37°C with 20 mg trypsin (Promega). Samples were subsequently cleaned with a cation exchange Cartridge Kit (Applied Biosystems).

Strong cation exchange chromatography

10 mM KH2PO4 (pH 2.7), 25% ACN buffer was added to the sample, which was then injected onto a Polysulphoethyl A strong cation exchange chromatography (SCX) column (PolyLC, Columbia, MD). The flow rate was set to 0.5 ml min-1. After equilibration, a 0-35% gradient of 10 mM KH2PO4, 25% ACN, 0.5 M KCl was added for 30 min. Each SCX fraction was then moved to autosampler vials (Dionex/LC Packings, Amsterdam).

One-dimensional reversed-phase chromatography with online mass spectrometry

A hybrid Quadruple-TOF LC-MS/MS mass spectrometer (QStar Pulsar I, MDS Sciex) with a nanoelectrospray ionization source (Proxeon, Odense, Denmark) was used to complete the liquid chromatography-mass spectrometry/mass spectrometry (LC-MS/MS) analyses. A C18AQ Nano LC and a Zorbax C18 guard column (Agilent Technologies) performed the chromatographic separation. The ANALYST QS v. 1.1 software service pack (ABI MDS SCIEX, Concord, Canada) gathered the data.

Mass spectrometry data analyses

the information dependent acquisition file was viewed using ANALYST v. 1.1 software, and the peak lists were assembled with a built-in MASCOT script (1.6b16 ABI--Matrix Science Limited). Spectra with less than 10 peaks were discarded. MASCOT v. 2.0 (Matrix Science Limited) was used to analyze the data. Spectrometry data were searched against a database of amino acid sequences from the CIPRO ( C iona intestinalis Protein) database

Comparison of dN/dSvalues between candidate GRPs and control proteins

Of proteins identified by proteomic analysis, we selected 39 proteins (30 candidate GRPs and 9 control proteins) for further analysis using proteomic and bioinformatic approaches (see Results section). Genomic sequences for the genes encoding the 39 proteins were located by performing tblastn searches to the Ciona intestinalis Ensembl genome server Primers were developed in coding regions, and the partial or entire coding regions of all 39 proteins were sequenced from cDNA of two Type A and two Type B individuals, all from allopatric populations. Testis tissues from these four individuals were collected in 2008 and immediately placed in RNAlater (Ambion). The tissue/RNAlater was frozen at -80°C within seven days of collection.

Total RNA was extracted from testis tissue with the RNAdvance Kit (Agencourt) and was used to synthesize single-stranded cDNA using SuperScript III reverse transcriptase (Invitrogen) and an oligo (dT) primer. A 5-fold dilution of the single-stranded cDNA was then PCR-amplified with TRsa and TS-PCR primers. The resulting PCR product was diluted 50-fold and used as the template for amplification of the coding regions for the 30 candidate GRPs proteins and 9 control proteins. The amplified coding regions were incubated with 0.25 μl each of Exonuclease I and Shrimp Antarctic Phosphatase at 37°C for 30 min, followed by 90°C for 10 min. The products were purified using CleanSeq beads (Agencourt), and the purified products were sequenced with a Big Dye Terminator Cycle sequencing kit and an Automated 3730 DNA Analyzer (Applied Biosystems). Sequences were edited, trimmed and aligned with Aligner (CodonCode Corporation, Dedham, MA). Primers and cycling conditions are available from the authors.

Once the sequences were obtained, the codeml program in PAML 4.4 [54] was used to obtain pairwise dN/dS values for each Type A vs. Type B combination (Type A Individual #1 vs. Type B Individual #1, Type A Individual #1 vs. Type B Individual #2, Type A Individual #2 vs. Type B Individual #1, Type A Individual #2 vs. Type B Individual #2). The average dN/dS value for all four combinations was calculated for each gene, and the dN/dS values of the putative GRPs and control proteins were found to be distributed non-normally using the Shapiro-Wilk test in R (version 2.10.1). The Shapiro-Wilk test is the most robust test of non-normality for small to medium sample sizes [55, 56]. The dN/dS values of the candidate GRPs and control proteins were therefore compared using a one-tailed Mann-Whitney U test in R (version 2.10.1).

Sympatric vs. allopatric divergence analyses

Samples were collected in 2005-2009 from allopatric and sympatric populations of Type A and Type B. The allopatric population of Type A was from Half Moon Bay, CA while sympatric populations of Type A were from Perros-Guirec, France and Concarneau, France. Allopatric populations of Type B were collected from Woods Hole, MA and Gosport, England.

To obtain genomic DNA, ovaries were dissected from freshly-collected individuals, cut into several pieces, immediately preserved in DMSO (dimethyl sulfoxide), and ultimately (within 10 d) stored at -80°C until needed. Total DNA was extracted from the ovaries using the Qiagen DNeasy® Tissue Kit (Qiagen Corporation, Santa Clarita, CA).

At least 10 individuals from each of 4 populations (allopatric Type A and B, sympatric Type A and B) were sequenced for the genes encoding three candidate GRPs: CIPRO37.40.1, CIPRO60.5.1, CIPRO100.7.1. The criteria for selecting these three candidate GRPs for the sympatric/allopatric comparison, from the 30 candidate GRPs used in the comparison of dN/dS values between GRPs and control proteins, are discussed in the Results section.

The same individuals were also sequenced for two control proteins: CIPRO53.35.1 and mitochondrial cytochrome oxidase I (mtCOI). Sequences from all five loci were deposited in GenBank (Accession Numbers HQ872081-HQ872454). A signature of enhanced selection in sympatry vs. allopatry could be due to selective processes or demographic processes (e.g. recent population growth). But demographic processes would affect all genes, not just candidate GRPs. So proteins not involved in the fertilization process (control proteins) were also subjected to divergence analyses. CIPRO53.35.1 was selected as a control protein because it contains a domain similar to a ribosomal L32 protein domain and is expressed in many different tissues; whereas mtCOI is an enzyme in the electron transport chain of the mitochondria. For all nuclear genes the PCR product for each locus was cloned using the pGEM®-T kit and up to eight clones were PCR-amplified and sequenced as described above to obtain both alleles. When only one allele was found in eight clones, the individual was assumed to be homozygous for that allele.

We used omegaMap version 0.5 [57] to determine whether omega values were significantly different for sympatric Type A vs. allopatric Type A populations and for allopatric Type A vs. allopatric Type B populations. We chose omegaMap over PAML for three reasons. First, omegaMap is designed to calculate omega values in the presence of recombination, and every one of the five genes (those encoding CIPRO37.40.1, CIPRO60.5.1, CIPRO100.7.1, CIPRO53.35.1, and mtCOI) has at least one population where intragenic recombination is present. Second, unlike PAML, omegaMap takes a population genetics approach, so we can use all the alleles we sequenced from each population in our calculation of omegaMap. Lastly, the Bayesian inference implemented in omegaMap provides a perfect framework for testing whether omega values are significantly different in sympatry vs. allopatry.

Using omegaMap, we calculated omega (dN/dS) values for sympatric Type A, sympatric Type B, allopatric Type A and allopatric Type B populations separately for each gene. We chose 250,000 iterations for each run, with thinning set to 1,000. We used an improper inverse distribution to specify priors for all parameters. Initial parameter values were μ = 0.1, κ = 3.0, ω = 0.5, ρ = 0.1. A constant model was used, so that all sites are assumed to have the same omega value. The number of iterations discarded as burnin varied across runs, but was determined by plotting the traces of μ and κ; iterations affected by the starting value of the parameter were discarded. Two independent runs were conducted for each population/gene (e.g. CIPRO37.40.1_Sympatric A Run 1 and 2). These two runs were combined in all cases, after it was determined that the mean and 95% credible interval for each parameter in the two runs matched closely.

After the omega values for each population/gene were estimated, we permuted the posterior distribution of CIPRO37.40.1 sympatric Type A omega values and the posterior distribution of CIPRO37.40.1 allopatric Type A omega values (for example), obtained a distribution of the differences between these values, and calculated the mean and 95% credible interval for this distribution. If the 95% credible interval contained zero, the estimated omega values were not significantly different. All calculations were done in R (version 2.10.1).

Sympatric vs. allopatric polymorphism analyses

Both alleles from at least 10 individuals from each of four populations (allopatric Type A and B, sympatric Type A and B) were sequenced for three candidate gamete recognition genes and two control genes as described in the "Sympatric vs. allopatric divergence analyses" section.

For each population and each gene, the summary statistics θ and π were calculated in DnaSP 5.10.1 [58]. We also employed the following tests: McDonald-Kreitman [59], Tajima's D [60], Fu and Li's D* and F* [61], and Fay and Wu's H [62] in DnaSP. Statistical significance of D, D*, F* and H were determined using 1,000 coalescent simulations in DnaSP. Estimates of per gene recombination for each population were made in DnaSP and were then imported into the simulations. 95% confidence intervals for D, D*, F* and H statistics were also recorded; sympatric and allopatric populations were determined to be significantly different for each statistic if the confidence intervals of the sympatric population did not contain the mean of the allopatric population.


  1. 1.

    Brown WL, Wilson EO: Character displacement. Syst Zool. 1956, 5: 49-64. 10.2307/2411924.

    Article  Google Scholar 

  2. 2.

    Howard DJ: Reinforcement: the origin, dynamics, and fate of an evolutionary hypothesis. Hybrid Zones and the Evolutionary Process. Edited by: Harrison RG. 1993, Oxford University Press, New York, 46-69.

    Google Scholar 

  3. 3.

    Rice WR, Hostert EE: Laboratory experiments on speciation: What have we learned in 40 years?. Evolution. 1993, 47: 1637-1653. 10.2307/2410209.

    Article  Google Scholar 

  4. 4.

    Dettman JR, Jacobson DJ, Turner E, Pringle A, Taylor JW: Reproductive isolation and phylogenetic divergence in Neurospora: Comparing methods of species recognition in a model eukaryote. Evolution. 2003, 57: 2721-2741.

    Article  PubMed  Google Scholar 

  5. 5.

    Le Gac M, Giraud T: Existence of a pattern of reproductive character displacement in Homobasidiomycota but not in Ascomycota. J Evol Biol. 2008, 21: 761-772. 10.1111/j.1420-9101.2008.01511.x.

    CAS  Article  PubMed  Google Scholar 

  6. 6.

    Kay K, Schemske DW: Natural selection reinforces speciation in a radiation of neotropical forest plants. Evolution. 2008, 62: 2628-2642. 10.1111/j.1558-5646.2008.00463.x.

    Article  PubMed  Google Scholar 

  7. 7.

    Cooley JR, Marshall DC, Hill KBR, Simon C: Reconstructing asymmetrical reproductive character displacement in a periodical cicada contact zone. J Evol Biol. 2006, 19: 855-868. 10.1111/j.1420-9101.2005.01056.x.

    Article  PubMed  Google Scholar 

  8. 8.

    Nosil PB, Crespi BJ, Sandoval CP: Reproductive isolation driven by the combined effects of ecological adaptation and reinforcement. Proc R Soc Lond B. 2003, 270: 1911-1918. 10.1098/rspb.2003.2457.

    CAS  Article  Google Scholar 

  9. 9.

    Kameda Y, Kawakita A, Kato M: Reproductive character displacement in genital morphology in Satsuma land snails. Am Nat. 2009, 173: 689-697. 10.1086/597607.

    Article  PubMed  Google Scholar 

  10. 10.

    Albert AYK, Schluter D: Reproductive character displacement of male stickleback mate preference: reinforcement or direct selection?. Evolution. 2004, 58: 1099-1107.

    CAS  Article  PubMed  Google Scholar 

  11. 11.

    Gabor CR, Ryan MJ: Geographical variation in reproductive character displacement in mate choice by sailfin mollies. Proc R Soc Lond B. 2001, 268: 1063-1070. 10.1098/rspb.2001.1626.

    CAS  Article  Google Scholar 

  12. 12.

    Lemmon EM: Diversification of conspecific signals in sympatry: geographic overlap drives multidimensional reproductive character displacement in frogs. Evolution. 2009, 63: 1155-1170. 10.1111/j.1558-5646.2009.00650.x.

    Article  PubMed  Google Scholar 

  13. 13.

    Coyne JA, Orr HA: "Patterns of speciation in Drosophila revisited". Evolution. 1997, 51: 295-303. 10.2307/2410984.

    Article  Google Scholar 

  14. 14.

    Coyne JA, Orr HA: Patterns of speciation in Drosophila. Evolution. 1989, 43: 362-381. 10.2307/2409213.

    Article  Google Scholar 

  15. 15.

    Dobzhansky T: Genetics and the origin of species. 1937, New York: Columbia University Press

    Google Scholar 

  16. 16.

    Higgie M, Blows M: The evolution of reproductive character displacement conflicts with how sexual selection operates within a species. Evolution. 2008, 62: 1192-1203. 10.1111/j.1558-5646.2008.00357.x.

    CAS  Article  PubMed  Google Scholar 

  17. 17.

    Noor MAF: Speciation driven by natural selection in Drosophila. Nature. 1995, 375: 674-675. 10.1038/375674a0.

    CAS  Article  PubMed  Google Scholar 

  18. 18.

    Coyne JA, Orr HA: Speciation. 2004, Sunderland: Sinauer

    Google Scholar 

  19. 19.

    Noor MAF: Reinforcement and other consequences of sympatry. Heredity. 1999, 83: 503-508. 10.1038/sj.hdy.6886320.

    Article  PubMed  Google Scholar 

  20. 20.

    Rundle HD, Schluter D: Reinforcement of stickleback mating preferences: Sympatry breeds contempt. Evolution. 1998, 52: 200-208. 10.2307/2410935.

    Article  Google Scholar 

  21. 21.

    Saetre GP, Moum T, Bures S, Kral M, Adamjan M, Moreno J: A sexually selected character displacement in flycatchers reinforces premating isolation. Nature. 1997, 387: 589-592. 10.1038/42451.

    CAS  Article  Google Scholar 

  22. 22.

    Levitan DR, Fukami H, Jara J, Kline D, McGovern TM, McGhee KE, Swanson CA, Knowlton N: Mechanisms of reproductive isolation among sympatric broadcast-spawning corals of the Montastraea annularis species complex. Evolution. 2004, 58: 308-323.

    Article  PubMed  Google Scholar 

  23. 23.

    Geyer LB, Palumbi SR: Reproductive character displacement and the genetics of gamete recognition in tropical sea urchins. Evolution. 2003, 57: 1049-1060.

    Article  PubMed  Google Scholar 

  24. 24.

    Springer SA, Crespi BJ: Adaptive gamete-recognition divergence in a hybridizing Mytilus population. Evolution. 2007, 61: 772-783. 10.1111/j.1558-5646.2007.00073.x.

    CAS  Article  PubMed  Google Scholar 

  25. 25.

    Swanson WJ, Vacquier VD: Reproductive protein evolution. Ann Rev Ecol Syst. 2002, 33: 161-179. 10.1146/annurev.ecolsys.33.010802.150439.

    Article  Google Scholar 

  26. 26.

    Nydam ML, Harrison RG: Introgression despite substantial divergence in a broadcast spawning marine invertebrate. Evolution. 2010,

    Google Scholar 

  27. 27.

    Nydam ML, Harrison RG: Genealogical relationships within and among shallow-water Ciona species (Ascidiacea). Mar Biol. 2007, 151: 1839-1847. 10.1007/s00227-007-0617-0.

    Article  Google Scholar 

  28. 28.

    Caputi L, Andreakis N, Mastrototaro F, Cirino P, Vassillo M, Sordino P: Cryptic speciation in a model invertebrate chordate. Proc Natl Acad Sci USA. 2007, 104: 9364-9369. 10.1073/pnas.0610158104.

    Article  PubMed  PubMed Central  Google Scholar 

  29. 29.

    Van Name W: The north and south American ascidians. B Am Mus Nat Hist. 1945, 84: 1-476.

    Google Scholar 

  30. 30.

    Kott P: The ascidians of Australia. Aus J Mar Fresh Res. 1952, 3: 206-233.

    Google Scholar 

  31. 31.

    Linne C: Vindobonae: Typis Ioannis Thomae, 1767-1770. Systema Naturae. 1767

    Google Scholar 

  32. 32.

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

    Google Scholar 

  33. 33.

    Inaba K, Padma P, Satouh Y, Shin-I T, Kohara Y, Satoh N, Satou Y: EST analysis of gene expression in testis of the ascidian Ciona intestinalis. Mol Repro Dev. 2002, 62: 431-445. 10.1002/mrd.10131.

    Article  Google Scholar 

  34. 34.

    Pond SK, Muse SV: Site-to-site variation of synonymous substitution rates. Mol Biol Evol. 2005, 22: 2375-2385. 10.1093/molbev/msi232.

    CAS  Article  PubMed  Google Scholar 

  35. 35.

    Swanson WJ, Wong A, Wolfner MF, Aquadro CF: Evolutionary expressed sequence tag analysis of Drosophila female reproductive tracts identifies genes subject to positive selection. Genetics. 2004, 168: 1457-1465. 10.1534/genetics.104.030478.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  36. 36.

    Clark NL, Gasper J, Sekino M, Springer SA, Aquadro CF, Swanson WJ: Coevolution of interacting fertilization proteins. PLOS Genetics. 2009, e:1000570-

    Google Scholar 

  37. 37.

    Walters JR, Harrison RG: Combined EST and proteomic analysis identifies rapidly evolving seminal fluid proteins in Heliconius butterflies. Mol Biol Evol. 2010

    Google Scholar 

  38. 38.

    Andres JA, Maroja LS, Harrison RG: Searching for candidate speciation genes using a proteomic approach: seminal proteins in field crickets. Proc R Soc Lond B. 2008, 275: 1975-1983. 10.1098/rspb.2008.0423.

    CAS  Article  Google Scholar 

  39. 39.

    Torgerson DG, Kulathinal RJ, Singh RS: Mammalian sperm proteins are rapidly evolving: evidence of positive selection in functionally diverse genes. Mol Biol Evol. 2002, 19: 1973-1980.

    CAS  Article  PubMed  Google Scholar 

  40. 40.

    Wyckoff GJ, Wang W, Wu CI: Rapid evolution of male reproductive genes in the descent of man. Nature. 2000, 403: 304-309. 10.1038/35002070.

    CAS  Article  PubMed  Google Scholar 

  41. 41.

    Lambert CC, Goudeau H, Franchet C, Lambert G, Goudeau M: Ascidian eggs block polyspermy by two independent mechanisms: one at the plasma membrane, the other involving the follicle cells. Mol Reprod Dev. 1997, 48: 137-143. 10.1002/(SICI)1098-2795(199709)48:1<137::AID-MRD16>3.0.CO;2-Y.

    CAS  Article  PubMed  Google Scholar 

  42. 42.

    Bigelow RS: Hybrid zones and reproductive isolation. Evolution. 1965, 19: 449-458. 10.2307/2406242.

    Article  Google Scholar 

  43. 43.

    Nydam ML, Harrison RG: Polymorphism and divergence within the ascidian genus Ciona. Mol Phyl Evol. 2010, 56: 718-726. 10.1016/j.ympev.2010.03.042.

    Article  Google Scholar 

  44. 44.

    Satoh N: Developmental Biology of Ascidians. 1994, Cambridge: Cambridge University Press, 1

    Google Scholar 

  45. 45.

    Sawada H, Yokosawa H, Lambert CC, eds: The Biology of Ascidians. 2000, New York: Springer, 1

    Google Scholar 

  46. 46.

    Hoshi M, De Santis R, Pinto R, Cotelli F, Rosati F: Sperm glycosidases as mediators of sperm-egg interaction in the ascidians. Zool Sci. 1985, 2: 65-69.

    CAS  Google Scholar 

  47. 47.

    Marino R, De Santis R, Hirohashi N, Hoshi M, Pinto MR, Usui N: Purification and characterization of a vitelline coat lysin from Ciona intestinalis spermatozoa. Mol Repro Dev. 1992, 32: 383-388. 10.1002/mrd.1080320412.

    CAS  Article  Google Scholar 

  48. 48.

    Suzuki M, Nishikawa T, Bird A: Genomic approaches reveal unexpected difference within Ciona intestinalis. J Mol Evol. 2005, 61: 627-635. 10.1007/s00239-005-0009-3.

    CAS  Article  PubMed  Google Scholar 

  49. 49.

    Boursot P, Auffray JC, Britton-Davidian J, Bonhomme F: The evolution of house mice. Ann Rev Ecol Syst. 1993, 24: 119-152. 10.1146/

    Article  Google Scholar 

  50. 50.

    Vacquier VD, Moy GW: Isolation of bindin: the protein responsible for adhesion of sperm to sea urchin eggs. Proc Natl Acad Sci. 1977, 74: 2456-2460. 10.1073/pnas.74.6.2456.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  51. 51.

    Lessios HA: Reproductive isolation between species of sea urchins. Bull Mar Sci. 2007, 81: 191-208.

    Google Scholar 

  52. 52.

    Lee Y H, Ota T, Vacquier VD: Positive selection is a general phenomenon in the evolution of abalone sperm lysin. Mol Biol Evol. 1995, 12: 231-238.

    PubMed  Google Scholar 

  53. 53.

    Yang Z, Swanson WJ, Vacquier VD: Maximum-Likelihood analysis of molecular adaptation in abalone sperm lysin reveals variable selective pressures among lineages and sites. Mol Biol Evol. 2000, 17: 1446-1455.

    CAS  Article  PubMed  Google Scholar 

  54. 54.

    Yang Z: PAML 4: Phylogenetic analysis by maximum likelihood. Mol Biol Evol. 2007, 24: 1586-1591. 10.1093/molbev/msm088.

    CAS  Article  PubMed  Google Scholar 

  55. 55.

    Conover WJ: Practical Nonparametric Statistics. 1999, New York; Wiley Publishing, 3

    Google Scholar 

  56. 56.

    Shapiro SS, Wilk MB: An analysis of variance test for normality. Biometrika. 1965, 52: 591-599.

    Article  Google Scholar 

  57. 57.

    Wilson DJ, McVean G: Estimating diversifying selection and functional constraint in the presence of recombination. Genetics. 2006, 172: 1411-1425. 10.1534/genetics.105.044917.

    CAS  Article  PubMed  PubMed Central  Google Scholar 

  58. 58.

    Rozas J, Sanchez-DelBarrio JC, Messeguer X, Rozas R: DnaSP, DNA polymorphism analyses by the coalescent and other methods. Bioinformatics. 2003, 19: 2496-2497. 10.1093/bioinformatics/btg359.

    CAS  Article  PubMed  Google Scholar 

  59. 59.

    McDonald JH, Kreitman M: Adaptive protein evolution at the Adh locus in Drosophila. Nature. 1991, 351: 652-654. 10.1038/351652a0.

    CAS  Article  PubMed  Google Scholar 

  60. 60.

    Tajima F: Statistical method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics. 1989, 123: 585-595.

    CAS  PubMed  PubMed Central  Google Scholar 

  61. 61.

    Fu YX, Li WH: Statistical tests of neutrality of mutations. Genetics. 1993, 133: 693-709.

    CAS  PubMed  PubMed Central  Google Scholar 

  62. 62.

    Fay JC, Wu CI: Hitchhiking under positive Darwinian selection. Genetics. 2000, 155: 14005-1413.

    Google Scholar 

Download references


We would like to thank Erin Newman-Smith at the Ascidian Stock Center (supported by grant # R24GM075049 from the National Institutes of Health) for supplying animals, the Smith Laboratory at UCSB and John Bishop at Plymouth MBA for providing logistical assistance, Dr. Kazuo Inaba for sharing information about the CIPRO database, Derek Smith and the staff of the University of Victoria Genome BC Proteomics Centre for generating the proteomics data, members of the Harrison Laboratory and Jose Andres for advice, and two anonymous reviewers for critiquing the manuscript. omegaMap was accessed through the Computational Biology Service Unit's Web Computing Resources at Cornell University. This work was funded by NSF Doctoral Dissertation Grant DEB-0907862 to M.L.N. and R.G.H. and NSF Grant DEB-0639904 to R.G.H.

Author information



Corresponding author

Correspondence to Marie L Nydam.

Additional information

Authors' contributions

MLN collected the samples, generated the data, performed the analyses and wrote the manuscript. RGH provided funding, and assisted with the design of the study and the editing of the manuscript. Both authors have read and approved the final manuscript.

Authors’ original submitted files for images

Below are the links to the authors’ original submitted files for images.

Authors’ original file for figure 1

Rights and permissions

This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Reprints and Permissions

About this article

Cite this article

Nydam, M.L., Harrison, R.G. Reproductive protein evolution in two cryptic species of marine chordate. BMC Evol Biol 11, 18 (2011).

Download citation


  • Control Protein
  • Sympatric Population
  • Allopatric Population
  • Sperm Protein
  • Fertilization Process