Principal Component Analysis and the Water Quality Index—A Powerful Tool for Surface Water Quality Assessment: A Case Study on Struma River Catchment, Bulgaria

Abstract

The water quality assessment of the surface water bodies (SWBs) is one of the major tasks of environmental authorities dealing with water management. The present study proposes a water quality assessment scheme for the investigation of the surface waters’ physicochemical status changes and the identification of significant anthropogenic pressures. It is designed to extract valuable knowledge from the Water Frame Directive (WFD) mandatory monitoring datasets. The water quality assessment scheme is based on the Canadian Council of Ministers of the Environment water quality index (CCME-WQI), trend analysis of estimated WQI values, and Principal Component Analysis (PCA) using calculated excursions during the determination of WQI values. The combination of the abovementioned techniques preserves their benefits and additionally provides important information for water management by revealing the latent factors controlling water quality, taking into account the type of the SWB. The results enable the identification of the anthropogenic impact on SWBs and the type of the corresponding anthropogenic pressure, prioritization and monitoring restoration measures, and optimization of conducted monitoring programs to reflect significant anthropogenic pressures. The proposed simple and reliable assessment scheme is flexible to introducing additional water quality indicators (hydrological, biological, specific pollutants, etc.), which could lead to a more comprehensive surface water quality assessment.

1. Introduction

Everyday human life necessitates the use of one of the Earth’s most vital resources—water. This use often leads to adverse impacts on the flora and fauna of the planet since about 80% of the used water is returned to the environment untreated, leading to decreased water quality. As a result, a worldwide increase in freshwater scarcity is observed. To reflect this problem, an adaptation of the relevant legislation has been initiated at the start of the new millennium—WHO Guidelines [1], Millennium Development Goals [2] and Directives [3] with the sole goal of improving the management and utilization of natural resources on the one hand, and to reduce the human pollution, on the other.

In addition to the anthropogenic influence, deterioration of the physical, chemical, or biological water quality parameters in surface waters is also due to natural processes such as the wearing of rocks, run-off, and ion exchange, which often makes the water resource unsuitable for its intended use. Therefore, to understand the changes in water quality, it is mandatory to monitor the possible increasing or decreasing trends in water quality to successfully perform water management. To achieve this goal, the responsible stakeholders should execute appropriate monitoring programs, including all required parameters to figure out the natural and anthropogenic processes controlling water quality [4]. The monitoring conducted by the environmental authorities produces large and complex data sets including water quality parameters with different sampling frequencies. Furthermore, it involves the determination of many water quality indicators at sampling locations positioned in different types of surface water bodies (SWBs) during different seasons [5,6].

The water quality index (WQI) approach is a preferred tool for water resource management as it estimates surface water quality in a single value. This value estimates the water quality at a specific place and at a specific time, taking into account water quality indicators’ adopted threshold values [7,8]. This way of presentation is convenient for the assessment of temporal surface water quality changes by combining WQI with the Mann–Kendall test [9,10,11] or with other techniques such as Spearman correlation tests [12] and rescaled range analysis [13]. The main disadvantage of such an approach is that the calculation of the integral WQI is performed by dimensionality reduction in the monitoring data set without taking into account its structure.

The WQI concept was introduced long ago in Germany [14], and since then, numerous WQIs have been developed [15,16], such as the National Sanitation Foundation WQI (NSFWQI) [17] and Oregon WQI (OWQI) [18]. They are applied in many regions of the world [19,20]. Different modifications of WQI have been introduced by the researchers, and comparisons of the different modifications exist in the scientific literature [21,22]. The Canadian Council of Ministers of the Environment (CCME) WQI is currently one of the most widely used, not only in all the provinces of Canada but in many countries worldwide [23,24,25,26], as it is one of the main recommended indices for assessing the quality of drinking water by the United Nations Environment Programme Global Environment Monitoring System (GEMS)/Water Programme [27]. Although there are many WQIs available worldwide, none is universally accepted due to the difference in climate and land use. The main advantages of CCME-WQI facilitating its application are its simplicity, flexibility to input water quality indicators and thresholds, and tolerance to missing data [7]. Thus, the adoption of CCME-WQI to different national environmental legislations could be easily performed.

To reveal the relationships within the data set and extract knowledge concerning surface water quality, a variety of multivariate statistical methods are added to the WQI. The Cluster Analysis (CA) is used to detect the groups of similarities between sampling locations and/or water quality parameters [12]. The Principal Component Analysis (PCA) is used to reveal the “hidden” factors controlling surface water quality [28,29,30,31,32]. Numerous groups of studies use both CA and PCA as additions to the WQI to introduce the specific outputs of both techniques to water quality assessment [11,33,34,35,36,37,38,39]. To investigate the temporal behavior of water quality indicators at regions with pronounced seasonality, Discriminant Analysis (DA) is used [13,40,41,42,43,44]. In some studies, multivariate methods for the selection of important water quality indicators, such as redundancy discriminant analysis [45] and the Boruta algorithm [10], are used. Where the monitoring scheme includes apportioning studies for enough sampling locations and specific for the investigated region water quality indicators [40,44], water quality risk assessment by prediction models is carried out [10,46,47]. In the case of high-density monitoring networks, GIS techniques are applied for the visualization of spatial patterns [48,49].

The above-mentioned multivariate statistical methods reveal the relationships between sampling points and/or water quality indicators in the original monitoring data but do not take into consideration the specificity of the investigated catchments and environmental legislation. None of these multivariate approaches added to the WQI in the proposed schemes for water quality assessment takes into account the adopted environmental thresholds, which could bias the estimation for catchments including different surface water bodies (SWBs) with different adopted thresholds. Additionally, the multivariate analysis of original monitoring data could not give a real insight into the “structure” of composite WQI, which could reflect different water quality scenarios [6].

The aim of this study is to develop a simple and reliable scheme for water quality assessment of a catchment that includes different types of SWBs. The case study of the Struma River catchment in Bulgaria is explored. The assessment is performed by the following steps: (i) calculation of WQI (using CCME-WQI) for each sampling location and estimation of temporal changes in water quality by the Mann–Kendall test; (ii) application of PCA using calculated excursions during WQI calculation to reveal latent factors controlling water quality, and (iii) spatial analysis based on factor score plots obtained by the PCA.

2. Materials and Methods

2.1. Struma River Catchment

The Struma River is one of the longest rivers in Bulgaria and is second only to the Maritsa River in terms of catchment area (10,797 km²). The catchment basin has a predominantly mountainous and hilly terrain. The river springs at a height of 2180 m at the southern slope of Vitosha and leaves Bulgarian territory at the village of Kulata (62 m above sea level). The river catchment has an oblong shape and asymmetry. It is located in Southwestern Bulgaria and includes parts of Kraishteto, Vitosha, Osogovo-Belasishka mountain range, Rila and Pirin, as well as many valleys closed between them (Figure 1).

The Struma River basin is divided into three parts [50]. The upper course covers the western and southwestern slopes of Vitosha, the Pernik, Radomir, and Kyustendil valleys. The middle course of the river starts from the southern periphery of the Kyustendil Valley to the southern border of the country. The lower course of the river is entirely on the territory of Greece.

The Struma River receives most of its right tributaries in Kraishteto (Konska, Svetla, Treklyanska, and Dragovishtitsa), and larger right tributaries in the middle course are Eleshnitsa, Logodashka, Lebnitsa, and Strumeshnitsa. The left tributaries of the river are significantly more abundant because they are formed in the high Rila and Pirin mountains.

The discharge of the river varies significantly along the course. At Pernik, the amount of water is around 2.2 m³/s, while at the border, it reaches 76.1 m³/s. The high water is observed in the period March–May, and the low water in July–October.

In 2021, about 412,000 people lived in the river catchment area, and compared to 2011, the population decreased by about 50,000. The urban population predominated (73 to 27%). The main part lives in the large regional cities—Pernik, Blagoevgrad, and Kyustendil, as well as in several municipal centers—Dupnitsa, Sandanski, Petrich, and Radomir.

Blagoevgrad district has the smallest share of the population connected to wastewater treatment plants (WWTPs)—32.7%. The Kyustendil region has 62.1%, and the Pernik region—73.8% [51]. There are 5 WWTPs downstream, which discharge directly into the Struma River or in its tributaries—WWTP Batanovtsi, WWTP Radomir, WWTP Dupnitsa, WWTP Kyustendil, and WWTP Blagoevgrad. The agricultural land, including arable land and pastures, covers 29.2% of the area of the Struma catchment.

2.2. Sampling Locations

The codes and descriptions of the 30 monitoring locations are presented in

Table 1. Sampling locations and water bodies’ types.

Code	Description	Type *
A1	Studena Dam (surface sample)	L3
A2	Struma River near Pernik town, after the Church quarter, before the mouth of the Rudarska River	R5
A2-1	Konska River before flowing into the Struma River	R5
A3	Struma River on the bridge near Batanovtsi town, after the Wastewater Treatment Plant	R5
A3-1	The Arkata River before flowing into the Struma River	R13
A4	Struma River at the bridge near Priboy village, after the confluence of the Arkata River, before the Pchelina Dam	R5
A5	Pchelina Dam (surface sample)	L13
A6	Struma River after the Pchelina dam, near Zabliano village	R5
A7	Struma River near Razdavitsa village	R5
B0-1	Dragovishtitsa River near the border	R3
B0-2	Sovolyanska Bistrica River before its mouth	R5
B0-3	Banshtitsa River after Kyustendil town, before flowing into the Struma River	R5
B1	Struma River near the village of Nevestino	R5
B1-1	Eleshnitsa river before its mouth, near Chetirtsi village	R3
B1-2	Dyakovo Dam (surface sample)	L13
B1-3	Dupnishka Bistrica River after the village of Bistrica near the town of Dupnitsa before its mouth	R5
B1-4	Razmetanitsa River before its mouth	R13
B1-5	German River before flowing into the Struma River, bridge Boboshevo town	R5
B2	Struma River before Blagoevgrad town	R5
B2-1	Blagoevgradska Bistrica River after Blagoevgrad town before flowing into the Struma River	R5
B3	Struma River after Blagoevgrad town	R5
B4	Struma River after Simitli town, road. bridge on E79 in the city of Orlovets	R5
B4-1	Brezhanska River before its mouth, Poleto village	R5
B5	Struma River before Kresna town	R5
B5-1	Sandanska Bistrica River before its mouth	R3
B6	Struma River after the confluence of the Sandanska Bistrica River, after Sandanski town	R5
B6-1	Lebnitsa River before its mouth	R3
B6-2U	Strumeshnitsa River near the border (bridge to the village of Gabrene)	R5
B6-2D	Strumeshnitsa River before its mouth, bridge to Mitinovo village	R5
B7	Struma River near the border with Greece (bridge to Topolnitsa village)	R5

Notes: * Types of water bodies: R3—mountain river types; R5—semi-mountainous river types; R13—plain river types; L3, L13—types of “lakes” with oligotrophic conditions.

. The letter in the code represents the river course (A—upper, B—middle). The sampling stations located on the Struma River are numbered from the beginning of the course (A1, B1, etc.), while for tributaries, an additional number is added to specify the number of inlets between two neighboring locations on the main river. (A 2-1, B 1-1, etc.). On the transborder river Strumeshnitsa, two sampling locations (B 6-2U and B 6-2U) are sampled. Additionally, the type of SWB where the corresponding sampling point is located is given.

2.3. Data Acquisition

In this study, data obtained from the mandatory monitoring program of the Ministry of Environment and Waters for the Struma River catchment, conducted at the sampling locations described in Table 1 between 2010 and 2021, is used. The mandatory monitoring data for surface water quality include the following physicochemical parameters: dissolved oxygen (DissO₂), pH, electrical conductivity (EC), ammonia (NH₄⁺), nitrates (NO₃⁻), nitrites (NO₂⁻), total nitrogen (N), phosphates (PO₄³⁻), total phosphorus (P), and biochemical oxygen demand after 5 days (BOD₅). The selected physicochemical parameters are analyzed at least 4 times a year [52] and are an important integral part of the ecological status assessment of SWBs. All analyses are performed in accredited laboratories every month according to national legislation requirements [52].

2.4. Water Quality Index

The water quality integral indicator applied in this study is the composite WQI developed by the CCME [53]. Since 2006 it has been recommended by the United Nations Environment Programme (UNEP) for application in assessing surface water quality [54]. The CCME-WQI is based on three factors characterizing the anthropogenic impact on the water quality:

• F₁ (Scope) represents the percentage of water quality indicators not meeting the regulatory guideline values (“failed variables”) over the total number of variables included in the water quality assessment;

$$F_1=\left(\frac{Number of failed variables}{Total number of variables}\right)\times 100$$

(1)

• F₂ (Frequency) represents the percentage of measurements in which a water quality indicator exceeds the guideline values (“failed tests”) over the total number of tests (measurements);

$$F_2=\left(\frac{Number of failed tests}{Total number of tests}\right)\times 100$$

(2)

• F₃ (Amplitude) represents the extent of deviation of the “failed tests” values relative to the corresponding guideline values. The amplitude is calculated utilizing a three-step algorithm, at the beginning of which an assessment is made of the magnitude of the deviations (excursion) of the so-called “bad samples” relative to the corresponding maximum allowable concentrations:

$${excursion}_i=\left(\frac{{Failed Test Value}_i}{{Objective}_j}\right)-1$$

(3)

or

$${excursion}_i=\left(\frac{{Objective}_j}{{Failed Test Value}_i}\right)-1,$$

(4)

where Failed Test Value_i is the value of the “bad” sample, and Objective_j is the reference value of the maximum allowable concentration for the corresponding quality indicator. Equation (4) is applied to calculate the amplitude of indicators whose reference values of maximum allowable concentration decrease when moving to a worse category (e.g., dissolved oxygen content). The normalized sum of deviations NSE is then calculated:

$$NSE=\sum _{i=1}^n\frac{{excursion}_i}{Number of tests}.$$

(5)

Finally, the amplitude (F₃) is calculated using the formula:

$$F_3=\left(\frac{NSE}{0.01NSE+0.01}\right).$$

(6)

The calculation of WQI is performed by aggregation of the obtained factors as follows:

$$WQI=100-\left(\frac{\sqrt{F_1^2+F_2^2+F_3^2}}{1.732}\right).$$

(7)

The equation above gives a value between 0 and 100, with 0 being the “worst” and 100 indicating the “best” water quality. Within this range, the original CCME-WQI scale [53] ranks water quality into 5 categories: poor, marginal, fair, good, and excellent (

Table 2. Modified WQI categorization scheme for Bulgarian legislation compliance.

Status [52]	WQI	Notes	Original CCME-WQI Scale [53]
Very good	80–100	Clean and conditionally clean waters—the water quality meets the reference values for “Good” status	Excellent (95–100) Good (80–94)
Good	65–79	Weakly polluted waters—the water quality randomly deteriorates from the reference values for “Good” status	Fair (65–79)
Moderate	0–64	Polluted waters—the water quality does not meet the reference values for “Good” status	Marginal (45–64) Poor (0–44)

). According to the Bulgarian legislation [52], harmonized with the European legislation through the WFD, the characterization of the river waters status based on physicochemical parameters is performed in 3 categories—moderate, good, and very good (Table 2) [6]. Currently, there is no national WQI in Bulgaria. Testing of various WQIs began in 2007, and the CCME-WQI was found to be the most flexible to be applied according to the requirements of national legislation [55,56].

The classical CCME-WQI approach [53] was harmonized with the adoption of the Bulgarian legislation reference ranges [52] using threshold values for “good” status in the subsequent calculations (see

Table 3. National Bulgarian guidelines for characterizing the ecological status of surface waters based on physicochemical indicators for quality [52].

Type *	Status	DissO2, mg/L	pH	EC µS/cm	NH4+ mg/L	NO3− mg/L	NO2− mg/L	N mg/L	PO43− mg/L	P mg/L	BOD mg/L O2
R3	Very good	10.5–8.0	–	650	<0.04	<0.2	<0.01	<0.2	<0.01	<0.012	<1
	Good	8.0–6.0	6.5–8.5	750	0.04–0.4	0.2–0.5	0.01–0.025	0.2–0.8	0.01–0.02	0.012–0.03	1–2.5
	Moderate	<6.0	–	>750	>0.4	>0.5	>0.025	>0.8	>0.02	>0.03	>2.5
R5	Very good	10.5–8.0	–	700	<0.04	<0.5	<0.01	<0.5	<0.02	<0.025	<1.2
	Good	8.0–6.0	6.5–8.5	750	0.04–0.4	0.5–1.5	0.01–0.03	0.5–1.5	0.02–0.04	0.025–0.075	1.2–3
	Moderate	<6.0	–	>750	>0.4	>1.5	>0.03	>1.5	>0.04	>0.075	>3
R13	Very good	9.0–7.0	–	700	<0.10	<0.7	<0.03	<0.7	<0.07	<0.15	<2
	Good	7.0–6.0	6.5–8.5	750	0.10–0.3	0.7–2	0.03–0.06	0.7–2.5	0.07–0.15	0.15–0.3	2–4
	Moderate	<6.0	–	>750	>0.3	>2	>0.06	>2.5	>0.15	>0.3	>4
L3, L13	Very good	10.5–8.0	–	650	<0.03	<0.2	<0.01	<0.2	0.007–0.0125	<0.0125	<1
	Good	8.0–6.0	6.5–8.7	750	0.03–0.08	0.2–0.5	0.01–0.025	0.2–0.8	0.0125–0.04	0.0125–0.04	1–2.5
	Moderate	<6.0	–	>750	>0.08	>0.5	>0.025	>0.8	>0.04	>0.04	>2.5

Note: * Types of water bodies as in Table 1.

). For most of the investigated physicochemical indicators (except for pH and DissO₂), the threshold value is the lower limit of the “good’ status range presented in Table 3. “Good” status based on the respective parameter is achieved when the measured value is lower than the threshold value (for example, less than 750 µS/cm for EC). When the measured value of the indicator is higher than this threshold, the status is classified as “moderate”. The excursions of these water quality parameters are calculated using Equation (3). For DissO₂, the opposite approach was adopted; a “good” status is achieved when being higher than the threshold value for the respective range (6 mg/L), and Equation (4) is used for calculations of the excursions. The “good” status is only achieved for pH when the measured value falls within the reference range (Table 3). The proposed harmonization was already presented and validated for another Bulgarian river catchment [6]. Following the aforementioned scheme, annual WQIs are calculated for each sampling point during the monitoring period.

2.5. Statistical Data Treatment

2.5.1. Mann–Kendall Test

The Mann–Kendall test is a common and practical nonparametric method to determine the monotonic trend in different environmental studies producing time series [57]. Its main advantage is that there is no requirement for a specific distribution of the data, and missing data is allowed [58]. In this study, the variance correction version of the Mann–Kendall test is performed [59]. The ‘modifiedmk’ R package (version 1.6) is used for all trend analysis calculations.

2.5.2. Principal Component Analysis (PCA)

A well-documented statistical approach in surface water studies, the Principal Component Analysis (PCA) [60,61,62], is used for the analysis of the obtained monitoring data. The data interpretation is achieved through the latent factors (principal components) formed by original variables (the water quality indicators) by looking at the major sources of variance in the data [63]. The latent factors form a matrix composed of factor loadings (weights of the original variables) and factor scores (projections of the sampling locations on the principal component axis). Thus, data can be interpreted using a fewer number of principal components than the number of original variables extracting useful information about hidden relationships between original variables and revealing the impact of latent factors on different sampling locations. All calculations for PCA were performed using PLS Toolbox 9.0 (Eigenvector Research, Manson, WA, USA) of MATLAB R2021a (The Mathworks, Natick, MA, USA).

3. Results

3.1. Statistical Analysis

The basic statistics for the 10 water quality indicators is presented in

Table 4. Basic statistics of the water quality indicators.

Indicator	Unit	n	Mean	Median	Minimum	Maximum	Standard Deviation
DissO2	mg/L	1300	8.30	8.30	2.39	15.99	2.12
pH	–	1301	8.08	8.12	6.24	11.39	0.45
EC	μS/cm	1300	467.16	424.00	25.00	1883.00	270.74
NH4+	mg/L	1284	0.45	0.10	0.00	16.80	1.12
NO3−	mg/L	1281	1.15	0.85	0.00	10.35	1.19
NO2−	mg/L	1272	0.05	0.02	0.00	2.70	0.10
N	mg/L	1272	2.20	1.56	0.00	39.00	2.36
PO43−	mg/L	1285	0.14	0.08	0.00	4.90	0.24
P	mg/L	1277	0.22	0.12	0.00	5.30	0.36
BOD5	mgO2/L	1298	3.50	2.49	0.25	79.00	4.65

. The comparison of the results obtained in this study (2010–2021) with previously published data for the period 1989–1998 outlines a decrease in measured water quality indicator values except for pH and DissO₂ [64]. The comparison between monitoring results and guideline regulation values reveals three groups of indicators responsible for the moderate status of the surface waters in the Struma River catchment. The first group consists of water quality indicators with more than 50% of exceedings: P (69.4%), PO₄³⁻ (68.4%), and N (51.2%). They are followed by NO₂⁻ (38.1%), BOD₅ (36.6%), NO₃⁻ (25.8%), and NH₄⁺ (23.7%). The pH (12.5%), DissO₂ (12.1%), and EC (11.1%) have the smallest contribution to the moderate water status.

3.2. WQI

The results of the calculated water quality indices for all 30 monitoring stations along the Struma River catchment are presented in Figure 2. They demonstrate that in more than half of the points, the surface water status is “moderate”, in about 30%—“good”, and in only 16%—“very good”. In most parts of the catchment, the water quality does not meet the goals set in the Struma River Basin Management Plan to improve the status to “good”. The proportion of “moderate” status in the upper course (A) of the river is particularly high (almost ¾).

The annual trend of WQI for the investigated sampling locations is presented in Figure 3. The performed Mann–Kendall test outlines 15 out of 30 sampling locations with a significant rate of change in the WQI, while only B 1-1 (the Eleshnitsa river) has a decreasing trend during the investigated period.

3.3. PCA

The relationships between water quality indicators and/or sampling locations were examined by PCA. The input data set used for PCA consists of variables—water quality indicators and objects—and sampling locations. The water quality indicators are presented by the excursions, calculated according to Equations (3) and (4). Only for pH was the original data used. Prior to the analysis, the input data underwent autoscaling pretreatment. The PCA reveals that the first three principal components explain almost 70% of dataset variation. The number of latent factors is determined based on their eigenvalues and the internal model validation error. In the formation of the first principal component, explaining 42.68% of data variance, the nutrients (N, P), PO₄³⁻, NH₄⁺, and BOD₅ have a significant impact (Figure 4). This latent factor could be conditionally named “urban waste waters”, but to a smaller extent, could also be indicative of agricultural activity, which is not recognized as driving pressure for the Struma River catchment. The second principal component (16.11%) is formed predominantly by EC and NO₃⁻ and reflects mainly the industrial wastewaters in the Struma River catchment. The last principal component (11.03%) is positively correlated with pH and negatively correlated with DissO₂. This latent factor could be related to the surface water pollution response at the investigated sampling locations. A high pH may inhibit the photosynthesis of algae, which will reduce the DissO₂ [65].

The factor score plots resemble the impact of the abovementioned latent factors at the sampling locations in the Struma River catchment. The PC1 score plot outlines the more pronounced impact of urban wastewater on sampling locations in the upper course (A) located on the Struma River—A3 (the Struma River at Batanovtsi town), A4 (the Struma River before the Pchelina Dam) and A5 (the Pchelina Dam), and on tributaries from the middle course (B) B 0-3 (the Banshtitsa River), B 2-1 (the Blagoevgradska Bistrica River), B 5-1 (the Sandanska Bistrica River), and B 6-2D (the Strumeshnitsa River before flowing in the Struma River) (Figure 5).

The influence of urban wastewater is presented in Figure 6. The most influenced are sampling stations located on the upper course of the Struma River—A3, A4, and A5 followed by the tributary Razmetanitsa River in the middle course (B1-4).

The impact of surface water pollution response does not discriminate sampling locations, but locations influenced by urban wastewaters: A5 (the Pchelina Dam) and B1-4 (the Razmetanitsa River), which possess more scattered factor scores (Figure 7).

4. Discussion

The Struma River catchment, one of the longest in Bulgaria, has a significant impact on water quality, and as such, it has been explored in this study. Our results show that there is a significant difference between surface water quality in the upper and the middle course of the catchment due to the specifically located urban and industrial pressures.

The upstream sampling location of the Studena dam (A1) in the upper course is used as a drinking water supply source for the town of Pernik and has had a “good” or “very good” status for the entire period of investigation. The PCA analysis does not reveal any anthropogenic pressures, and the trend analysis results show a small increasing trend in WQI of 0.22 units per year. The following sampling locations situated along the Struma River (A2–A6) have moderate physicochemical status during the investigated period. The A2 (Struma River at Pernik) and A3 (Struma River at Batanovtsi town) are impacted by the urban and industrial wastewaters from the Pernik region, while the wastewater from the Radomir region affected the downstream sampling stations A4 (the Struma River before the Pchelina Dam), A5 (the Pchelina Dam), and A6 (the Struma River at Zabliano village). The Pchelina Dam was built in 1975 as a precipitator of the polluted waters from the communal and industrial pollution of the Struma River waters from the Pernik and Radomir regions. The complex nature of the anthropogenic pressure is confirmed by the previous studies on surface waters [66] and sediments of the Pchelina Dam [67]. The high variation of PC3 factor scores for A5 is potentially caused by the changes in the surface water quality and operation of the small hydroelectric station located near the Pchelina Dam. The WQIs of the sampling locations downstream of the Struma River are improved after passing through the reservoir—a conformation of the study in the period between 1989 and 1998 [68]. The sampling location of the Struma River at Razdavitsa village (A7) has had a good status for 6 out of 12 investigated years. The significant water quality improvement trend (0.43 WQI units per year) for A3 could be related to the improvement of the WWTP of Pernik and the related sewage network. The tributary Konska River (A2-1) has a pattern similar to the sampling location A2 (Struma near Pernik town), and its “moderate” status is due to urban wastewaters and small industries from Breznik town and surrounding villages. The tributary Arkata River (A3-1) has dominantly “good” and “very” good status with a significant positive water quality trend. The “moderate” status during 2017 and 2018 is caused by the small industries located in the tributary catchment (see Figure 7).

Based on the obtained results for the upper course of the Struma River catchment, the main mitigation measures for improving surface water quality should be directed to the more efficient treatment of industrial wastewaters from the regions of Pernik and Radomir, accompanied by the improvement of their sewage networks.

The sampling locations of the Struma River in the middle course (B1–B7) have predominantly “good” physicochemical status with increasing downstream water quality tendency, similar to the tendency of sampling stations in the upper course, located after the Pchelina Dam. The PCA results do not reveal any significant urban or industrial impact on the abovementioned sampling locations. There is a significant improvement in the water quality trend for the Struma River at Nevestino village (B1), while the most pronounced within the entire catchment of the Struma River is at Kresna town (B5) with 0.61 WQI units per year. The first two tributaries in the middle course—the Dragovishtitsa River (B0-1) and the Sovolyanska Bistrica River (B0-2)—have significant positive water quality trends. It should be mentioned that the Dragovishtitsa River is a transboundary river, and the WQI categories of B0-1, located close to the border, reflect the anthropogenic pressures on the territory of the Republic of Serbia. The Banshtitsa River (B0-3) has a “moderate” status during the whole investigated period, which is caused by the urban wastewater of Kyustendil town. Extreme exceedings (over 25 times) for PO₄³⁻ and P are observed. This is in agreement with the results previously reported [69], where the water quality of the Banshtitsa River is classified as a WWTP effluent. The tributaries that flow into the Struma River between sampling locations B1 and B2 could be divided into three groups. The Eleshnitsa River (B1-1) and the Dyakovo Dam (B1-2) have “good” and “very good” statuses without evidence of any anthropogenic impact from the PCA analysis. The Dupnishka Bistrica River (B1-3) and the Razmetanitsa River (B1-4) have varying “moderate” to “good” water quality statuses during the period of investigation. The reasons for the estimated “moderate” status are urban wastewater of Dupnitsa town and industrial wastewater of the thermal power plant Bobov dol for B1-3 and B1-4, respectively. The variation in the surface water quality of the Razmetanitsa River is supported by the scattered surface water pollution response factor scores (PC3). The German River (B1-5) has a “moderate” status during the whole period of observation. The PCA results do not detect upstream urban and industrial pressures at sampling points B1-3 and B1-4. It seems that surface water processes that take place between sampling locations B1-4 and B1-5 “change” the relationships between the physicochemical parameters from that included in PC1 and PC2. This assumption is supported by the registered extreme exceedings (over 25 times) for PO₄³⁻ and P in the German River and from the very high measured concentration of total organic carbon during 2021–9.48 mg/L. The introduction of additional water quality indicators specific to upstream pressures could overcome their “masking”. The next two downstream tributaries, influenced by urban wastewater, are the Blagoevgradska Bistrica River (B2-1) and the Sandanska Bistritsa River (B5-1). Both sampling locations have a significant increase in the WQI. This improvement is due to the putting into the operation of WWTP Blagoevgrad in 2012 and the improvement of the sewage network of Sandanski town for B2-1 and B5-1, respectively. The last two tributaries sampling locations where PCA results indicate urban wastewater impacts are B6-2U and B6-2D, located on the transboundary Strumeshnitsa River. The “moderate” status of the surface waters coming from the Republic of North Macedonia (B6-2U), which is characterized by PO₄³⁻ and P exceedings, does not undergo significant change at B6-2D.

In contrast with the upper course, the main anthropogenic pressures in the middle course of the Struma River have urban origins. Since the introduced mitigation measures for the treatment and collection of urban waters from Blagoevgrad and Sandanski town led to the improvement of surface water quality, the reconstruction of the sewage network of Kyustendil town is necessary to be planned in further river catchment management.

In general, along the course of the main river, the anthropogenically impacted tributaries do not lead to the deterioration of the water quality of the Struma River. The PC1 and PC2 factor scores have similar values for the first Studena Dam (A1) and the last Struma River at Topolnitsa village (B7) sampling points in the Struma River catchment.

5. Conclusions

The present study introduces an applicable and reliable scheme for water quality assessment of catchments, including different types of SWBs using WFD mandatory monitoring datasets. The assessment is performed by combining CCME-WQI, the Mann–Kendall test, and PCA, exploring in this way their advantages and complementarity.

The proposed approach provides important knowledge to the environmental authorities for (i) identification of the anthropogenic impact on SWBs and the type of the corresponding anthropogenic pressure, (ii) prioritization and monitoring of mitigation measures, and (iii) optimization of conducted monitoring programs to reflect significant anthropogenic pressures. The obtained knowledge could be used for future river catchment management. Additionally, the proposed water quality assessment scheme could be applied to monitor the effectiveness of the implemented mitigation measures.

The assessment scheme is flexible to introducing additional water quality indicators (priority pollutants, hydrological, and biological), which could lead to a more comprehensive surface water quality assessment. The opportunity for introducing special pollutants for a particular river basin opens the door for the proposed approach to be implemented not only in another river catchment but to be used on the river basin management level.