Analyses on the Temporal and Spatial Characteristics of Water Quality in a Seagoing River Using Multivariate Statistical Techniques: A Case Study in the Duliujian River, China

In the Duliujian River, 12 water environmental parameters corresponding to 45 sampling sites were analyzed over four seasons. With a statistics test (Spearman correlation coefficient) and multivariate statistical methods, including cluster analysis (CA) and principal components analysis (PCA), the river water quality temporal and spatial patterns were analyzed to evaluate the pollution status and identify the potential pollution sources along the river. CA and PCA results on spatial scale revealed that the upstream was slightly polluted by domestic sewage, while the upper-middle reach was highly polluted due to the sewage from feed mills, furniture and pharmaceutical factories. The middle-lower reach, moderately polluted by sewage from textile, pharmaceutical, petroleum and oil refinery factories as well as fisheries and livestock activities, demonstrated the water purification role of wetland reserves. Seawater intrusion caused serious water pollution in the estuary. Through temporal CA, the four seasons were grouped into three clusters consistent with the hydrological mean, high and low flow periods. The temporal PCA results suggested that nutrient control was the primary task in mean flow period and the monitoring of effluents from feed mills, petrochemical and pharmaceutical factories is more important in the high flow period, while the wastewater from domestic and livestock should be monitored carefully in low flow periods. The results may provide some guidance or inspiration for environmental management.


Introduction
River ecosystems have received more attention in recent years because they not only provide water resources for socioeconomic development, but they also play an important role in ecological integrity. The surface water quality in rivers plays an important role in the productivity and life of surrounding human societies. However, water quality deterioration in rivers has been extremely serious in developing countries due to the intensive anthropogenic activities along the river and lack of water treatment in the entire basin [1][2][3]. Nowadays, abundant literature focusing on the water quality in inland rivers or coastal seas can be easily found, while there is little reported information on seagoing rivers due to the complexity of their ecosystems and diversity of pollution sources [2,4]. Due to the transition from fresh water to sea water, seagoing rivers not only have the pollution As the main water source supplying two important coastal wetland ecological environment conservation areas, namely the Beidagang Wetland Nature Reserve and the Tuanbo Birds and Nature reserve, the Duliujian River plays an important role in ecological environment protection. However, the Duliujian River also receives a large number of point and non-point sources pollutants throughout the year since that the up-middle reaches goes through the urban population living quarters and the downstream reach crosses some industrial parks that contain many textile, chemical, pharmaceutical, petroleum chemical, oil refinery, iron-steel and electroplating factories [20,28]. Furthermore, the wastewater from the aquaculture along the river and the surrounding areas also aggravates the pollution status of the Duliujian River.

Sampling and Analysis
In order to monitor the temporal and spatial changes of the water quality in the Duliujian River, 15 sampling cross-sections (S1-S15) were set along the whole Duliujian River considering natural conditions and human activities ( Figure 1). Further, three sampling sites on each sampling cross-section were determined at 1/4, 1/2 and 3/4 across the river width. However, it should be stated that the sites S1.1-S2.3 on S1 and S2 in spring (May 2017) were not available due to dry up. Therefore, seasonal water samples from 45 sampling sites were collected in winter (December 2016, based on S1-S15), summer (July 2017, based on S1-S15) and autumn (September 2017, based on S1-S15), and 39 in spring (May 2017, based on S3-S15), respectively. During each sampling period, the work was carried out from S1 to S15 in sequence. In addition, weather forecasts for at least 15 days were checked before each sampling, and then seasonal sampling was completed in the consecutive days that were most suitable for sampling, in order to reduce the measurement errors caused by weather. At each sampling site, some water quality parameters were measured in situ like water temperature (WT), water depth (WD), and others, and at the same time, 1.5 L water was collected 1.0 m (1/2 of WD when WD < 1 m) below the water surface using a polymethylmethacrylate sampler and stored in a polyethylene plastic bottle. These sampling bottles were precleaned with detergent, rinsed three times with deionized water and dried under 30 • C in an oven for 24 h. Then, the water samples were placed in insulated boxes containing ice bags and transported to the laboratory for further analyses. The preservation and handling of the samples were all performed according to the Water Quality Sampling-Technical Regulation of the Preservation and Handing of Samples (HJ493-2009) norms [33]. Due to the normally closed sluices along the river, the average flow velocity was very slow. As such, in comparison with the nearest distance between two adjacent sections (3.0 km), the effects of pollutant transfer with water flow should be very slight and therefore ignored. reserve, the Duliujian River plays an important role in ecological environment protection. However, the Duliujian River also receives a large number of point and non-point sources pollutants throughout the year since that the up-middle reaches goes through the urban population living quarters and the downstream reach crosses some industrial parks that contain many textile, chemical, pharmaceutical, petroleum chemical, oil refinery, iron-steel and electroplating factories [20,28]. Furthermore, the wastewater from the aquaculture along the river and the surrounding areas also aggravates the pollution status of the Duliujian River.

Sampling and Analysis
In order to monitor the temporal and spatial changes of the water quality in the Duliujian River, 15 sampling cross-sections (S1-S15) were set along the whole Duliujian River considering natural conditions and human activities ( Figure 1). Further, three sampling sites on each sampling crosssection were determined at 1/4, 1/2 and 3/4 across the river width. However, it should be stated that the sites S1.1-S2.3 on S1 and S2 in spring (May 2017) were not available due to dry up. Therefore, seasonal water samples from 45 sampling sites were collected in winter (December 2016, based on S1-S15), summer (July 2017, based on S1-S15) and autumn (September 2017, based on S1-S15), and 39 in spring (May 2017, based on S3-S15), respectively. During each sampling period, the work was carried out from S1 to S15 in sequence. In addition, weather forecasts for at least 15 days were checked before each sampling, and then seasonal sampling was completed in the consecutive days that were most suitable for sampling, in order to reduce the measurement errors caused by weather. At each sampling site, some water quality parameters were measured in situ like water temperature (WT), water depth (WD), and others, and at the same time, 1.5 L water was collected 1.0 m (1/2 of WD when WD < 1 m) below the water surface using a polymethylmethacrylate sampler and stored in a polyethylene plastic bottle. These sampling bottles were precleaned with detergent, rinsed three times with deionized water and dried under 30 °C in an oven for 24 h. Then, the water samples were placed in insulated boxes containing ice bags and transported to the laboratory for further analyses. The preservation and handling of the samples were all performed according to the Water Quality Sampling-Technical Regulation of the Preservation and Handing of Samples (HJ493-2009) norms [33]. Due to the normally closed sluices along the river, the average flow velocity was very slow. As such, in comparison with the nearest distance between two adjacent sections (3.0 km), the effects of pollutant transfer with water flow should be very slight and therefore ignored. In Table 1, the monitored water quality parameters and their abbreviations, analytical methods and units are summarized. All parameters were measured or analyzed in accordance with relevant In Table 1, the monitored water quality parameters and their abbreviations, analytical methods and units are summarized. All parameters were measured or analyzed in accordance with relevant standard methods [34]. After filtering through 0.7 µm glass fiber filter and extraction for 12 h using 90% (v/v) acetone by a freeze-thaw method, the absorbance of the water samples was measured using a UV spectrophotometer, by which, the concentration of Chla could be calculated according to a relevant standard [35]. NH 4 + -N, NO 3 − -N, NO 2 − -N and PO 4 3− -P were determined from adequate samples filtered through 0.45 µm cellulose ester membranes. DIN was employed to represent the concentration of inorganic nitrogen and it was calculated using Equation (1): In order to ensure the quality of the experimental data, quality control procedures were performed throughout the sampling and analysis processes. In addition, two replicates of the parameters tested in laboratory for each sampling site were measured to minimize the experimental error and enhance the reliability of the measured results. Key statistics of these water quality parameters measured in this study are listed in Table 2.

Data Processing and Statistical Analyses
Prior to statistical analyses, the missing values were filled with the mean values of each quarter sequence. Then, logarithmic transformations were performed to standardize the environment variables so that the data could meet the assumptions of normality and homoscedasticity for statistical analyses (except for pH). Logarithmic transformation is a kind of standardization which can eliminate the influence of different units and scales of the measurements, make the data dimensionless and tends to reduce the effect of outliers [15,36]. The logarithmic transformation equation can be expressed as: S ij = Log 10 (x ij + 1), I = 1, . . . , 12, j = 1, . . . , 45, where x ij is the i-th variable of the j-th sample; S ij is the standardized value of x ij . Whether the data obeyed a normal distribution could be determined by the standardized skewness and standardized kurtosis. The environmental variables, whose values of standardized skewness and standardized kurtosis were outside the range of −2 to +2, were considered to correspond to a non-normal distribution [23]. With this method, only 12 parameters with proper distribution were selected although more than 20 water environmental parameters were available. For example, chemical oxygen demand (COD) was not chosen for further analyses due to its low confidence in the log-normal distribution test. Spearman correlations (a statistical test) were employed to reveal the relationships among the water quality parameters. A correlation coefficient near −1 or 1 means a strong negative or positive linear relationship between two variables, while one closes to 0 demonstrates a weak relationship. In this study, the relationships with correlation coefficients greater than 0.5 at the significant level of p < 0.01 were chosen for further detailed analysis.
Cluster analysis (CA) is a multivariate statistical technique that can investigate the similarities and dissimilarities between sampling sites, which groups the relevant objects into different aggregations based on their interdependent variables or characteristics [29,37]. Hierarchical agglomerative cluster analysis (HACA), the most common approach used in CA, was chosen to provide intuitive similarity relationships between temporal and spatial patterns [2,23]. HACA analysis was employed on the standardized dataset by the Ward's method using squared Euclidean distances without any prior knowledge of the number of clusters [22,23]. The cases in each group clustered by this method have the smallest sum of squares deviation, and the clustering processes are presented in a dendrogram. The categories were based on the linkage distance D link /D max × 100, where D link is the distance between two clusters or objects, and D max is the maximum D link .
One-way analysis of variance (ANOVA) was performed to compare the differences of water quality variables in different seasons or sampling sites (p < 0.05; least-significance difference, LSD of assuming homogeneity of variance; Tamhane's T2 of non-homogeneity of variance). ANOVA results were employed to show the reliability of CA results.
Principal component analysis (PCA) is another multivariate statistical technique that can reduce dimensionality of the original data set to extract the information about correlation among variables analyzed in the water samples with a minimum loss of the original information [30,38]. PCA was carried out to explore the main parameters that can explain the data variance in the groups divided by CA. The Kaiser-Meyer-Olkin (KMO) and Barlett's tests were performed to test the suitability of the data. All the data calculations and statistical analyses in this study were performed by IBM SPSS 19 (SPSS Inc., Chicago, IL, USA) and Statistica 10 (Statsoft Inc, Tulsa, OK, USA).

Physicochemical Characteristics and Correlation Analyses
In order to analyze the spatial variations and the correlations of these water quality parameters monitored in this research, the values of these 12 water quality parameters were plotted in Figure 2 corresponding to four different seasons. The Spearman correlation coefficients between these 12 physicochemical parameters were also calculated and are shown in Table 3.
As can be seen in Figure 2, WT, with a range of 1.9 • C to 32.1 • C across the study period, varied due to the season. It was also noticeable that WT had significant positive correlations with EC and TDS, which may rely on the enhancement of irons activity due to increasing WT. In addition, seawater intrusion due to the opening of estuarine sluices in summer and autumn (with higher WT) may also be responsible to the much higher EC (with an annual value of 20841.1 µs/cm) and TDS (with an annual value of 13,765.4 mg/L) than that of fresh water (TDS of drinking water is no more than 1000 mg/L [39]). As such, the correlations of WT with EC and TDS observed in this research were different from some other studies due to the effects of seawater intrusion [31,40]. than 1000 mg/L [39]). As such, the correlations of WT with EC and TDS observed in this research were different from some other studies due to the effects of seawater intrusion [31,40].

Figure 2.
Spatial and temporal variability of water quality in the Duliujian River. Note: The major ticks of the x-axis are the third site at each section, and the minor ticks before each major ticks are the first and second site. i.e., the major tick number 1 means S1.3, two minor ticks before number 1 means S1.1, S1.2. The major tick labels are also corresponding to their section number.
With an average of 8.7, pH varied from 7.8 to 9.9 and indicated the slightly alkaline characteristic of water. The effects of pH on PO4 3− -P adsorption should be responsible for strong positive correlations between TP and PO4 3− -P. Under alkaline conditions, phosphorus was mainly released by ion exchange, by which, metal cations combined with phosphate are replaced by OH − therefore leading to a higher dissolved PO4 3− -P concentration [41,42]. Meanwhile, the negative charges on the surface sediment increased gradually with the increase of pH value [20], and the increased repulsion between the more negatively charged phosphate species and negatively charged surface sediment resulted in the lower adsorption of PO4 3− -P [43]. PO4 3− -P, therefore, was released into water and raised the TP concentration. The correlations between pH and TN, DIN (Table 3) indicated that the increase of pH might have a negative impact on biochemical processes related to dissolved organic nitrogen (DON) such as phytoplankton secretion, nitrogen fixation and absorption of aquatic plants [5,6,44]. Spatial and temporal variability of water quality in the Duliujian River. Note: The major ticks of the x-axis are the third site at each section, and the minor ticks before each major ticks are the first and second site. i.e., the major tick number 1 means S1.3, two minor ticks before number 1 means S1.1, S1.2. The major tick labels are also corresponding to their section number.
With an average of 8.7, pH varied from 7.8 to 9.9 and indicated the slightly alkaline characteristic of water. The effects of pH on PO 4 3− -P adsorption should be responsible for strong positive correlations between TP and PO 4 3− -P. Under alkaline conditions, phosphorus was mainly released by ion exchange, by which, metal cations combined with phosphate are replaced by OH − therefore leading to a higher dissolved PO 4 3− -P concentration [41,42]. Meanwhile, the negative charges on the surface sediment increased gradually with the increase of pH value [20], and the increased repulsion between the more negatively charged phosphate species and negatively charged surface sediment resulted in the lower adsorption of PO 4 3− -P [43]. PO 4 3− -P, therefore, was released into water and raised the TP concentration. The correlations between pH and TN, DIN (Table 3) indicated that the increase of pH might have a negative impact on biochemical processes related to dissolved organic nitrogen (DON) such as phytoplankton secretion, nitrogen fixation and absorption of aquatic plants [5,6,44]. The relationship between EC and TDS is well known and they had similar spatial trends in this study as the values increased along the river with significant increases at the last three sites in summer, autumn and winter, while a sharp growth from S6.3 (EC: 6940.0 µs/cm, TDS: 4251.0 mg/L) to S7.1 (EC: 42,128.0 µs/cm, TDS: 27,618.5 mg/L) was seen in spring. Besides, the correlation between DIN and EC might be attributed to microbial activity. SD and WD also had similar spatial variation tendencies except for spring when a rubber dam under construction led to artificial effects on WD (S6.3, S7.1). The minimum values of SD and WD were 14.2 cm (winter, S4.1) and 0.7 m (spring, S15.3), respectively. Specifically, there were no strong correlations between SD/WD and other factors.
DO, an important index in water quality, exhibited a higher concentration in winter and a lower concentration in summer with a minimum value of 0.05 mg/L at S7.1. The physical processes of water should be responsible to the moderate negative correlations between DO and EC, TDS. Furthermore, the changes of aquatic plants and plankton and therefore their production ability of oxygen may also be another important factor influencing the DO concentration.
Chla displayed a large temporal variation with a lowest concentration of 2.2 µg/L at S13.2 in spring. As an important index representing the phytoplankton biomass, Chla did not exhibit a significant correlation with TN indicating that phytoplankton were affected by other limiting environmental factors rather than TN, which required more further study.
TP and PO 4 3− -P showed significant downtrend along the river across the seasons, varying from 0.05 to 1.12 mg/L and 0.00 to 0.81 mg/L, respectively. TP and PO 4 3− -P were associated well with WT, which indicated that these correlations might be caused by sediment release or aquatic life absorption, or the frequent anthropogenic activities such as agricultural land fertilized or aquaculture waste water discharge in spring and summer.

Spatial Sites Grouping and Source Identification
CA was used to categorize the spatial similarity groups of all the monitored sites according to their water quality similarities. EC and PO 4 3− -P, highly correlated with TDS and TP, respectively, were not used in CA to avoid repeated counting of the contribution rates of similar water quality parameters. As can be seen in Figure 3, the 45 sampling sites were categorized into four clusters at (D link /D max ) × 100 < 32, and these sampling sites showed a reasonable consistency in their cross sections. Besides, ANOVA was performed to demonstrate the significance and reliability of the grouping result (Table 4). It was illustrated that the results of CA were effective for the significant differences of each water quality parameter among each cluster at p-value < 0.05 level. As shown in Table 5 and Figure 4, PCA was employed to find out the potential factors corresponding to these four clusters. In Table 5, it can be seen that four, two, four and four principal components for the water quality parameters were extracted in cluster 1, cluster 2, cluster 3 and cluster 4, respectively. In order to identify the main pollution source of each river reach, the spatial CA and PCA results can be comprehensively analyzed. These clusters were determined by judging their water quality, which was primarily influenced by land use. Sites within each cluster have similar features and pollution source types. Cluster 1 contained S1.1-S2.3 (S1-S2) which were located in the upstream part of the Duliujian River. A large number of villages and fewer factories in this area indicated that these sampling sites were polluted by domestic waste water and agricultural activities. Cluster 1 was perceived as a relatively low pollution region because, as Figure 2 shows, the concentrations of nutrient like TN, TP of the S1.1-S2.3 were relatively lower in the sampling sites, especially in comparison with the following three sections. PCA shows that PC1 including WT, DO, EC, TDS, WD, TP, PO 4 3− -P and DIN, which explained 44.75% of the total variance, was the dominant component and attributed to the soil erosion, domestic sewage and agricultural activities. Soil erosion may increase the ion concentration leading to a significant effect on DO, EC and TDS. The high loadings of TP, PO 4 3− -P and DIN might be attributed to domestic sewage and surface runoff from farmland containing a large amount of nutrients [45]. PC2 explaining 27.57% of the total variance was attributed to the physicochemical indexes related to phytoplankton growth for the moderate and high loadings of TDS, pH, WD, TN and Chla. Cluster 2 contained S15.1-S15.3 (S15) located on the last cross-section at the estuary. It mainly contributed to the exceptional EC and TDS due to the sea water intrusion (PC1, 51.44%) [17]. Because of the groundwater overpumping, land subsidence and natural geological conditions, the hydrodynamic balance between seawater and fresh water was destroyed and resulted in sea water intrusion. River and soil were therefore degenerated as a result of serious sea water intrusion.
Cluster 3 included S7.1-S14.3 (S7-S14) located on the middle and lower reaches of the investigated river. PC1, which accounts for 38.34% of the total variance, was attributed to the frequent intrusions of effluents from fisheries and livestock activities and mixing of seawater and river water. Due to the thriving fisheries [23], excessive fish feed and large amounts of fish waste containing a large number of nutrients were directly dumped into the river. As for livestock, sewage from washing farming sites and untreated poultry manure and urine were dissolved in rainfall runoff and then discharged into the river. Hence, DIN, as the dominant component of fish feed and poultry manure, was noticeable in PC1. Besides, from the textile, pharmaceutical, petroleum chemical and oil refinery factories, etc. along the river reach, a large number of industrial effluents containing acidic pollutants and phosphorus were discharged into the river, and became significant point contamination sources. As the PC2 shows, 28.10% of the total variance was explained by this component with moderate and high loadings of WT, pH, TN, TP and PO 4 3− -P. Although there are a lot of point and non-point source pollutions, this reach was grouped as moderate pollution, this must significantly rely on that all these sampling sites were located in the Beidagang wetland nature reserve where the self-purification ability and assimilative capacity were reasonably strong.    Cluster 4 comprising S3.1-S6.3 (S3-S6) was perceived as highly polluted because of the high nutrient concentration (PC1, 41.20% and PC2, 26.79%). In this region, there are many small factories producing feed, furniture and pharmaceuticals. High loading values of TP and PO4 3− -P in this reach should be attributed to the effluents from feed mills which contained a lot of nutrient, and high loading values of DO was attributed to the effluents from furniture and pharmaceutical factories containing rich organic substance (as seen in Table 6, high emission of COD), which may consume a large amount of oxygen and resulted in low-oxygen environment [1] and therefore were the mainly point pollution sources in this region. The relatively intensive residential population surrounding the up-middle reaches of the investigated river also made that the domestic wastewater became an important point pollution source. Besides, livestock activities were a notable contributor to the nonpoint pollution sources since there are many poultry farms. PC2 with high loading values of TN, Chla and DIN was attributed to point pollution sources and nutrients related to phytoplankton. As shown in Table 6, this cluster has a high emission of NH4 + -N which increased the TN and provided important nutrients for phytoplankton growth.  should be attributed to the effluents from feed mills which contained a lot of nutrient, and high loading values of DO was attributed to the effluents from furniture and pharmaceutical factories containing rich organic substance (as seen in Table 6, high emission of COD), which may consume a large amount of oxygen and resulted in low-oxygen environment [1] and therefore were the mainly point pollution sources in this region. The relatively intensive residential population surrounding the up-middle reaches of the investigated river also made that the domestic wastewater became an important point pollution source. Besides, livestock activities were a notable contributor to the non-point pollution sources since there are many poultry farms. PC2 with high loading values of TN, Chla and DIN was attributed to point pollution sources and nutrients related to phytoplankton. As shown in Table 6, this cluster has a high emission of NH 4 + -N which increased the TN and provided important nutrients for phytoplankton growth.

Temporal Periodic Grouping and Source Identification
The temporal similarity grouping results are shown in Figure 5. The 45 sampling sites in four seasons were grouped into three clusters at (D link /D max ) × 100 < 60. EC and PO 4 3− -P were not includedin CA for the same reason as 3.2. The ANOVA was performed to demonstrate the significance and reliability of the grouping result (Table 7). All parameters contributed to the CA results since the differences of these water quality parameters among three clusters were significant at the p-value < 0.05 level. Hence, the grouping results of CA were proved to be effective. It was noteworthy that the temporal samples were grouped into three different clusters rather than four (seasons). These clusters can be explained by the hydrological characteristics of mean, high and low flow periods which were divided by the mean precipitation in 1980-2010 [46], consistent with [28]. In the Duliujian River watershed, spring and autumn with the seasonal precipitations of 0.17 × 10 8 m 3 , 0.21 × 10 8 m 3 were classified as mean flow period, while summer (0.82 × 10 8 m 3 ) and winter (0.03 × 10 8 m 3 ) were the high and low flow periods, respectively. Sampling data, anthropogenic activities and natural conditions may provide some explanations for the grouping results. In addition, PCA was employed to find out the main factors for these three clusters. PCA extracted four, three and four principal components for the water quality parameters in the mean, high and low flow periods, respectively. Further, these PCs with eigenvalues larger than 1 explained 82.40%, 79.13% and 82.49% of the total variance according to the relevant data sets (Table 8, Figure 6). All the significances of Bartlett's sphericity test were below 0.001 and the corresponding results of Kaiser-Meyer-Olkin (KMO) were 0.53, 0.60 and 0.51, which indicated the reliability of the PCA results. In this study, the loadings were classified as 'strong', 'moderate', and 'weak' according to the loading values of >0.70, 0.70-0.50, and 0.50-0.30, respectively [47,48]. In detail, Cluster 1 consisted of almost all the sampling sites in May 2017 which corresponding to the mean flow period. Except that, some of the results that not divided into this cluster were mainly affected by the construction of water conservancy project. During the sampling period, a rubber dam was being built between the cross S6 (S6.1, S6.2 and S6.3) and S7 (S7.1, S7.2 and S7.3). Hence, the samplings S3.1-S6.3 were taken from a lake-like area with no water supply and discharge, which resulted in the fact that the S3.1-S6.3 with relatively lower EC, TDS and higher nutrient were grouped into Cluster 2. Besides, S7.1-S15.3 showed higher EC and TDS due to the seawater intrusion in comparison with the other three seasons. In PCA results for mean flow period (corresponding to Cluster 1), the first principal component (PC1) explained 38.96% of the total variance and exhibited strong positive loadings on EC, TDS and TN, while a strong negative loading on PO 4 3− -P and moderate negative loadings on DO, pH, TP and Chla. PC1 was mainly attributed to natural processes and nutrient pollution. Seawater intrusion brought a large number of ions such as Na+, Mg 2+ , Ca + as well as Cl − into the river, which significantly increased the values of EC and TDS [17]. The positive loading with pH was attributed to the seawater intrusion exacerbated the soil erosion, which improved the dissolved CO 3 2− and HCO 3− from rock dissolution [28,40]. Besides, nutrients like TN, TP, and PO 4 3− -P were mainly from domestic wastewater, agricultural and livestock activities [20]. PC2, which explained 24.69% of the total variance, had high negative loadings on WT and moderate negative loadings on SD and TN DIN, while high positive loadings on DIN. This component may represent the metabolism of microorganisms since DIN was associated with nitrification and denitrification closely [7]. PC3, accounting for 11.66% of the total variance showed moderate negative loading on DO and positive loading on SD. This component was regarded as the oxidation process since the effluents from furniture factories and pharmaceutical factories contained a lot of organic substance, whose oxidation process consumed a large amount of DO, and reduce the SD. PC4 was interpreted as the physical process component and may be associated with the seasonal effect of the flow. Cluster 2 consisted of almost all the sampling sites in December 2016, which corresponds to the low flow period. In this cluster, DO was relatively high and the content of TN and TP were moderate. PC1 for this cluster had strong positive loadings on EC and TDS, while strong negative loadings on DO, TN, PO 4 3− -P, and DIN and moderate loadings on WT, pH, TP and Chla. This component was attributed to natural processes (seawater intrusion), metabolic process of microorganisms (anaerobic fermentation processes [36], nitrification and denitrification) and anthropogenic activities (effluents from factories, wastewater from livestock and domestic). PC2 had strong loading on Chla, moderate loadings on WT, EC, TDS SD, and TP. This component was considered as phytoplankton growth processes. Compared to the high flow period, EC and TDS became important factors, and TP became the limited factor for phytoplankton. PC3, explaining 10.87% of the total variance with the negative loadings on DO and pH, representing the oxidation of wastewater. Cluster 3 included almost all the sampling sites in July 2017 and September 2017 corresponding to the high flow period. It can be seen from the PCA results for this cluster (corresponding to high flow period), PC1 accounting for 31.53% of the total variance, had strong negative loadings on TN, TP, PO 4 3− -P and DIN, moderate loadings on WT and WD. This component was interpreted as due to anthropogenic processes [19]. Effluents from factories, waste water from domestic, fisheries and livestock activities containing a lot of nutrients were probably the mainly reason. In this period, the average NH 4 + -N and TP emission concentrations of the 24 swage outlets along the river were 2.49 mg/L and 0.47 mg/L according to the data from Tianjin Eco-environmental monitoring center, which demonstrated high nutrients input. PC2 had strong positive loadings on EC, TDS, while moderate negative loadings on TP and PO 4 3− -P, representing the natural processes such as seawater intrusion and sediment absorption. PC3 explaining 16.40% of the total variance, had strong loadings on SD and Chla, moderate loading on TN, reflecting the phytoplankton growth process. Chl a was an index represented phytoplankton biomass, and SD and TN were the important factors for plankton growth [36]. PC4, accounted for 13.84% of the total variance, had strong positive loading on pH, moderate loadings on WT, DO and WD, reflecting the oxidation processes. Effluents from feed mills, petrochemical and pharmaceutical factories may contain carbohydrates which hydrolyze to produce acids or organic acids, thus decreasing the pH values.
and sediment absorption. PC3 explaining 16.40% of the total variance, had strong loadings on SD and Chla, moderate loading on TN, reflecting the phytoplankton growth process. Chl a was an index represented phytoplankton biomass, and SD and TN were the important factors for plankton growth [36]. PC4, accounted for 13.84% of the total variance, had strong positive loading on pH, moderate loadings on WT, DO and WD, reflecting the oxidation processes. Effluents from feed mills, petrochemical and pharmaceutical factories may contain carbohydrates which hydrolyze to produce acids or organic acids, thus decreasing the pH values.      Overall, PCA indicated the temporal variations of water pollutants. TN, TP and PO 4 3− -P, were always in PC1 which indicated that nitrogen and phosphorus were the main basic factors in all periods, should therefore more attention be paid to nutrient control and effluents from sewage outlets. For the high flow period, microbes and phytoplankton metabolized and reproduced quickly with the high nutrient concentration. According to the algal measurements in the study period, the density of Microcystis, which can reflect the condition of algal bloom, was much higher in summer (3.69 × 10 6 ind./L) and autumn (57.79 × 10 6 ind./L) than the spring (0.03 × 10 6 ind./L) and winter (0.17 × 10 6 ind./L). More attention should be paid to the potential risk of algal blooms. Significant loadings of EC and TDS in PC1 except for the high flow period indicated that a natural process (seawater intrusion) was the main basic factors influencing water quality in mean and low flow periods, and large flows can reduce the impacts of seawater intrusion.

Conclusions
In this study, the spatial and temporal variations in water quality of Duliujian River were identified using statistical techniques (Spearman correlation, CA, ANOVA and PCA). This study revealed that the water quality parameters varied reasonably with season and spatial position.
The results of spatial CA and PCA revealed that the most polluted reach was the up-middle reaches, including S3.1-S6.3. This reach was mainly polluted by nutrients from feed mills, furniture and pharmaceutical factories as well as livestock activities. The middle and lower reaches of the river, S7.1-S14.3, were moderately polluted. Though it was contaminated by complex pollutants from up-middle reaches as well as surrounding factories, it showed good self-purification ability since the river crosses some wetlands. The upstream of the river including S1.1-S2.3 was lightly polluted, since it was only polluted by domestic wastewater and agricultural activities. The estuary of the river, S15.1-S15.3, was influenced significantly by seawater intrusion.
The results of temporal CA showed that the temporal samples were grouped according to their hydrological characteristics into mean, low and high flow periods rather than the four seasons. Temporal PCA indicated that nitrogen and phosphorus were mainly the basic factors in all periods, and nutrient control and effluents from sewage outlets should therefore be paid more attention. Specially to the high flow period, more attention should be paid to the potential risk of algal blooms since that phytoplankton metabolized nutrients and reproduced quickly with high concentration nutrient. Significant loadings of EC and TDS in PC1 indicated that seawater intrusion was the mainly basic factors influencing water quality in mean and low flow periods.
Proper treatment for domestic and livestock wastes as well as polluted surface runoff should be undertaken to reduce the nutrient concentration. Management needs to be carried out to minimize effluent pollution from factories. These results can help managers gain insights into the temporal and spatial variations in water quality of a seagoing river ecosystem and consequently make better ecological management decisions.