{'fraction_non_alphanumeric': 0.03281318106529509, 'fraction_numerical': 0.012187254816493766, 'mean_word_length': 4.416052880075543, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 28, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Morphological classifiers

This work proposes a new type of classifier called Morphological Classifier (MC). MCs aggregate concepts from mathematical morphology and supervised learning. The outcomes of this aggregation are classifiers that may preserve shape characteristics of classes, subject to the choice of a stopping criterion and structuring element. MCs are fundamentally based on set theory, and their classification model can be a mathematical set itself. Two types of morphological classifiers are proposed in the current work, namely, Morphological k-NN (MkNN) and Morphological Dilation Classifier (MDC), which demonstrate the feasibility of the approach. This work provides evidence regarding the advantages of MCs, e.g., very fast classification times as well as competitive accuracy rates. The performance of MkNN and MDC was tested using p -dimensional datasets. MCs tied or outperformed 14 well established classifiers in 5 out of 8 datasets. In all occasions, the obtained accuracies were higher than the average accuracy obtained with all classifiers. Moreover, the proposed implementations utilize the power of the Graphics Processing Units (GPUs) to speed up processing.


Introduction
Supervised learning consists of analysing labelled instances in order to generate a classification model capable of producing labels for unlabelled instances. Supervised learning is well established and vastly employed in fields such as computer science and engineering, and in a great number of sub areas, such as robotics [1] , computer vision [2,3] , data mining, knowledge discovery [4] , health care [5][6][7][8] and many others [9] .
Classification algorithms are built upon disctinct concepts and theories. Some are based on decision-trees, others on neural networks, Bayes theorem, mathematical functions, rules, regressions, etc. The performance of each classifier varies according to the used dataset. According to the no free lunch theorem, there is no one model that works best for every optimization problem. The assumptions of a great model built for a specific problem may not hold for another problem [10] . Usually, multiple models are evaluated in order to determine the best for a particular problem.
In this work, we propose the use of mathematical morphology (MM) in the context of classification. Mathematical Morphology [11,12] is commonly employed to process and analyse geometric structures in images and is based on set theory. We introduce the theory of morphological classifiers and propose two paradigms for this new classification approach, one inspired mostly by k-Nearest Neighbours (k-NN) [13] and morphology and the other by morphological dilations [14] .
The work is organized as follows. In the next section, we perform a state-of-the-art review in supervised learning and mathematical morphology. In Sections 3 and 4 , we define the main proposal of this work and present the MkNN and MDC algorithms. Section 5 is concerned with experiments and performance comparison to other state-of-the-art approaches. Finally, in the last section, we provide the conclusions and propose avenues for future work.

Literature review
A brief introduction to supervised learning, followed by an overview of the various types of classification algorithms is introduced in this section. Next, we introduce the essential operations of Mathematical Morphology and provide examples of its application.

Classification algorithms
Lazy algorithms [15] generate classification models in real time. That is, the algorithm generates a model right before classifying  unlabelled instances. Their main weakness stems from that fact, i.e., generating new models when unlabelled instances are encountered requires rework. However, this results into a significant strength, i.e., there is no requirement to store a model in memory.
k-NN is one of the most popular lazy algorithms and probably the simplest one. It requires the k nearest neighbours of an unlabelled instance to be examined in order to determine its class. The predominant class among the k neighbours of an unlabelled instance is chosen as the target class. k-NN is a robust algorithm in terms of accuracy, yet it does not deal well with noise in some specific cases and for some values of k . One such occasion is when a small number of instances of one class is surrounded by several instances of other classes. A simple k-NN with k = 3 is depicted in Fig. 1 . The unlabelled instance is checked against the 3 nearest labelled instances and is classified as the most frequently encountered class that, in this case, is the square class.
Decision tree based classification relies on decision tree models. Essentially, these algorithms generate one or more trees that, in turn, contain several nodes. Each node split is usually associated to a formulated Boolean criterion, which depends on one or more attributes of the dataset. For instance, a i > s , where a i is the attribute at index i and s is a constant. The split governs the flow of the decision, and the leaves of the trees represent the various classes. A popular decision tree algorithm is C4.5 [16,17] (J48 in Weka [18] ). Random Forests [19] is also frequently employed in classification. It generates several decision trees that, when combined together, compose an ensemble learning method [20] . Fig. 2 illustrates a possible decision tree based classification model.  Another type of classificier is neural networks. Neural networks are often categorized as function-based classifiers. They associate input patterns to output patterns by using the representation of intermediate layers, activation functions and weights that connect the various nodes, mimicking the organization of the brain. Two popular neural networks are the Radial Basis Function Network (RBF Network) [21,22] and the Multilayer Perceptron (MLP) [23,24] . Fig. 3 shows a feedforward neural network with n input units, m hidden units and two outputs. The outputs of the hidden units are the weighted sums of the input attributes after they pass through the activation functions of these units. Finally, the outputs of this neural network are calculated as the weight sums of the input from the hidden units after they pass through the activation function of the output units. Support Vector Machines (SVM), on the other hand, are examples of function-based classifiers. In summary, SVM traces hyperplanes that segregate the data, and uses the created clusters or segregations to predict classes. In its basic form, SVM works in a straightforward fashion. Fig. 4 illustrates a possible solution, shown by the dashed line, which separates the data optimally. Optimality in this case refers to higher accuracy and margin, i.e., gap size g . Since this problem is 2-dimensional, the hyper-plane is a simple line. Algorithms that derive from SVM such as Sequential Minimal Optimization (SMO) [25] also bear the same optimization objective.   6. A decision table that predicts the class of unlabelled instances based on the rules r 0 , r 1 and r 2 . If r 0 and r 1 are true, the instance is probably a circle. Some of the previously described algorithms can be thought of as probabilistic models. That is, some decision tree algorithms, for instance, output a probability on the correctness of classification of a certain node, or the respective path of the tree that a given instance followed. These probabilities can be computed in various ways and therefore, these algorithms could be considered probabilistic. However, some classifiers are based on conditional probabilities, more specifically, on Bayes Theorem [26] . These are explicitly categorized as probabilistic. Here, we refer to them as Bayesbased algorithms, with the popular algorithms being the Naive Bayes [27,28] and Bayes Net [6,29] . In Bayes theorem, P ( E i ) denotes the probability of a given event E i being true. The event could be, for instance, the probability of an attribute being equal, greater, lower or different than a certain value. P ( E a | E b ) denotes the conditional probability, which is the probability of E a being true given that E b is true. A simple Bayes net is shown in Fig. 5 . The joint probabilities at the right of the figure indicate the probability of the analysed instance belonging to a certain class, given the described events. The exclamation marks indicate the NOT operator.
Lazy, decision tree, function, and Bayes classifiers are probably the most commonly used classifers in the state-of-the-art. However, there is another category of classifiers that often produces very good results, yet it is often disregarded. These are rule-based classifiers such as Decision Table [30] and Conjuctive Rule [18] . A decision table can be extraordinarily simple. Such table can be constructed using the default rule mapping to the majority class. In general, decision tables are relatively easy to understand as generated rules are often simple and constructed using perfectly understandable logical constraints. The construction process involves heuristical search of feature subsets for rule finding, as well as cross-validation [30] . Fig. 6 shows a decision table model. The first variables r 0 , r 1 and r 2 represent rules based on training instances. In the second  column, if rules r 0 and r 1 are true, then the table classifies the instance as a circle. In the fifth column, on the other hand, if an unlabelled instance has rules r 0 and r 1 equal to false and r 2 true, then it is classified as a square. The rule itself can be a simple comparison split of an attribute (e.g., if an attribute value is higher than a constant), as in decision tree models, or more complex associations.

Mathematical morphology (MM)
Mathematical morphology [12,31] is concerned with geometric structures present in images. The morphological analysis of images is performed using set theory [32] . The structuring element that is used in most mathematical morphology operations is defined by a mathematical set of elements that represents its pixels. Dilations and erosions are the two principal and fundamental operations in MM. A dilation of a binary image A = { a 1 , a 2 , . . . , a n } , containing n foreground pixels, by a structuring element B , is given by [14] : where A b represents the elements of A translated by b . In other words, dilation is given by the union of A + b for every element b in B . Each element represents the position of a pixel. Fig. 7 illustrates a dilation of a binary image A by a cross-shaped structuring The dilation operation sets the intensity of any background pixels via the superposition of the structuring element on each of the foreground pixels of the binary image to the foreground value. In this specific case, since the spatial extent of the structuring element is fairly limited, it spreads the foreground pixels that already exist to their neighbouring pixels, thus filling the gaps of the original image.
Erosion is the opposite operation to dilation, and is defined as [14] : Fig. 8 shows the erosion of an image A by the same crossshaped structuring element B . Erosion erases pixels of an image respecting the structuring element. In this specific case, since the spatial extent of the structuring element is fairly limited, it erodes the pixels near the original ones, so that holes or noise in the image become more evident.
Mathematical morphology involves a wide variety of operations such as opening, closing, top-hat, reconstruction and others. These operations are derived from the essential operations of dilation and erosion. For instance, opening is a combination of erosions and dilations, while closing is the opposite, a combination of dilations and erosions, respectively. Top-hat is the result of subtracting the original image from its opening [33] . Morphological reconstruction, on the other hand, reconstructs the shapes of the image using seed points, dilations and a mask [34] .
In general, mathematical morphology is often employed in the literature for geometric feature filtering, e.g., to segment structures [35] , extracting and enhancing blood vessels [36,37] , etc. Liu et al., [38] use grey-level mathematical morphology operations and a genetic algorithm to create a thresholding method for images. Furthermore, it is also used to remove noise from images [39,40] . The noise removal potential is evidenced in Fig. 7 , where a great deal of noise from the original image is removed with a simple dilation. The closing operation is also used to remove noise from images [41] .
The closest that mathematical morphology comes to classification in the existing literature is by providing image descriptors for further classification [42,43] . The power of mathematical morphology has not yet been explored in the context of classification models. In this work, we provide a theoretical foundation for mathematical morphology-based classifiers, and propose two such models.

Theoretical framework
The essence of the approach is to expand pattern sets of labelled instances such that the expansion covers the entire feature space while preserving the morphology of the original clusters of labelled instances. We commence a process that subsets are expanded using set operations from their original state i = 0 , to their final state, i = I l , for each class l . For this, we focus on data-driven classification functions, based on the p -dimensional shape of the original class clusters. The expansion in a specific direction may be terminated when a suitable termination criterion is reached. This criterion could be, for instance, the classification error P err > t err , where t err is a predefined misclassification threshold.
We consider a class of labelled patterns represented by a number of disjoint clusters as a multiset. A multiset is a generalization of the set concept, which allows the presence of repeated instances. For example, {1, 2, 2} is a different multiset in comparison to {1, 2} although they are the same set. Thus, given a training dataset X = { x 1 , . . . , x n } that is a multiset containing instances such that x ∈ C , where C ⊂ R p , p being the number of attributes of the instances, with corresponding labels c(x) ∈ { 1 , . . . , L } , an unlabelled instance x is assigned a label l ∈ { 1 , . . . , L } , such that: where I l corresponds to the final iteration for class l and A I l l is the expanded set partition of instances belonging to class l , which consists of J l disjoint subsets a I l l, j such that: In what follows, we describe the methodology for obtaining the expanded set partition A I l l . Let us assume a set of pattern points A 0 l belonging to class l , with index 0 referring to the original set partition of training samples:  The original set partition of samples belonging to class l , A 0 l , is represented by the union of the disjoint subsets a 0 l, j . An example of this in two dimensions is shown in Fig. 9 . Here, the class set partition A 0 l is the union of { a 0 l, 1 , a 0 l, 2 , a 0 l, 3 , a 0 l, 4 } .
Next, we consider the expanded set partition A i +1 l , corresponding to the application of a set operator f at iteration i : . . . . . .
where f ( A ) is a function that expands set A . Therefore, given I l iterations of the expanded class set partitions, we have L equivalence classes as follows: A candidate set morphology-preserving operation f based on mathematical morphology is the dilation [14] : by a suitable structuring element B. The result of this for iteration 1 is shown in Fig. 10 , assuming that iteration 0 is shown in Fig. 9 , where in this specific case B = In a dilator classifier, all A 0 l can be expanded simultaneously, where intersections of labels are erased when they occur. The result of the dilation of L classes until they cover the entire domain of C is a classification model.
Operation f may be any other operation apart from strict dilations. For instance, it could be a combination of morphological operations such as dilations and erosions, a pseudo-random or random function operating on the set partition, etc.
In the case of p -dimensional set partitions A 0 l , the patterns are represented in a p -dimensional grid or hypercube. Each original set partition A 0 l is represented by a unique class label l , and an appropriate set operator, such as a dilation, shown in Eq. 9 , is applied to produce the expanded set partition A i l at iteration i . The set expansion process is controlled by a convergence function, which monitors a performance measure, e.g., classification error. The evolution of a class set partition at iteration i is terminated, according to the convergence function, and the membership of the class set partition is reset to that at iteration i − . . , L } , its accuracy rate is given by: where ˆ c (y) is the predicted class of instance y and | Y | stands for the cardinality of set Y . Theorem 3.1 illustrates a possible methodology for the expansion of the set partitions. The fundamental concept of this approach is that it is possible to exploit redundancy, i.e., the complement of the union of L − 1 set partitions representing one of the set partitions, which reduces memory and computational requirements.
However, with L set partitions, L options for the complement set partition are possible, which leads us to an optimization problem. Theorem 3.1 proves that the selection of the complement partition affects accuracy.

Theorem 3.1. Impact of complement set partition A q on classification performance:
Consider the L expanded set partitions A i l at iteration i , l ∈ There is a set partition A I l q , q ∈ { 1 , . . . , L } , at final iteration I l , that maximizes accuracy for a specific classification problem.
Proof. Eq. 13 introduces the use of the complement of a set partition in the expansion process, which reduces memory and processing time requirements. For instance, in binary classification, the partitions are A i 1 and A i 2 . The expansion is performed in one of the two partitions, while the other corresponds to the complement of the first partition.
Assuming that X l is a multiset of instances with labels equal to l , then the accuracy rate for the binary case is given by: The previous accuracy concerns to A is expanded, which provides the following accuracy rate: Thus, we want to prove that the choice of A i q has an impact on the accuracy. For this to be true, we assume that Eq. 18 is true and derive contradiction. If this does not hold for at least a single case, then the theorem is proved.
Let us assume a similar case for both set accuracy computations, where A I 1 1 and A I 2 2 , respectively, are expanded until domain C is fully covered, such that A is a contradiction as X 1 and X 2 can be different multisets. Therefore, the accuracy can be affected by the choice of the complement set partition.

Methods
In this section, we propose two classification algorithms based on morphology. The first algorithm consists of an adaptation of the k-NN algorithm to the proposed morphological framework. The second algorithm is a mathematical morphology classifier based on simple dilations.
Before introducing the proposed classifiers, we need to consider that the input dataset is preprocessed so that the classification space resembles a discrete grid in order to speed up computations and for ease of implementation. Given a dataset X , each position of the discrete grid G contains a number of repeated instances in X and is defined by: (20) where ν is a precision vector ν = ν 1 , . . . , ν p and round ( ν•x) is also a p -dimensional vector. ν•x is the Hadamard product of these two vectors: (21) where in this specific case x k stands for the k th attribute of x.
The grid is a subset of natural numbers representing C , which increases the processing performance of the algorithms and facilitates data compression.

Morphological k-NN (MkNN)
In this section, we present an algorithm similar to k-NN, however conforming to the theoretical definitions of Section 3 . This algorithm generates a classification model, which does not require re-training and thus, is no longer a lazy classifier. However, as in the traditional k-NN classifier, an unlabelled instance y has its surroundings examined until the algorithm checks at least k other instances. Finally, it classifies y as the majority class of the visited neighbouring instances.
Formally, given an unlabelled instance y, the neighbouring instances are visited so that at iteration i , labelled instances that are at a distance ≤ i from y are visited. If at any iteration i , the number of visited labelled instances is ≥ k , where k is an input parameter as in the k-NN algorithm, then the unlabelled instance y is said to be of the same class as the majority of the visited labelled neighbouring instances. The f operation in this case is given by: The algorithm proceeds by analysing and examining the nearest labelled instances. For a single unlabelled position y in grid G , the surrounding labelled instances x are examined, as shown in function (): where d is a distance function and mode () returns the most frequent label in the multiset. This process is repeated for every possible unlabelled position in the grid before the algorithm reaches convergence, i.e., it terminates when k nearest instances have been visited. The convergence function for class l is then given by: where | C | indicates the total size of the subset C or grid G , i.e., the total amount of elements in the grid. Finally, the labelling function of the MkNN algorithm conforms to Eq. 11 . The pseudo-code for constructing this classification model is shown in Algorithm 1 .

Morphological dilator classifier (MDC)
Morphological Dilator Classifiers inherit more of the properties of mathematical morphology in contrast to MkNN. This algorithm dilates positions in grid G until a termination criterion is reached.
The f operation for MDC is governed by a rule β that indicates the orientation in which the instances should be dilated. Function f is then given by: Data : k and γ are input parameters, where k stands for how many neighbouring instances should be visited, γ stands for a weight applied to the central grid position and σ represents the maximal number of iterations.
T 1 , . . . , T L are counters for the classes, X is the training dataset, A 0 l respects Eq. 6 and c(x) returns the class label of instance x . where the operator dil refers to a p -dimensional dilation of the expanded set partition A i l at iteration i by the p -dimensional structuring element b . Structuring element b is a set whose structure is defined by the choice of distance metric d , distance measure i and orientation factor β.
Examples of expansion rules for β, which alter the orientation of the structuring element and hence the expansion of the associated set partition for the case of two dimensions, are: where e 1 represents the x coordinate of e and e 2 its y coordinate.
Combinations of these rules are also possible. Different expansions and p -dimensional rules can also be applied. For the MDC implementation in this work, we consider these 4 different orientations and any combination of them, resulting in a total of 16 distinct growing patterns. At low level overview, the f function MDC is given by: Besides, the convergence function of MDC is given by: where P err is a classification error measure and t err is an error threshold, indicating the maximum accepted error prior to the classifier terminating training. σ forces expansions to terminate if i exceeds σ . ∨ represents the Boolean operation OR . The individual error P instErr (y) of an instance y is given by: where T l is a counter of neighbouring instances of y that have been engulfed by the expansion and belong to class l. τ is a scaling factor that, such as in GRASP [44] , allows the expansion to continue even if T l < T u . The total error of the set partition is the sum of P instErr such that: The pseudo-code for generating the MDC classification model is shown in Algorithm 2 .
Data : γ , τ , β and σ are input parameters. γ stands for a weight applied to the central instance, τ controls the convergence as shown in Eq. 30, β a orientation or expansion altering rule, as shown in Eq. 27, and σ is a threshold for the number of partition set expansions. T 1 , . . . , T L are counters for the classes 1 . . . L , X is the training dataset and A 0 l respects Eq. 6. q indicates the complement set partition.
Apply the complement operation to set partition A I q q as described in Equation 13;

Iterative 2D distance
The algorithms presented up to this point consider pdimensional multi-labelled data. In the experiments section, however, we decided to adopt processing p -dimensional datasets as combinations of 2-dimensional models. Several 2D "weak" learners can be created from p -dimensional data and combined to form a p -dimensional classification model. The Random Forest [19] algorithm utilizes a similar concept. It generates several trees with potential reductions on data dimensionality and combines them together to form a general classification model.
Moreover, 2D models are faster to compute and more efficient to store. Not every pairwise combination of attributes has to be explored. Even if they are, it is still more computationally efficient than strict p -dimensional dilations (proven in Theorem 4.1 ). Heuristics can be used to eliminate some of these combinations.
In this section, we introduce a 2-dimensional iterative distance that is efficient when computed using GPUs. The order in which the neighbours are visited in both MkNN and MDC algorithms respects this distance and is shown in Fig. 11 , for the upper-right quadrant, where y is the central instance. The 3 remaining quadrants follow the same ordering pattern.
The grey shaded pixels of Fig. 11 can be obtained using the sum of arithmetic progressions, as shown in Eq. 32 , where n is the index of the number in the shaded column or line. For instance, number 15 is in the 5 th line from bottom to top and also in the 5th column from left to right, thus, n would be equal to 5 in this case. a 1 is the initial value of the progression and r is its rate which, in this case, both are constants equal to 1.
Thereafter, we can state that the numbers to the right of the grey shaded ones are the increments by 1 of their left numbers. So, by fixing the y -axis, the x index can be summed to S n , according to n , to produce the distance value between two elements. This pattern is repeated with respect to the diagonal line of the table, but the axes must be inverted. Thus, by considering these aspects, the where,

Complexity analysis
Theorem 4.1 provides the worst-case complexity and associated analysis of the proposed algorithms assuming that they use the combination of 2-dimensional models to generate a p -dimensional classification model. The number of the 2-dimensional models is equal to the number of pairwise combinations of pattern dimensions. Proof. The neighbourhood of each instance in all implementations is checked respecting d t , which accesses each of the neighbouring elements once, respecting the order given by d t . Since classifiers are run in the GPU, each grid instance is theoretically processed in parallel.
The worst case scenario would be dilating an instance in the image at the top-left position v = (0 , 0) , and its neighbourhood would need to be grown until it covers the entire image. The position w represents the pixel at the right-bottom corner, such that , where x t , y t are the dimensions of the image or grid with respect to the x and y axes. The number of iterations until the image is fully covered is given by the distance function in Eq. 37 . (37) Assuming y t ≥ x t , then We can assume that x t ≥ 1 and y t ≥ 1 to avoid null or empty datasets. In this case, the modulus can be removed, since the result inside them would be positive regardless. | − y t + 1 | outputs values equal or greater than 0. Thus, by assuming that y t is at least 1, we modify our parenthesis to accommodate this assumption. This is rewritten as ( y t − 1 ): For y t < x t the process is exactly the same, and in this case we obtain: Therefore, the worst-case complexity for constructing 2D clas- In other words, the complexity of 2D implementations is O (y 2 The worst case complexity for the ensemble approach, where 2D dimensional models are combined considering every attribute in p -dimensions is: which is the cost for generating a 2D model multiplied by the number of combinations of every two attributes. Expanding the equation: and finally, we obtain: As previously discussed, combining 2D models is much better in practice than computing dilations for p -dimensions directly. If we apply the expansions in 3 dimensions, for instance, we get a cubic complexity of O (x 3 t + y 3 t + z 3 t ) . Usually, y t and x t , which represent the size of the grid in 2 dimensions, are substantially greater than the number of dimensions p . Let us suppose a hypothetical but frequent scenario whereas y t = x t = z t = sp, where s is a scaling factor that relates the size of the grid to the number of dimensions of the problem. The proposed ensemble approach gives p 2 ((sp) 2 + (sp) 2 ) = 2 s 2 p 4 (43) operations. In the case of expanding the set partitions in pdimensions directly, operations are required. It is clear that if s is not small, the number of operations in the second case (3 s 3 p 3 ) is much greater than in the first case (2 s 2 p 4 ). In fact, if s is greater than 2 p 3 , which is the case for every practical classification problem, the number of operations of the latter case is already greater than the former.
In any case, if the number of dimensions p of the problem is equal to 1, which is theoretically the best case, the complexity in the 2D combinations case is 2 s 2 and in the p -dimensional case 3 s 3 .
In the best possible case, it is already clear (2 s 2 < 3 s 3 ), i.e., it is more efficient to use the 2D combinations approach.

2D model combinations in practice
Let us suppose that an arbitrary input dataset has 4 dimensions, labelled as {d 0 , d 1 , d 2 , d 3 }, excluding the class label, such as the original UCI Iris dataset. In this case, we have a total of 6 combinations of 2 pairs of features (d 0 and d 1 , d 0 and d 2 , d 0 and d 3 , d 1 and d 2 , d 1 and d 3 , and d 2 and d 3 ). Each one of these combinations is a 2D dataset. 2-dimensional models are trained with each of these combinations, and hence each combination produces a respective classification or predictive model. In total, we have 6, 2-dimensional classification models. Now, let us suppose that we would like to classify instances using these generated predictive models. Each combination produces a vote for the label of an unlabelled instance. These votes can be combined using different paradigms. A simple approach would be to label unlabelled instances with the label that received the most votes. This is exactly what is done in the p -dimensional experiments in Section 5 .
We consider that some of these 2D models do not produce very accurate votes. Thus, just the votes of the best 2D models are considered, i.e., the ones that achieved the highest accuracies in training and testing with one of the 2D combinations (the training and validation is performed using 10-fold cross validation). The number of 2D models that should have their vote considered was chosen empirically for each experiment. We can fairly say that numbers close to 0 work far better than numbers close to p 2 , as determined empirically.
The number of combinations is related to the number of features in the dataset, which is not related at all to the number of instances. Iris has 4 features and 150 instances. It is possible to increase the number of instances indefinitely without altering the number of features and vice-versa. Thus, the potential increase in complexity using the 2D combinations would be a vote (which is O(1), since in 2D it is possible to use matrices, and a single position in the matrix is visited), multiplied by the number of combinations. The complexity for generating a 2D model is O (y 2 as previously highlighted. Thus, the complexity for generating 6 of these models is given by the multiplication of this value with the number of combinations, which is given by p 2 , where p is the number of dimensions (4 for the Iris dataset). The resultant complexity is O (p 2 (y 2 t + x 2 t )) , which is the worst-case complexity reported in Theorem 4.1 .

2D model combination theory
Multi-label p -dimensional classification is defined as the aggregation of 2 dimensional binary classification models. A multi-label classification, where class labels l belong to { 1 , . . . , L } , is decomposed in L binary problems. For binary problem l = 1 , the classification models predict if an unlabelled instance belongs to class 1 or not. For binary problem 2, the same follows, the associated models are trained to predict if an unlabelled instance belongs to class 2 or not. This is repeated for all L class labels.
The final class label is given by the aggregation of these L votes, respecting the accuracy of each binary problem, obtained using 10fold cross validation during the training phase of the algorithm. The votes of classes that obtain the highest accuracies are computed first. Let us suppose we have the following occasion: Accuracy : 90 . 4% ← label : 1 Accuracy : 70 . 1% ← label : 2 Accuracy : 92 . 6% ← label : 3 (45) At first, since class 3 has the highest accuracy (92.6%), its votes are considered first in detriment of the other classes. That is, if this classification model labels an unlabelled instance with label 3, then it is predicted to be label 3. If this is not the case (i.e., not label 3), then the instance is tested with the model that has the second highesr accuracy, i.e., 90.4%, and the same procedure is repeated. If it is labelled as 1, then it is said to be label 1. Otherwise, it is labelled as 2. For this procedure to work well, true positive rates are expected to be greater than true negative rates in the binary models. Otherwise, a contrary logic should be used instead. Logical combinations of these votes can be used as well.
For each binary labelling problem, the dataset is decomposed in pairwise combinations of features, as previously discussed. p 2 2dimensional classification models are produced, and each 2D combination provides a vote for the binary problem.
A total of p 2 or less votes are considered for each binary labelling procedure. Some votes can be disregarded according to the accuracy of each 2D model, since they may provide weak or bad votes. At first, all 2D models are ranked according to their accuracy obtained during the 10-fold cross validation process in the training phase. Next, a threshold parameter is set, which indicates how many models should be considered in the voting process. Let us consider the case where there is a total of 22, 2D models (20 combinations of pairwise features), and the threshold indicates that just the top 10 models are trustworthy, then the remaining 12 models are disregarded. In this case, for the given binary problem, the algorithm provides 10 votes, each one indicating whether the instance belongs to a class or not. The majority of the votes indicates the predicted label.

Classification model compression
The classification models produced by the proposed approach are, in theory, sets that contain p -dimensional instances. Storage of the model, in practice, can be performed using the said sets or using p -dimensional matrices.

Rectangular compression for sets
In the first case, the classification model is given by: . . . , h p ) 1 , . . . , (h 1 , . . . , h p ) n } (46) where p represents the dimension of the problem and n represents the total amount of instances in the classification model. For p dimensions, pnk bits are required to store the model as a set, as demonstrated below: nk bits . . . , h p ) 1 , . . . , (h 1 , . . . , h p ) n } : pnk bits (47) where k is the amount of bits designed for each variable that, in turn, represent the coordinates of the elements in the set. However, if the rectangles are stored in a compressed notation (e.g., (x : x + θ x , y : y + θ y ) for two dimensions), where θ j represents the size of the rectangle in dimension j, k additional bits are required to indicate this extended notation and a further pk additional bits are required to store all θ j . The required amount of bits is given by: where r represents the total amount of rectangles.
Thus, a single rectangle of size (θ 1 , . . . , θ p ) requires (2 pk + k ) bits. Following this reasoning, r(2 pk + k ) bits are required to store r rectangles. In contrast, storing the same r rectangles using Eq. 47 , assuming p = dimensional rectangles requires r(θ 1 × θ 2 × . . . × θ p ) elements, which is the area of the p -dimensional rectangles multiplied by their number. This results to a r(θ 1 × θ 2 × . . . × θ p ) k bits requirement. The equivalence of these two equations is given by: the compression approach is much more efficient, which is virtually always the case (very small rectangles being the only expection).

Partition trees for matrices
The second option is to consider images or matrices to store the classification models. In this case, another type of compression can be used. We propose the use of 2 p -trees to separate the matrix in p -dimensional quadrants, i.e., orthants. One of the positive aspects of working with matrices is that most implementions of morphological operators are based on matrices or, for the 2D case, images. Fig. 12 shows a 2 p -tree representing a 4x4 binary image. Fig. 13 shows the same 2 p -tree in a flat perspective, where the used orientations are the left-right and top-bottom.
The 2 p -tree representation consists of dividing the entire matrix in 2 p orthants, recursively, until the desired representation precision is reached. In Fig. 12 , the image is divided in 4 quadrants at first. The top-right and bottom-left quadrants contain uniform colors, black and white, respectively, and these nodes are promptly turned into leaves. The top-left and bottom-right quadrants are further divided in 4 quadrants each due to the fact that the pixel values in these quadrants are not uniform. The next level of the tree is the deepest level in this case, where the quadrants cannot be divided any more. The leaves of this tree contain values for the binary case (black and white), however it is possible to generalize the representation for grey-scale images, i.e., multi-label classification models, where the leaves of the tree store the respective grey-level (or label) value.
It is clear that when the tree branches halt prematurely, a storage compression is achieved. In Fig. 12 , two of the nodes of the second level are actually leaves, which were prematurely pruned due to the uniform pattern of the data in their respective quadrants. Thus, the 8 elements of the original image are compressed into 2. The worst-case complexity for reaching any element in the tree is O(log 2 p n ), which is very efficient. The image shown in Fig. 12 is written as: • "XX01X10111011" in depth-first order, which contains 13 values. 'X' represents a branch or parent node. The original image contains 16 values and the compressed version contains 13. In this case, a compression of 3 elements is obtained. A more expressive compression can be observed in large uniform classification models. The classification models obtained with the proposed algorithms, as shown in the results section, are highly prone to this type of compression, i.e., they are uniform in most regions.
For p -dimensional classification matrices, each element in the string representation must represent all possible classes in the problem, i.e., k must respect 2 k ≥ L , where L is the amount of possible distinct classes.

Run length encoding
Finally, a text-compression method, such as run-length encoding [2,7,45] can be applied to both approaches, i.e., the one based on sets and the one utilizing p -dimensional matrices, before storing the classification model. The set, as well as the 2 p -trees, can be written as a single string and therefore text-compression, in general, and run-length encoding, in particular, apply.
Run-length encoding consists of replacing characters repetitions with a unique occurrence of the character along with the the number of times that the character occurrence happens in sequence. The sequence "AAAABBC" can be written as "3A2BC", which requires less memory. Run-length encoding can also be applied along with some pattern recognition approaches, where patterns that occur most frequently are replaced by a predefined set of patterns that requires less memory. In summary, a number of approaches may be used to mitigate the memory burden of the proposed classification model, despite the fact that this is by and large resolved by using the combinations of 2D models instead of the direct pdimensional approach. An appropriate discretization/quantization of the classification space also mitigates the storage requirements.

Parameter selection
The performance of classifiers depends on the chosen set of input parameters. In order to fairly compare the performance of morphological classifiers to other classification algorithms and to visu-alize the obtained results, we performed two separate sets of experiments. At first, experiments with 2 dimensional datasets were performed. That is, two p -dimensional datasets were reduced to 2 dimensions, and the parameters of each classifier were tuned using an evolutionary algorithm (EA) [46] .
Parameters were optimized using two methods. The first approach consists of using EA for tuning of numerical parameters, while categorical attributes are manually configured. The second approach uses Auto-Weka [17] , which determines the best suited algorithm from the Weka [18] machine learning framework and an appropriate set of parameters, including a subset of the dataset attributes, if it leads to a performance improvement.

Evolutionary algorithm
The proposed EA consists of a combination of an elitist strategy, mutation and cross-over. The size of the population carried on to the next generation was set to 20. At first, 20 individuals are randomly generated. In each generation, a single individual is created; this receives randomly selected parameters of two randomly selected parents (cross-over). After that, mutation is performed on the selected parameters. This individual has its objective function computed to confirm whether it may enter in the top 20 members of the population. If this is true, then the individual with the smallest objective function value among the top 20 is removed. The objective function is defined as the accuracy of the classifier on 10-fold cross validation. Algorithm 3 illustrates the overall pro-Data : is the population at each generation, () is a function that generates a random number ∈ [0 , 1] , ε represents individuals and represents a maximal value of each parameter, chosen according to the classification algorithm. cess of the EA. Each EA run was limited a total of 12 hours. This was also the case for Auto-Weka. Since MkNN and MDC were faster than the algorithms in Weka, for these cases, the genetic EA was set to converge in 6 hours.

Auto-Weka
As opposed to EA, the results of Auto-Weka were not manually tweaked. Although Auto-Weka performs a 10-fold cross validation internally, it presents its indexes (accuracy, true positives, etc) tested on the entire training set. Therefore, the reported results are somewhat unfair with regards to the first methodology (EA plus manual tweaking), where the results are presented as tested on the 10 partitions of the 10-fold cross validation. The authors of Auto-Weka reported that this was defined due to how Weka works internally [17,47] . If the final results reported by Auto-Weka were derived from cross validation alone, the reported indexes would be lower in comparison.

Experimental results
Weka [18] , a Java based machine learning framework, was used to run and test the classification algorithms apart from MDC and MkNN. MkNN and MDC were implemented in CUDA/C/C++ and their source code is available at [48] . All the implementations proposed in this work conform to the use of the complement approach shown in Theorem 3.1 due to the previously described benefits.

2-dimensional experiments
This section addresses the details of datasets used in our experiments as well as the performance of the classifiers, visual results and comparisons. The parameters were selected, as described in Section 4.7 . We performed experiments with 2-dimensional subspaces of p -dimensional datasets. In this regard, predictive models can be visualized and analysed as an image. In the next section, experiments with p -dimensional datasets are performed.
The first implementation of the proposed classifiers considers input datasets as multisets, that is, repeated instances, i.e., instances where the feature values are identical, influence the set expansion. The second implementation excludes these repetitions. The subscript Rep in the name of the algorithms indicates the version that uses the multiset concept. As a remark, combining 2D models is more efficient than using p -dimensional dilations, as shown in Section 4.4 .

Used datasets
Two reduced versions of datasets from the UCI repository [49] were considered in this set of experiments. The first one is the widely known Iris dataset. The number of attributes of this dataset was reduced to 2, i.e., attributes sepallength and petallength were randomly selected and the remaining features were excluded, except for the class label. Instances whose classes were equal to Irissetosa were also removed from the dataset in order to transform this dataset into a binary classification problem.
The second dataset was Diabetes , and the selected attributes were plas and insu . The original dataset already consists of a binary classification problem and thus no instances were removed. Table 1 illustrates this information and the total number of instances of each dataset.

Iris dataset
The results obtained with the Iris dataset and each algorithm (optimized with EA along with manual adjustments) as well as the Auto-Weka approach, which is highlighted in italic -refer to Section 4.7 for more information, are shown in Table 2 .
The parameters of the morphological classifiers were also automatically selected using the EA approach ( Section 4.7 ). Among Table 2 Performance comparison on the Iris dataset.

Classifier
Acc (%) TP (%) TN (%) Multilayer Perceptron [50] 95 92 98 RBF Network [21] 94 94 94 Conjunctive Rule [50] 94 92 96 SPegasos [50] 93 92 94 k-NN (IBk) [13] 93 94 92 SVM [51] 93 92 94 SMO [25] 93 92 94 Auto-Weka [17] 93 92 94 REPTree [50] 92 88 96 Random Forest [19] 91 92 90 Decision Table [ Rep stands for the version of the algorithm that considers instance repetitions (multiset). The algorithms highlighted in bold are the ones proposed in this work. the classifiers in comparison, Auto-Weka selected Logistic Model Trees (LMT) [57] as the best suited approach. The results were obtained using 10-fold cross validation (except for Auto-Weka, refer to Section 4.7 ). Apart from morphological classifiers, Multilayer Perceptrons achieved the highest accuracy, i.e., 95%. Three implementations of the proposed algorithms also achieved 95% of accuracy with this dataset (MkNN, MDC and MDC Rep ). In most cases, MkNN and MDC produce high true positive rates, usually higher than the remaining evaluated classification algorithms while maintaining high accuracy (this is more evident with the MDC algorithm). Fig. 14 shows the utilised 2D Iris dataset after the preprocessing described in the beginning of Section 4 . The classes that are equal to 1 are shown as red circles while the ones that are equal to 2 are depicted as green squares. 10-fold cross validation separates the dataset in 10 partitions. Next, 10 iterations follow (the folds). In the first iteration, a classification model is created trained on 10 − 1 of the partitions and tested on the remaining 1. This process is repeated 9 more times while always switching the test partition to a partition that has not been previously used as test set. Fig. 15 shows one of the classification models of the 10 iterations (or folds) of the Iris dataset shown in Fig. 14    The dark grey area in Fig. 15 represents what does not belong to A I 1 1 , while the light grey represents instances that belong to A I 1 1 . Assuming that repeated instances are discarded and conceiving the fact that this dataset almost does not contain repeated instances, then the larger the area where the light gray intersects instances labelled as 1 (red circles) and excludes the ones labelled as 2 (green squares), the higher the accuracy of the model. Fig. 16 shows the 2D Diabetes dataset. Fig. 17 , on the other hand, shows the classification models of one fold of each proposed algorithm. In the case of MDC Rep , the orientation parameter β chosen by the EA was equal in all directions but the left. Therefore, the dilator clearly grows the structure in the remaining 3 directions. In the case of MDC, β restricted the growth downwards. Table 3 shows the obtained results on the 2D Diabetes dataset. In this case, MkNN Rep and SMO outperformed the remaining algorithms, and tied on the accuracy. However, MDC still produces a high true positive rate. However, the accuracy obtained with MDC was 71.61%, which is considerably less than the 74.86% obtained by SMO. Auto-Weka, in this case, also selected LMT as the best classifier.

P-Dimensional Experiments
In this subsection, experiments and comparisons are performed with p -dimensional binary and multi-label datasets. We considered 8 datasets from the UCI public repository. The proposed pdimensional approach uses a combination for 2-dimensional clas-  Rep stands for the version of the algorithm that considers repeated instances (multiset). The algorithms highlighted in bold are the ones proposed in this work. sification models, due to previously described advantages (see Theorem 4.1 ). Each pairwise combination of features generates a binary predictive model. The total number of models is given by p 2 , where p represents the number of features or dimensions of the dataset. An unlabelled instance is classified using the votes of all p 2 models, where each vote corresponds to a binary vote, belonging or not to a class.
Multi-label datasets were separated in L binary problems, where L is the number of labels of the dataset. That is, binary classifiers were trained for each of the L classes. Let us consider   Next, for the binary problem of class 1, the approach follows as previously described. p 2 binary classification models are generated for the class 1 problem. Each model provides a vote for an unlabelled p -dimensional instance. In total, we have p 2 or less votes (depending on the significance of the model/vote), which decide whether the instance belongs to class 1 or not. The procedure is repeated for the other two classes, i.e., class 2 and class 3. At the end, each binary problem (class 1, class 2 and class 3) provides a vote for each one of the class labels. The majority of votes decides the predicted label of the unlabelled instance. Section 4.5 describes this decomposition logic in more detail. Table 4 contains the accuracies obtained with each classifier for each dataset. Table 5 shows the training and testing times of MDC and MkNN. Training times include the training phase and the time required for building the 2D models for all pairwise feature combinations as well as the multi-label training, as applicable for the dataset, considering 9/10ths of the original dataset. Test times represent the time required for testing 1/10th of the dataset (one fold of the 10-fold cross validation). Although training phases require additional time, testing times are remarkably fast due to the nature of the classification model, which are matrices whose elements are accessed with O(1) complexity in the worst case.

Conclusion
This work proposes a new category of classification algorithms called Morphological Classifiers (MC). MCs are an aggregation of mathematical morphology and supervised learning. The classifiers rely on the fact that datasets may contain well defined shapes and structures, which naturally come under the remit of mathematical morphology and could be a positive aspect to explore. In addition to establishing the theoretical framework of MCs, we propose two GPU implementation approaches.
This work is the first to propose such a type of classifier, which was shown to outperform some well established classifiers in the performed experiments. Experiments presented herein focused on proposing and evaluating 2-dimensional models and further aggregating them to form p -dimensional classification models. MDC and MkNN classifiers tied or outperformed 14 well established classifiers in 5 out of 8 p -dimensional UCI datasets. In all of these cases, the proposed classifiers obtained accuracies that were higher than the average accuracy across all classifiers.
Although the proposed theoretical framework and algorithms readily work for p -dimensional datasets, we proved that it is faster, in terms of required steps in the worst-case scenario, to generate 2-dimensional classification models and to combine them, rather than working directly on p -dimensional grids and creating p -dimensional classification models using p -dimensional dilations. These 2-dimensional models are combined using a voting scheme, as an ensemble classifier.
The complexity for checking whether an instance belongs to a certain class using matrices is O p 2 in the worst case, where p stands for the number of features or dimensions of the dataset. This is a very fast instance classification that does not relate to the number of instances in the dataset. The classification time is potentially faster than decision-tree based classification times, which are already fairly fast in comparison to other classification paradigms.
Drawbacks include the fact that training morphological classifiers may be slow, particularly when dilations are performed in p -dimensions directly, for instance. However, due to the nature of the operations performed, such classifiers are amenable to parallelization. Moreover, in a fashion similar to the training of deep learning networks, training times may not be considered a major disadvantage, as the actual classification process is instantaneous and classification models demonstrate high accuracy rates. As a final remark, the code for the proposed classifers in CUDA/C/C++ is available as open-source at [48] along with the utlized datasets.

Future work and discussion
Other opportunities to explore in regards to morphological classifiers include introducing more rules for the expansion process of MDC. The current MDC proposal is fairly naive, and this simplicity certainly introduces more classification errors. Growing in a fractal fashion was also evaluated and produced acceptable results. However, the fractal approach was too simple and did not outperform most classifiers. However, the idea of exploring fractal information in datasets and using it in classification problems has not been fully explored. The morphological classifiers framework allows for this possibility.
Another direction for further research is to explore algorithms that are dilated in p -dimensions using p -dimensional morphological operations such as p -dimensional dilations, conditional thickening, etc. We performed a p -dimensional dilation experiment with the Iris dataset, and the dilator algorithm obtained the same accuracy as the one obtained with the 2-dimensional models combination approach, shown in Table 4 . Direct p -dimensional approaches are possible but require more memory and processing power. Thus, future research may focus on reducing this burden by using approximations or real-time functions that generate the classification model in classification time, hence enabling p -dimensional morphological approaches to be performed in real time for a broad class of datasets.
Furthermore, data compression techniques can be applied to reduce the burden of memory requirements. We intend to deeply explore the application of such techniques and to propose efficient ways of storing classification models for big datasets in future work. It is necessary to highlight, however, that the 2D combination approach works well even with reasonably large datasets, as the memory and processing burdens are mitigated. Testing and training times using the combination approach can be found in Table 5 . Morphological functions could also be used. In this case, it is not necessary to store a discretized classification space. These functions could define whether an instance belongs to a certain area, and hence, to a class. The key question here is whether it is viable to merge different functions and to construct such models during the classification process. This approach will also be further explored in future research.