Environmental Factors and Genetic Diversity as Drivers of Early Gonadal Maturation: A Gonadosomatic Index Based Investigation on Indian Shad, Tenualosa ilisha Population of Bangladesh

In recent years, attaining gonadal maturation in smaller Hilsa (Tenualosa ilisha) has become a burning issue for Hilsa fishery of Bangladesh. Causes of early maturation are not yet clearly understood. Along with environmental parameters, genetic differentiation within the population was hypothesized as the main driver, and therefore, assessing the correlation between gonadosomatic index (GSI) and environmental factors and analyzing genetic diversity were set as objectives of the present study. To address these complex issues, six diverse habitats across Bangladesh were chosen for Hilsa sample collection. For GSI, gonad was dissected from fresh fish and preserved in Bouin’s fluid for histological observation. Water quality parameters such as temperature, dissolved oxygen, pH, and salinity were also assessed. 35 fish from each habitat were used to extract and amplify DNA through the PCR technique, and genetic diversity was examined. Further, to draw a firm conclusion, the phylogenetic tree of the Hilsa population was developed by the unweighted pair-group method of arithmetic mean method based on the Cyt b gene of mitochondrial DNA. Results of GSI studies revealed that peak spawning months of T. ilisha were in October and February, where October showed the highest values in all six habitats. Histological examination showed different stages of gonadal development in different sizes and ages of Hilsa. Among all sampling sites, no statistical difference was observed for GSI value; however, smaller sized and aged Hilsa being ripped were evident in Gaglajur Haor and Kali River. Among the observed water quality parameters, temperature correlated with GSI strongly. Increased GSI was observed with temperature augmentation from downstream to upper stream, irrespective of body size and age. A perplex correlation between dissolved oxygen of observed habitats and GSI was executed. Other physico-chemical parameters viz. pH and salinity exhibited weak and moderate positive association with the GSI, respectively. Haplotype diversity of mitochondrial DNA divided the Hilsa population into three possible sub-populations, where the most distant group (Gaglajur Haor and Kali River) was subjected to early gonadal maturity. Results of this study make clear conclusions regarding the role of environmental and genetic factors on early gonadal maturations, pointing fingers at the curse of climate change and anthropogenic stressors for the migration of the Hilsa fishery of Bangladesh.

In recent years, attaining gonadal maturation in smaller Hilsa (Tenualosa ilisha) has become a burning issue for Hilsa fishery of Bangladesh. Causes of early maturation are not yet clearly understood. Along with environmental parameters, genetic differentiation within the population was hypothesized as the main driver, and therefore, assessing the correlation between gonadosomatic index (GSI) and environmental factors and analyzing genetic diversity were set as objectives of the present study. To address these complex issues, six diverse habitats across Bangladesh were chosen for Hilsa sample collection. For GSI, gonad was dissected from fresh fish and preserved in Bouin's fluid for histological observation. Water quality parameters such as temperature, dissolved oxygen, pH, and salinity were also assessed. 35 fish from each habitat were used to extract and amplify DNA through the PCR technique, and genetic diversity was examined. Further, to draw a firm conclusion, the phylogenetic tree of the Hilsa population was developed by the unweighted pair-group method of arithmetic mean method based on the Cyt b gene of mitochondrial DNA. Results of GSI studies revealed that peak spawning months of T. ilisha were in October and February, where October showed the highest values in all six habitats. Histological examination showed different stages of gonadal development in different sizes and ages of Hilsa. Among all sampling sites, no statistical difference was observed for GSI value; however, smaller sized and aged Hilsa being ripped were evident in Gaglajur Haor and Kali River. Among the observed water quality parameters, temperature correlated with GSI strongly. Increased GSI was observed with temperature augmentation from downstream to upper stream, irrespective of body size and age. A perplex correlation between dissolved oxygen of observed habitats and GSI was executed. Other physico-chemical parameters viz. pH

INTRODUCTION
Bangladesh has one of the most extensive and active deltas, fed by three mighty rivers: the Padma, the Meghna, and the Jamuna. This contributes to a high potential for fresh and brackish-water fish, in addition to the vast marine resources (Hossain, 2001;Haldar and Islam, 2003). There are almost 293 fish species where 18 fish from Clupeiformes were reported from Bangladesh territory by Hossain et al. (2012), and Hilsa, Tenualosa ilisha (Hamilton, 1822), is the dominating species of the group. Hilsa is mainly an anadromous fish that migrates from saline water to freshwater for the purpose of breeding. During the breeding period, it comes to the river channels of the adjacent sea countries and returns to the sea after spawning. After hatching, the juvenile Hilsa remain in freshwater for a few months for their development and then return back to the sea. It owns the largestspecies fishery in Bangladesh, contributing 12% (0.53 million MT) of the total fish production (Department of Fisheries, 2019).
With a wide geographical distribution of Hilsa ranging from the Persian Gulf to the Arabian Sea, this transboundary fish takes shelter in the shallow water of the Bay of Bengal (Blaber, 2000;Arai and Amalina, 2014). For breeding purposes (with a further intention of nursing the juvenile), Hilsa follow their main migration route [Meghna river (lower part)-Meghna river (upper stream)-Padma (lower stream of the Ganges)-Jamuna (New Brahmaputra River)] to reach in upper stream rivers, which facilitate them with an uncompromising favorable environment in purpose of their breeding success. Apart from this main migration route, there is also evidence that Hilsa takes entry in freshwater through some other unpopular and secondary coastal migration route (Pashur River-Madhumati River-Gorai River-Padma River; reviewed in Sarker et al., 2021). With an unknown number of mysterious issues about the migration intention of Hilsa, the fish have also been found to use unfocused water bodies, including marshy wetland ecosystems known as Haor and hill stream rivers (e.g., Someshwari River) bordering the northeastern part of Bangladesh as their travel route (Sarker et al., 2021).
Though it can be found in almost every river, estuary, and marine environment in the country (Ahmed et al., 2008), the number of small sized Hilsa with matured gonad in natural water bodies are increasing alarmingly (Almukhtar et al., 2016). The reason behind this menacing problem observed in the Hilsa of Bangladesh is still unknown. Whether these small sized Hilsa are adults, or environmental factors (e.g., climate, hydrological, and ecological changes) or genetical factors have forced the fish to develop premature gonad in the juvenile stage is a big question (Almukhtar et al., 2016). Environmental factors are very much related to climate change, and Bangladesh, a low-lying floodplain country bordering the GBM basins (Ganges, Brahmaputra, and Meghna), is considered as one of the most vulnerable countries to the consequences of climate change (Agrawala et al., 2003;World Bank, 2010). Intensive rainfall during monsoons and other natural disasters (for example, cyclones, floods, storm surges, etc.) becomes a regular incidence in Bangladesh due to geographical south attachment with Bay of Bengal and north attachment with Himalayas (Ferdous and Baten, 2011). With a sharp increase of 0.103 • C temperature per decade of the country in the last four decades (Shahid, 2010), it has been projected that the temperature will continue to increase by 1 • C by 2030, 1.4 • C by 2050, and 2.4 • C by 2100 (IPCC, 2007). Rising temperature has also been altered by the seasonal variability in rains with their spatial distribution (Shahid and Khairulmaini, 2009) by augmenting air moisture-holding capacity.
Usually, the old literature indicates that Hilsa stock of the Indian subcontinent became mature between 16 and 17 cm in case of males and for females the range was 19-20 cm, when their age was 2+ or 3+ (Pillay, 1958). However, in an alteration of seasonal pattern with a combined effect of vigorous nourishment during their migratory journey, the age range has diverged for the sexual maturity of Hilsa, which now can be appear at different sizes and ages of fish (Zhang et al., 2009;Mohamed and Qasim, 2014). Furthermore, genetic factors, feed competition, or their persistence of stay in the river may stimulate early gonadal development. It may also be possible for some small Hilsa to reside in the upper stream instead of returning to the downstream and thus, attain gonadal maturity in freshwater (Ahmed et al., 2018). However, these hypotheses need intensive research for confirmation. Moreover, determining the gonadosomatic index (GSI), observation on gonadal maturation, genetic diversity of the selected stocks as well as assessment of the water quality parameters of the habitat of Hilsa are reasonably necessary for the validation of these hypotheses. Therefore, it was the center of the attention of the present study.
Reproduction is a continuous developmental process throughout ontogeny, requiring energetic, ecological, anatomical, biochemical, and endocrinological adaptations (Caputo et al., 2000). The GSI is a good indicator of reproductive activity being used to determine gonadal maturation stages (Hismayasari et al., 2015). It represents details regarding the reproductive status and determining reproduction time of fish (Mohan and Jhajhria, 2001;Shankar and Kulkarni, 2005). GSI also provides data on seasonal variations in multiple spawning and variation of food consumption (Pradhan, 1965) and facilitates the effective conservation of any species (Ghaedi et al., 2013;Rahman, 2014). Zhang et al. (2009) confirmed that the spawning activities of fish can be determined based on the monthly variations of GSI, and it is proportionate to each maturity stage. Bladon et al. (2019) in their study, confirmed that a significant drop of GSI from high to low is a signal of egg depositing activity in fish and that is why the GSI is a convenient method to estimate the spawning period in fish. Variations in the GSI are related to the changes in the water quality parameters and food abundances, reflecting the size variation of Hilsa (Hasan et al., 2016). In addition, water quality parameters control the maturation of gonads in fish (Hazard and Eddy, 1951). Therefore, fish nutritional status and gonad maturation are a result of their relationship with environmental factors (Ribeiro-Silva et al., 2018). Once again, histological study on the gonadal changes of fish is an essential tool for identifying reproduction level, which includes different developmental phases. Unver and Saraydın (2012) found a strong correlation between histological change in gonad and the sexual maturation of Hilsa. Therefore, a combined study of GSI and gonad histology can precisely determine the maturation timing of any fish species. Study of vitellogenic stage (VG) of oocyte and increased gonad weight (GW) and histology of gonad are the important aspects of reproductive biology that indicate the maturation level and can estimate approximate spawning time (Brown-Peterson et al., 2011). Knowledge on the fish size and age at first sexual maturity is essential to estimate the spawning stock's size (Neja, 1992). Hence, size and age at first sexual maturity are closely linked to lifetime fecundity of fish (Stearns, 1992;Bernardo, 1993). In this regard, Ahmed et al. (2020) initiated a good attempt to determine the size and age at first sexual maturity of T. ilisha by otolith examination and confirmed that Hilsa of 1.2 years with a total length (TL) of 26 cm is sexually matured, which is a clear indication of early gonadal maturation.
Water quality parameters such as water temperature, pH, and oxygen are the crucial drivers of the gonad development for several fish species in the tropical region (Bailly et al., 2008). For instance, the water temperature can trigger the reproductive processes of many fish species and, therefore, their gonadal development may be connected to rising sea water temperature (Mondal, 2012). Mylonas and Zohar (2007) reported that a change in water temperature had a high correlation with the GW of fish and GSI and sexual maturation. Seasonal variations of physico-chemical factors of water also directly influence the gonad maturation of fish since they affect the availability of food and trigger a number of neuro-endocrinological processes in fish (Araujo et al., 2019). Roomiani et al. (2014) again reconfirmed that the differences in the spawning maturation of fish in different areas may be due to genetical and environmental factors like temperature, pH, DO, etc. After maturation for spawning purpose Hilsa shad follow two types of migration pattern-(1) fluvial potamodromous (fish remain in freshwater habitats and breed therein) and (2) marine ecotype (fish migrate to salt water for breeding purposes) in Bangladesh (Raja, 1985). Based on these, different levels of genetic variations and diversities have been recognized by previous studies (Majumder and Alam, 2009;Asaduzzaman et al., 2019Asaduzzaman et al., , 2020Sarker et al., 2021). However, interpretation on whether these variations are caused by different migratory patterns or is the cause of reciprocal migration is not yet understood. The genetic architecture of vertebrates, including fish, is a ubiquitous and potent caliber that shapes ovarian maturation, leading to defining the next generation (Barson et al., 2015;Mobley et al., 2021). It is well established that the gonadal maturation of migratory Atlantic Salmon is strongly influenced by genetic factors affecting the span of different life stages (Healy et al., 2019). However, the relationship between genetic diversity and early gonadal maturity of Hilsa shad remains unboxed. Due to retaining maternal inheritance throughout the generation, preventing recombination, and showing a better alteration rate, mitochondrial DNA is considered more appropriate to measure genetic deviation than nuclear DNA (Majumder and Alam, 2009). To represent genetic differentiation among the six different habitats and find an underlying relationship with gonadal maturation, the Cyt b gene of mitochondrial DNA were analyzed to calculate haplotypes and genetic diversity.
In compliance with the above research, a better understanding of early gonadal development in T. ilisha during the natural reproductive cycle is crucial to determine the actual cause of early gonadal maturation of Hilsa in Bangladesh. Therefore, the present study was aimed to (1) determine the seasonal variation of GSI, (2) assess the gonad maturity of T. ilisha across six diverse habitats of Bangladesh through histological observation, (3) determine the variation in the water quality parameters in those habitats during the reproductive period, and (4) analyze the genetic diversity to identify whether the selected Hilsa stock maintained single or multiple sub-populations.

Sample Collection, Age, and Growth Determination of T. ilisha
Fresh Hilsa were collected randomly from six different locations of Bangladesh, namely the Meghna River Estuary, Chandpur; the Bay of Bengal, Cox's Bazar; the Kali River, Kishoreganj; the Tetulia River, Barisal; the Padma River, Mawa, Munshiganj; and the Gaglajur Haor, Mohanganj, Netrokona (Figure 1). Fish were collected with the help of local fishers by using fishing boats from January to December, 2018. All samples were carried to the Fisheries Biology and Genetics Laboratory of Bangladesh Agricultural University (BAU), Mymensingh. Thirty samples from each habitat were collected for analysis. TL were taken to the nearest mm using a measuring scale attached to a board. The individual whole fish was measured as body weight (BW) using an electronic digital balance (MONOBIOC, B2002-S). The age and growth pattern of different sized Hilsa were determined based on the lunar ring in the otolith and length-weight relationships, respectively (Ahmed et al., 2020). Thereafter, the samples which were under normal maturation were used for the present study. To visualize the rings in the otolith, they were ground and polished using a combination of the acid contact method (Sinha and Jones, 1967) and the wet-dry sandpaper method (Ahmed et al., 2020). The use of too much acid was avoided to hinder contact with acid to the edge of the otolith. During this grinding procedure, the otolith was sporadically examined under a microscope to check for the appearance of the rings. The growth pattern of Hilsa was determined by the length weight-relationship for each of the four habitats. The mean values of TL, weight, age, and their growth pattern were recorded for tabulation.

Determination of Gonadosomatic Index
The fish's abdomens were opened by cutting the abdomen with a sharp scissor, and the muscle, fat tissue, digestive organs, blood veins, etc. were appropriately removed. Both left and right gonads were measured together as GW. The collected gonads were kept in Bouins fixative about 8-10 h. The ratio of gonad and fixative was 1:10. After 8-10 h, gonads were rinsed with 10% formalin. At the end, these gonads were preserved with 10% formalin in small vials. The GSI of each fish specimen was calculated by the following formula (Howaida et al., 1998): Where, GSI = Gonadosomatic Index, GW = Gonad Weight, and BW = Body Weight.

Relationship Between Body Weight, Total Length, Age, and Gonadosomatic Index
For establishing the relationship of GSI with different body parameters viz: TL, BW, and age, all data were recorded with sophisticated instruments. TL and BW were measured as described in section "Sample Collection, Age and Growth Determination of T. ilisha" in this manuscript. Age of each fish were counted from their otolith. GSI were calculated by the method as described in section "Determination of Gonadosomatic Index." All data were processed in Microsoft excel to generate a comparative graph between GSI and TL, BW, age of Hilsa.

Histological Observation of the Ovary
The histological observation was performed with several steps such as dehydration, cleaning, infiltration, embedding, trimming, sectioning, staining, and mounting. The dehydration, clearing, and infiltration of the sample were done using an automatic tissue processor in a series of increasing alcohol solutions, according to Johnston et al. (2009). After that, two changes of xylene and three series of changes in paraffin wax were completed in the processor. Finally, samples were left overnight in the laboratory embedded in and covered by paraffin mass, which eventually solidified and supported the segment during sectioning. Paraffinembedded blocks were cut by a microtome knife at 4-5 µm size and left the sections into a water bath at a temperature of 37 • C as described by Ahammad et al. (2021). the side of the paraffin blocks were trimmed by a sharp knife to obtain a suitable cuboidal shape in which the undesirable wax (paraffin) layer of the embedded blocks was removed. Sectioning was done by using a manual microtome machine (Biobase, Bk-MT260M).
After obtaining sections in ribbons, it was correctly placed in a hot water bath having hot water of 50 • C. At the end, the slide was stained with xylene, hematoxylin followed by periodic acid Schiff reaction. All images from histological observation were processed following the method of Moallem et al. (2015), using a high magnification microscope (CX 41, Olympus) attached with a high-resolution camera.

Observation of the Water Quality Parameters
Water samples were collected in glass bottles throughout the year from six selected habitats of Hilsa. Samples were triplicated from each location for minimizing standard error. After sampling, the bottles were marked with the respective identifying numbers.
Water quality parameters such as temperature, dissolved oxygen, and pH were recorded immediately, and salinity was measured in the laboratory. Water temperature ( • C) from each habitat was recorded by using a Celsius thermometer (SMART SENSOR * AR867, China). Dissolved oxygen (mg/l) and pH were measured using a digital DO meter (HANNA, Model-HI 9146, Romania) and a direct reading digital pH meter (Milwaukee PH meter, Model-PH55/PH56, United States), respectively. Salinity was measured by a Hand Refractometer (ERMA Hand Refractometer, Model-MASTER-S10a, Mumbai, India).

Determination of Genetic Diversity of Selected Hilsa Populations
Genomic DNA was extracted from the fin samples of Hilsa (N = 210; 35 samples/populations) using GeneJET Genomic DNA Purification Kit (Thermo Fisher Scientific, United Kingdom). Next, DNA was measured by GENEVA NANO micro-volume spectrophotometer, and amplification was performed by PCR using the universal primer set for the Cyt b gene of the mitochondrial DNA. Then amplified PCR products were electrophoresed and visualized on 1.0% agarose gels containing 2 µl ethidium bromide. The PCR products were purified according to FavorPrep GEL/PCR Purification Kit (Favorgen Biotech Corp) protocol and sequenced bidirectionally using an ABI 3730 Genetic Analyzer following the manufacture's protocol. The raw nucleotide sequences of all the specimens were trimmed and aligned among themselves using the ClustalW software in MEGA 7.0 (Kumar et al., 2018). Haplotypes were calculated using the Kimura-2-Parameter model (Kimura, 1980). The complete aligned dataset was calculated for polymorphic sites (Ps), the number of haplotypes (H), haplotype diversity, and nucleotide diversity for each population using DnaSP 5.10 (Librado and Rozas, 2009). We estimated Tajima's D (Tajima, 1989) and Fu's Fs (Fu, 1997) to detect signatures of population bottlenecks or expansions and deviations from the pattern of polymorphism expected from the neutral model of evolution. The deviation from neutrality assessed by both D and Fs indices was determined in 1,000 simulations. The observed and expected distributions of the number of pair-wise mutational differences between haplotypes were calculated to their relative frequencies (mismatch distributions) under the sudden expansion model (Rogers and Harpending, 1992).

Statistical Analysis
Growth and GSI data of Hilsa were analyzed using MS Excel 2016. Comparison between TL, weight, age and GSI were established using MS Excel 2016. Significant difference of GSI between the habitats was done by one-way analysis of variance followed by duncan's multiple range test at 95% level of confidence. Scatter plot representing the water quality parameters against GSI of T. ilisha was attained from the linear model. The Spearman rank correlation test was executed to measure the linear or monotonic association (Al-jabery et al., 2020) between water quality parameters and GSI of T. ilisha. The statistical tests were two-sided, and statistical significance was considered at P < 0.05. Statistical analyses were performed in R (RStudio Team, 2018). The genetic distance values were used to construct the unweighted pair-group method of arithmetic mean (UPGMA) dendrogram using the software POPGENE version 1.32 and MEGA version 10.0.5 (2018).

Age, Size, and Growth Pattern
Mean TL was found highest in the case of Hilsa from the Tetulia River. However, the highest standard deviation (in the case of mean weight) indicates that the size variation was very high in this habitat, and the fish were not in uniform size. On the other hand, the TL and weight was lowest in the case of fish from the Gaglajur Haor. Mean age was highest in the fish from the Tetulia River, which was on average 41.8 months. The lengthweight relationship indicates a positive allometric growth for the Meghna River, the Bay of Bengal, the Tetulia River, and the Padma River. On the other hand, the Kali River and the Gaglajur Haor showed isometric and negative allometric growth, respectively ( Table 1).

Gonadosomatic Index During Spawning Season
In all the habitats, GSI values were found apparently similar (P > 0.05; Supplementary Table 1). Hilsa from the Meghna River, the Padma River, the Tetulia River, the Bay of Bengal, the Gaglajur Haor, and the Kali River showed GSI values ranged from 4. 9-11.6, 4-12.4, 4.9-13.3, 4.7-12.4, 4.5-11.2, and, 4.7-11.6, respectively, recorded from January to December. The highest GSI values were recorded in October, and the lowest in April. However, in the case of the Bay of Bengal the lowest GSI was recorded in June. Mean GSI of the Tetulia River showed higher values in August (8.4) and October (13.3) than fish from other habitats in the same months. Clearly, the inconsistencies for the Gaglajur Haor and the Kali River were visible because they both showed lower curves than the other four habitats.
It indicates that fish from these two habitats consist of more negligible differences between their body and GW (P > 0.05). The mean GSI value is higher in February and October in all the habitats indicating these 2 months as peak GSI and thereby the spawning time for Hilsa (Figure 2). Nevertheless, between these 2 months, GSI in October was higher than February, indicating October as the major spawning season for Hilsa in all six habitats.
Length-Weight, Age, and Gonadosomatic Index Relationship Body weight, TL, and age with respect to GSI showed a regular and uniform pattern in all Hilsa samples from six habitats (P > 0.05). However, in the case of Hilsa from the Kali River and Gaglajur Haor, during the time of peak spawning season, the BW was the lowest (Figure 3). On the other hand, the BW of Hilsa from other habitats showed the highest value in different seasons of the year. In relation to GSI, all graphs showed a proportional relationship with TL, BW, and age as these three parameters increased with GSI in all six samples. In all the collected samples, fish were older and held the highest GSI in October compared to other months of the year.

Ovarian Maturity
In all the six habitats, at the first stage, oocytes were spherical, oval, and multi-faceted, and the large nucleus occupies most of the cell. Prominent nucleoli were found sporadically in the nucleoplasm. The cytoplasm was strongly basophilic and stained with hematoxylin ( Figure 4A). Nucleolus and Yolk nucleus were clearly visible in the pronucleus stage ( Figure 4B). Small alveoli were formed around the nucleus, which increased in number and size during the cortical alveoli stage ( Figure 4C). In primary vitellogenesis stage, vitellin granules around the nucleus clearly appeared ( Figure 4D). Fat vacuoles were seen in large numbers around the pore. The number of granules increased and moved to the center of the cell in secondary vitellogenesis ( Figure 4E). Cell layers were thickened, and zona radiata was stained with eosin heavily (Figure 4F). In the tertiary vitellogenesis stage, oocytes were seen with an unaided eye (Figure 4G), and they were more spherical in shape ( Figure 4H).
In the comparison of ovarian development with age from the previous study by Ahmed et al. (2020), Hilsa from the Meghna River, the Bay of Bengal, the Padma River, and the Tetulia River showed regular ovarian development in respect to their age. No fish were found with a fully matured ovary at the age below 2 years. However, in the case of the Kali River and  the Gaglajur Haor, irregular ovarian development was recorded, and fully developed gonads were confirmed by identifying the tertiary vitellogenesis stage in some small size Hilsa with an age of 14 months.

Temperature
The uphill pattern of the regression line with R > 0.476 in habitats except Gaglajur Haor indicates a moderate positive relationship between temperature and GSI of T. ilisha. The relation between temperature and GSI of T. ilisha in Gaglajur Haor was found to have a strong positive relationship (R > 0.93; Figure 5). In all habitats, p-values were lower than 0.05, indicating that there is a very strong correlation between the temperature and gonadal maturity of T. ilisha.

Dissolved Oxygen
The downhill pattern of the regression line with R < −0.968 in four habitats from downstream except for the Meghna River indicates a negative relationship between dissolved oxygen and GSI of T. ilisha. The negative relationship found in the Padma River was significant (P < 0.05) and in the other two habitats, the Tetulia River, and the Bay of Bengal, the negative association was non-significant (P > 0.05; Figure 6). For the Meghna River and other two habitats from upper stream (the Kali River and the Gaglajur Haor), the uphill pattern of the regression line with a R < 0.476 indicates a moderate positive relationship between dissolved oxygen and GSI of T. ilisha. In the case of a positive relationship, the Meghna River was found to have a significant impression (P < 0.05), where in both habitats from the upper stream, a positive but non-significant relationship was found (P > 0.05; Figure 6).

pH
Among the four downstream habitats, the uphill pattern of the regression line with R < 0.476 in the Padma River, the Meghna River, and the Bay of Bengal apprise the moderate positive relation within pH and GSI of T. ilisha. The relationship between pH and GSI for the Tetulia River and other observed habitats downstream exhibited very weak positive relationships for having R∼0.07. And these positive relationships for all four habitats downstream were found to have non-significant impressions (P > 0.05). The relationship between GSI of T. ilisha and pH for Kali River was significantly moderate negative (R = −0.67, P < 0.05). In another habitat, Gaglajur Haor, from the upper stream, the relationship was very weakly negative (R∼ −0.13) and non-significant (P ≥ 0.05; Figure 7).

Salinity
The relationship between GSI of T. ilisha and salinity in two saline habitats (Tetulia River and Bay of Bengal) from downstream was found to be moderately positive (R < 0.476) and non-significant (P > 0.05; Figure 8).

Spearman Rank Test
Data regarding the GSI of T. ilisha and water quality parameters of all observed habitats were included in the Spearman rank correlation test to obtain their specific correlation coefficient (SCC) values. In the case of the Spearman rank test, the habitats were categorized into upper stream (Kali River and Gaglajur Haor) and downstream (Padma River, Meghna River, Tetulia River, and Bay of Bengal) and their average values were considered for both GSI and all water quality parameters. Among the average water quality parameters of the upper stream, only the average temperature showed a strong correlation (SCC: 0.615), with a very high statistical significance (P < 0.05) to GSI of upper stream T. ilisha. For the same pair of variables (average GSI and average temperature), the same strong (SCC: 0.581, P < 0.05) correlations were observed in T. ilisha downstream ( Table 2).
In the same dataset, the average temperature was found to have strong but negative (SCC: −0.628, P < 0.05) correlation with moderate dissolved oxygen content for upper stream. Moreover, the average temperature downstream had a strong correlation (SCC: 0.693, P < 0.05) with average pH ( Table 2).

Genetic Diversity of Selected Populations of Hilsa Shad Using mtDNA Cyt b Gene Sequence
The genetic diversity of six different populations of Hilsa shad, T. ilisha was presented in Table 3. However, the highest number of haplotypes were observed in the Bay of Bengal population (24) followed by the populations of Tetulia River (18), Meghna River (16), Padma River (14), Gaglajur Haor (11) Table 3).
From the pairwise mismatch distributions of mtDNA haplotypes, we observed smooth and unimodal (bell-shaped) distribution of all Hilsa populations (Bay of Bengal, Tetulia River, Meghna River, Padma River, Gaglajur Haor, and Kali River; Figure 9). No mismatch distribution was found in the wild population.
The UPGMA dendrogram showed that the six Hilsa shad populations and wild T. ilisha species had been grouped into two major clusters, as shown in the dendrogram (Figure 10). The wild T. ilisha was in one cluster and separated by the highest genetic distance of D = 0.80. In contrast, the remaining cluster for four rivers (i.e., Tetulia River, Meghna River, Padma River, and Kali River) including the Bay of Bengal and Gaglajur Haor populations was subdivided into three clusters and separated by the genetic distance of D = 0.15.

DISCUSSION
Gonadosomatic index is an indicator of the state of gonadal development of fish and can be estimated during spawning season (Duarte et al., 2007;Oso et al., 2011;Amtyaz et al., 2013). Narejo et al. (2008) studied GSI for female T. ilisha from River Indus, Pakistan and recorded GSI data separately for summer (from April to September) and winter (from October to March). They observed that in the females, GSI index was higher in July and February, which is inconsistent with the results obtained in the present study. This might be due to the different geographical habitats where the studies were conducted. Flura et al. (2015) conducted a study on GSI of T. ilisha from the Meghna River, Chandpur, Bangladesh where they found the peak GSI value for females in June, October, and February. The highest mean GSI value was found in October; hence it determined October as the Pearson correlation coefficient is reported as "R" while probability value is reported as "P." main spawning season of T. ilisha. In the present study, from October to February is reported as peak spawning period for Hilsa, which is supportive of the findings of Flura et al. (2015). For all six habitats, GSI values started to increase gradually from August and reached their highest peak in October and then decreased, which indicates the spent stage of the fish. Thereafter, the GSI indices increased again in the month of February, which is determined as the winter spawning period. This report of one major spawning period in October and one minor spawning period in February corroborated with the previous studies conducted by Flura et al. (2015) and Bladon et al. (2019) in the case of Hilsa. However, it is evident from the present study that no early or late peak spawning season has been occurring in these six habitats. In all sampling sites, the maturation of the ovary was found linear with the age and size of the Hilsa. But in the case of Kali River and Gaglajur Haor, smaller sized and aged Hilsa were found with matured ovaries. In general, Hilsa gets its first sexual maturity in 31-33 cm of TL (Hossain et al., 2018). But this length of fish does not necessarily represent the age of fish in some stock (Ahmed et al., 2020). In a previous work with the same fish sample, Ahmed et al. (2020) conducted a study from the same six locations using the same fish sample as in the present study. In that study, the author reported age of Hilsa in respect to their length and weight. In the present study, histological observation was carried out using the same ovary sample from the fish used by Ahmed et al. (2020) to determine age through otolith examination where some Hilsa were found with matured ovary Pearson correlation coefficient is reported as "R" while probability value is reported as "P." before reaching their usual age of sexual maturity. Except for the Kali River and Gaglajur Haor, Hilsa showed no incongruity in the case of ovarian development in respect to their size and age. Mean GSI of Kali River and Gaglajur Haor were higher than that of the rest four sampling site, indicating that ovaries of kali and Gaglajur Haor sized larger when compared to their respective body size. Data showing that smaller sized and aged fish bear larger GSI, elucidates the evidence of early gonadal maturity in these two habitats. However, ripended fish ranging smaller than 25 cm (declared as Jhatka-juvenile Hilsa by Department of Fisheries under Ministry of Fisheries and Livestock, Bangladesh) provokes alarm from considering sustainable fisheries management.
A few fish were found to have developed gonad at small size but were mature in age. This irregularity in ovarian FIGURE 8 | Scatter plot and linear regression for the year-round means of Gonadosomatic index (GSI) of Tenualosa ilisha against mean salinity for two diverse habitats across Bangladesh: (A) Tetulia river, (B) Bay of Bengal. For other four habitats (Kali river, Gaglajur Haor, Padma river, and Meghna river) scatter plot and linear regression was avoided for their observed data recorded of zero salinity. Pearson correlation coefficient is reported as "R" while probability value is reported as "P." development may occur due to their separate environments because these two-sampling sites are situated in the upper stream of the country, which may cause difference in availability of nutrients. Some environmental factors also may drive this ovarian development in Hilsa. According to Nair and Nancollas (1958) and Mondal (2012), sea water temperature rise is an important factor that may push early gonadal maturity in Hilsa.
As the temperature of the Bay of Bengal is rising day by day (Kay et al., 2018), early maturity of Hilsa in this habitat is possible. This matured fish may migrate to the upper streams, the Kali River and the Gaglajur Haor. According to Hasan et al. (2016), anomaly of current velocity in the upper stream areas also influence sexual maturity in Hilsa. These variations could be attributable to genetic makeup, such as inbreeding, genetic depression, or changes to the expression levels of genes related to gonad development of Hilsa (Ahmed et al., 2020).
Variation in gonad development of fish may also be possible due to climate change, food availability, chemicals, agricultural fertilizers, industrials and urban waste water, various endocrine disrubters or a combination of these factors (Yoneda and Wright, 2005;Caballero-Gallardo et al., 2016;Servili et al., 2020). It can be assumed that the Hilsa in the Kali River and Gaglajur Haor could not migrate to the sea after spawning and were somehow trapped in these water bodies. However, this phenomenon needs further investigation for confirmation by transcriptomic analysis of gonad, kidney, brain, liver, and muscle samples of different sized Hilsa. Environmental temperature variation of aquatic habitats over the year is fundamentally important to fish for their ectothermic physiological process to which they are ambient (Neuheimer and Taggart, 2007;Buisson et al., 2008;Durance and Ormerod, 2009). Successful reproduction of T. ilisha is influenced by different extrinsic factors such as temperature, pH, salinity, and other properties of water (Jahan et al., 2017) and limits to attain successful artificial reproduction. Non-uniform water temperature across different habitats throughout the year leaves a deep impression on the early reproductive process of T. ilisha as observed in this study (Supplementary Tables 2-7). This alteration in temperature may influence the gamete maturation, ovulation, and spawning of T. ilisha like other fish (Hilder and Pankhurst, 2003;Pankhurst and King, 2010). In the case of female T. ilisha, the increased or substandard temperature of their ambient environment possibly results in conformational change of cytoskeleton dynamics like Danio rerio (Nonnis et al., 2021) and this cytoskeleton participates in oocyte maturation following ooplasm segregation (Ivanenkov et al., 1987). The maturation accelerated by an aggregation of microtubules in high temperature and halted by breakdown of microtubules in decreased temperature (Lessman, 1987;Nogales et al., 1995). These interfere in nucleus migration for active participation of actin microfilaments in cytoplasmic polarization (Becker and Hart, 1996), formation of the first polar body (Maro and Verlhac, 2002), fertilization and embryo morphogenesis (Klymkowsky and Karnovsky, 1994) entangled with final oocyte maturation. In contrast to the existing relationship between temperature and gonadal maturity found in this study, Jonsson et al. (1991) proved no influence of differing temperature on the maturity of Salmo trutta. Meanwhile, geomorphological and ecological changes in the habitats altering the reproduction pattern of Hilsa. Moreover, temperature augmentation of sea water have initiated the maturation of Hilsa (Mondal, 2012). In previous years it took about minimum 2 years to reach gonadal maturation for Hilsa, however, in recent years they are maturing within 1 year (Mondal, 2012). This might be due to the factors that are adversely influencing the upstream migration are siltation in habitats, newly constructed manmade barriers on the rivers and change in temperature variation over the years between downstream and upper stream (Mondal, 2012). This improves our understanding of the largest GSI in the smallest fish of the youngest group. The present study reveals weak and moderate positive relationships of GSI with pH and salinity, respectively. Hossain et al. (2016) reported the 6-7.5, the optimum pH for different stages of T. ilisha; though in this study, the pH ranged from 6.4 to 9.65 in different habitats. In the case of Padma River, Meghna River, and Bay of Bengal, the pH was very much close to the mentioned optimum range by Hossain et al. (2016; Supplementary Tables 2-7). For other habitats this relationship ranged between weak positive to moderate negative possibly for the reason that high pH driven waterbodies acidification differentially affect physiology, reproduction, biochemistry, and survival of organisms (McClellan-Green et al., 2007;Ellis et al., 2014;Lane et al., 2015) of those habitats. Especially the female T. ilisha possibly becomes vulnerable to stress from alteration of their resource energy location and high production cost of eggs compared to sperm like of other fish (Marčeta et al., 2020). Moderate positive relationships between GSI and salinity of two saline habitats (Tetulia River and Bay of Bengal) of this study suggest that gonadal development of T. ilisha was promoted by salinity. In Eriocheir sinensis (Long et al., 2017) and Crassostrea gasar (Paixão et al., 2013), gonadal development was stimulated by their subservient habitats salinity. There are a number of potential explanations that could have contributed to this relationship; (1) the resemblance of osmolality of brackish water and T. ilisha hemolymph spare less energy consumption on osmoregulation (Dutta et al., 2019) make the way to use nutrient/energy for gonadal development; (2) gonadal development of T. ilisha is substantially fueled by calcium ion enriched brackish water. This fact has been substantiated from poor gonadal development of Crangon crangon (Gelin et al., 2001), and Grapsid crabs (Bas and Spivak, 2000) in low salinity.
We observed that freshwater populations of Hilsa especially Kali River, Gaglajur Haor, Meghna River, Padma River, and Tetulia River showed the lower number of haplotype, Ps, haplotype, and nucleotide diversity compared to the Bay of Bengal population (see Table 3 and Figures 9, 10), indicating low genetic diversity existing in the Hilsa stock of freshwater habitats. Majumder and Alam (2009) also observed similar results, reporting the lowest genetic variability occurred in the Meghna River compared to the Bay of Bengal population of Hilsa shad in Bangladesh, while using D-Loop region of mitochondrial DNA. Their results suggested that Hilsa shad have at least two, one marine and the other riverine-estuarine, or possibly three (riverine, estuarine, and marine), stocks of Hilsa in Bangladesh. Dahle et al. (1997) revealed three sub-populations of Hilsa shad in Bangladesh waters, such as Freshwater, Brackish water, salt water (Bay of Bengal), when using RAPD markers. Preliminary  genetic differentiation is more evident in distant groups (Gaglajur Haor, Kali River, Padma River, Meghna River, and Tetulia River) when we consider the Bay of Bengal as natal habitat and can opine that migration patterns of Hilsa might influence the genetic variation. Hilsa from the Gaglajur Haor and Kali River might not follow the complete migration to the marine environment due to environmental and/or anthropogenic stressors and return to their habitat from middle streams. On the other hand, Hilsa near the coastal region (lower part of Meghna, Padma, and Tetulia rivers) might have easier access to marine habitats during the juvenile stage and return during the breeding season. As a result, while Hilsa from Padma, Meghna, and Tetulia rivers showed close genetic diversity to the Bay of Bengal, fish from Gaglajur Haor, and Kali River showed abnormal diversity, getting off from wild population. With this anomaly, the sub-population of Hilsa as distinct and the distant groups (Kali River, Gaglajur Haor, upper part of Meghna and Padma Rivers) were found to be subjected to early gonadal maturation. Unrooted phylogenetic tree based on Cyt b gene of mitochondrial DNA clearly indicates that early gonadal maturation in Hilsa shad of Bangladesh might be controlled by genetic factors. However, this incidence needs more clear investigation of upregulatory or downregulatory genes involved in early gonad maturity of Hilsa shad through transcriptomic analysis.

CONCLUSION
For stock assessment of any fish, the evidence-based calculation is necessary. From the aforementioned results, the obvious and caliber dominance of environmental factors, particularly temperature and genetic diversity on the early gonadal maturity of Hilsa shad can be concluded. While GSI is strongly correlated with body size and age of Hilsa shad, early gonadal maturation evident in distant population groups (Gaglajur Haor and Kali River) is linked to the haplotype diversity. Premature gonad ripening of Hilsa, a menacing problem in Bangladesh, is clearly unpacked in the present research. It will support the researchers in developing further genetic research and assist policymakers in developing more applicable strategies to strengthen production and conserve future stocks.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author/s.

ETHICS STATEMENT
The animal study was reviewed and approved by the Animal Welfare and Bangladesh Agricultural University Ethics Committee (Ref. no. BAURES/ESRC/Fish/10).

AUTHOR CONTRIBUTIONS
AA: conceptualization, methodology, investigation, writingoriginal draft, visualization, and fund acquisition. NH: formal analysis, data curation, writing-original draft, and visualization. MH: conceptualization, validation, writing-review, and editing. AB: methodology, data arrangement, and writing-original draft. MAh: sample collection, methodology, validation, and data analysis. MAl: sample collection, writing-review, and editing and data analysis. MAs: writing-review, and editing and data analysis. MB: funding acquisition and writing-review, and editing. YM: funding acquisition and writing-review, and editing. All authors contributed to the article and approved the submitted version.

ACKNOWLEDGMENTS
Authors would like to acknowledge Bangladesh Fisheries Research Institute (BFRI) for the financial support to carry out the research.