Cytochrome b gene based population study of Tenualosa ilisha (Hamilton) in the Brahmaputra river system of India

Tenualosa ilisha, an anadromous fish, inhabited in the lower region of the estuaries and the foreshore areas migrating to river for spawning. The Hilsa population entered the Brahmaputra River system from Bay of Bengal. Various studies have shown the Hilsa population differs significantly from one another found in different rivers. In the present communication attempt has been made to investigate population of Hilsa shad in the Brahmaputra River using cytochrome b gene. Analysis of 425 bp of Cytochrome b gene fragment revealed 4 haplotypes with two clusters in dendogram. Existence of three polymorphic sites in the haplotypes sequences obtained from two groups of Hilsa shad showed low genetic differentiation in the Brahmaputra River.


Introduction
Hilsa shad, Tenualosa ilisha (Hamilton) is a popular fish found in India, Bangladesh, Pakistan, Myanmar, Indonesia, Sumatra, Kuwait, Iraq and Iran. The fish lives in the sea, but during breeding period they migrate to freshwater river. Hilsa fishing supports livelihood in Bangladesh and India, and to maintain sustainable Hilsa production, control of wanton killing of the breeding population is needed (Sunny et al., 2019). In India they migrate to the Brahmaputra and Ganga River system for spawning (Pillay and Rosa, 1963) and due to intense fishing and exploitation of brood fish during their upstream breeding migration its number is declining (Bhaumik, 2017). In Brahmaputra River Hilsa fishing is mainly done during the monsoon season, but during the later part of winter also, the mature fish come to river for spawning. The breeding season of Hilsa fish is mainly from May to August with a shorter season from January to February and sometimes it stretches into March. Different studies on population of Hilsa have been done by various authors encompassing morphometric and genetic aspects and found contradictory results. The Hilsa population was recognized as three different ecotypes by Majumdar (1939), three stocks by  and existence of different morphotypes by Ghosh et al, (1968). The molecular studies based on allozyme markers by Dhale et al. (1997) revealed genetic differences of Hilsa shad population between freshwater and brackish waters in Bangladesh, but similar studies by Lal et al. (2004) report the existence of a single panmictic population in India sampled from five rivers including Brahmaputra. Differences in Hilsa population in freshwater and brackish waters in Bangladesh were reported by Hussain et al. (1998) while working with RAPD marker. Salini et al. (2004) reported presence of genetically distinct Hilsa population in Bangladesh, Kuwait and Indonesia. The genetic study conducted by Rahman and Naevdal (2000) differentiated the marine group from freshwater and estuarine localities of the Hilsa population in Bangladesh. Similarly, Behera et al. (2015) found five haplotypes of in the Hilsa population in Bhagirathi and Tapti River by cytochrome b gene analysis. Verma et al. (2016) reported the presence of genetic differentiation among the Hilsa populations with high haplotype diversity in the Arabian Sea and Bay of Bengal from the sequencing of mitochondrial control region. In another study, SNP (single nucleotide polymorphism) analysis showed the presence of three genetically diverged local populations in different environmental habitats (Asaduzzaman et al., 2019;2020). The ecology and adaptations of Hilsa were studied through genomic data by Mohindra et al. (2019) to know the adaptation mechanisms in different salinity levels during migration. The present study of Hilsa population has aimed to find out any genetic variability in the Brahmaputra River system because the past studies unable to describe the population genetic structure. Therefore, partial sequencing of cytochrome b gene of mitochondrial genome was done to find out the presence of genetically distinct Hilsa population in Brahmaputra River as it has been regarded as an important tool of genetic variability studies and analyzing the population structure and phylogenetic studies in fish and other animal species (Irwin et al., 1991;Li et al., 2009;He et al., 2011;Behera et al., 2015).

Sample collection
A total of 80 samples of T. ilisha of varying size and weight were collected during the spawning period from different landing sites of the Brahmaputra River system in Dhubri and Goalpara district of Assam during the period of July 2017 to November 2017. From Goalpara district 27 (GH001-GH027) samples collected from four different locations, namely Madartary ghat, Krishnai, Joleswar ghat and fish Market of Goalpara town; and 53 (DH001-DH053) samples were collected from Netai Dhubuni Ghat, Dhubri New ghat, Jogomaya Ghat and Dhubri fish Market of Dhubri town. After taking weight and length of each fish, tissue samples of Hilsa were dissected out immediately and immersed in pure grade of ethanol in the sample bottle. After transporting to laboratory the vials were kept at -80° Celsius until DNA extraction.

DNA extraction
Total genomic DNA was isolated from muscle tissue by QIAGEN-DNeasy Blood and Tissue kit by manufacturer protocol. For each sample approximately 20 mg of tissues was cut into small pieces and genomic DNA was extracted by proteinase K digestion. The extracted DNA was stored at -20°C for PCR amplification. The extraction of DNA was followed by PCR amplification.

Primers
The partial cytochrome b gene sequences of mitochondrial DNA were isolated via PCR. The primers used in the reaction was taken from Irwin et al. (1991) and the primers were L14724 (5'-CGAAGCTTGATATGAAAAACCATCGTTG-3') and H15149 (5'-AAACTGCAGCCCCTCGAATGATATTTGTCCTC A-3').

PCR amplification
The PCR reaction was done in15µl reaction containing purified water 3.75µl, taq buffer 1.5µl, MgCl2 1.5µl, dNTPs 1.5µl, primers 1.5µl, BSA 1.5µl, taq DNA polymerase enzyme 0.75µl and 3µl of extracted DNA. The PCR cycling profile was initial denaturation at 94°C for 5 minutes followed by 35 cycles comprising subsequent denaturation at 94°C for 30 sec, annealing at 53.4°C for 30sec, extension at 72°C for 50 sec, incubation at 72°C for 10 minutes and finally at 4°C forever. 76 samples result positive PCR amplification. The amplified PCR products were then resolved on 2% agarose gel in 1XTBE buffer containing 2μl of ethidium bromide. After electrophoresis the PCR products were purified by adding 1μl of Exo-SAP (Exonuclease-I and Shrimp Alkaline Phosphatase in buffer) per 10μl reaction and set in PCR machine as 37° C for 20 minutes, 80° C for 15 minutes and finally at 4°C. After this, the purified PCR products were sequenced in the ABI Prism BigDye Terminator Cycle Sequencing kit (Applied Biosystems). Asian J Agric & Biol. 2021(2).

Data analysis
Amplified cytochrome b regions were sequenced in both forward and reversed directions to check the validity of nitrogenous bases. Alignment of sequences was done automatically using CLUSTALW program (Thompson et al., 1994), inbuilt in the genetic analysis package MEGA 7 . The sequences were checked with the Chromas software (Technelysium, Australia) and then visually refined. The sequenced data were analyzed with the help of BLAST (www.ncbi.nlm.nih.gov.) for homology search. Sequence variation sites, nucleotide diversity and haplotype diversity within population were estimated using DnaSP v5 (Librado and Rozas, 2009). After finding of haplotypes, a Neighbor-joining tree was constructed for all haplotypes with the Gadusia chapra sequence (Accession number: AP011603) as an outgroup based on the Kimura's 2 parameter model using MEGA 7.0. By applying this same model genetic distance among the haplotypes were also calculated. To assess confidence in the branching of dendrogram bootstrap analysis (1000 data sets) was performed. For a visualization of the relationship among the haplotypes Median-joining (MJ) network from all Hilsa sequences was constructed using the Network 5.0 (Bandelt et al., 1999). To explain the genetic structure of the Hilsa population collected from Goalpara and Dhubri, the allele frequencies (Fst) and AMOVA analysis were done in Arlequin 3.5 (Excoffier and Lischer, 2010).

Results
The Nucleotide positions of the sequenced segment of DNA were assigned from 14348 th to 14773 th position according to the complete mitochondrial genome reference sequence of the GenBank accession no. AP011610 and NCBI Blast result. From the sequenced cytochrome b region of T. ilisha samples, four haplotypes were obtained. The haplotype sequences have three variable (polymorphic) sites; out of these two are singleton variable sites and one parsimony informative sites. The average nucleotide diversity of the haplotypes is 0.002496 ± 0.001852 and the average nucleotide composition of all haplotype sequences are 25.92% Adenine, 27.53% Thymine, 17.49% Guanine and 29.06% Cytosine. The mean number of pairwise difference is 1.060702 ± 0.710370. The haplotype diversity is 0.533 and the average number of nucleotide differences is 1.667, indicating the low genetic diversity in Hilsa population in the Brahmaputra River. The similarity or percentages of the identity of haplotypes obtained were compared ( Table 1) with some of the Genbank sequences of the respective part of the DNA sequence. This implies that the haplotypes H_1 and H_4 has 100% similarity with few Genbank sequences. The Neighbor-Joining tree (Fig.1) has shown medium bootstrap values among the 4 haplotypes where Gadusia chapra was used as an outgroup.  (Fig.2) showed the positions of haplotypes. The MJ network showed the interrelation of haplotypes and their mutation sites. A single mutation site was observed by the network for all four haplotypes. The four haplotypes obtained in this study comprising of different frequencies (Fig. 3). At both the locations similar haplotypes were obtained except the H_3 which was only recorded in Dhubri. The haplotype 4 (H_4) comprises the largest number of sequences (46.05%) and percentage of occurrences of other three haplotypes namely H_1, H_2 and H_3 are 19.74%, 28.95% and 5.26% respectively. While conducting genetic structure analysis of Hilsa population, the Fst value (0.03604) showed that there is a little difference between the groups (Wright, 1965). In the AMOVA analyses the genetic variation within the population was found 68.94%, whereas between the groups was 3.70%. The mean pairwise difference of Dhubri and Goalpara groups are 0.6556 and 0.6400 respectively.

Discussion
A detailed knowledge of the genetic makeup is of primary importance for management as well as survival of a species. The future of a population is dependent upon its genetic diversity. The genetically less diverged population utilized limited available habitat, thereby create competition among themselves (Johnson et al., 2016). The present analysis of the 76 cytochrome b sequences comprising two groups showed that the there exists genetic variation with one transition and two transversions in the Hilsa population in the Brahmaputra River. The population differentiation between the two groups was less significant, although H_3 was only found in one group. A number of similar studies suggested the existence of different races of Hilsa shad in India and Bangladesh. Brahmane et al. (2013), studied cytochrome b gene of Hilsa population from Ganga and Hooghly Rivers and found low genetic diversity with single population. In contrast Behera et al. (2015) found five cytochrome b gene haplotypes of Hilsa in Bhagirathi and Tapti River where the differentiation of population was established in a single nucleotide substitution in the 307bp length of cytochrome b gene.
A genetically variable population of Hilsa in Bangladesh was also obtained by Rahman and Naevdal (2000) using   Verma et al. (2016) while working on the mitochondrial control region found broad structuring of Hilsa population between the Arabian Sea and Bay of Bengal but did not detect structure among riverine populations. However the recent SNP study based on NextRAD sequencing suggested that existence of local adaptive riverine populations called genetically structured ecotypes in differing environmental habitats (Asaduzzaman et al., 2019). In addition to this, two genetically structured clusters, one from marine and estuarine and other from riverine collection sites were described in where the riverine population was further subdivided into turbid freshwater and the clear freshwater ecotypes (Asaduzzaman et al., 2020). In the current study, 4 haplotypes of Hilsa were obtained, but the differentiation between haplotypesa is nominal. However, it is also important to note that a 425bp nucleotide study is not sufficiently reasonable to make a conclusion over a species. But it can be ascertained that there exists a genetically variable Hilsa population in the River Brahmaputra. Owing to this low genetic differentiation of the Hilsa population, the haplotypes obtained in this study cannot be regarded as a separate race or sub population and thereby it is difficult to demonstrate the Hilsa population as distinct one and separate from each other Rivers. The low differentiation observed might be due to inbreeding population structure and other evolutionary forces such as migration, selection, genetic drift and geographic barriers. Therefore, it can be said that Hilsa population of Brahmaputra River is genetically polymorphic though it did not reveal much genetic variation in the population. A combined study of physical, ecological, morphometric, and genetic study with more molecular marker may solve the question of existence of ecotypes and riverine population differentiation of Hilsa shad.

Conclusion
The present study has revealed significant information on the genetic structure of Hilsa population in Brahmaputra River and can be concluded that a genetically variable Hilsa population is available in the Brahmaputra River.