SciELO - Scientific Electronic Library Online

vol.37Type specimens of Limnophorini (Diptera: Muscidae) deposited in the Museum für Naturkunde, Humboldt-Universität zu Berlin (Berlin, Germany)Occurrence and ecological implication of a tropical anguillid eel, Anguilla marmorata, in Brunei Darussalam, Borneo Island author indexsubject indexarticles search
Home Pagealphabetic serial listing  

Services on Demand




Related links


Zoologia (Curitiba)

Print version ISSN 1984-4670On-line version ISSN 1984-4689

Zoologia (Curitiba) vol.37  Curitiba  2020  Epub Sep 07, 2020 


Genetic analysis of whole mitochondrial genome of Lateolabrax maculatus (Perciformes: Moronidae) indicates the presence of two populations along the Chinese coast

1Fujian Key Laboratory of Genetics and Breeding of Marine Organisms, Xiamen University. 361000 Xiamen, China.

2Shenzhen Research Institute, Xiamen University. 518000 Shenzhen, China.

3State Key Laboratory of Large Yellow Croaker Breeding, Ningde Fufa Fisheries Company Limited. 352130 Ningde, China.

4South China Sea Fisheries Research Institute, Chinese Academy of Fishery Sciences. 510300 Guangzhou, China.


The whole mitochondrial genome of Lateolabrax maculatus (Cuvier, 1828) was used to investigate the reasons for the observed patterns of genetic differentiation among 12 populations in northern and southern China. The haplotype diversity and nucleotide diversity of L. maculatus were 0.998 and 0.00169, respectively. Pairwise FST values between populations ranged from 0.001 to 0.429, correlating positively with geographic distance. Genetic structure analysis and haplotype network analysis indicated that these populations were split into two groups, in agreement with geographic segregation and environment. Tajima’s D values, Fu’s Fs tests and Bayesian skyline plot (BSP) indicated that a demographic expansion event may have occurred in the history of L. maculatus. Through selection pressure analysis, we found evidence of significant negative selection at the ATP6, ND3, Cytb, COX3, COX2 and COX1 genes. In our hypotheses, this study implied that demographic events and selection of local environmental conditions, including temperature, are responsible for population divergence. These findings are a step forward toward the understanding of the genetic basis of differentiation and adaptation, as well as conservation of L. maculatus.

KEY WORDS: Adaptation; genetic differentiation; selection


The Chinese sea bass, Lateolabrax maculatus (Cuvier, 1828), is an endemic species in Northeast Asia, mainly distributed along the coasts of China, Japan and Korea (Niu et al. 2015). The species is associated with reefs in shallow oceanic waters that widely vary in temperature and salinity. This fish spawns on rocky reefs and juveniles may enter rivers. Because of its delicious meat, fast growth rate and short breeding period, L. maculatus is a valuable recreational and aquaculture fishery species. It is commonly reared in cages and ponds in China, with yields exceeding 10×104 tons from 2010 to 2014. In recent years, cage-reared L. maculatus have been imported from China into Korea and Japan for human consumption (An et al. 2014). Germplasm degradation and a reduction in disease resistance of cultured L. maculatus restricts its sustainable development (Wang et al. 2017), and calls for an assessment of the population structure of this fish (Barhoumi et al. 2017). The results of such study can be used to prevent overfishing and ultimately declines in the numbers of L. maculatus.

An earlier study using isozyme and RAPD analyses had revealed differences among L. maculatus Chines populations from north to south (Li et al. 2005). In addition, preliminary phylogenetic analysis demonstrated that populations of L. maculatus clustered into different clades, and suggested low genetic diversity within populations (Zhao et al. 2018). Genetic diversity and population structure were investigated based on a single mitochondrial gene (COI) in five geographical populations (Wang et al. 2017). The maternally inherited mitochondrial genome, with not only low recombination frequency but also high mutation rate, is a feasible marker to evaluate the genetic relationships among individuals or populations (Dong et al. 2015, Stoljarova et al. 2016).

Revealing the phenotypic and genetic differentiation of a species in distinct populations is a striking topic in evolutionary biology. Convincing evidence suggests that natural selection usually results in adapted phenotypes and genetic divergence (local adaptation) among populations inhabiting diverse environments (Hélène and Luca 2011, Hansen 2010, Kawecki and Ebert 2010). Several studies have ascertained selection in the mitochondrial genes of marine fish. The Cytochrome c oxidase (COX2) was selected to analyze tuna and billfish populations because of the high aerobic performance required of them in natural conditions (Dalziel et al. 2006). To analyze populations of Pacific salmon species, the ND2 and ND5 genes were used since they are under positive selection, illustrating the considerable role of the ND gene complex in the evolution (Garvin et al. 2011).

The whole mitochondrial genome was used to analyze the population structure and selection in the Atlantic herring and Japanese Swellshark. In contrast, a systematic study of the whole mitochondrial genome of L. maculatus has not yet been performed. A more comprehensive sequencing effort involving the mitochondrial genome may increase the resolution of not only population structure but also intraspecific selection dynamics of marine fish.

In this study, we sequenced a large number (n = 85) of whole mitochondrial genome sequences of L. maculatus from 12 locations, covering the entire Chinese coast (Fig. 1). The population structure, genetic diversity and selection pressure were analyzed applying the whole mitochondrial genome. We also compared 12 protein-coding genes in the mitochondria between northern and southern populations and to identify potential genes implicated in divergence and environmental adaptation. The results will be critical and helpful to resource conservation and fisheries management of this target species for stock enhancement in the China Sea (Chen et al. 2019).

Figure 1 Locations of 12 sampling sites of Lateolabrax maculatus. 85 individuals of L. maculatus were collected from 12 geographic locations including Tianjin (TJ), Yantai (YT), Wendeng (WD), Lianyungang (LY), Zhoushan (ZS), Wenzhou (WZ), Shantou (ST), Shenzhen (SZ), Zhanjiang (ZJ), Haikang (HK), Tieshan (TS), Fangcheng (FC) along Chinese coastal waters. The average sea surface temperature (SST) of 1981-2010 at 12 sampling sites were retrieved from NOAA Optimum Interpolation (OI) Sea Surface Temperature (SST) V2 ( 


Sampling and sequencing

A total of 85 individuals of L. maculatus were collected from 12 geographic locations including Tianjin (TJ), Yantai (YT), Wendeng (WD), Lianyungang (LY), Zhoushan (ZS), Wenzhou (WZ), Shantou (ST), Shenzhen (SZ), Zhanjiang (ZJ), Haikang (HK), Tieshan (TS), Fangcheng (FC) along Chinese coastal waters (Fig. 1). The sample numbers of every location are shown in Table 1. After eugenol anesthesia, white muscles were collected and used for DNA extraction (Zhou et al. 2019). Total DNA extracted from muscle was used to construct a library. Sequencing libraries with 250 bp insert size were constructed according to the manufacturer’s instructions (Illumina, San Diego, CA, USA) and paired-end sequencing was performed using the Illumina HiSeq2500 platform with a read length of 2 × 150 bp. All reads containing adaptor sequences were discarded before the removal of uncertain or low-quality bases using SolexaQA++ (version v. (Cox et al. 2010). The retained and trimmed reads (clean reads) were used for variants calling.

Table 1 Genetic diversity parameters for 12 populations of Lateolabrax maculatus based on whole-mitochondrial sequences. 

Locality Abbreviation Numbers Number of haplotypes (h) Sample number of variable sites (S) Haplotype diversity (Hd) Nucleotide diversity (π) Tajima’s D Fu’s Fs
Fangchenggang FC 9 9 37 1.000 0.00068 -0.847 -0.784
Tieshangang TS 10 10 42 1.000 0.00058 -1.711* -3.595
Haikang HK 10 9 38 0.978 0.00062 -1.142 -0.250*
Zhanjiang ZJ 4 4 14 1.000 0.00044 -0.403 -0.037
Shenzhen SZ 5 5 50 1.000 0.00155 0.407 0.814
Shantou ST 4 4 12 1.000 0.00036 -0.841 -0.288
Wenzhou WZ 4 4 44 1.000 0.00140 -0.361 1.273
Zhoushan ZS 5 5 55 1.000 0.00151 -0.401 0.778
Lianyungang LY 4 4 80 1.000 0.00252 -0.432 1.898
Yantai YT 10 10 97 1.000 0.00159 -1.136 0.347
Tianjin TJ 10 9 129 0.978 0.00205 -1.291* 0.799
Wendeng WD 10 10 140 1.000 0.00252 -1.049 -0.544
All 85 78 243 0.998 0.00169 -1.451* -26.785

South population including FC, TS, HK, ZJ, SZ, ST. North population including WZ, ZS, LY, YT, TJ, WD. *Significant at level (p < 0.05).

This study was approved by the Animal Care and Use committee at College of Ocean and Earth Sciences, Xiamen University. All the methods used in this study were carried out in accordance with approved guidelines.

Assembling mitochondrial genome assembly and annotation

Producing mitochondrial genome assembly from whole genome sequencing (WGS) data is an accurate and efficient approach. Data from the whole genome sequencing of L. japonicus were downloaded from NCBI SRA database (SRR6040926 and SRR6040927) for assembly. Data sets SRR6040926 and SRR6040927 were sequenced on the Illumina HiSeq2500 platform and consisted of paired end reads with a read length of 150bp. Mitochondrial genome sequence assembly based on short pair-end reads was performed using NOVOPlasty (version 2.7.2) (Dierckxsens et al. 2016) by taking L. japonicus mitochondrial genome (downloaded from NCBI, accession number NC_042503) (Lavoué et al. 2014) as a reference sequence. Subsequently, the mitochondrial reference genome of L. maculatus was annotated by DOGMA ( (Wyman et al. 2004). Finally, the newly assembly genome was annotated by MITOS (version 999) (Bernt et al. 2013), a web server for the annotation of metazoan mitochondrial genomes.

Variants calling

The clean reads of 85 individuals were aligned to the reference mitochondrial genome of L. maculatus using BWA (version 0.7.17-r1188) (Li and Durbin 2009) with default parameters. Then GATK (version (Mckenna et al. 2010) was applied for SNP and indel discovery and genotyping across all samples simultaneously using standard hard filtering parameters according to GATK Best Practices recommendations (Der Auwera et al. 2013). Further filtering was performed using vcftools (version 0.1.15) (Danecek et al. 2011) to remove all SNPs with minor allele counts less than 2, genotype missing counts greater than 2, or multiple alleles.

Genetic diversity, population structure and selection analysis

The number of haplotypes (h), haplotype diversity (Hd), number of variable sites (S) and nucleotide diversity (π) were calculated using DnaSP 5.0 (Librado and Rozas 2009). Fu’s Fs statistics (Fu 1997) and Tajima’s D (Tajima 1989) with 1000 permutation were used to test for neutrality by software Arlequin 3.0 (Excoffier et al. 2005) in historical demographic evaluation. The changes in effective population size were investigated by Bayesian skyline plot (BSP) method (Drummond et al. 2005) in BEAST 1.10 (Suchard et al. 2018) to reconstruct the demographic dynamics. We selected the HKY nucleotide substitution model, a strict clock (Soares et al. 2012), and a substitution rate of 2×10-9 (Burridge et al. 2008). MCMC analyses were run for 2 × 108 generations when necessary to achieve adequate effective sample size (i.e., ≥ 200), sampling every 1,000 generations and discarding the first 10% as burn-in. BEAST outputs were visualized with the Tracer v.1.7.1 ( To explore molecular variance, the analysis of molecular variance components (AMOVA) and fixation index (FST) were carried out within and among different populations. Isolation by distance was assessed by correlating pairwise genetic distance (FST) between populations with geographical distance (calculated from longitude and latitude) using the Mantel tests (100,000 permutations) implemented by the software Arlequin 3.0 (Excoffier et al. 2005). Pairwise Kimura 2-parameter distances were calculated using Mega 5.0 (Tamura et al. 2011). Except for ND6 gene, which has a different codon usage biases (Hasegawa et al. 1998), the Ka/Ks values of the other 12 protein-coding genes in the mitochondria were evaluated independently utilizing DnaSP 5.0 (Librado and Rozas 2009). From each haplotype block generated by based “Confidence Intervals” method (Gabriel et al. 2002), one tag SNP were selected using Haploview (version 4.2) (Barrett et al. 2005). Then software FRAPPE version 1.1 (Tang et al. 2005) was used to estimate the genetic ancestry of each sample based on the genotypes of tag SNPs. The maximum iteration of EM (expectation maximization) was set as 10,000; convergence threshold was set as 10,000, and the number of populations (K) was set as two and three for each calculation.

Phylogenetic trees were constructed using Bayesian inference methods to explore the phylogenetic relationship among different individuals from 12 geographical locations. Lateolabrax japonicus (GenBank: AP006789.1) was used as an outgroup. Phylogenetic analysis was carried out with MrBayes3.1.2 (Huelsenbeck and Ronquist 2001). The program Modeltest version 3.7 (Posada and Crandall 1998) was used to choose an appropriate substitution model of sequence evolution. The (TrN+I) model was selected as the best fit model. Four independent Markov chains were run for 10,000,000 generations by sampling one tree per 1000 generations and allowing adequate time for convergence. The first 2500 trees (25%) were discarded as part of a burn-in procedure (that was determined by the likelihood of being stationary), the remaining 7500 sampling trees were used to construct a 50% majority rule consensus tree. Two independent runs were used to provide additional confirmation of the convergence of the Bayesian posterior probabilities (BPP) distribution. A median-joining network was constructed using Network 4.51 (Bandelt et al. 1999) to infer relationships among 78 haplotypes of L. maculatus.


Mitochondrial genome assembling and variants calling

The total length of the mitochondrial genome of L. maculatus was 16,601 bp, which is close to L. japonicus (16,593 bp). The mtDNA comprised 13 protein-coding genes (PCGs) with an average length of 878.46 bp as well as two rRNA genes and 22 tRNA genes (Fig. 2). About 1032 G clean bases were obtained by next generation sequencing for 85 individuals from 12 geographic locations along Chinese coast. After aligning all clean reads to the newly assembled mitochondrial genome of L. maculatus, the average mapping ratio was 0.2% while the mapping depth was approximately 1422 (Supplementary Table S1). We obtained 249 high-confidence single nucleotide polymorphism (SNP) sites and 24 indels (Fig. 2). At last, combining all SNPs, indels and the reference genome of L. maculatus, the complete mitochondrial sequences of 85 individuals were generated. The genotype data was given in Supplementary Table S2. The mtDNA sequences of 85 individuals are reliable and used for further analysis.

Figure 2 The structure of Lateolabrax maculatus mitochondrial genome. The total length of mitochondrial genome of L. maculatus was 16,601 bp comprising 13 protein-coding genes (PCGs), 2 rRNA genes and 22 tRNA genes. 249 high-confidence single nucleotide polymorphism (SNP) sites and 24 indels was identified in 85 individuals. 

Genetic diversity of L. maculatus

A total of 85 whole mitochondrial sequences were obtained in 12 populations (FC, TS, HK, ZJ, SZ, ST, WZ, ZS, LY, YT, TJ, WD), including 78 haplotypes (Table 1). Among them, 72 haplotypes were ‘specific haplotypes’ while the remaining six were shared by more than two individuals. The haplotype diversity (Hd) and nucleotide diversity (π) of all samples were 0.998 and 0.00169, respectively. For different populations, Hd ranged from 0.978 (HK, TJ) to 1.000 (TS, HK, ZJ, SZ, ST, WZ, ZS, LY, YT, WD), and π ranged from 0.00036 (ST) to 0.00252 (WD, LY). Haplotype diversity (Hd) and nucleotide diversity (π) for all populations are shown in Table 1. LY population had the highest haplotype diversity and nucleotide diversity (Hd = 1.000, π = 0.00252) (Table 1).

Population genetic structure and historic demography of L. maculatus

To evaluate the degree of differentiation among populations, pairwise FST values and genetic distances were calculated based on mitochondrial sequences. According to SST in the sampling places, we can divide populations into two groups, South (FC, TS, HK, ZJ, SZ, ST) and North (WZ, ZS, LY, YT, TJ, WD) (Fig. 1). The pair wise FST values ranged from 0.001 (SZ-ZS) to 0.429 (TS-YT). Moreover, most FST values between northern and southern populations were larger than 0.05, which indicated significant genetic differentiation between them (Table 2). The pairwise Kimura 2-parameter distances ranged from 0.019 (HK-TS) to 0.038 (YT-FC, YT-HK, YT-ST, YT-TS) (Table 2). In our study, the test of isolation by distance (IBD) showed a significant correlation (R² = 0.6909, p = 1.196E-17) between genetic and geographical distances among populations (Fig. 3). We used hierarchical AMOVA grouping the populations according to their hypothesis (north and south). All populations were divided into two groups, north (WZ, ZS, LY, YT, TJ, WD) and south (FC, TS, HK, ZJ, SZ, ST), according to SST in the sampling places. AMOVA analysis indicates that 24.66% of the genetic variation is between northern and southern groups (p < 0.05). The genetic variance within populations is 72.39% and among populations is 2.95% (Table 3).

Figure 3 Plot of pairwise estimates of genetic (FST) and geographical distance between populations of Lateolabrax maculatus

Table 2 Pairwise Kimura 2-parameter distance (above diagonal) and fixation index (FST) (below diagonal) between localities. 

FC - 0.011 0.011 0.012 0.025 0.014 0.019 0.021 0.032 0.038 0.035 0.037
TS -0.024 - 0.010 0.011 0.025 0.013 0.019 0.020 0.031 0.038 0.035 0.037
HK -0.024 -0.027 - 0.012 0.025 0.013 0.019 0.020 0.031 0.038 0.035 0.037
ZJ -0.027 0.004 0.035 - 0.025 0.014 0.018 0.020 0.031 0.037 0.034 0.035
SZ 0.282* 0.331* 0.326* 0.225 - 0.026 0.025 0.026 0.033 0.030 0.031 0.033
ST 0.021 0.031 0.042 -0.016 0.188 - 0.021 0.021 0.155 0.038 0.035 0.037
WZ 0.061 0.120* 0.124 -0.030 -0.044 -0.005 - 0.024 0.032 0.035 0.033 0.036
ZS 0.113* 0.145* 0.154* 0.025 0.001 0.018 -0.095 - 0.031 0.034 0.033 0.035
LY 0.271* 0.312* 0.311* 0.174 0.009 0.009* -0.003 -0.022 - 0.037 0.036 0.036
YT 0.397* 0.429* 0.427* 0.325* -0.003 0.309* 0.122* 0.112* 0.032 - 0.035 0.037
TJ 0.299* 0.331* 0.332* 0.210* -0.042 0.196 0.026 0.032 -0.037 -0.022 - 0.037
WD 0.285* 0.315 0.314* 0.194* -0.019 0.181 0.034 0.036 -0.078 0.005 -0.053 -

*Significant at level (p < 0.05).

Table 3 Analysis of molecular variance (Amova) analysis in 12 populations of Lateolabrax maculatus from north and south coast of China. 

Source of variation d.f. Sum of squares Variance components Percentage of variation F-statistics
Among groups 1 208.762 4.55 24.66 FSC= 0.039*
Among populations within groups 10 170.593 0.544 2.95 FST= 0.276*
Within populations 73 962.05 13.361 72.39 FCT= 0.246*
Total 84 1341.405 18.457 100

*Significant at level (p < 0.05).

A clustering algorithm implemented in Structure (Pritchard et al. 2000) was used to infer the relationships among populations. To reduce the redundant information of highly linked SNPs, we identified 119 haplotype blocks from all SNPs and selected one tag SNP from each block. The number of populations, K, was set to two and three to reveal any hierarchical population structure. The most likely population structure was observed for K = 3, and discriminated the two dominating regions, i.e. southern and northern populations (Fig. 4). The phylogenetic analysis (Bayesian tree) based on whole-mitochondrial sequences was in accordance with their geographic distribution, of which WD, TJ, YT, LY, ZS and WZ clustered into group “North”; FC, TS, HK ZJ, and ST clustered into group “South” (Fig. 5). It is noteworthy that population SZ is distributed in the southern and northern regions. To understand the phylogenetic and geographical relationships among all identified haplotypes, a network, based on the median-joining method, was constructed (Fig. 5). Consequently, a slightly star-like topology with a great number of unique haplotypes (91.8%) was obtained. Haplogroups South and North included most haplotypes from southern and northern populations, respectively, except for SZ, ZS and WZ populations, which had individuals in both haplogroups (Fig. 6).

Figure 4 Admixture analysis among all populations derived from 85 whole-mitochondrial sequences. The K value was set 2 and 3. 

Figure 5 Bayesian tree constructed based on 86 whole-mitochondrial sequences of Lateolabrax maculatus. Each line represents one individual in the population. The reseda area and orange area represent north population and south population, respectively. Lateolabrax japonicus was used as outgroup. 

Figure 6 The median-joining network constructed based on 78 haplotypes of Lateolabrax maculatus. Each cycle represents a haplotype, the area of the circle is proportional to the frequency of haplotype. Different geographical populations were shown in different colors. 

The neutrality tests of Tajima’s D, Fu’s Fs were carried out to detect the historic demography of L. maculatus. For all samples, the neutrality test obtained negative Tajima’s D value (-1.451). Almost Tajima’s D values for 12 populations were negative, and the values for TS and YT were significant (Table 1). In Fu’s Fs tests, a very large negative Fu’s Fs index (-26.785) was presented when all samples were pooled together. These values signified an excess of low frequency polymorphisms, indicating population size expansion and/or negative selection acting in the population. Meanwhile, the BSP analysis showed the species underwent population expansions (Fig. 7).

Figure 7 The changed trend of effective population numbers with the time based on Bayesian skyline plot method. X-axis is the timescale before present, and Y-axis is the estimated effective population size. Solid curves indicate median effective population size; the shaded range indicates 95% highest posterior density intervals. 

Mitochondrial genomic selective pressure analysis

To investigate the evolutionary rate differences between southern and northern clusters, the number of nonsynonymous substitutions per nonsynonymous site (Ka) and the number of synonymous substitutions per synonymous site (Ks) were calculated for 12 mitochondrial PCGs. The ratio of Ka/Ks is normally deemed to be an indicator of selective pressure and evolutionary relationship among species at the molecular level (Li et al. 2015, Shen et al. 2007). According to previous studies (Yang and Nielsen 2000), the Ka/Ks value > 1, Ka/Ks = 1, and Ka/Ks < 1 suggest positive (adaptive) selection, neutral evolution and negative (purifying) selection, respectively. In our study, the Ka/Ks values of 12 PCGs varied approximately from 0.0177 (ATP6) to 0.1065 (ATP8), which indicated negative (purifying) selection (Fig. 8). Our result of the Ka/Ks ratio implied that the abundant genes evolved under negative selection, which means natural selection against deleterious mutations with negative selective factors (Yang and Bielawski 2000). The Ka/Ks values were higher in ATP8 and ND5 genes among all populations, while the values were lower in ATP6, ND3, Cytb, COX3, COX2 and COX1 genes. From this, ATP8 and ND5 genes may be under the less selective pressure while ATP6, ND3, Cytb, COX3, COX2 and COX1 genes might have undergone adaptive evolution to cope with different environment. But there were no significant genetic differences in certain mitochondrial genes between the southern and northern populations.

Figure 8 The intensity of purifying selection of 12 mitochondrial genes of Lateolabrax maculatus. The different colors represent different geographical populations. 


Genetic variance in mtDNA between south and north populations of L. maculatus

The genetic diversity of a population is correlated with the survival status of the species. A reduction or loss of genetic diversity is a great threat to wild populations living in a changeable environment. Haplotype diversity (Hd) and nucleotide diversity (π) are two important parameters to reflect the genetic diversity of a population (Liu et al. 2010). Both northern and southern populations of L. maculatus have high haplotype diversity. However, the nucleotide diversity of the southern population is lower than of the northern population. The FST value is an effective index when assessing gene flow and genetic diversity between populations. The normal value, FST of 0.05, was defined as inappreciable genetic differentiation, while greater than 0.25 means high genetic differentiation within a population (Holsinger and Weir 2009). Based on the pairwise FST values of populations from the same major geographic location, there is no significant genetic variation (South: 0.004-0.122; North: 0.003-0.078) except SZ (Table 2). The FST value of SZ was closer to the northern populations than to the southern populations. This probably demonstrated that gene flow was prominent among populations along the coastline without physical barriers. Considering the aquaculture of L. maculatus in Shenzhen using specimens introduced from the northern regions, the SZ population may have spread from the northern population.

By contrast, comparing the southern and northern populations of L. maculatus from different habitats revealed significant differences (0.003-0.429). The pairwise genetic distances between populations ranged from 0.010 to 0.038, which revealed a significant genetic divergence between geographically distant populations. Similar patterns were observed in the Fst value. In our research, a significant relationship between pairwise genetic distance (FST) and geographic distance (p < 0.05) was identified (Fig. 3). The greater the geographic distance between two populations, the higher the Fst value. When we used hierarchical AMOVA analysis between northern and southern groups, there was a 24.66% genetic variation. Adaptation to local conditions leads to a reduction in the within-population diversity as well as an increase in the among-population diversity (Storz 2010). Admixture analysis (Fig. 4), Phylogenetic trees (Fig. 5) and Median-joining network (Fig. 6) presented consistent results that individuals from different populations clustered with their geographic clade. The analysis indicated high level of genetic divergence between southern and northern geographic populations of L. maculatus, which was mainly in agreement with previous population genetic studies of L. maculatus based on markers from nuclear (Shao et al. 2009, Xu et al. 2017).

Population expansion of L. maculatus in history

The Tajima D-test (Tajima 1989) and FS test of Fu (Fu 1997) are used to examine whether the hypothesis of neutrality holds. Significant negative D and FS statistics can be interpreted as symbols of population expansion. Tajima’s D statistics is more sensitive than Fu’s Fs tests on population expansion (Fu 1997). A negative and significant Tajima’s D value implies that the sequence contains more nucleotide changes than the neutral evolution model, which may suggest a population expansion event in history and explain the wide distribution of several shared haplotypes (Zhang et al. 2017). Tajima’s D values for 12 populations of L. maculatus were mostly negative (Table 1). When all samples were pooled together, a highly negative Fs index (-26.785) and a significant negative D value (-1.451) were obtained, which indicated that a demographic expansion event in the past. The BSP analysis (Fig. 7) also supported the hypothesis of population expansion (Slatkin and Hudson 1991).

Negative selection in mitochondrial protein-coding genes

The effects of IBD in the genetic differentiation of populations was originally illustrated in 1948 (Wright 1943). Since then, IBD has been observed and reported in many population genetic studies, for example, the population genetic studies Nile Tilapia throughout west Africa and the coral-specialist yellowbar angelfish from the Northwestern Indian Ocean, the genetic distances between populations were clearly correlated to their geographic distance (Lind et al. 2019, Torquato et al. 2019). In this study, the IBD test performed by Mantel tests showed a significant correlation between genetic and geographical distances among L. maculatus populations. A high correlation value (R2 = 0.6909) suggested that the gene flow in L. maculatus could be the impact of geographical distances, which limited the dispersion of wild populations. Apart from this, the biological characteristics (e.g., dispersal ability and adaptive capacity), ecological factors (e.g., complicated topography and climatic environment) and anthropogenic factors (e.g., aquaculture and agriculture cultivation) might account for the genetic isolation (Liu et al. 2019). For L. maculatus, one of the apparent selective factors is sea surface temperature (SST), which can influence enzyme kinetics, physiology and ecological fitness (Powers et al. 1991). In winter, the SST is 10 °C in the Bohai Sea Gulf and 25 °C in the Southern Chinese Sea demonstrated a steep temperature gradient from north to south along the coast of China (Zhao et al. 2018). The lower water temperature in the Bohai Sea and Yellow Sea might limit the dispersal capabilities of L. maculatus (Wang et al. 2017). SST may be the principal environmental factor to adapt for specific population of L. maculatus (Strohm et al. 2015, Sun et al. 2011). Moreover, potential hydrological barriers could also result in a geographical segregation and lead to genetic divergence (Ma et al. 2010).

The rejection of the historical and equilibrium neutral models implicates natural selection. The southern population may display significant genetic differences in certain protein coding genes compared with the northern population, which is related to their respective adaptations to specific environments. To accurately detect the genomic footprints left by natural selection, we performed selection pressure analysis across the mitochondrial genome. The Ka/Ks ratios of 12 protein-coding genes were calculated but failed to detect significant evidence of selection between northern and southern populations. This could be a function of the fact that these mitochondrial genes are particularly conserved (James Bruce et al. 2008).

For comparison, the Cytb gene of wild boars was significantly different between the Siberian and the Vietnamese clades (Chen et al. 2018). Since Siberian wild boars can survive the cold subarctic climate, it is plausible to think that genes related to energy production might contribute to their local adaptation. Mitochondria are the powerhouses of cells, playing essential roles in both thermogenesis and thermoregulation (Gambert and Ricquier 2007, Yoshihiko et al. 2008). However, the molecular mechanism of temperature adaptation may be not necessarily in the mitochondrion. For example, a comparison of the whole genomic sequences of Leuciscus waleckii inhabiting extremely alkaline water and freshwater revealed significant genetic differences in some genes and genome regions that have undergone positive selection in the alkaline water population (Xu et al. 2017).

Correlations between genetic diversity, differentiation and environmental factors were also of interest due to the strong environmental gradients found on the coast of China, as well as the potential adaptive value of mitochondrial variability. However, there is very little published information on the adaptive roles of individual mitochondrial genes in fish. Further genetic and evolutionary studies utilizing whole genome approaches are necessary for better understanding the local adaptation of L. maculatus in the coast of China.

The whole mitogenomes of L. maculatus from different latitudes have been sequenced and used for population genetic analysis. Apparent genetic differences exist between the northern and southern populations. Population expansion in L. maculatus was estimated to infer their phylogenetic status, selection pressures and historical demography. No significantly different selection pressures were identified in the mitochondrial PCGs between the northern and the southern population. Studies on the population genetics of L. maculatus is significant for the effective management of fishery genetic resources and insights into population structure. Our results afford a crucial step for future analyses of the genetics of adaptation and conservation efforts in the L. maculatus species.


We acknowledge financial support from the Knowledge Innovation Program of Shenzhen City (Fundamental Research, Free Exploration, #JCYJ20170818142601870), Guangdong Basic and Applied Basic Research Fund (Guangdong Natural Science Fund, #2019A1515010731), and Fundamental Research Funds for the Central Universities, Xiamen University (#20720190108, 20720190099 and 20720190102).


An HS, Lee JW, Kim HY, Kim JB, Chang DS, Park JY, Myeong JI, An CM (2014) Population genetic structure of the sea bass (Lateolabrax japonicus) in Korea based on multiplex PCR assays with 12 polymorphic microsatellite markers. Genes & Genomics 36: 247-259. ]

Bandelt H, Forster P, Rohl A (1999) Median-joining networks for inferring intraspecific phylogenies. Molecular Biology and Evolution 16: 37-48. ]

Barhoumi M, Khoufi W, Kalai S, Ouerhani A, Essayed S, Zaier G, Jaziri H, Meriem SB, Fehribedoui R (2017) The use of Fourier analysis as a tool for Oblada melanura (Linnaeus, 1758) stock unit separation in the south central Mediterranean Sea. Journal of the Marine Biological Association of the United Kingdom 98: 1725-1732. ]

Barrett JC, Fry B, Maller J, Daly MJ (2005) Haploview: analysis and visualization of LD and haplotype maps. Bioinformatics 21: 263-265. ]

Bernt M, Donath A, Juhling F, Externbrink F, Florentz C, Fritzsch G, Putz J, Middendorf M, Stadler PF (2013) MITOS: improved de novo metazoan mitochondrial genome annotation. Molecular Phylogenetics and Evolution 69: 313-319. ]

Burridge CP, Craw D, Fletcher D, Waters JM (2008) Geological Dates and Molecular Rates: Fish DNA Sheds Light on Time Dependency. Molecular Biology and Evolution 25: 624-633. ]

Chen B, Li Y, Peng W, Zhou Z, Shi Y, Pu F, Luo X, Chen L, Xu P (2019) Chromosome-Level Assembly of the Chinese Seabass (Lateolabrax maculatus) Genome. Frontiers in Genetics 10: 275. ]

Chen J, Pan N, Thi TNT, Kamaldinov EV, Zhao S (2018) Selective constraints in cold-region wild boars may defuse the effects of small effective population size on molecular evolution of mitogenomes. Ecology & Evolution 8: 8102-8114. ]

Cox MP, Peterson DA, Biggs PJ (2010) SolexaQA: At-a-glance quality assessment of Illumina second-generation sequencing data. BMC Bioinformatics 11: 485-485. ]

Dalziel AC, Moyes CD, Fredriksson E, Lougheed SC (2006) Molecular evolution of cytochrome c oxidase in high-performance fish (teleostei: Scombroidei). Journal of Molecular Evolution 62: 319-331. ]

Danecek P, Auton A, Abecasis G, Albers CA, Banks E, DePristo MA, Handsaker RE, Lunter G, Marth GT, Sherry ST (2011) The variant call format and VCFtools. Bioinformatics 27: 2156-2158. ]

Der Auwera GAV, Carneiro MO, Hartl C, Poplin R, Angel GD, Levymoonshine A, Jordan T, Shakir K, Roazen D, Thibault J (2013) From FastQ data to high confidence variant calls: the Genome Analysis Toolkit best practices pipeline. Current protocols in human genetics 43: 1-33. ]

Dierckxsens N, Mardulyn P, Smits G (2016) NOVOPlasty: de novo assembly of organelle genomes from whole genome data. Nucleic Acids Research 45: e18-e18. ]

Dong C, Xu J, Wang B, Feng J, Jeney Z, Sun X, Xu P (2015) Phylogeny and Evolution of Multiple Common Carp (Cyprinus carpio L.) Populations Clarified by Phylogenetic Analysis Based on Complete Mitochondrial Genomes. Marine Biotechnology 17: 565-575. ]

Drummond AJ, Rambaut A, Shapiro B, Pybus OG (2005) Bayesian Coalescent Inference of Past Population Dynamics from Molecular Sequences. Molecular Biology and Evolution 22: 1185-1192. ]

Excoffier L, Laval G, Schneider SW (2005) Arlequin (version 3.0): An integrated software package for population genetics data analysis. Evolutionary Bioinformatics 1: 47-50. ]

Fu YX (1997) Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection. Genetics 147: 915-925. ]

Gabriel S, Schaffner SF, Nguyen HA, Moore J, Roy J, Blumenstiel B, Higgins JM, Defelice M, Lochner AL, Faggart M (2002) The structure of haplotype blocks in the human genome. Science 296: 2225-2229. ]

Gambert S, Ricquier D (2007) Mitochondrial thermogenesis and obesity. Current Opinion in Clinical Nutrition and Metabolic Care 10: 664-670. ]

Garvin MR, Bielawski JP, Gharrett AJ (2011) Positive Darwinian selection in the piston that powers proton pumps in complex I of the mitochondria of Pacific salmon. Plos One 6: e24127. ]

Hélène C, Luca F (2011) Evidence for morphological and adaptive genetic divergence between lake and stream habitats in European minnows (Phoxinus phoxinus, Cyprinidae). Molecular Ecology 20: 4490-4502. ]

Hansen MM (2010) Estimating the long-term effects of stocking domesticated trout into wild brown trout (Salmo trutta) populations: an approach using microsatellite DNA analysis of historical and contemporary samples. Molecular Ecology 11: 1003-1015. ]

Hasegawa M, Cao Y, Yang Z (1998) Preponderance of slightly deleterious polymorphism in mitochondrial DNA: nonsynonymous/synonymous rate ratio is much higher within species than between species. Molecular Biology & Evolution 15: 1499-1505. ]

Holsinger KE, Weir BS (2009) Genetics in geographically structured populations: defining, estimating and interpreting FST. Nature Reviews Genetics 10: 639-650. ]

Huelsenbeck JP, Ronquist F (2001) MRBAYES: Bayesian inference of phylogenetic trees. Bioinformatics 17: 754-755. ]

James Bruce S, Christoph F, Elson JL, Nils G, Ran L (2008) Purifying selection of mtDNA and its implications for understanding evolution and mitochondrial disease. Nature Reviews Genetics 9: 657-662. ]

Kawecki TJ, Ebert D (2010) Conceptual issues in local adaptation. Ecology Letters 7: 1225-1241. ]

Lavoué S, Nakayama K, Jerry DR, Yamanoue Y, Yagishita N, Suzuki N, Nishida M, Miya M (2014) Mitogenomic phylogeny of the Percichthyidae and Centrarchiformes (Percomorphaceae): comparison with recent nuclear gene-based studies and simultaneous analysis. Gene 549: 46-57. ]

Li H, Durbin R (2009) Fast and accurate short read alignment with Burrows-Wheeler transform. Bioinformatics 25: 1754-1760. ]

Li M, Zhao M, Zhong A, Xue L (2005) The analysis of genetic variation of Lateolabrax japonicus from Rizhao of Shandong and Xiamen of Fujian by isozyme and RAPD methods. Acta Oceanologica Sinica 27: 119-123. [ Links ]

Li X, Huang Y, Lei F (2015) Comparative mitochondrial genomics and phylogenetic relationships of the Crossoptilon species (Phasianidae, Galliformes). BMC Genomics 16: 42. ]

Librado P, Rozas J (2009) DnaSP v5: a software for comprehensive analysis of DNA polymorphism data. Bioinformatics 25: 1451-1452. ]

Lind CE, Agyakwah SK, Attipoe FY, Nugent C, Crooijmans RPMA, Toguyeni A (2019) Genetic diversity of Nile tilapia (Oreochromis niloticus) throughout West Africa. Scientific Reports 9: 1-12. ]

Liu Q, Liu L, Song N, Wang X, Gao T (2019) Genetic structure in the marbled rockfish (Sebastiscus marmoratus) across most of the distribution in the northwestern Pacific. Journal of Applied Ichthyology 35: 1249-1259. ]

Liu S, Wang K, Yuan X, Wang D, Yue X (2010) Genetic diversity and differentiation between populations of Glyptothorax zanaensis in the middle and lower reaches of the Nujiang River. Hereditas 32: 254-263. ]

Ma C, Cheng Q, Zhang Q, Zhuang P, Zhao Y (2010) Genetic variation of Coilia ectenes (Clupeiformes: Engraulidae) revealed by the complete cytochrome b sequences of mitochondrial DNA. Journal of Experimental Marine Biology & Ecology 385: 14-19. ]

Mckenna A, Hanna MG, Banks E, Sivachenko A, Cibulskis K, Kernytsky AM, Garimella K, Altshuler DS, Gabriel S, Daly M (2010) The genome analysis toolkit: a mapreduce framework for analyzing next-generation DNA sequencing data. Genome Research 20: 1297-1303. ]

Niu S, Liu Y, Qin C, Wang X, Wu R (2015) The complete mitochondrial genome and phylogenetic analysis of Lateolabrax maculatus (Perciformes, Moronidae). Mitochondrial Dna 28: 173-175. ]

Posada D, Crandall KA (1998) MODELTEST: testing the model of DNA substitution. Bioinformatics 14: 817-818. ]

Powers DA, Lauerman T, Crawford D, Dimichele L (1991) Genetic mechanisms for adapting to a changing environment. Annual Review of Genetics 25: 629-659. ]

Pritchard JK, Stephens M, Donnelly P (2000) Inference of population structure using multilocus genotype data. Genetics 155: 945-959. [ Links ]

Shao C, Chen S, Xu G, Liao X, Tian Y (2009) Eighteen novel microsatellite markers for the Chinese sea perch, Lateolabrax maculatus. Conservation Genetics 10: 623-625. ]

Shen X, Ren J, Cui Z, Sha Z, Wang B, Xiang J, Liu B (2007) The complete mitochondrial genomes of two common shrimps (Litopenaeus vannamei and Fenneropenaeus chinensis) and their phylogenomic considerations. Gene 403: 98-109. ]

Slatkin M, Hudson R (1991) Pairwise comparisons of mitochondrial DNA sequences in stable and exponentially growing populations. Genetics 129: 555-562. ]

Soares P, Alshamali F, Pereira JB, Fernandes V, Silva N, Afonso C, Costa MD, Musilova E, Macaulay V, Richards MB (2012) The Expansion of mtDNA Haplogroup L3 within and out of Africa. Molecular Biology and Evolution 29: 915-927. ]

Stoljarova M, King JL, Takahashi M, Aaspõllu A, Budowle B (2016) Whole mitochondrial genome genetic diversity in an Estonian population sample. International Journal of Legal Medicine 130: 67-71. ]

Storz JF (2010) Using genome scans of DNA polymorphism to infer adaptive population divergence. Molecular Ecology 14: 671-688. ]

Strohm JHT, Gwiazdowski RA, Hanner R (2015) Fast fish face fewer mitochondrial mutations: Patterns of dN/dS across fish mitogenomes. Gene 572: 27-34. ]

Suchard MA, Lemey P, Baele G, Ayres DL, Drummond AJ, Rambaut A (2018) Bayesian phylogenetic and phylodynamic data integration using BEAST 1.10. Virus Evolution 4: vey016. ]

Sun Y, Shen Y, Irwin DM, Zhang Y (2011) Evaluating the roles of energetic functional constraints on teleost mitochondrial-encoded protein evolution. Molecular Biology and Evolution 28: 39-44. ]

Tajima F (1989) Statistical method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics 123: 585-595. [ Links ]

Tamura K, Peterson DS, Peterson N, Stecher G, Nei M, Kumar S (2011) MEGA5: Molecular evolutionary genetics analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods. Molecular Biology and Evolution 28: 2731-2739. ]

Tang H, Peng J, Wang P, Risch N (2005) Estimation of Individual Admixture: Analytical and Study Design Considerations. Genetic Epidemiology 28: 289-301. ]

Torquato F, Range P, Benhamadou R, Sigsgaard EE, Thomsen PF, Riera R, Berumen ML, Burt JA, Feary DA, Marshell A (2019) Consequences of marine barriers for genetic diversity of the coral-specialist yellowbar angelfish from the Northwestern Indian Ocean. Ecology and Evolution 9: 11215-11226. ]

Wang W, Ma C, Chen W, Zhang H, Kang W, Ni Y, Ma L (2017) Population genetic diversity of Chinese sea bass (Lateolabrax maculatus) from southeast coastal regions of China based on mitochondrial COI gene sequences. Biochemical Systematics & Ecology 71: 114-120. ]

Wright S (1943) Isolation by Distance. Genetics 28: 114-138. ]

Wyman SK, Jansen RK, Boore JL (2004) Automatic annotation of organellar genomes with DOGMA. Bioinformatics 20: 3252-3255. ]

Xu J, Li JT, Jiang Y, Peng W, Yao Z, Chen B, Jiang L, Feng J, Ji P, Liu G (2017) Genomic Basis of Adaptive Evolution: The Survival of Amur Ide (Leuciscus waleckii) in an Extremely Alkaline Environment. Molecular Biology & Evolution 34: 145. ]

Yang Z, Bielawski JP (2000) Statistical methods for detecting molecular adaptation. Trends in Ecology & Evolution 15: 496-503. ]

Yang Z, Nielsen R (2000) Estimating synonymous and nonsynonymous substitution rates under realistic evolutionary models. Molecular Biology and Evolution 17: 32-43. ]

Yoshihiko O, Yoshiaki K, Yukie A, Takanori I, Miyuki M, Yuka I, Megumi I, Kazushige M, Yusuke K, Hiroyuki K (2008) Functional coexpression of the mitochondrial alternative oxidase and uncoupling protein underlies thermoregulation in the thermogenic florets of skunk cabbage. Plant Physiology 146: 636-645. ]

Zhang R, Leng P, Wang X, Zhang Z (2017) Molecular analysis and genetic diversity of Aedes albopictus (Diptera, Culicidae) from China. Mitochondrial Dna A Dna Mapp Seq Anal 29: 1-6. ]

Zhao Y, Peng W, Guo H, Chen B, Zhou Z, Xu J, Zhang D, Xu P (2018) Population genomics reveals genetic divergence and adaptive differentiation of Chinese sea bass (Lateolabrax maculatus). Marine Biotechnology 20: 45-59. ]

Zhou Z, Han K, Wu Y, Bai H, Ke Q, Pu F, Wang Y, Xu P (2019) Genome-wide association study of growth and body-shape-related traits in large yellow croaker (Larimichthys crocea) using ddRAD sequencing. Marine Biotechnology 21: 655-670. ]

Publication Notes

Available online: August 25, 2020

Zoobank Register:

Publisher: © 2020 Sociedade Brasileira de Zoologia. Published by Pensoft Publishers at

Received: December 03, 2019; Accepted: March 31, 2020; Published: August 25, 2020

Corresponding author: Peng Xu (

Editorial responsibility:

Walter A.P. Boeger

Author Contributions:

JG and BC contributed equally to this work; PX conceived and supervised the study; JG and BC analyzed the data and drafted the manuscript; DZ and QK collected the samples; ZZ and YS construct the ddRAD library; PX and BL revised the manuscript. All authors have read and approved the manuscript.

Competing Interests:

The authors have declared that no competing interests exist. The authors alone are responsible for the content and writing of the article.

Creative Commons License This is an open-access article distributed under the terms of the Creative Commons Attribution License