Compiled data were tested for normality before performing statistical analyses. Statistical analyses comprised Kruskal-Wallis test (H), Principal Component Analysis (PCA) and Non-Metric Multi Dimensional Scaling (NMDS).  All tests were performed using the R-software (R Development Core Team, 2006).Box plots are used throughout this report to visually summarize data. These graphs depict the following statistical measures: median, upper and lower quartiles. The box plot is interpreted as follows: The box itself contains the middle 50 percent of the data. The upper edge of the box indicates the 75th percentile of the data set, and the lower edge of the box indicates the 25th percentile. The line in the box indicates the median value of the data. If the median line within the box is not equidistant from the edges of the box, then the data are skewed. “Gridding” is the operation of spatial interpolation of scattered 2D data points onto a regular grid. Gridding allows the production of a map showing a continuous spatial estimate. The spatial coverage of the map is generated automatically as a square covering the data points. Non-metric multidimensional scaling is based on Bray-Curtis distance matrix was performed for classifying the sites across river basins. In NMDS, data points are placed in 2 or 3 dimensional coordinates system preserving ranked differences.

Data were tested for normality before performing statistical analyses. Statistical analyses comprised Kruskal-Wallis test (H), Principal Component Analysis (PCA) and Non-Metric Multi Dimensional Scaling (NMDS).  All tests were performed using the R-software (R Development Core Team, 2006). The non-parametric Kruskal-Wallis test was used to assess whether species richness, species diversity and turnover across water quality regimes were significantly different. Temporal variation in diatom assemblages in each site was analyzed by NMDS using absolute abundance data. NMDS is an ordination method well suited to data that are non-normal or are arbitrary or discontinuous and for ecological data containing numerous zero values (Minchin, 1987; McCune & Grace, 2002). Results were visualized showing the most similar samples closer together in ordination space (Gotelli & Ellison, 2004). A final stress value, typically between 0 and 15, was evaluated as a measure of fitted distances against the ordination distance, providing an estimation of the goodness-of-fit in multivariate space. Changes in species composition or percentage turnover (T) were used to indicate community persistence. T was calculated as T=(G + L)/(S1 + S2) times 100 where G and L are the number of taxa gained and lost between months respectively, and S1 and S2 are the number of taxa present in successive sampling months (Diamond & May, 1977; Brewin et al. 2000; Soininen & Eloranta, 2004). The relationship between local population persistence, local abundance in terms of relative abundances, and regional occupancy were examined using correlation analysis (Soininen & Heino, 2005). For the species distribution model the species were classified as core species as species that occurred in over 90% of sites, and satellite species as species that occurred in fewer than 10% of sites (McGeoch & Gaston, 2002, Soininen & Heino, 2005). Local occupancy of each diatom species was calculated by their percentage of occurrence at each site across the seasons. Seasonal diatom community was related to the water quality parameters using multiple linear regressions. Finally, water quality variables were used in PCA to elucidate the spatial water quality variation.