Genetic structure of the snakehead murrel, Channa striata (channidae) based on the cytochrome c oxidase subunit I gene: Influence of historical and geomorphological factors

Nucleotide sequences of a partial cytochrome c oxidase subunit I gene were used to assess the manner in which historical processes and geomorphological effects may have influenced genetic structuring and phylogeographic patterns in Channa striata. Assaying was based on individuals from twelve populations in four river systems, which were separated into two regions, the eastern and western, of the biodiversely rich state of Perak in central Peninsular Malaysia. In 238 specimens, a total of 368-bp sequences with ten polymorphic sites and eleven unique haplotypes were detected. Data on all the twelve populations revealed incomplete divergence due to past historical coalescence and the short period of separation. Nevertheless, SAMOVA and FST revealed geographical structuring existed to a certain extent in both regions. For the eastern region, the data also showed that the upstream populations were genetically significantly different compared to the mid- and downstream ones. It is inferred that physical barriers and historical processes played a dominant role in structuring the genetic dispersal of the species. A further inference is that the Grik, Tanjung Rambutan and Sungkai are potential candidates for conservation and aquaculture programmes since they contained most of the total diversity in this area.


Introduction
Channa striata, locally known as haruan or snakehead murrel, an eminent tropical freshwater fish widely used for medicinal and pharmaceutical purposes (Mat Jais et al., 1994;Michelle et al., 2004), is also an important food source in the Asia-Pacific region (Froese and Pauly, 2008;Hossain et al., 2008). This carnivorous airbreather species is encountered in rivers, swamps, ponds, canals, drains, reservoirs, rice fields, small streams, mining pools, roadside ditches and lakes, across southern Asia, southern China, Indochina and the Sunda Islands (Mohsin and Ambak, 1983;Lee and Ng, 1994;Hossain et al., 2008). In Malaysia, and due to its abundance in nature, it is normally marketed alive, fresh from the catch. Aquaculture itself is only significant in certain neighbouring countries, such as Thailand, Pakistan, Taiwan, the Philippines, Vietnam, Cambodia and India (Wee, 1982;Hossain et al., 2008). The aquacultural potential of the species in Malaysia, has, as yet, not been fully exploited despite its many advantageous characteristics, notably high market price, air breathing ability, hardiness and high tolerance to adverse environmental conditions (Samantaray and Mohanty, 1997;Ali, 1999;Froese and Pauly, 2008). As is frequently observed in many important food-fish, overharvesting and other anthropogenic factors have resulted in severe damage to its natural habitat, with a subsequent decline in the indigenous stock of the species (Nagarajan et al., 2006;Hossain et al., 2008). Furthermore, unsystematic hatchery programmes can also lead to inbreeding depression, with a possible reduction in fecundity, adaptation ability and survival rate (Beaumont and Hoare, 2003;Sun et al., 2004). Hence, in order to effectively conserve and manage the fish, vital information on relevant population genetics is required, specifically through assessment of its genetic diversity and structuring for potential brood-stock identification. Previous studies have focused on reproductive biology (Ali, 1999), medical and pharmaceutical properties (Baie and Sheikh, 2000;Michelle et al., 2004), biochemical composition (Zuraini et al., 2006;Zakaria et al., 2007), ecology (Lee and Ng, 1994;Amilhat and Lorenzen, 2005), breeding (Haniffa et al., 2000), diet (Roshada, 1994;Arul, 2008) and morphological characters (Chandra and Banerjee, 2004), with only limited available information on population genetics of the species, especially in Malaysia. Ambak et al. (2006) and Mat Jais et al. (2009), when examining the genetic structure of snakeheads from Peninsular Malaysia based on RAPD and RFLP data, reported a positive corelation between hydro-geographic factors (mainly due to the Titiwangsa Mountain Range) and population genetic differentiation, respectively. Hara et al. (1998) also reported similar geographic structuring of C. striata in Thailand, based on allozymic data.
A typical example of the common aspect of natural landscapes of Peninsular Malaysia (Mohsin and Ambak, 1983), the state of Perak, situated in the central-western region, is overspread with pronounced geographical features, complex natural ecological heterogeneity and dense short river tributaries connecting and running into the Strait of Malacca, thereby sustaining rich biodiversity. Little is known regarding the influence of Pleistocene climate fluctuations and geomorphological processes involving drainage re-arrangements, ecological changes and natural or anthropogenic physical barriers, that have lead to the current biotic composition of Malaysia in general and Perak State specifically. The Chenderoh dam, and the Bintang and Keledang mountain ranges, as well as the Kerian and Perak rivers and their tributaries, are among the significant geographical features dominating the local landscape. These are probably the factors mainly involved in shaping the genetic pattern of regional biotic systems. Moreover, deglaciation during the late Pleistocene is believed to have periodically separated the two main rivers in the region, the Kerian and Perak, as well as their tributaries. According to the prevailing hypothesis, during the Pleistocene these two rivers coalesced, thus constituting a single course that ran north to the Andaman Sea (Voris, 2000). Therefore, Pleistocene climate fluctuations possibly played a significant role in genetic distribution, through changes in the physical land surface, habitats and the natural range of many regional species.
In this study, the mitochondrial cytochrome c oxidase subunit I (COI) gene was partially sequenced to examine the genetic diversity and structure of C. striata populations within Perak state in the central Peninsular Malaysia and particularly to infer the mechanisms or forces most likely to have been involved in shaping these populations, thus providing critical genetic information for brood-stock management and species conservation.

Sample acquisition
A total of 238 individuals of C. striata, representing twelve populations throughout Perak state, were sampled from the wild during the period 2007 to 2009. These populations were classified into four categories, according to their placing along the respective river system. The Grik, Tanjung Rambutan and Sungkai populations, through being the closest to their origins, were then considered as up-stream populations, whereas, the Tapah, Kubu Gajah, Lenggong and Ulu Kinta were clustered as upper-mid stream, the Kuala Kangsar and Gopeng as mid-stream, and finally, the Parit Buntar, Kampung Gajah and Teluk Intan as downstream populations. These populations were further divided into two major regions, the eastern and western, effectively separated by the Bintang mountain range. Only one single river system, the Kerian river, was investigated in the western region. Since the eastern covered a larger area, with several distinct tributaries, this was further divided into three riverine systems, namely the Perak, Kinta and Batang Padang -Sungkai. Two further significant geographical landmarks in this area are (1) the Keledang mountain range, which partially separates the Perak and Kinta rivers along the upper-mid stream, and (2) the Chenderoh dam, built across Perak river. Details of sampling locations, and region and sample sizes, are shown in Table 1 and Figure 1.

DNA analysis
An AquaGenomic DNA isolation Kit (BioSyntech, Salt Lake City, Utah, USA) was used for total DNA extraction from fins and muscles, according to manufacturer's in- Jamsari et al. 153  structions. DNA integrity and quantities were assessed on 0.8% (w/v) agarose gel, with a spectrophotometer (U-1900 UV/VIS spectrophotometer 200V -Hitachi, Tokyo, Japan). The DNA was then stored at -20°C until use. A segment of the COI mtDNA gene was PCR-amplified using the primers L6154 (5'-AYC ARC AYY TRT TYT GRT TCT-3') and H6556 (5'-TGR AAR TGI GCI ACW ACR TA-3') (Teletchea et al., 2006). Amplification was with an MJ PTC-200 thermal cycler (MJ Research, Waltham, MA, USA) at a total volume of 25 mL containing 1X PCR buffer, 4.2 mM MgCl 2 , 0.2 mM dNTPs, 0.6 mM of each primer, 0.08 U of Taq polymerase and 50 to 100 ng of template DNA. The reaction programme was carried out initially at 95°C for 5 min followed by 30 cycles with the following profile: 94°C for 60 s, 50°C for 60 s, 70°C for 60 s and finally 5 min of final extension at 72°C. DNA amplification products were separated in 1.5% (w/v) agarose gels at 100 V with 0.5X Tris-borate-EDTA (TBE) buffer, stained with ethidium bromide and visualized under UV illumination. PCR products were purified using a QIAquick PCR purification kit (Qiagen, Valencia, CA, USA), and sequenced using Big Dye Terminator v3.1 and an ABI3730XL Genetic Analyzer (Applied Biosystems, Foster City, CA, USA).

Data analysis
Initial editing of ambiguous bases was undertaken with MEGA 4.0 software (Tamura et al., 2007). The edited sequences were then aligned by using Clustal W 1.6 implemented in the same software. The alignments thus obtained were further visually cross-checked. Amino acid sequence translation (vertebrate mitochondrial code) was applied to evaluate the accuracy of COI sequences, and then translated back for subsequent analysis. Collapse 1.2 (Posada, 2006) was used to determine identical haplotypes in the aligned matrix. All haplotypes were deposited in GenBank (accession numbers GQ244413 to GQ244422 and GQ334376). Molecular diversity indices [number of haplotypes, polymorphic sites, transitions, transversions, and haplotype (h) and nucleotide (p) diversities], besides a population comparisons by pairwise F-statistics, were calculated by using Arlequin version 3.1 software (Schneider et al., 2000), in order to reveal the level of genetic variation and population structure. To correct for multiple comparisons, the sequential Bonferroni correction was applied. Spatial analysis of molecular variance (SAMOVA) was performed using SAMOVA v.1.0 (Dupanloup et al., 2002), indicated the amount of genetic variation, as well as geographically homogeneous population clusters. Finally, the CONTRIB version 1.02 software (Petit et al., 1998) was used to evaluate the contribution of each population to total diversity, as measured by allelic richness (CTR).

Results
Genetic variability, haplotype distribution and haplotype network relationships A total of 368 bp of unambiguous COI sequence alignments, with ten (2.7%) variable nucleotide positions and six (1.6%) parsimoniously informative sites, were obtained. Nine of the base substitutions were transitions, whereas only one was a transversion. All the variable sites occurred at the third codon position. Altogether, eleven unique putative haplotypes were identified from all the 238 individuals sampled. The mean base composition (%) was 24.7A: 20.2G: 29.9T: 25.2C. In the eastern region (Perak River, Batang Padang -Sungkai river and Kinta river systems), haplotype 01 was predominant and unique, whereas haplotype 07 was common but not exclusive to the western Kerian river system. Most of the haplotypes were shared by multiple populations, with only three haplotypes being population specific. Nine of the 11 haplotypes (except 04 and 07) were unique to either the eastern or western regions. Haplotypes 01, 07 and 04 were the most common, contributing with 59.2%, 21% and 9.7%, respectively, to total haplotype occurrence. 154 Genetic structure of Channa striata The genetic variability parameters, viz., nucleotide and haplotype diversities, in populations, river systems and the total population, are presented in Table 2. Nucleotide and haplotype diversities ranged from 0.0007 to 0.0067 and 0.1000 to 0.7879, respectively. As a whole, relatively low mean population nucleotide (0.0048) and haplotype (0.5964) diversities were observed. Although much higher genetic variation was observed in upstream populations, this was only moderate downstream ones, reaching the lowest in upper-middle and mid-stream populations. In general, genetic diversity was slightly higher in upper-mid stream than in mid-stream populations. In the upstream populations, diversity was high (h: 0.7126 to 0.7879 and p: 0.0052 to 0.0067), thereby contributing significantly to total genetic diversity (Figure 2), whereas in the remainder, this was low to moderate (h: 0.1000 to 0.3333 and p: 0.0013 to 0.0030). Low genetic variability (h: 0.2541 to 0.4914 and p: 0.0025 to 0.0035) was also observed when all the populations in the same river system were treated as a single entity.

Population differentiation (F ST ) and the spatial analysis of molecular variation (SAMOVA)
Overlapping haplotype distribution, haplotype sharing and low nucleotide mutation among all the populations, were apparent from the data. Nevertheless, substantial genetic differentiation was observed among the western and eastern lineages (F ST ranging from 0.4535 to 0.8526 and SAMOVA with F CT : 0.6210) (Tables 3 and 4). However, with the exception of the Grik, Tanjung Rambutan and Sungkai populations, no significant genetic differentiation among populations within river systems and regions themselves was detected (Table 3). In the Perak river, genetic differentiation between the upstream Grik population and the downstream population Kampung Gajah, was insignificant. Interestingly, no significant differentiation (F ST : 0.0000 to 0.0244) was observed when assessing pairwise comparison of the Grik, Sungkai and Tanjung Rambutan populations, although the three belonged to different river systems, i.e., they were genetically similar to each other but not to other members within their own river system and region. Furthermore, in the eastern region, where all the various populations, through occupying one and the same river system, were treated as a single group, low and insignificant mutual population differentiation (F ST ranging from 0 to 0.0122) was recorded. Limited statistical analysis could be applied to the western region, consisting of only one river system (the Kerian river) comprising two populations. The Parit Buntar, Kubu Gajah, Grik and Sungkai populations were the main contributors to population differentiation among all those studied (Figure 2). SAMOVA analysis indicated that 62.10% of the total variation existed between the two regions, 4.84% among populations within the region itself, and 33.06% within the populations, with significant support at all hierarchical levels. Furthermore, the SAMOVA tests also showed that a large proportion of the molecular variance was attributable to genetic differences among populations within the total (F ST : 0.6694) and to differences between the eastern and western regions (F CT : 0.6210). However, inspection of the data revealed that the significant (although low) differences among populations within regions (F SC : 0.1277) were contributed by the three genetically significantly different upstream populations, Tanjung Rambutan, Grik and Sungkai.

Genetic variation
Despite their abundance and wide distribution, as well as several biologically advantageous attributes, genetic variation was relatively low (mean h: 0.5964, p: 0.0048) in the C. striata populations investigated. The three upstream populations near the riverheads, upstream Grik, Tanjung Rambutan and Sungkai, revealed high genetic diversity, whereas in the three downstream populations this was moderate, and very low in the six uppermiddle and mid-stream ones. From this it can be inferred that upstream populations are colonizing centers or refugia. These refugia harbour high variability (with h: 0.7126 to 0.7879 and p: 0.0052 to 0.0067), thereby implying a large initial effective population size. Pleistocene glaciation periods gave origin to substantial habitat re-organization, subsequently leading to the displacement of populations into glacial refugia, thereby giving rise to high intraspecific diversity, both through secondary contact between differentiated assortments, as also with the local gene pool, as appraised by Grant and Bowen (1998) and Petit et al. (2003) in marine fishes and European trees shrubs, respectively. Nevertheless, no refugium was detected, neither in the western region nor in the Batang Padang river in the eastern, probably due to limited sampling effort, although it is very likely such a parallel pattern does exist. Interestingly, all putative refugia were located in the upstream reaches, thereby indicating altitudinal shifts in refugee movements during the interglacial period. Jamsari et al. 155 The lower genetic variability recorded in the upper-middle, mid-and downstream populations may reflect their re-colonisation via successive founder events, each with a relatively low Ne value, and by either a single or only few organelle-lineage colonisation (Grant and Bowen, 1998;Wang et al., 2000), in a passive downstream dispersal-migration route. This downstream migration pattern had already been documented in a study by Halls et al. (1998) of C. striata in Bangladesh. Alternatively, the low variability observed in these three populations could also be the result of a historical bottleneck event that may have almost eradicated these populations. Notwithstanding, genetic variation in three downstream populations, viz., the Parit Buntar, Kampung Gajah and Teluk Intan, was higher than in the upper-middle and mid ones. Presumably their downstream location was more propitious for the overlapping of colonisation routes, with concomitant admixture of haplotypes or alleles from several separate refugia, or other connected river systems. This had already been observed in other studies of the same scenario (Nguyen et al., 2006).

Phylogeographic relationships
Based on F ST and SAMOVA analyses, three genetically different groups were identified: (1) the western region, (2) an upstream population in the eastern region, and (3) a mid downstream population, also in the eastern region. These populations were, however, incompletely divergent, evident by overlapping haplotype distribution, haplotype sharing and low nucleotide mutation. The lack of genetic divergence within geographical regions (except in refugia), incompatible with the low migratory behavior of C. striata (Halls et al., 1998;Amilhat and Lorenzon, 2005), implied factors other than free gene flow, possible alternatives being either recent population expansion with insufficient time for coalescence (cf. Wang et al., 2000), or interconnection of the areas studied. Mid and lower stream topology may not have been a sufficient hindrance against gene flow, when compared to the situation upstream. Furthermore, this region is also liable to high flooding (Mohsin and Ambak, 1983;Liu and Chan, 2003), thereby facilitating lateral dispersal to flood plains, watersheds or between Jamsari et al. 157 (Hurwood and Hughes, 1998;Wang et al., 2000).
Two topographically isolated stocks of C. striata were revealed. Geographical separation by the Bintang mountain range, and changes in sea-water level during the Pleistocene, constitute two significant elements which have probably influenced the genetic distribution pattern of C. striata between the western and eastern regions of this highland. This mountain range, composed of several prominent peaks up to 1800 m high, only permitted genetic exchange between the systems at the ancient confluence, now submerged into the bed of the Strait of Malacca. Periodic changes in sea level witnessed the disappearance of this historical connection, with subsequently isolated populations forming distinct evolutionary units. Overlapping haplotype distribution and sharing, as well as low nucleotide mutation, provide ample evidence of this connection. It is apparent that, following construction of the Chenderoh dam across the Perak river in the 1930's, insufficient time has passed for any marked differentiation to have occurred among the isolated populations (cf. Hurwood and Hughes, 1998;Sun et al., 2004). Nevertheless, this event could potentially lead to substantial genetic restructuring over time, especially with the complete isolation of the genetically richer upstream sites from downstream populations. Similarly, the Keledang Mountain Range has not contributed to any obvious genetic divergence between the Perak and Kinta Rivers.
In theory, population segregation will normally contribute to genetic differentiation among isolated populations (cf. Hurwood and Hughes, 1998;Wang et al., 2000). Thus, it was unusual to discover that the three presumptive refugia in different rivers were genetically very closely related. Due to the small-sized study area, one possible factor could be the absence of the suitable ecological heterogeneity for stimulating overall local adaptation. Furthermore, all the populations were either of a common origin, or originated from a single colonization unit, seeing that during the glaciation period, each refugium maintained or accumulated all (or almost all) the haplotypes either from the original population or from other sites. Thirdly, the cytochrome oxidase subunit 1 (COI) was incapable of detecting population variability in this species, probably due to its relatively low mutation rate (Hebert et al., 2003;Hellberg, 2006), but also to the short time-span since population contraction into refugia during the last Pleistocene, insufficient for observable genetic differentiation to have occurred.
Overall, physical barriers to gene flow are thought to have played an important role in establishing evolutionary units between the two separate regions, whereas within each region, it appears that historical processes leading to environment disturbance and subsequent habitat re-organisation, have played a more predominant role in structuring the genetic dispersal of the species.

Conservation and management implications
The conservation of species is a source of great concern, especially regarding those endangered and/or of high economic value. This study provided useful and critical genetic information for planning management, conservation and ranching guidelines for C. striata in the biodiversityrich Perak State. Nevertheless, based on limited genetic information, all inferences should be reconfirmed using faster evolving markers, such as microsatellites. Based on F ST and SAMOVA analyses, two genetically identified stocks, the eastern and the western, were detected, both requiring separate monitoring and management (Moritz, 1994). Translocation between the two is not recommended, to thus avoid genetic disruption (Jørstad and Farestveit, 1999), adverse competition (Minckley, 1995), or the introduction of diseases (Mahidol et al., 2007). Furthermore, in the case of restoration involving local stock enhancement, special precautions must be taken, as this species reportedly exerts a negative ecological impact on aquatic communities, mainly due to its carnivorous behaviour (Cagauan, 2007;Froese and Pauly, 2008). With respect to the eastern region, if the upstream population as a source of gene pool for downstream populations could be verified, stock transfer from the former could be a beneficial option. Otherwise, they should be managed separately. Furthermore, as the more elevated upstream areas harbor a large proportion of the total genetic variation in this region, maximum priority for conservation should be given to this area. A strict regulation of harvest and fishery management should be implemented to protect genetically depauperate downstream populations. Habitat protection from agricultural and industrial activities, including the development of surrounding areas, as well as the construction of river-dams which prevent population connectivity, must be carefully regulated, seeing that the latter brings about habitat fragmentation, thereby potentially causing further losses in genetic diversity or an increase in inbreeding depression. As an alternative, the introduction and improvement of C. striata aquaculture is a fine option for ensuring the maintenance of this species. Nevertheless, hatchery facilities and operation, as well as brood-stock selection, must be systematically controlled, so as to prevent fugitive genetic contamination into the wild and other adverse genetic risks (Mahidol et al., 2007). Based on its high genetic variability and positive contribution towards genetic variation in this region, Grik, Sungkai and Tanjung Rambutan are the most appropriate candidate population of wild species to receive priority as a baseline stock for selective breeding.