Among other things, it has been used in the past to study the distribution patterns of phyto- and bacterio-plankton in different seasons [48], as well as to identify the most important habitats of endangered and endemic fish inhabiting coral reefs [49]. Self-organizing maps were presented in 1982 by [PERSON] in [50]. These are a type of neural network used to create a transformation that preserves the features of objects in a multidimensional space into a map. These maps are usually represented as a two-dimensional grid on which observations from the multidimensional space are marked. As this is intended to be a transformation that preserves the characteristics of objects, the data that are similar in the original space are represented close together in a self-organizing map [51]. The Kohonen map consists of \(p\cdot q\) neurons, where the numbers \(p\) and \(q\) are the dimensions of this map. These numbers are established in advance at the beginning of the algorithm. Each neuron is an \(n\)-element vector (where \(n\) is the number of explanatory variables in the original set). Usually, neurons in Kohonen maps have a rectangular, hexagonal, or circular topology. Training a self-organizing map involves reducing the direct neighborhood of each neuron while reducing the weights assigned to individual neurons, which are randomly assigned at the start of the algorithm [52]. Kohonen networks have been used to visualize, cluster, and analyze multivariate environmental engineering data, among others, in [53; 54; 55]. There are many algorithms for hierarchical clustering, the main differences between them are the choice of the metric in which the distances between the observed objects are calculated and the choice of the method for determining clusters. The most commonly used metric is Euclidean, which was also used in this work; the distance between two vectors \(x\), \(y\) of length \(n\) is determined using the formula:  
\[d(x,y)=\sqrt{\sum_{i=1}^{n}(x_{i}-y_{i})^{2}}. \tag{6}\]  
There are many methods for determining clusters, for example, [PERSON]'s method presented in [56], as well as clustering by single [57] and full linkage methods [58]. In this paper, the full linkage method, which is an agglomerative method, was used. It relies on the fact that at the beginning of the algorithm, each observation is in a separate cluster. The distances between clusters \(A\) and \(B\) are calculated as follows:  
\[D(A,B)=\max_{x\in A,\ y\in B}d(x,y). \tag{7}\]  
Thus, joining of two clusters depends on the farthest pair of observations. In this method of object clustering, as the cluster grows in size, it becomes increasingly difficult to attach more objects to the cluster, because the increase in cluster size automatically distances it from other observations and clusters considered in the analysis [59]. The hierarchical clustering method was used in [41] to classify the water quality in the Fuji River and in [60] to assess the sources of heavy metal contamination in water. In an MLP (multilayer perceptron) neural network, the results are obtained not only through a set of \(n\) explanatory variables derived from the input data, but also through a set of auxiliary variables, which are latent units. Each of them is a linear combination of the original variables, usually transformed with a hidden \(g\) function, which is not linear in nature:  
\[h_{k}(x)=g\left(\beta_{0k}+\sum_{i=1}^{n}x_{i}\beta_{jk}\right). \tag{8}\]  
The coefficient \(\beta_{jk}\) denotes the effect of the \(j\)-th variable on the formation of the \(k\)-th hidden unit. Assuming that there are \(H\) hidden units in the network, the result of the neural network will be a linear combination of the following form:  
\[f(x)=\gamma_{0}+\sum_{i=1}^{H}\gamma_{k}h_{k}. \tag{9}\]  
Optimization of the coefficients of such a network involves minimizing the sum of squares of the errors with an additional correction attaching a predetermined regularization parameter \(\lambda\) to control the overfitting of the model to the training data.