Introduction

Sundarbans, the largest deltaic stretch of mangrove forest in the world, is formed at the estuarine phase of the Ganges—Brahmaputra river system (Chaudhuri et al. 2012). The area of Sundarbans is characterized by a low flat alluvial plain covered with mangrove swamps and marshes, and it is crisscrossed by a large number of tidal rivers, estuaries, creeks and saline water courses (Bagchi 1972). The biodiversity of Indian Sundarbans is extremely rich and also consists of extremely rich coastal, estuarine and marine fisheries. Taking this into consideration, it has been declared as world heritage site by United Nations Educational, Scientific and Cultural Organization (1987) which shall go a long way to protect and conserve biodiversity (Sreelekshmi et al. 2020). Recently, in February 2019, the Indian Sundarbans has also been designated as ‘Wetland of International Importance’ (Ramsar site) (Sreelekshmi et al. 2020). Gopal and Chauhan (2006) stated the extensive spatiotemporal changes in hydrological regimes, topography and texture of the substratum, salinity and their interactions cause habitat heterogeneity in the mangrove ecosystems which affects the aquatic vegetation and fauna directly, and also causing damage especially to the seabirds. The author also stated that there is continuous threat to biodiversity due to the burgeoning human population, pollution, growing industrialization, oil spills, oil exploration, global climate change and sea level rise.

Tidal creeks are extremely abundant in Sundarbans estuarine systems which serve as feeder rivers, providing food control, storm water drainage and habitat to wildlife. It is ecologically important for the transfer of materials from terrestrial to marine biome, habitat to nursery ground for commercial fisheries (Mallin and Lewitus 2004; Saliu and Ekpo 2006). There has been a considerable loss of creeks due to overuse, pollution, diversion, land filling, etc., and Indian Sundarbans is also not an exception thereof. The nature and distribution of biodiversity are primarily ruled by variations in the physicochemical properties of water (Murugan and Ayyakkannu 1991). It is also responsible for supporting aquatic life and maintaining the aquatic environment and is likely to have an impact on the growth of the organisms as well as the productivity of the ecosystems (Das 2012). The physicochemical properties of the estuarine environment are highly dynamic due to the influx of fresh water. Surface waters are the most vulnerable and are exposed to rapid pollution due to their easy accessibility both the natural and anthropogenic causes (Singh et al. 2004; Barakat et al. 2016). These creeks are also acting as the probable channels for carrying pollutants (municipal and industrial wastes run off from agricultural land) to the estuarine core. Barletta et al. (2019) opined that the estuaries are often found to have the hard-hit zones with respect to the heavy metals and pesticides, and usually it has a high residence time for pollutants. In addition, the seasonal variation of precipitation, surface runoff, is a seasonal phenomenon which also supports to constitute the polluting source, largely affected by climate. Since estuaries and the supported creeks are highly fecund ecosystems and provides congenial home for the inherent biota, it is necessary to prevent and control the sources of pollution and to have reliable information on quality of water for effective management. Therefore, it is imperative to understand the changes in water quality variables from time to time for the quality assessment of the ecosystem (Chang 2008). The gradual increase in human activities for developing coastal areas also have resulted in loss of ecological values of the tidal creeks (Vernberg and Vernberg 2001). The large variability in hydrological parameters in coastal water influences, topographic heterogeneity and their interactions, has resulted in a biodiversity explosion in the Sundarbans mangrove eco-region (Ansari et al. 2017). Although such creek ecosystems are highly vulnerable to environmental changes and anthropogenic interferences, yet their ecological significance is underestimated as very less studies on these systems have been so far conducted relative to larger, better-known estuarine systems (Mallin and Lewitus 2004). Therefore, the present study evaluated the water quality status of four creeks (Costala, Gangasagar, Chemaguri and Mooriganga) in the Sundarban eco-region. Reports are available on fluctuations of physicochemical parameters in various parts of Sundarbans (Sarkar et al. 1986; Mitra et al. 2009; Manna et al. 2010; Das 2012; Chaudhury et al. 2012; Rakshit et al. 2015); however, the studies on spatiotemporal variations of water quality in tidal creeks, at Sagar Island are still scarce. The island had been affected by natural coastal processes, storm waves and artificial constructions such as seawalls and jetties, etc. (Mondal et al. 2017). The effects of climate change are also well documented in and around Sagar Island (Mitra et al. 2009).

The multivariate statistical approaches have become popular in recent times due to their ability to treat large volume of spatial and temporal data from variety of monitoring sites (Barakat et al. 2016). Different multivariate statistical techniques such as cluster analysis, principal component analysis, factor analysis and discriminant analysis are used to interpret the complex data matrices for a better understanding of the water quality and ecological status of the studied ecosystems (Vega et al. 1998; Noorie et al. 2010; Wang et al. 2014). It also allows the identification of possible factors that have an impact on water systems and provides an essential tool for reliable management of water resources (Wunderlin et al. 2001; Reghunath et al. 2002; Simeonova et al. 2003). Helena et al. (2000) stated that multivariate statistical techniques facilitate to characterize and evaluation of surface water quality, which will be useful further in verifying spatial and temporal patterns linked to seasonality. Wang et al. (2007) evaluated the spatiotemporal variations of water quality of 19 rivers by applying the cluster analysis and factor analysis. Similarly, Kumarasamy et al. (2014) and Khan et al. (2016) applied cluster analysis and factor analysis to investigate the variability of hydrochemistry from Tamrapani river, Southern India and Ramganga river (Ganga basin), respectively. The correlation analysis, principal component analysis and canonical correlation components were used by Noori et al. (2010) to study the seasonal variation of water variables, relationship between physical and chemical parameters and selection of principal and non-principal monitoring stations in the Karoon river basin, Iran. Sharma et al. (2015) also applied former multivariate techniques to identify the potential sources of pollution and clustering of monitoring stations of Ganga and Yamuna River in Uttarakhand, India. Furthermore, usefulness of cluster analysis, principal component analysis and factor analysis for interpreting the complex data sets, identifying the attributes responsible for pollution as well as water quality management previously reported by several workers (Shreshtha and Kazma 2007; Noori et al. 2010, 2012; Wang et al. 2012, 2014; López‐López et al. 2014; Bostanmaneshrad et al. 2018). In the recent times, several authors have evaluated the underlying relationship between macro-scale parameters (land use, population density, geology, erosion, etc.) with the micro-scale attributes, i.e., water quality parameters (Bostanmaneshrad et al. 2018; Chimwanza et al. 2014; Liu et al. 2016). Bostanmaneshrad et al. (2018) quantified the significant relationship between land use and phosphorus, total solids and turbidity, erosion levels and electrical conductivity, and erosion and total solids in their study from Siminehrood river basin. All the previous studies showed that the multivariate statistical methods are important tools to determine the relationships between the water quality parameters, and identify particular attributes responsible for the relationship. With this view, the present study is an attempt to understand the present ecological health and status of tidal creeks in the mangrove ecosystems of Sagar Island by approaching multivariate techniques. The basic objective of the study is to evaluate the spatiotemporal variability in physicochemical parameters of surface water samples collected from the selected tidal creeks.

Study area

Sagar Island (between latitude 21º 37′ and 21º 52ʹ N and longitude 88º 03′ and 88º 11′E), part of the Hooghly-Matlah estuary is subjected to intense tidal and wave actions (Majumdar et al. 2002). The island has unique mangrove vegetation and surrounded by the river Hooghly on the north and the west, Mooriganga on the eastern side, while the southern part faces the Bay of Bengal. The island is inundated by various natural tidal creeks and artificial canals which are a source of saline surface water in the island (Chakraborty 1995). The island experiences a semidiurnal tidal regime with high tide zone ranges from 5 to 6 m (Paul and Bandhopadhay 1987). The island is also highly vulnerable to environmental degradations, burgeoning population and erosion which accounted for approximately 30 km2 net land loss in the last three decades (Hazra et al. 2002). Among the tidal creek system at Sagar, Chemaguri creek is the most prominent one with planted mangroves along its length (Mitra et al. 2009). The present study was carried out from September 2015 to August 2016 in 4 (four) major creeks viz. Costala (21º 48ʹ 829ʺ N 88º 09ʹ 93ʺ E), Gangasagar (21º 41ʹ 222ʺ N 88º 02ʹ 80ʺ E), Chemaguri (21º 40ʹ 745ʺ N 88º 07ʹ 695ʺ E) and Mooriganga (21º 48ʹ 828ʺ N 88º 09ʹ 93ʺ E) (herein referred as S1, S2, S3 and S4, respectively) at Sagar Island, Indian Sundarbans. The sampling sites were designed to cover a wide range of determinants of key sites, which reasonably represent the water quality of the creeks in the Island. The sites were selected mainly based on their mangrove vegetation covers and also ranging from region with high riverine influence (Costala and Mooriganga) to brackish water regions and mangrove vegetation (Chemaguri and Gangasagar). The geographical locations of the sampling sites and study area in Sagar Island, Sundarbans are shown in Fig. 1.

Fig. 1
figure 1

The map showing the geographical locations of the sampling sites and study area in Sagar Island, Sundarban

Full size image

Sampling procedure

Surface water samples were collected in pre-cleaned acid washed polyethylene bottles (1.0 L) during low tide before noon from the selected locations. Water temperature (WT) was measured on field by using a mercury thermometer (0–100 °C); pH with a digital pH meter (Model No. 101 E), and turbidity (Nephelometric turbidity unit) was measured by a turbidity meter (Model No. EI 331E). The in situ estimation of dissolved oxygen (DO) and salinity (ppt) was done by the titrimetric method (APHA 2005). Collected water samples were brought to the laboratory in cold conditions for nutrient analysis. Nutrients such as nitrate, phosphate and silicate were analyzed following standard methods described in APHA (2005). For analysis of chemical oxygen demand (COD) and biochemical oxygen demand (BOD), water samples were collected separately in 500 ml polyethylene bottles. COD was measured by the open reflux method, and BOD was estimated by a 3-days BOD incubation (BOD3) (APHA 2005). All the methods were standardized as per ambient conditions and blank measurements were taken into consideration for the estimation procedures. HPLC grade water was used throughout the laboratory works. Since determination of BOD (BOD3 and/or BOD5) is a long-time estimation process, several attempts have been made to shorten the process by applying alternate analytical methods (Noori et al. 2011, 2013). The analytical model such as reduced-order neural network model (RONNM) and reduced-order adaptive neuro-fuzzy inference system (ROANFIS) model fed with input vectors [latitude (N), longitude (E), suspended solids, water temperature, river discharge, DO, electrical conductivity, nitrate and total phosphorus] is developed which indicated the best fitted online predictions of BOD5 (Noori et al. 2011, 2013).

Data treatment

Statistical computations such as one-way ANOVA, Pearson correlation, clustering dendrogram and factor analysis were performed with the help of software applications (SPSS v. 23 and XLSTAT v. 2018). The environmental data were normalized by transforming log(x + 1) prior to analysis. The water quality parameters of different seasons were subjected to one-way analysis of variance (ANOVA) and Duncan's multiple range tests (DMRT) post-hoc test using SPSS v.23. The post-hoc test was carried out to find out the significant variations of measured water variables between seasons. Pearson correlation (2-tailed) was conducted to comprehend the interdependence among the physicochemical parameters. This minimizes the effect of between stations correlations and between sampling campaigns relationships (Barakat et al. 2016). Further, the correlation coefficient − 1 or 1 indicated the strongly negative or positive relationship between the two variables, whereas its value closer to 0 defined no linear relationship between them at 0.05 significance level (Kumar et al. 2006; Barakat et al. 2016). Hierarchical agglomerative cluster analysis was performed to understand a measure of similarity using Euclidean distances (Zhao et al. 2011; Singh et al. 2004; Sharma 1996). The hierarchical cluster analysis is the most common approach where the clusters are formed based on the similar nature attributes and forming higher clusters step-by-step (Willet 1987). Generally, the Euclidean distance portrayed the proximity between two samples and ‘distance’ which is represented by the ‘difference’ between the obtained values from both the samples (Otto 1998), and it is illustrated through a dendogram (Forina et al. 2002). Further, non-metric dimensional scaling (NMDS) was also performed to visualize the dimensions to understand the similarity or dissimilarities among data sets. Likewise, factor analysis was performed (XLSTAT v.18) on the water quality data to find out the factors related to influence each of them. The first step of factor analysis is to standardize the raw data and make them dimensionless. The correlation coefficient matrix, eigen values and eigen vectors were determined to yield the covariance matrix (Wang et al. 2014). The factors with eigen values that exceed 0.5 were considered in this study, whereas only three factors obtained that exceed one. Varimax rotation is an important second step in factor analysis used for maximizing the sum of the variance of the squared loadings, where ‘loadings’ means correlations between variables and factors. The varimax rotation on the factor loading infers the rotated principal factors (Singh et al. 2004, 2012; Noori et al. 2010; Ghahremanzadeh et al. 2018).

Results and discussion

The variation of water quality parameters in different seasons (pre-monsoon, monsoon and post-monsoon) is documented in Table 1. One-way ANOVA showed that, out of ten parameters, five parameters (WT, salinity, DO, turbidity and COD) varied significantly between the seasons (p ≤ 0.05). The water remained neutral to alkaline throughout the study period well within the predicted range between 7.36 and 8.16. The result is in agreement with the findings of Manna et al. (2010) and Choudhury et al. (2012) from Sundarban waters. The rate of photosynthesis process, which removes CO2 from water through bicarbonate degradation, freshwater influx and decomposition of organic matter, could have contributed to the fluctuation of pH (Rahman et al. 2013). The alkaline nature of water (pH > 7.0) with low variations between the sampling stations, suggested that the water mass remained well buffered throughout the study period, and it indicated the presence of biodegradable organic matter in the water column. Surface water temperature varied significantly (p ≤ 0.05) from post-monsoon to pre-monsoon. Water temperature of the tidal creeks showed a forecasted range of variability with lowest in post-monsoon in the month of January (19.34 °C), and had a maximum in pre-monsoon in the month of June (31.13 °C). In tidal influenced estuarine environment of the Sundarbans, temperature values undergo a wide diurnal and seasonal variation. A gradual increase in salinity was observed from post-monsoon to pre-monsoon, recorded maximum in May (28.08 ppt) and lowest in October (0.79 ppt). The distinct temporal variation of salinity in the present study was similar to Saravanakumar et al. (2008); Choudhury et al. (2012) and Gogoi et al. (2019). Saravankumar et al. (2008) stated that low salinity level in brackish water habitats such as backwaters, estuaries and mangrove waters is due to influx of freshwater from land runoff caused by precipitations or by tidal variations, which supported our present findings. Manna et al. (2010) also reported similar findings with the lower salinity during winter months (November-February; avg. 16.7 PSU) as compared to summer months (March to May; 23.5 PSU) from Bara Herobhanga Khal (Jharkhali waters), Sundarbans. Moderate concentrations of dissolved oxygen were recorded across the stations with its maximum in October at the station Chemaguri (7.8 mgl−1) and minimum in December at Mooriganga (3.29 mgl−1). Dissolved oxygen and water temperature were influenced to a greater extent by atmospheric conditions and sampling times (Debels et al. 2005); as a result, its concentrations vary during various moments of the day and under different weather conditions. Lower concentrations of DO during post-monsoon may be due to increased salinity, temperature and biological activity (Levinton 2001). Our observation is in conformity with Sarkar and Bhattacharya (2010) where the authors stated that dissolved oxygen concentration in the surface water of Sundarbans mangrove wetland ranged from 5.18 to 6.49 mgl−1 with an average value 4.0 mgl−1, indicating that surface waters are moderately oxygenated. Turbidity was recorded in the expected level with its peak during monsoon which might be due to riverine influx bearing excess suspended load coupled with high monsoon precipitations. The steady expansion of waterway transportation/vessel trafficking was also one of the reasons for the enhancement of turbidity at the stations under study. The mean value of turbidity was highest at Chemaguri in the month of November and the lowest at Gangasagar in December. The observed range of COD was 8.0–40.0 mgl−1 across the stations with an average value 19.60 ± 7.0 mgl−1, which exceeds the permissible range of World Health Organization, i.e., < 10 mgl−1. The precipitation causes dilution during monsoon which could be the reason for lower concentration of COD during monsoon. The moderate concentration of COD in the present study was indicating toward pollution of the aquatic ecosystem due to natural and anthropogenic stresses. Bhattacharya et al. (2015) reported high range of COD values (101.28–111.31 mgl−1) in their studies from Jambu Island, Sundarbans and stated that it is derived from the primary production of the dense mangrove forests. High COD values in a tropical coastal wetland in Southern Mexico were also previously reported by Hernández-Romero et al. (2004), and it was associated with mangrove-enriched organic matters which supported the present findings. The BOD3 ranges between 1.62 and 2.23 mgl−1 with an average 1.96 mgl−1, where maximum and minimum were recorded at Gangasagar and Mooriganga, respectively. An average value of 2.04 mgl−1 BOD was reported previously by Rahman et al. (2015) from Sundarban waters, Bangladesh which is in line with our present observations.

Table 1 Seasonal variation of water quality parameters in the study sites of Sundarbans
Full size table

Non-uniform seasonal trend for nutrients was observed throughout the study period. In general, nutrients are sufficient to support primary productivity in an estuarine environment. Tidal fluctuation is the main driving force for the distribution of nutrients in Indian Sundarbans, where its dispersion primarily depends on the tidal motions of the water body (Chaudhuri et al. 2012). No significant (p ≥ 0.05) variations of nutrients (phosphate, nitrate and silicate) were observed across the seasons. Phosphate concentration was estimated to be highest at station Mooriganga (11.7 mgl−1) in the month of October while the lowest at Gangasagar (0.4 mgl−1) in March. Lower concentration of phosphate during pre-monsoon season might be due to its utilization by photoautotrophs, and buffering actions of sediments in various environmental conditions (Rajasegar 2003). Silicate and nitrate showed higher concentration in monsoon as compared to post-monsoon and pre-monsoon during the study period. Higher values of nitrate and silicate in monsoon may be attributed to freshwater influx and land runoff (Sathpathy et al. 2009; Ramakrishnan et al. 1999). Similarly, Karuppasamy and Perumal (2000) stated that higher magnitude of nitrate during monsoon might be caused by the organic matter being enriched by the monsoonal flow and decomposition of terrestrial runoff from the catchment areas. Sarkar and Bhattacharya (2010) also opined that land-based nutrients, especially from the agricultural fields in proximity, contributed to a greater extent to the enrichment of nutrient concentrations during monsoon. The inverse trend was observed with regard to phosphate concentrations that recorded the maximum during post-monsoon. Silicate values were found to be relatively higher as compared to the nitrate and phosphate in the present study, which is in line with Vajravelu et al. (2018). Maximum concentration of nitrate and silicate was observed in the month of October (3.20 mgl−1; at Chemaguri) and August (31.7 mgl−1 at Costala), and lowest in the month of April (0.62 mgl−1 at Gangasagar) and March (4.24 mgl−1 at Chemaguri), respectively. The high silicate and phosphate ratio (Si:P) in the present study is corroborated with the findings of Chaudhuri et al. (2012) and Sarkar and Bhattacharya (2010). However, the Si:P ratio was somewhat lower than ‘modified’ Redfield ratio (15:1) (Brzezinski 1958), and our present observations in seasonal estimates of silicate indicated broad seasonal fluctuations in the creek environment. Similarly, Choudhury and Bhadury (2015) reported that seasonal estimates of N:P ratio remained below the Redfield ratio (16:1), indicated nutrient limited (nitrogen) environment in Sagar island which is in line with our present observation from Sagar Island. Arumugum et al. (2016) also showed low N:P ratio than the Redfield ratio in Muthupet mangrove waters depicted a low bio-availability of nitrogen for primary production which confirms to our estimates. Redfield (1958) approached to a conclusion that atomic ratios of elements in the biochemical cycle of plankton were statistically uniform and follow the stoichiometric ratio of C:N:P = 106:16:1 and this variation is mainly due to synthesis or decomposition of organic matter. The estimated ratio of nitrate, phosphate and silicate is shown in Table 2.

Table 2 Nutrient parameters: nitrate, phosphate and silicate (N:P:Si) ratio in the study sites of Sundarbans
Full size table

Intra-relationship of various environmental variables (Karl Pearson’s correlation coefficient) showed significant positive correlation of water temperature with salinity (r = 0.42; p ≤ 0.01), DO (r = 0.55; p ≤ 0.01), silicate (r = 0.40; p ≤ 0.01); COD (r = 0.38, p ≤ 0.01), and it was negatively correlated with phosphate (r = -0.11). According to Osibanjo et al. (2011) and Şener et al. (2017), one very common but essential parameter, pH governs various other parameters like alkalinity, solubility and also the hardness. The pH had significant positive correlation with DO (r = 0.31; p ≤ 0.05). Similar observation was also reported by Gogoi et al. (2019) in their study from Kailash Khal, Sundarbans; depicted that pH had a positive correlation with water temperature and also salinity with the water temperature and DO. The strong positive correlation was observed between salinity and COD (r = 0.39 p ≤ 0.01), and nutrients such as nitrate (r = 0.35 p ≤ 0.01) and silicate (r = 0.41 p ≤ 0.01). Gogoi et al. (2020) also found similar observation in their study from Sundarban waters, where salinity had a positive correlation with the nutrient parameters (phosphate, silicate and sulfate).

Multivariate approach

Cluster analysis was performed season-wise to understand the percentage similarity among various stations. The hierarchical group average showed the variations in water parameters by clustering similar nature of samples/sampling stations. It formed two clusters with the maximum similarity between the stations Chemaguri and Costala during pre-monsoon, while showed maximum dissimilarity between the station Chemaguri during monsoon and Gangasagar during pre-monsoon (Fig. 2). The freshwater influx with its high allochthonous inputs during the rainy season also could be a reason for variations of water attributes during monsoon. Further, NMDS also showed the similar nature of ordinance pattern among samples/sampling stations that were observed in cluster analysis (Fig. 3). The stress value was found less than 0.14 which is good ordinance pattern that distance among items/samples are perfect. In this study, water parameters such as dissolved oxygen and turbidity were found to be strong parameters influencing their variables in the system. The other important parameters were found to be COD, BOD3 and phosphate concentrations during the study period (Table 3).

Fig. 2
figure 2

Dendogram showing the similarity of stations in seasons

Full size image
Fig. 3
figure 3

NMDS showing the environmental variables recorded in various stations and seasons (1: pre-monsoon; 2: monsoon; 3: post-monsoon)

Full size image
Table 3 Intra-relationship of environmental variables in Sundarbans
Full size table

Results of the factor analysis showed, four factors were extracted which explained 53.21% of the total variance in the observed data (Table 4 and Fig. 4). The eigen value and cumulative percentage variability of the extracted factor from the data set are shown in the scree plot (Fig. 5). Factor 1 accounts for 22.62% variance where the parameters (loadings) such as pH, WT, salinity, DO, silicate and COD were positively correlated, while the attributes viz., turbidity, phosphate and nitrate were negatively correlated. Similar observations also reported by Sarkar and Bhattacharya (2010), where they stated that, the factors appear (negative correlation of nitrate, phosphate and BOD) to be originated from the combined effect from the anthropogenic activities. Factor 2 explained the 15.88% total variance, and it was negatively loaded with variables; DO, turbidity and nitrate. The causes may be due to gradual increase in water transportation/fishing vessels and agriculture expansions. Likewise, factor 4 represents the 6.68% total variance. The factor loading revealed that parameters such as phosphate, COD and BOD3 were negatively correlated with each other. Since these causes of these parameters (phosphate and COD) are primarily from the source of agricultural land and industrial activities. BOD represents the high organic load, which enriching the growth of microorganisms in the mangrove patches. From this analysis, it is reflected that the major causes of the water quality variations/degradation seemed to be more usage of nitrogenous and phosphate-rich fertilizers, discharge of untreated industrial effluents, increasing waterway transportation which accelerates the oil spills in river/estuary and also the other non-point sources. With the best uses of Sundarban waters in different aspects such as fishing, irrigation, ferry service, boating and tourism, some water quality parameters are still well within the permissible limits, where it is seemed to be suitable for bathing, aquatic biodiversity, fisheries and recreation, etc., according to the criterion placed from the various best designated use prescribed by the World Health Organization, US Environmental Protection Agency and Pollution Control Board, India (Table 5).

Table 4 Factor analysis of the water quality parameters after varimax rotation
Full size table
Fig. 4
figure 4

Scree plot of the factor analysis after varimax rotation

Full size image
Fig. 5
figure 5

Factor loading of water quality parameters after varimax rotation

Full size image
Table 5 Water quality standards for aquatic biodiversity propagation
Full size table

Conclusion

The present study reflected the changes in some important hydrological factors round the year which hinted toward the dynamic nature of the estuarine ecosystem due to climate change and human-induced impacts. Salinity and turbidity showed the maximum annual range of variation followed by dissolved oxygen in the present study. Other parameters did not show any significant annual range of variation in the present study, though phosphate and COD values were higher in post-monsoon and pre-monsoon seasons, respectively. As there were no significant industrial activities in the near vicinity; therefore, the higher COD might be due to the presence of mangroves at sampling sites. Water variables such as dissolved oxygen, water temperature, salinity, turbidity and nitrate were the effective factors influencing their variables in the creeks, which was evident from factor analysis. The higher magnitude of COD in this study doesn’t seem to indicate pollution in the study area. Higher turbidity seemed to be due to monsoonal runoff and ferrying activities in the creeks. Long-term monitoring approach is required to look for time-scale changes to comprehend the climatic emergencies and anthropogenic stress. Also, studies on anthropogenic interferences such as pollution, coastal development, destructive fishing, ferrying, overexploitation, etc., resulted in unique information that may provide mitigating measures to conserve biodiversity. Water quality criteria for wildlife propagation and fisheries need to be addressed for the sustenance of endemic biota inhabiting in the Sundarban Biosphere Reserve, which is also opined by many previous workers. Although the creek waters (Sundarbans) are impacted by various human-induced activities (such as fishing, irrigation, ferry service, boating and tourism), some water quality parameters are still well within the permissible limits, where it is suitable for bathing, aquatic biodiversity, fisheries and recreation purpose.