Water Quality of Inflows to the Everglades National Park over Three Decades (1985–2014) Analyzed by Multivariate Statistical Methods

The Everglades, a vast subtropical wetland, dominates the landscape of south Florida and is widely recognized as an ecosystem of great ecological importance. Data from seven inflow sites to the Everglades National Park (ENP) were analyzed over three decades (1985–2014) for temporal trends by the STL (integrated seasonal-trend decomposition using LOESS) method. A cluster analysis (CA) and principal component analysis (PCA) were applied for the evaluation of spatial variation. The results indicate that the water quality change trend is closely associated with rainfall. Increasing rainfall results in increasing flow and thus, decreasing concentrations of nitrogen and phosphorus. Based on 10 variables, the seven sampling stations were classified by CA into four distinct clusters: A, B, C, and D. The PCA analysis indicated that total nitrogen (TN) and total phosphorus (TP) are the main pollution factors, especially TN. The results suggest that non-point sources are the main pollution sources and best management practices (BMPs) effectively reduce organic nitrogen. However, TN and TP control is still the focus of future work in this area. Increasing the transfer water quantity can improve the water quality temporarily and planting submersed macrophytes can absorb nitrogen and phosphorus and increase the dissolved oxygen (DO) concentration in water, continuously improving the water quality.


Introduction
Surface water, especially from rivers and lakes, is closely related to human activities, so the water quality cannot be ignored. Congress created the Everglades National Park (ENP) in 1947 to protect and preserve portions of the south Florida ecosystem. Because the Everglades is located in the upstream region of the ENP, the area is subject to upstream water management practices. In order to protect and monitor the ENP waters, Water Conservation Areas (WCAs) and Stormwater Treatment Areas (STAs) were constructed to improve the water quality in the north of the ENP.
Some researchers have focused on water quality improvement in WCAs and STAs. For example, Entry [1] studied the impact of STAs and best management practices (BMPs) on nutrient concentration in the Florida Everglades using Kendall tau or Tobit trend analyses for the entire year. The trend analyses suggested that agricultural BMPs have led to a large initial reduction in nutrient concentrations in the ENP. Chen et al. [2] examined annual phosphorus (P) removal in the six Everglades STAs during 1995-2011. STA performance, in terms of outflow total phosphorus (TP) concentration, was shown to depend on three primary variables: the hydraulic loading rate (HLR), the inflow TP concentration, and the loading rate (PLR). However, the mechanism through which these modifications have affected

Study Area and Data Source
In this paper, selected chemical and physical water quality parameters were examined according to a monthly time series at seven stations at L29 and L31W canals in the ENP in the last three decades   (Figure 1; The relative position and functions of L29, L67C and L67A canals refer to http://www.evergladeshub.com/). Among the seven stations, S12D and S332 sites were cancelled after 2008 and 2006, respectively. Water quality parameters and hydrogeologic data were downloaded from the South Florida Water Management District's database (SFWMD: http://www.sfwmd.gov/portal/ page/portal/pg_grp_Sfwmd_era/pg_sfwmd_era_dbhydrobrowser). The monitoring sites are named S333, S332, S18C, S12A, S12B, S12C, and S12D based on their structures. Sample handling and chemical analyses were done in accordance with published protocols [18]. The detection limits and analytical methods for selected parameters are listed in Table 1. The statistical evaluations included dissolved oxygen (DO), pH, turbidity (TURB), total nitrogen (TN = TKN + NO X -N), ammonia nitrogen (NH 4 -N), total Kjeldahl nitrogen (TKN), total phosphorus (TP), ortho-phosphorus (PO 4 -P), and rainfall. structures. Sample handling and chemical analyses were done in accordance with published protocols [18]. The detection limits and analytical methods for selected parameters are listed in Table 1.

Seasonal Trend Decomposition Using the Loess (STL) Method
The time series of indicators were analyzed by seasonal trend decomposition using the Loess (STL) method. LOESS is considered to be a smoother and a nonparametric statistical method. In STL, observed time series are decomposed into the trend component, seasonal component, and residual component [7]. The trend component with low frequency can be considered to be a change tendency. The seasonal component with high frequency can be viewed as variation with stable seasonal disturbance. The residual component with random disturbance can be regarded as irregular variation. Therefore, removing the seasonal and residual components to obtain the trend component is very useful for understanding the change tendency of variables [19].

Seasonal Trend Decomposition Using the Loess (STL) Method
The time series of indicators were analyzed by seasonal trend decomposition using the Loess (STL) method. LOESS is considered to be a smoother and a nonparametric statistical method. In STL, observed time series are decomposed into the trend component, seasonal component, and residual component [7]. The trend component with low frequency can be considered to be a change tendency. The seasonal component with high frequency can be viewed as variation with stable seasonal disturbance. The residual component with random disturbance can be regarded as irregular variation. Therefore, removing the seasonal and residual components to obtain the trend component is very useful for understanding the change tendency of variables [19].

Cluster Analysis (CA) and Principal Component Analysis (PCA)
In this study, hierarchical cluster analysis (CA) and principal component analysis (PCA) were performed with the commercial statistics software package SPSS version 19.0 for Windows. CA is an extremely useful tool which can classify a set of observations into several mutually exclusive unknown groups based on a combination of internal variables. CA allowed the grouping of river water samples based on their similarities in chemical composition [20]. PCA is a very powerful technique for reducing the dimensionality of a data set with many inter-related variables [12]. An eigenvector is weighted linear combination of the original variables, which is equal to coefficients multiplied by the original correlated variables to obtain new uncorrelated (orthogonal) principal components. PCA allows easier interpretation of a given multidimensional system through reducing the number of correlated variables to a smaller set of orthogonal factors [21].

Basic Statistics of the Monthly Measured Data at Each Station
The basic statistics of the monthly measured data on water quality are summarized in Table 2. DO was measured from 8:00 a.m. to 3:30 p.m. The maximum value of DO was 15.4 mg L −1 at S12A, and the minimum value was 0.1 mg L −1 at S332. The maximum values of DO at all seven stations exceeded 8 mg L −1 which shows very high primary productivity in the water. However, the minimum values of DO were only 0.1-0.6 mg L −1 which shows that the water was sometimes under the anaerobic condition. The magnitudes of the mean values at each station were in the following order: S12B > S18C > S12A > S12C >S12D > S333 > S332. The pH value was relatively stable across all stations in the range of 7. 19-7.45. Turbidity values were relatively low, and the mean value was below 2 NTU, except at S332. The maximum concentration of NH 4 -N was 1.36 mg L −1 at S12B; the minimum value was below 0.005 mg L −1 , and the mean value at all stations was below 0.1 mg L −1 , except at S332. The mean TKN and TN values were followed the same order of magnitude: S333 > S12D >S12C > S12B > S12A > S332 > S18C. The PO 4 -P and TP concentrations were very low at all stations and higher at the stations on L29 than in L31W and C111. The water quality standard in the Everglades is 0.01 mg L −1 for TP (Florida Administrative Code, FAC 62-302.540) and 1.27 mg L −1 for TN [22]. The mean values of TN at S12C, S12D and S333 were 1.29, 1.37, and 1.44 mg L −1 , respectively, which are higher than the standard value, and TP exceeded this criterion at all stations except S332 and S18C. The mean values of TP and PO 4 -P at both S332 and S18C in C111 were much lower than those at S333 in L31W. Table 2. Basic statistics of the monthly measured data at seven stations Parameters S12A S12B S12C S12D S333 S332 S18C

DO, pH, and Turbidity
The trend changes in DO, pH, and turbidity at S12A, S12B, S12C, and S12D were similar ( Figure 2), as were the DO trend changes at S332, S333, and S18C. There were upwards trends in DO from 1985 to 1990 and 2005 to 2010 and downwards trends in DO from 1995 to 2003 and after 2007 for S12A, S12B and S12C. The DO values were higher in 1990 and 1996 and lower in 1993 and 2003 at S12D. There was a very clear downward trend in DO values from 1990 to 1995 at all stations, especially at S332. On the whole, the pH values showed similar upwards trends at all stations with the lowest point reached in 1992, except for at S12A. A downwards trend in turbidity was found at S18C, although there were some fluctuations. All stations showed an obvious downwards trend in turbidity from 1990 to 1998 and an upwards trend at S12A, S12B, S12C, and S12D after 2005. Almost all stations had their minimum values at around 1995 and 2005.

DO, pH, and Turbidity
The trend changes in DO, pH, and turbidity at S12A, S12B, S12C, and S12D were similar ( Figure 2), as were the DO trend changes at S332, S333, and S18C. There were upwards trends in DO from 1985 to 1990 and 2005 to 2010 and downwards trends in DO from 1995 to 2003 and after 2007 for S12A, S12B and S12C. The DO values were higher in 1990 and 1996 and lower in 1993 and 2003 at S12D. There was a very clear downward trend in DO values from 1990 to 1995 at all stations, especially at S332. On the whole, the pH values showed similar upwards trends at all stations with the lowest point reached in 1992, except for at S12A. A downwards trend in turbidity was found at S18C, although there were some fluctuations. All stations showed an obvious downwards trend in turbidity from 1990 to 1998 and an upwards trend at S12A, S12B, S12C, and S12D after 2005. Almost all stations had their minimum values at around 1995 and 2005.   Figure 3 shows the change trends in NH 4 -N, TKN, and TN at the seven stations. There were similar change trends in NH 4 -N, TKN, and TN at S12A, S12B, S12C, S12D, and S333 in L29. There was no NH 4 -N change trend at S12A, S12B, S12C, or S12D, while there was an obvious downward change at S332 and S18C after 1999. The NH 4 -N concentration was obviously higher at S332 than at other stations. TKN and TN change trends were similar at S332 and S18C. There were downwards trends in TKN and TN at all stations over the three decades, but changes were slow after 1995. TN was lower from 1995 to 2005. There was almost no change in NH 4 -N at S12A, S12B, S12C, S12D, and S333, but TKN and TN decreased, indicating that organic nitrogen decreased at all stations Figure 3 shows the change trends in NH4-N, TKN, and TN at the seven stations. There were similar change trends in NH4-N, TKN, and TN at S12A, S12B, S12C, S12D, and S333 in L29. There was no NH4-N change trend at S12A, S12B, S12C, or S12D, while there was an obvious downward change at S332 and S18C after 1999. The NH4-N concentration was obviously higher at S332 than at other stations. TKN and TN change trends were similar at S332 and S18C. There were downwards trends in TKN and TN at all stations over the three decades, but changes were slow after 1995. TN was lower from 1995 to 2005. There was almost no change in NH4-N at S12A, S12B, S12C, S12D, and S333, but TKN and TN decreased, indicating that organic nitrogen decreased at all stations

PO4-P and TP
There was a downwards Loess trend in PO4-P for all stations (Figure 4). The PO4-P concentration decreased drastically before 1990, and gently after that at S12A. The decrease was linear at S12B. There was almost no change in PO4-P before 2008, a decrease from 2008 to 2010, and almost no change in recent years at S12C. There were two and three small peaks-1993, 2000 at S12D and 1993, 2000, 2005 at S333. The PO4-P concentration did not change at S332. There was a distinct peak in PO4-P at S18C in 1995 which decreased after 2008. TP fluctuated; there was an upwards trend in 1985-1990, 1995-2000, and 2005-2008 and a downwards trend in 1990-1995 and 2000-2005 and after 2008 for S12 A to S12D. There was an upwards trend in TP from 1985 to 1990 and from 1995 to 2000, and a downwards trend from 1990 to 1995 and after 2000 at S333 and S332. TP decreased slowly from 1985 to 2005 at S18C.

PO 4 -P and TP
There was a downwards Loess trend in PO 4 -P for all stations (Figure 4). The PO 4 -P concentration decreased drastically before 1990, and gently after that at S12A. The decrease was linear at S12B. There was almost no change in PO 4 -P before 2008, a decrease from 2008 to 2010, and almost no change in recent years at S12C. There were two and three small peaks-1993, 2000 at S12D and 1993, 2000, 2005 at S333. The PO 4 -P concentration did not change at S332. There was a distinct peak in PO 4 -P at S18C in 1995 which decreased after 2008. TP fluctuated; there was an upwards trend in 1985-1990, 1995-2000, and 2005-2008 and a downwards trend in 1990-1995 and 2000-2005 and after 2008 for S12 A to S12D. There was an upwards trend in TP from 1985 to 1990 and from 1995 to 2000, and a downwards trend from 1990 to 1995 and after 2000 at S333 and S332. TP decreased slowly from 1985 to 2005 at S18C.  Figure 3 shows the change trends in NH4-N, TKN, and TN at the seven stations. There were similar change trends in NH4-N, TKN, and TN at S12A, S12B, S12C, S12D, and S333 in L29. There was no NH4-N change trend at S12A, S12B, S12C, or S12D, while there was an obvious downward change at S332 and S18C after 1999. The NH4-N concentration was obviously higher at S332 than at other stations. TKN and TN change trends were similar at S332 and S18C. There were downwards trends in TKN and TN at all stations over the three decades, but changes were slow after 1995. TN was lower from 1995 to 2005. There was almost no change in NH4-N at S12A, S12B, S12C, S12D, and S333, but TKN and TN decreased, indicating that organic nitrogen decreased at all stations

PO4-P and TP
There was a downwards Loess trend in PO4-P for all stations (Figure 4). The PO4-P concentration decreased drastically before 1990, and gently after that at S12A. The decrease was linear at S12B. There was almost no change in PO4-P before 2008, a decrease from 2008 to 2010, and almost no change in recent years at S12C. There were two and three small peaks-1993, 2000 at S12D and 1993, 2000, 2005 at S333. The PO4-P concentration did not change at S332. There was a distinct peak in PO4-P at S18C in 1995 which decreased after 2008. TP fluctuated; there was an upwards trend in 1985-1990, 1995-2000, and 2005-2008 and a downwards trend in 1990-1995 and 2000-2005 and after 2008 for S12 A to S12D. There was an upwards trend in TP from 1985 to 1990 and from 1995 to 2000, and a downwards trend from 1990 to 1995 and after 2000 at S333 and S332. TP decreased slowly from 1985 to 2005 at S18C.

Cluster Analysis
Based on the 10 variables, the seven sampling stations were classified by CA into four distinct clusters: A, B, C, and D ( Figure 5). Cluster A was formed with S12A and S12B, Cluster B included S12C, S12D, and S333, while S332 and S18C were ascribed to Clusters C and D. The five stations located on the northern boundary were divided into two groups. S332 and S18C were divided separately. The seven stations can be reduced to four monitoring stations with the same monitoring function. Thus, more useful information could be obtained from as few stations as possible.

Cluster Analysis
Based on the 10 variables, the seven sampling stations were classified by CA into four distinct clusters: A, B, C, and D ( Figure 5). Cluster A was formed with S12A and S12B, Cluster B included S12C, S12D, and S333, while S332 and S18C were ascribed to Clusters C and D. The five stations located on the northern boundary were divided into two groups. S332 and S18C were divided separately. The seven stations can be reduced to four monitoring stations with the same monitoring function. Thus, more useful information could be obtained from as few stations as possible.

Principal Component Analysis
To identify the main pollution factors that influenced the identified group, S12A, S12D, S332, and S18C were selected as the typical stations for the PCA analysis ( Figure 6). The four factors in Cluster A (S12A) explained 75.24% of the total variance. The first factor explained 33.91% of the total variance with heavy loading on TURB, TKN, TN and TP of 0.792, 0.759, 0.768, and 0.806, respectively. Factor 2 was dominated by DO and pH, accounting for 14.93% of the total variance. Factor 3 was moderately loaded by NO2, NO3, and PO4-P, explaining 13.40% of the total variance. For Cluster A, the turbidity induced by total nitrogen and phosphorus and other non-point sources was concluded to be main pollution parameter. For Cluster B (S12D), the four factors explained 69.95% of the total variance. The first factor explained 29.48% of the total variance; TN had strong loading and TP, TKN, NO2, and NO3 had moderate loading. Factor 2 was dominated by TURB, which accounted for 15.28% of the total variance. Factor 3 was moderately loaded by PO4-P, which accounted for 13.62% of the total variance. Nitrogen and phosphorus were identified as the main pollution parameters in Cluster B. For Cluster C (S332), the four factors explained 76.08% of the total variance. The first factor explained 23.63% of the total variance with TN and TKN loaded strongly. For Factor 2, NH4 was strongly loaded and DO was moderately loaded, accounting for 17.46% of the total variance. Factor 3 accounted for 13.24% of the total variance with moderately loaded pH and NO2. Nitrogen was the most dominant factor and phosphorus was the second most dominant factor in Cluster C (S332). For Cluster D (S18C), the four factors explained 72.74% of the total variance. The first factor explained 25.40% of the total variance with TN and TKN being strongly loaded and NH4-N being moderately loaded. Factor 2 was dominated by DO, pH, and NO2 which accounted for 19.82% of the total variance. Factor 3 was moderately loaded by PO4-P and TP which accounted for 15.69% of the total variance. The analysis indicated that TN and TP are the main pollution factors. S12D S333 S12C S12A S12B S332 S18C Figure 5. Dendrogram showing the clustering of sampling sites according to the water quality parameters of selected sites.

Principal Component Analysis
To identify the main pollution factors that influenced the identified group, S12A, S12D, S332, and S18C were selected as the typical stations for the PCA analysis ( Figure 6). The four factors in Cluster A (S12A) explained 75.24% of the total variance. The first factor explained 33.91% of the total variance with heavy loading on TURB, TKN, TN and TP of 0.792, 0.759, 0.768, and 0.806, respectively. Factor 2 was dominated by DO and pH, accounting for 14.93% of the total variance. Factor 3 was moderately loaded by NO 2 , NO 3 , and PO 4 -P, explaining 13.40% of the total variance. For Cluster A, the turbidity induced by total nitrogen and phosphorus and other non-point sources was concluded to be main pollution parameter. For Cluster B (S12D), the four factors explained 69.95% of the total variance. The first factor explained 29.48% of the total variance; TN had strong loading and TP, TKN, NO 2 , and NO 3 had moderate loading. Factor 2 was dominated by TURB, which accounted for 15.28% of the total variance. Factor 3 was moderately loaded by PO 4 -P, which accounted for 13.62% of the total variance. Nitrogen and phosphorus were identified as the main pollution parameters in Cluster B. For Cluster C (S332), the four factors explained 76.08% of the total variance. The first factor explained 23.63% of the total variance with TN and TKN loaded strongly. For Factor 2, NH 4 was strongly loaded and DO was moderately loaded, accounting for 17.46% of the total variance. Factor 3 accounted for 13.24% of the total variance with moderately loaded pH and NO 2 . Nitrogen was the most dominant factor and phosphorus was the second most dominant factor in Cluster C (S332). For Cluster D (S18C), the four factors explained 72.74% of the total variance. The first factor explained 25.40% of the total variance with TN and TKN being strongly loaded and NH 4 -N being moderately loaded. Factor 2 was dominated by DO, pH, and NO 2 which accounted for 19.82% of the total variance. Factor 3 was moderately loaded by PO 4 -P and TP which accounted for 15.69% of the total variance. The analysis indicated that TN and TP are the main pollution factors.  Factor loadings for the first factor (Factor 1) and the second factor (Factor 2). (A-D) indicate S12A, S12D, S332, and S18C respectively.

Linkage of Water Quality Change with Precipitation
There are only two seasons (dry season and wet season) in the study area. The water quality change trend was closely associated with rainfall, for example, S12D and S332. The monthly rainfall trend is illustrated in Figure 7. Increasing rainfall resulted in increasing flow which resulted in a decreasing concentration of TP. Figure 7 shows that the rainfall was low in 1990 and 2000, and the corresponding TP concentration peaked in those two years ( Figure 4). In 1995, the rainfall was high and the TP concentration values at all stations were the lowest. The TP concentration was negatively correlated with rainfall. Linear regression between the annual average TP concentration and annual precipitation at S12D (from 1985 to 2014) and S332 (from 1985 to 2005) is illustrated in Figure 8. A significant negative correlation (p = 0.01) was observed at S12D. In the wet season, the chemistry and bioavailability of N and P in stormwater are subjected to further settling and biological degradation during the course of flow conveyance from agricultural fields to the basin outlet through the drainage canal network [23]. In this case, the catchment can serve as a functional physicochemical filter that regulates the fate and transport of nutrients [24]. Thus, the concentration of TP was lower in high flow years than low flow years. The concentration of N and P is also connected with the hydraulic residence time and land types and so on. In the area with high agricultural non-point sources, the TP concentration had a significant increasing trend with agricultural runoff in the wet season [25]. In this study area, the TP concentration declined in the wet season, which showed that the phosphorus content in the soil was low. There was no obvious correlation between rainfall and the PO4-P concentration. The Pearson correlation coefficients of water quality variables and monthly rainfall at S12D and S332 are shown in Table 3. DO and rainfall showed a significant negative correlation at the two stations (p = 0.01).  Figure 6. Factor loadings for the first factor (Factor 1) and the second factor (Factor 2). (A-D) indicate S12A, S12D, S332, and S18C respectively.

Linkage of Water Quality Change with Precipitation
There are only two seasons (dry season and wet season) in the study area. The water quality change trend was closely associated with rainfall, for example, S12D and S332. The monthly rainfall trend is illustrated in Figure 7. Increasing rainfall resulted in increasing flow which resulted in a decreasing concentration of TP. Figure 7 shows that the rainfall was low in 1990 and 2000, and the corresponding TP concentration peaked in those two years ( Figure 4). In 1995, the rainfall was high and the TP concentration values at all stations were the lowest. The TP concentration was negatively correlated with rainfall. Linear regression between the annual average TP concentration and annual precipitation at S12D (from 1985 to 2014) and S332 (from 1985 to 2005) is illustrated in Figure 8. A significant negative correlation (p = 0.01) was observed at S12D. In the wet season, the chemistry and bioavailability of N and P in stormwater are subjected to further settling and biological degradation during the course of flow conveyance from agricultural fields to the basin outlet through the drainage canal network [23]. In this case, the catchment can serve as a functional physicochemical filter that regulates the fate and transport of nutrients [24]. Thus, the concentration of TP was lower in high flow years than low flow years. The concentration of N and P is also connected with the hydraulic residence time and land types and so on. In the area with high agricultural non-point sources, the TP concentration had a significant increasing trend with agricultural runoff in the wet season [25]. In this study area, the TP concentration declined in the wet season, which showed that the phosphorus content in the soil was low. There was no obvious correlation between rainfall and the PO 4 -P concentration. The Pearson correlation coefficients of water quality variables and monthly rainfall at S12D and S332 are shown in Table 3. DO and rainfall showed a significant negative correlation at the two stations (p = 0.01).

Correlation of Water Quality and Management
Water flows to the S12 structures through S333 and then enters the ENP. Water entering the L29 originates from WCA-3. There are many natural variations in geology, hydrology, and vegetation, and differences in water management and land uses in the two different regions which explains the differences in observed values between the northern and eastern structure groups. For example, water flowing into WCA-3 through the STA system resulted in improved water quality. The water inflowing from the eastern structures had a considerable residence time in the canals before entering the ENP, and land uses adjacent to these canals include agricultural and urban uses, which were managed by best management practices (BMPs) [15]. A possible reason for TN concentration not changing after a change in rainfall may be the change in management that occurred in the early 1990s. Non-point source pollution was effectively controlled after that. The temporal changes in nutrients are complex and affected by many factors (discharge, transport, and transformation processes), so it is not possible to completely explain them [26]. The TN concentration in water decreased sharply from 1985 to 1995 and decreased slowly after 1995 (Figure 3). TKN changed . Rainfall trend at S12D and S332.

Figure 8.
Linear regression between the annual average TP concentration and annual precipitation at S12D and S332.

Correlation of Water Quality and Management
Water flows to the S12 structures through S333 and then enters the ENP. Water entering the L29 originates from WCA-3. There are many natural variations in geology, hydrology, and vegetation, and differences in water management and land uses in the two different regions which explains the differences in observed values between the northern and eastern structure groups. For example, water flowing into WCA-3 through the STA system resulted in improved water quality. The water inflowing from the eastern structures had a considerable residence time in the canals before entering the ENP, and land uses adjacent to these canals include agricultural and urban uses, which were managed by best management practices (BMPs) [15]. A possible reason for TN concentration not changing after a change in rainfall may be the change in management that occurred in the early 1990s. Non-point source pollution was effectively controlled after that. The temporal changes in nutrients are complex and affected by many factors (discharge, transport, and transformation processes), so it is not possible to completely explain them [26]. The TN concentration in water decreased sharply from 1985 to 1995 and decreased slowly after 1995 (Figure 3). TKN changed Figure 8. Linear regression between the annual average TP concentration and annual precipitation at S12D and S332.

Correlation of Water Quality and Management
Water flows to the S12 structures through S333 and then enters the ENP. Water entering the L29 originates from WCA-3. There are many natural variations in geology, hydrology, and vegetation, and differences in water management and land uses in the two different regions which explains the differences in observed values between the northern and eastern structure groups. For example, water flowing into WCA-3 through the STA system resulted in improved water quality. The water inflowing from the eastern structures had a considerable residence time in the canals before entering the ENP, and land uses adjacent to these canals include agricultural and urban uses, which were managed by best management practices (BMPs) [15]. A possible reason for TN concentration not changing after a change in rainfall may be the change in management that occurred in the early 1990s. Non-point source pollution was effectively controlled after that. The temporal changes in nutrients are complex and affected by many factors (discharge, transport, and transformation processes), so it is not possible to completely explain them [26]. The TN concentration in water decreased sharply from 1985 to 1995 and decreased slowly after 1995 (Figure 3). TKN changed similarly to TN, which indicated that BMPs effectively reduced organic nitrogen. However, according to the PCA, nitrogen pollution control still cannot be ignored.

Management Implications
From the trend analysis above, we inferred that organic nitrogen decreased at all stations, but the PCA analysis showed that nitrogen in water should still be considered to be a principal factor. The TP concentration was higher in the dry season than in the wet season. So, increasing the transfer water quantity is necessary to maintain water quality. Secondly, submersed macrophytes can absorb nitrogen and phosphorus in the water and sediment and increase the DO in the water. Thirdly, although the nitrogen and phosphorus are lower and the turbidity is relatively high in the eastern structure group, it is necessary to strengthen the role of the constructed ecological slope in the riparian buffer zones. The ecological slope protects the soil slope from collapsing, and vegetation on the slope can intercept the suspended solids in runoff and then reduce the turbidity in the water when it rains. Finally, to streamline administration, the seven stations can be reduced to four monitoring stations, and management practices should still focus on controlling nitrogen and phosphorus content in the future.

Conclusions
This paper demonstrated the temporal changes in water quality from seven monitoring stations at inflows to the ENP from 1985 to 2014 with multivariate statistical analysis. The STL method is perfectly adaptable to the water quality time series trend, and CA and PCA are useful tools for understanding the complex nature of water quality issues by identifying groupings in the set of data. We conclude that (i) the DO concentration range was 0.1-15.4 mg L −1 with a very clear downwards trend from 1990 to 1995 at the seven stations, and a correlation with rainfall. The pH values showed similar slight upwards trends at all stations, and the mean values were in the range of 7.19-7.45. All stations showed an obvious downwards trend in turbidity from 1990 to 1998, and there was an upwards trend at S12A, S12B, S12C, and S12D after 2005. (ii) The mean values of TN at S12C, S12D, and S333 were higher than the standard value, and TP exceeded this criterion at all stations except at S332 and S18C. The concentrations of TN at S12C, S12D, and S333 were relatively higher than at other stations. There was almost no change trend in NH 4 -N, but downwards trends in TKN and TN indicated that organic nitrogen decreased at all stations. The mean values of TP and PO 4 -P at both S332 and S18C in C111 were much lower than those at S333 in L31W. There was a downwards trend in PO 4 -P at all stations. The TP concentration was negatively correlated with rainfall. The increase in the TP concentration from 1985 to 1990 and from 1995 to 2000 may be attributed to the decline in precipitation. (iii) The cluster analysis suggested that seven inflow sites could be classified into four distinct clusters. This finding is useful for designing water quality monitoring programs that maximize the amount of variability captured in as few stations as possible. The principal component analysis showed that TP, and especially TN, are the main pollution factors. Therefore, management practices should still focus on controlling nitrogen and phosphorus contents in the future.