Genetic connectivity of the endangered brown sea cucumber Isostichopus fuscus in the northern Gulf of California revealed by novel microsatellite markers

José Miguel Ochoa-Chávez a, Miguel Ángel Del Río-Portilla b, Luis Eduardo Calderón-Aguilera c, Axayácatl Rocha-Olivares a, *

HTML PDF XML

a Department of Biological Oceanography, Molecular Ecology Laboratory, Centro de Investigación Científica y de Educación Superior de Ensenada, Carretera Ensenada-Tijuana No. 3918, 22860 Ensenada, Baja California, Mexico
b Department of Aquaculture, Centro de Investigación Científica y de Educación Superior de Ensenada, Carretera Ensenada-Tijuana No. 3918, 22860 Ensenada, Baja California, Mexico
c Department of Marine Ecology, Coastal Zone Ecology and Fisheries Laboratory, Centro de Investigación Científica y de Educación Superior de Ensenada, Carretera Ensenada-Tijuana No. 3918, 22860 Ensenada, Baja California, Mexico

* Corresponding author: arocha@cicese.mx (A. Rocha-Olivares)

Abstract

Isostichopus fuscus is an economically important sea cucumber that has been highly exploited along its distribution range in the eastern Pacific. The significant population decline is responsible for its listing as endangered in the IUCN Red List. Despite its importance for management and conservation, information about its population genetic structure is unavailable, largely due to the lack of suitable genetic markers. Here we develop species-specific microsatellite markers and use them to assess the genetic connectivity between populations in the Gulf of California. Next generation sequencing (Illumina) was used to shotgun-sequence the genome of 2 sea cucumbers. From these data, we identified and characterized 19 polymorphic microsatellite loci; which were tested in organisms from Bahía de los Ángeles, on the western shore of the Gulf of California. The number of alleles ranged from 5 to 22, observed heterozygosity from 0.35 to 1, and 4 loci deviated from Hardy-Weinberg equilibrium. We determined high levels of genetic connectivity between this locality and San Felipe, in the upper gulf (Amova φst = 0.002; p > 0.05) with a subset of 8 markers. The newly designed microsatellites are suitable for multiplexing panels and will be useful for the future genetic assessment of this important tropical sea cucumber.

Keywords:

Holothurid; Fisheries; Conservation genetics; Next generation sequencing; Genome sequencing

(http://creativecommons.org/licenses/by-nc-nd/4.0/).

Conectividad genética del pepino de mar café Isostichopus fuscus en peligro de extinción en el norte del golfo de California usando nuevos marcadores microsatelitales

Resumen

Isostichopus fuscus es un pepino de mar de importancia económica que ha sido altamente explotado a lo largo de su distribución en el Pacífico oriental. La reducción significativa en sus poblaciones es responsable de su inclusión en la Lista Roja de la UICN como especie amenazada de extinción. A pesar de la importancia que tiene para su manejo y conservación, no hay información sobre su estructura genética poblacional, en gran parte debido a la falta de marcadores genéticos adecuados. Para este trabajo desarrollamos marcadores microsatelitales específicos para la especie y los usamos para evaluar la conectividad genética entre poblaciones del golfo de California. Usamos secuenciación de siguiente generación (Illumina) para secuenciar por “shotgun” el genoma de 2 pepinos de mar. A partir de esos datos, identificamos y caracterizamos 19 loci microsatelitales polimórficos que fueron evaluados en organismos de Bahía de los Ángeles, en la costa occidental del golfo de California. El número de alelos varió entre 5 y 22, la heterocigosidad observada entre 0.35 y 1, y 4 loci se desviaron del equilibrio de Hardy-Weinberg. Usando 8 de estos marcadores se determinaron altos niveles de conectividad genética entre Bahía de los Ángeles y San Felipe, en el alto golfo (Amova φst = 0.002; p > 0.05). Los nuevos loci microsatelitales son adecuados para ser usados en páneles de multiplexación y serán útiles para la futura evaluación genética de este importante pepino del mar tropical.

(http://creativecommons.org/licenses/by-nc-nd/4.0/).

Palabras clave:

Holotúrido; Pesquerías; Genética de la conservación; Secuenciación de siguiente generación; Secuenciación genómica

Introduction

Isostichopus fuscus Ludwig 1875, is a commercially important species due to its high demand and price in Asian markets, where it is considered a delicacy. It is distributed in the eastern central Pacific coast from the Gulf of California, Mexico, to Ecuador including the Galapagos Islands (Maluf, 1988). Its shallow habitat and sedentary habits make this sea cucumber easily accessible to fishing efforts. Populations of I. fuscus have severely declined due to overfishing in México (Ibarra & Soberón, 2002) and Ecuador, including Galápagos (Toral-Granda, 2008). Consequently, I. fuscus has been included as endangered in the IUCN Red List. Despite its importance, the genetic status of its populations has been little studied, with only a handful of samples having been analyzed throughout its distribution (Lohr, 2003), and in Mexico it remains unknown. Here, we present the first I. fuscus species-specific polymorphic microsatellite loci identified from genomic data produced with next generation sequencing (NGA, Illumina), and use them to test genetic structure in the northern Gulf of California. These markers will prove useful for a variety of genetic studies.

Materials and methods

DNA extraction from the tentacles of 2 organisms collected in La Paz, Mexico (24°44’48.58” N, 110°40’36.27” W) was performed with lithium chloride protocol (Hong et al., 1995). De novo assembly of 84,843,094 reads obtained from NGS massive parallel sequencing (HiSeq 2500 Illumina) resulted in a total of 581,221 contigs averaging 898 base pairs (bp) in length. Contigs were processed with MSATCOMMANDER (Faircloth, 2008) searching for potential di- to hexa-nucleotide microsatellite loci with ≥ 6 repeats and enough flanking regions for primer design. A total of 11,224 microsatellite loci were obtained, and unique primers were designed for 4,920 of them. Subsequently, di-, tri- and tetra-nucleotide microsatellite loci were selected based on coverage (number of reads stacked per residue in contig) ranging from 10 to 100 and the number of motif repeats ≥ 10, to increase chances of high polymorphism in loci and their usefulness for population genetics studies (Gardner et al., 2011); producing 29 candidate microsatellites (15 di, 8 tri and 6 tetra-nucleotide) loci to be tested for amplification success and polymorphism.

DNA was extracted from tentacles and tube feet of 22 organisms from Bahía de los Ángeles (BLA) (24°44’48.58” N, 110°40’36.27” W), Mexico and from 19 organisms from San Felipe (SF) (29°48’18” N, 114°22’19” W), Mexico using PureLink ™ Genomic DNA Mini Kit (Invitrogen). The M13 primer sequence (TGTAAAACGACGGCCAGT) was attached to the 5’ end of each forward primer in order to use fluorescent dyes (FAM, NED, VIC, PET) following Schuelke (2000). PCR conditions consisted of 1 min initial denaturation at 94 °C, 30 cycles of 30 s denaturation at 94 °C, 30 s annealing at 63 °C (for all loci) and 30 s elongation at 72 °C, with 15 min elongation at 72 °C at the end of the last cycle. A second PCR was performed replacing the forward primer with FAM, NED, VIC or PET marked M13 primers, using the same conditions but an annealing temperature of 53 °C. Amplification products were subject to fragment analysis in an ABI Prism 3100 genetic analyzer (Applied Biosystems) and genotyped using GeneMarker v. 2.4.0 (Holland & Parson, 2011).

We used organisms from BLA to fully characterize microsatellite polymorphism in all loci, where the number of alleles per locus was evaluated with MStools (Park, 2001). Observed (H_O), and expected heterozygosities (H_E), Hardy-Weinberg equilibrium (HWE) and linkage disequilibrium (LD) were assessed using ARLEQUIN v. 3.5 (Excoffier et al., 2005) and heterozygote deficiency (Fis) with GENEPOP v. 4.4 (Rousset, 2008). Micro-checker v. 2.2.3 (Van Oosterhout et al., 2004) was used to infer the presence of scoring errors, large allele dropout and null alleles.

Finally a subset of 8 microsatellites (Ifus2-01, Ifus2-06, Ifus2-11, Ifus2-14, Ifus3-06, Ifus4-01, Ifus4-03 and Ifus4-06) was used to assess genetic diversity (N_A, H_O, H_E and HWE) in SF, and its genetic differentiation (φst) with BLA with an analysis of molecular variance (Amova) using ARLEQUIN v. 3.5 (Excoffier et al., 2005).

Table 1

Primers and polymorphism characterization of 19 novel microsatellite loci for Isostichopus fuscus in samples from Bahia de Los Angeles.

Locus	Primer sequences (5′ – 3′)	Repeat motif	Size range (bp)	Dye	N	N_A	H_O	H_E	P_HW	Fis	Null
Ifus2-01	F: AGCAGTTTCGGATGGTTGAG R: GTGACCGCGCCTTAATAACC	(AC)12	139-153	FAM	20	7	0.70	0.67	0.90	-0.05	-0.04
Ifus2-02	F: GCGCTTTCTAACACCTCGAC R: AACCTTTGTATTCACGCGGG	(AC)11	444-460	NED	19	7	0.68	0.77	0.01	0.11	0.02
Ifus2-03	F: AACCATCGTGATTCACCGTG R: TATGTCAAAGCCTCCCTCCC	(AC)11	440-468	VIC	20	10	0.95	0.86	0.02	-0.11	-0.08
Ifus2-04	F: AAGCGTCCATCATTGCTACC R: CGATTGTGTGGTACGGAAGC	(AC)10	290-304	VIC	22	6	0.50	0.56	0.24	0.11	0.02
Ifus2-05	F: TCTGATGGAGCTGGTCTAGG R: CGGTGGTTCTGAATCTGCAC	(AG)14	448-474	VIC	22	9	0.77	0.87	<0.001*	0.11	0.04
Ifus2-06	F: TGGTTCCCAGAAGATGAGTGG R: TGACCAAGCCAGCCTATGAC	(AG)13	374-386	NED	20	6	0.35	0.64	<0.001*	0.46	0.21*
Ifus2-07	F: ACTGACTTACGTACATGCCAG R: ACAGAGCACGCGTAATGAGG	(AG)10	441-453	FAM	22	7	0.91	0.81	0.06	-0.13	-0.09
Ifus2-10	F: TGCCAACTAATTTCACGGGC R: GGGAGGGCTGAGATTTGTAAG	(AT)10	392-402	PET	22	5	0.55	0.61	0.39	0.11	0.05

Locus	Primer sequences (5′ – 3′)	Repeat motif	Size range (bp)	Dye	N	N_A	H_O	H_E	P_HW	Fis	Null
Ifus2-11	F: TCACAATGATGCATGCGGTG R: TCATTACAACACTGGCAACCG	(AG)10	286-296	PET	20	6	0.40	0.57	0.06	0.31	0.15
Ifus2-14	F: TTTACAGCAATCGTCACGCC R: CACATTCATCTGGGCACGTG	(AC)12	360-386	VIC	20	10	0.80	0.84	0.88	0.05	0.01
Ifus3-01	F: AGTGCAGTGTACTCTCAGGC R: TGGGTTTCCTGCTAAGCTTC	(AGG)11	238-262	NED	20	9	0.65	0.88	0.004	0.26	0.12*
Ifus3-02	F: TCTTCGTCCTCTCGCTCAAC R: TGACGATGATGAAGGTGGCC	(ATC)21	361-412	FAM	17	13	0.47	0.94	<0.001*	0.50	0.24*
Ifus3-04	F: AATGGGTGTATTTGCCGCTG R: TCGCTTGCACCATCATTGTC	(ATC)20	295-352	FAM	21	15	0.62	0.92	0.003	0.33	0.16*
Ifus3-06	F: TGTGTTTGGCCGAGTTCAAC R: AGCAGTAAGGTTGGTCTGAGG	(ATC)11	170-233	NED	19	10	0.68	0.82	0.08	0.17	0.07
Ifus4-01	F: CGACGTCACTATGCCACAAG R: CACACTGCTCTGACGTTCAC	(ACAG)11	353-389	NED	20	9	0.60	0.85	0.04	0.30	0.14*
Ifus4-02	F: GGGAGGCCAATCAGATGTTG R: CATGTTTCAAGATGGCAACCC	(ACAG)18	403-447	VIC	15	12	0.73	0.91	0.23	0.20	0.08
Ifus4-03	F: ACAAAGCAGGCATCACAAAC R: AGGACCAATGAAACCACTGTG	(ACAT)12	273-389	VIC	20	22	1.00	0.97	1.00	-0.03	-0.03
Ifus4-04	F: ACTGCATATAGAACGCGTGC R: TGAGGTAAATTTGGGCCAAGG	(AGAT)18	396-544	VIC	12	15	0.42	0.97	<0.001*	0.58	0.27*
Ifus4-06	F: ACGAGTTTACAGAGGTGCCC R: CTGTCAATCCATCTGTCTGCC	(AGAT)18	186-278	FAM	20	16	0.90	0.91	0.12	0.02	0.00

Dye, fluorescent dye of M13 primer; N, number of individuals with successful amplification; N_A, number of alleles; H_O, observed heterozygosity; H_E, expected heterozygosity; P_HW, P-value for Hardy-Weinberg equilibrium (* denotes significant deviation after Bonferroni correction, p < 0.0026); F_IS, heterozygote deficiency; Null, null allele frequency (* inference of null alleles). GenBank accession numbers KX525270 to KX525288.

Results

Optimization was possible for nineteen (di-, tri-, and tetra-nucleotide repeats) microsatellite loci (Table 1). High polymorphism was evident in all loci, showing a number of alleles from 5 to 22, with an average of 10.2. Heterozygosities ranged from 0.35 to 1 (H_O) and from 0.56 to 0.97 (H_E). After Bonferroni correction, 4 loci (Ifus2-05, Ifus2-06, Ifus3-02 and Ifus4-04) significantly deviated from HWE, and the pair Ifus2-01- Ifus2-03 showed to be potentially linked. Scoring errors and large allele dropout were undetected; however, evidence of null alleles was found in 6 loci.

The 8 loci assessed in samples from SF showed levels of variation similar to BLA. The number of alleles ranged from 2 to 19 (9.5 average). Heterozygosities varied from 0.23 to 0.83 (H_O) and from 0.32 to 0.95 (H_E). Loci Ifus2-14 and Ifus4-01 deviated from HWE, showing heterozygote deficiency associated to null alleles and couples Ifus2-01 – Ifus2-11, Ifus3-06 – Ifus4-03, Ifus2-01 – Ifus2-06 and Ifus3-06 – Ifus4-01 presented potential LD. The Amova results showed no differentiation between BLA and SF localities with only 0.2 % of genetic variance related to interpopulation level (Table 2).

Table 2

Amova results of population differentiation (φst), using 8 novel microsatellite loci, between brown sea cucumber samples from Bahía de Los Angeles and San Felipe localities in the northern Gulf of California.

Source of variation	Variance component	Percentage of variation	φst	p-value
Among populations	0.01	0.2	0.002	0.47
Within populations	2.47	99.8

Discussion

Optimization was possible for 19 (di-, tri-, and tetra-nucleotide repeats) microsatellite loci and involved designing multiplexing panels of loci using the same annealing temperature, producing different amplicon sizes, as well as combining a variety of fluorescent dyes (Table 1).

The high genetic diversity and the absence of genetic differentiation between BLA and SF suggest the existence of a large panmictic population in the region and the presence of high connectivity in the recent evolutionary history. Despite the decrease of population density of I. fuscus in this area in the past 2 decades, genetic connectivity could have been sustained by the short distance between locations (≈150 km) and the highly dynamic circulation in the Ballenas Channel and the northern Gulf of California (Lavín & Marinone, 2003).

A broader evaluation considering a larger study area and several molecular markers is needed to determine the genetic status of I. fuscus in the Gulf of California. The new microsatellite markers reported herein will be useful for future population genetics analyses of this endangered and economically important sea cucumber along its distribution.

Acknowledgments

This research was financially supported by “Fondo sectorial de investigación en materia agrícola, pecuaria, acuacultura, agrobiotecnología y recursos fitogenéticos” (Sagarpa-Conacyt 0223255). The first author benefited from a graduate fellowship from Conacyt to support his M. Sc. program in Marine Ecology at Centro de Investigación Científica y de Educación Superior de Ensenada (CICESE). We thank the research staff from the Coastal Zone Ecology and Fisheries Laboratory at CICESE for fieldwork and sampling and Hector Reyes Bonilla (UABCS) and Dinorah Herrero Perezrul for providing specimens for genomic library construction.

References

Excoffier, L., Laval, G., & Schneider, S. (2005). Arlequin ver. 3.0: An integrated software package for population genetics data analysis. Evolutionary Bioinformatics Online, 1, 47–50.

Faircloth, B. C. (2008). MSATCOMMANDER: detection of microsatellite repeat arrays and automated, locus-specific primer design. Molecular Ecology Resources, 8, 92–94.

Gardner, M. G., Fitch, A. J., Bertozzi, T., & Lowe, A. J. (2011). Rise of the machines – recommendations for ecologists when using next generation sequencing for microsatellite development. Molecular Ecology Resources, 11, 1093–1101.

Holland, M. M., & Parson, W. (2011). GeneMarker ® HID: a reliable software tool for the analysis of forensic STR data. Journal of Forensic Sciences, 56, 29–35.

Hong, Y. K., Kim, S. D., Polnefuller, M., & Gibor, A. O. (1995). DNA extraction conditions from Porphyra perforata using LiCl. Journal of Applied Phycology, 7, 101–107.

Ibarra, A. A., & Soberón, G. R. (2002). Economic reasons, ecological actions and social consequences in the Mexican sea cucumber fishery. SPC Beche-de-Mer Information Bulletin, 17, 33–36.

Lavín, M. F., & Marinone, S. G. (2003). An overview of the physical oceanography of the Gulf of California. En O. Velasco-Fuentes, J. Sheinbaum, & J. Ochoa (Eds.), Nonlinear processes in geophysical fluid dynamics (pp. 173–204). Dordrecht, Netherlands: Springer.

Lohr, H. R. (2003). Genetic variation among geographically isolated populations of the commercially important sea cucumber, Isostichopus fuscus, in the Eastern Pacific (M. S. Thesis). San Francisco State University, CA, USA.

Maluf, L. Y. (1988). Composition and distribution of the central eastern Pacific echinoderms. Los Angeles, CA: Natural History Museum of Los Angeles County.

Park, S. (2001). MStools v 3 (Excel spreadsheet toolkit for data conversion). Smurfit Institute of Genetics. Trinity College, Dublin.

Rousset, F. (2008). GENEPOP ‘ 007: a complete re-implementation of the GENEPOP software for Windows and Linux. Molecular Ecology Resources, 8, 103–106.

Schuelke, M. (2000). An economic method for the fluorescent labeling of PCR fragments. Nature Biotechnology, 18, 233–234.

Toral-Granda, V. (2008). Population status, fisheries and trade of sea cucumbers in Latin America and the Caribbean. In V. Toral-Granda, A. Lovatelli, & M. Vasconcellos (Eds.), Sea cucumbers. A global review of fisheries and trade. FAO Fisheries and Aquaculture Technical Paper (pp. 213–229). Rome: FAO.

Van Oosterhout, C., Hutchinson, W. F., Wills, D. P. M., & Shipley, P. (2004). MICRO-CHECKER: software for identifying and correcting genotyping errors in microsatellite data. Molecular Ecology Notes, 4, 535–538.