A Metaheuristic Approach to the C1S Problem

Given a binary matrix, finding the maximum set of columns such that the resulting submatrix has the Consecutive Ones Property (C1P) is called the Consecutive Ones Submatrix (C1S) problem. There are solution approaches for it, but there is also a room for improvement. Moreover, most of the studies of the problem use exact solution methods. We propose an evolutionary approach to solve the problem. We also suggest a related problem to C1S, which is the Consecutive Blocks Minimization (CBM). The algorithm is then performed on real-world and randomly generated matrices of the set covering type.


Introduction
The problem of consecutive ones submatrix on a binary matrix has been known since the 1950s. It was suggested by Fulkerson and Gross [1] and described as follows: Let be an incidence matrix, and then can we rearrange its columns so that each row has a single block of ones? 1.1 Basic concepts of C1P and C1S Definition 1.1. A block of 1's (block of 0's) in a binary matrix A is any maximal set of consecutive one entries (zero entries) appearing in the same row [2]. A binary matrix A is said to have the Consecutive Ones Property (C1P) if the ones in every row appear consecutively [3].

An illustration of C1S and current solution approaches 2.1 The seriation problem
The incidence matrix A is an input to this problem. It is a binary matrix with m rows and columns (grave sites and object types, respectively); the entry of A has a value of one if grave has any object of type and zero otherwise. If the rows (graves) were in the suitable time successive order, 1 A problem belongs to the class NP if a solution to it can be generated and verified quickly, i.e. in polynomial time. It is believed that NP-hard problems are intractable, i.e. that there is no efficient algorithm to solve them, [7]. , then each column (object) would indicate a time period during which that object was common. Assume that each object represents the time interval from the appearance of the object to the time when it becomes uncommon. The chronological ordering of the graves is reached by getting rows permutation of that transforms it into Petrie form. This will minimize the total temporal range of an object summed over all the varieties. This range for an artifact type is the span from the first to the last appearance of a "1" in the column of corresponding to that artifact type.
Kendall [19] introduced a solution by using a Similarity Matrix of size for graves, which is nonnegative and symmetric. He found out that a permutation P can convert matrix A into a Petrie form (i. e., ), as it can transform the grave matrix Q into Robinson form (Robinson Matrix 1 , i. ). The final matrix has the C1P where each column displays the correct order of the artifact type from its appearance to it disappears.

Polynomial-time local-improvement heuristic for CBM
The problem of CBM-decision is NP-complete even if it is limited to binary matrices containing two ones in each row [13,23]. However, Haddadi introduced a polynomial-time heuristic which gives a permutation such that the optimum solution and the consecutive blocks do not differ by more than 50%.
Haddadi and others found a polynomial-time local-improvement heuristic for the CBM. They introduced two ( ) size local neighborhood search, where the blocks number of a neighbor is provided in ( ) operations [2,3].

A GA-based solution approach
The GA has been used in related cases. Here, we develop an implementation to deal with C1S. GA is a meta-heuristic inspired by the process of natural selection that belongs to the larger class of evolutionary algorithms (EA). It was proposed by Holland [24]. Genetic algorithms are commonly used to generate high-quality solutions to optimization and search problems by relying on bio-inspired operators, such as crossover, mutation, and reproduction [25].

Solution representation
The basic is to keep permuting the columns of the given binary matrix to get as many consecutive ones as possible in every row. An integer is used here to represent a solution or chromosome of size n, an array of the column indices of the matrix. The position of each gene in the chromosome corresponds to a column in this matrix. Figure-1 illustrates this representation.

Genetic operators
After randomly generating the initial population, we breed successive generations of offspring. This can be achieved by performing genetic operators which are crossover, mutation, and reproduction. Crossover operator. We use the 'Order Crossover' of Gen and Cheng [26] to breed a new child. This operation begins by picking a substring (subchromosome) from one of the two parents randomly, then copying it into its corresponding location in the child. The genes in the second parent that appeared in the substring are deleted to prevent repeating the gene. In the end, the genes are placed into the unfixed location of the child from left to right, relying on the order of the sequence (Figure-2).
Mutation operator. It is used to allow the genetic algorithm to search in the global solution space and prevent trapping in local optima, via changing one or more genes. Here, two genes from an individual are randomly chosen and then swapped. Figure-3 illustrates the mutation operator. 1 Robinson matrix is a matrix with entries which do not increase as one progresses along a row beyond the main diagonal and do not decrease as one continues to progress along that row towards the main diagonal. Reproduction. Without any modification, some fit chromosomes are copied via this operator into the next generation.

Fitness functions
The matrix does not include enough information to decide whether it has the property or not. However, we will use the entries of the matrix to build a fitness function. Many fitness functions are tested and the most reliable one with respect to the number of C1S columns and run time is selected. Fitness Function FF1. It computes the number of 1's in each row using a simple formula. Summing the number of 1's for all the rows will give the fitness value of the matrix. Suppose that we have the following rows with four ones each row 1 = (1 0 1 1 1 0), row 2 = (1 0 1 0 1 1). The first row has two blocks of 1's, k 1 , and k 2 , where k 1 is one element and k 2 is three 1's. We deal with every block as a sequence of 1's. The associated series is defined as the ordered formal sum ∑ where is the length of the block. Counting the elements for the two blocks will be as ∑ ∑ ( ) Applying the formula for the two rows gives fitness values of 7 and 5, respectively, where ∑ ∑ and ∑ ∑ ∑ . Applying this formula for a matrix gives different fitness values for every row. The fitness value for the whole matrix is found by accumulating the fitness values of all the rows, as shown in the next formula The larger the fitness, the better the matrix with respect to consecutive ones blocks. With different permutations of columns in each generation and computing the fitness function of the matrix, the consecutive solutions are improved every time by pushing the ones together to form less blocks. In the final generation, it is expected that a matrix with best consecutive ones in every row is obtained. Fitness Function FF2. This function relies on counting the ones in the different regions of the matrix. That is, by finding the label of the connected components of the matrix, where each maximal connected region is assigned as a unique label [27,28]. The matrix below has three regions of ones. The first connected components in the matrix has elements (label 1), the second has three elements (label 2), and the last has one element (label 3), as shown in the region matrix. We also deal with every region as a sequence of 1's, like the fitness FF1. The fitness values of the regions are 15, 6, and 1, respectively, which means that the fitness function for the matrix is 22 This function pushes the ones in different separated regions, which may produce a variety of fitness values more than those of FF1. It accumulates the ones, not only in the rows, but also in the columns. This helps to accumulate the ones in large blocks. Function FF1 may give similar fitness values for many rows, and consequently similar fitness values for different matrices. This causes that the worse matrix with less columns having the C1P property may be selected for the next generation. With FF2, there is a higher chance to distinguish between matrices. It may give better results than the FF1 function. Although it gives good results, unfortunately it costs more time than FF1. Thus, it is not the desired function for the proposed algorithm. Fitness Function FF3. This function counts the zeros, instead of the ones, in every row and then for the whole matrix. The same formula that is used for finding the maximal sum of 1's in FF1 is used here. Maximizing the fitness value leads to aggregating the zeros in blocks, which results in minimizing the number of gaps between the ones. This fitness almost gives the same result but with longer time, because the matrix is sparse so the calculations cost longer time than FF1. As we seek better results with shorter time, this fitness function is also not suitable for our implementation. Fitness Function FF4. The fitness function that Gargano and Lurie used for solving the Petrie Seriation Problem [22] is used to find the consecutive ones in the columns by permutation of the rows. Although it gives good results for small matrices, they did not mention whether it is suitable for large matrices with different densities. Thus, we decided to test it. Their fitness is evaluated with the so called Petrie Range Index (PRI) of the permuted binary matrix. We formulate the fitness for our matrix by choosing the first column and the last column , where in each row there is at least one element. The PRI for that row is . The fitness for the matrix is the sum of these values over all the rows. The best solution is reached through minimizing the PRI. It is noticeable that if the matrix has the same number of ones in each column, then the fitness value for the matrix will be the same for any permutation. As a result, good matrices may not cross to the next generation, which causes not having very good results. This approach is more suitable for matrices with high density with different number of ones in each column. It generates less columns with C1P but in shorter time as compared to the other functions. It is formulated as follow ∑ Here again, FF4 is not suitable for our implementation of GA. Fitness Function FF5. This function is meant to solve the CBM problem which is equivalent to C1S. The aim is to permutate the matrix columns to minimize the number of blocks of consecutive ones. Building the function relies on counting the blocks of the ones in each row, then for the whole matrix. The same two rows, and that are used in FF1 function for illustrations, are used here. The first row has two blocks of ones and the second one has three blocks . The formula for the matrix is ∑ ∑

Abo-Alsabeh and Salhi
Iraqi Journal of Science, 2021, Vol. 62, No. 1, pp: 218-227 223 Applying the formula for the two rows gives values of 2 and 3, respectively. The fitness value for the whole matrix can be obtained by summing up the fitness values of all the rows. This fitness also pushes the 1's together to make blocks. The smaller the fitness value, the lesser number of blocks the matrix will have, thus producing a larger C1S submatrix.

Stopping criterion
A common criterion used to stop the proposed genetic algorithm is the maximum number of generations.

Selection procedure
It usually implements a roulette wheel which is biased towards fit individuals. This means that good individuals are likely to be parents.
The proposed GA algorithm is shown in a pseudocode as shown in Algorithm 1.

Algorithm 1: Algorithm C1P
1 Input a binary matrix A; 2 Input the rate of crossover and the rate of mutation; 3 Generate a random population of permutations of columns of matrix A; 4 Evaluate the fitness of the permutations according to some fitness function; 5 Rank individuals according to their fitness; 6 Select parents from the population according to some selection procedure; 7 Generate a new population by applying the following operators: crossover, mutation, and reproduction; 8 Compute the fitness of the individuals of the new population; 9 Until (The sopping criteria are satisfied) Repeat from 5. 10 Return Permuted matrix.

Testing the quality of the fitness functions
To compare the quality of the fitness functions that are stated above, we implement them separately in Algorithm1. They are applied to the same matrices with a fixed number of generations and initial population. We run the GA on 5 different matrices for each size. The size of the population and the generation are fixed to 40 for the first two matrices and 100 for the rest. Table-1 records the number of columns (Nbcols) with C1P rounded to the nearest integer. We refer to the matrices size as (Mat.) and to the density as (Dens.). The table shows that FF1 and FF5 produce better results in a shorter time. Although both of them give almost the same results, the latter is superior in CPU time.

Computational experience 3.1 Results of Algorithm 1 on the CBM problem
This algorithm is implemented for matrices of the set covering problem with different densities. We use these matrices to serve as a test bed to run the GA implementation. The results shown in Table-2 are obtained by applying Algorithm 1 to ten different matrices of each size. We run the GA 10 times on every matrix, each time with a random initial population. In fact, the results of the algorithm rely on the arrangement of the binary matrix and, therefore, it may produce better results when it is repeated with different initial populations. The run time depends slightly on the density of

Abo-Alsabeh and Salhi
Iraqi Journal of Science, 2021, Vol. 62, No. 1, pp: 218-227 224 the matrix and the number of rows, but more on the number of columns. We set the number of generations and population size to 100. The heuristic gives the results of the number of blocks and columns. The first five columns have the initial information of the matrices with the number of generations (Gen) in column 3. The rest of the table shows the number of final blocks with the number of improved blocks, the number of columns (Nbcols) with C1P, and the time. For small matrices, the number of columns with C1P is good, but for matrices with size the results are not up to expectation. From our experience, the results can be improved by the following: 1. Increasing the number of generations and the initial population size. 2. For large matrices, it is better to have an initial population which covers more than half the columns of the matrix to obtain about the same number of columns having the C1P. This reduces computing time. This can be achieved by reducing the chromosome that presents the matrix to the half. 3. Dividing the matrix into submatrices and applying Algorithm 1, separately, then applying it for the whole matrix put together. This however requires more time. 4. The initial population can be seeded by applying Algorithm 1 many times, then making the final population from the best of each run and use it as the initial population. This also computationally expensive.
Overall, this algorithm is performed to different instances of real-world and randomly generated matrices. The nonsymmetric matrices that are created from the set covering problem are not checked for the consecutive ones property, so the optimal solutions are not known. Consequently, we cannot discuss the quality of our results. The rest of matrices, Real-world data, ( -) whose sizes are given in Table-3, arise from the problem of stop location, posted by a German railway company [29]. The problem is formulated as a Set Covering Problem (SCP) [29,30]. These cases produce matrices with 0, 1 entries, that are assumed to have almost C1P. Other small size matrices of types B and C, which are randomly created by Ruf and Schöbel [29], are both sparse and almost have the C1P with density values of 3% and 5%, respectively (Table-3). Finally, the algorithm is implemented as follows: The results of Table-4 show that the columns of the B and C matrices from Ruf and Schöbel could not satisfy the C1P. Also, applying this procedure alone is not enough to improve the number of consecutive blocks.
The results on the real-world data of Ruf and Schöbel are illustrated in Table-5. The first matrix R 1km almost satisfies the CBM, since the number of final blocks is near to , and thus optimal. With respect to the rest of matrices of this type, although the final numbers of blocks and columns differ from the optimal values, they are close to the lower bound , but not to the upper . However, it still produces large submatrices with C1P. From the results shown in the three tables, we can say that: 1. Performing Algorithm 1 improves the number of blocks but not that of the columns with C1P.

2.
Originally, any matrix should have C1S submatrix of at least two columns. The result from Algorithm 1 does not rely on the size of the C1S only, but also on the structure of the matrix. 3. Minimizing the number of blocks does not imply finding a large C1S submatrix. We can say that improving the C1S can improve the CBM, but the converse is not always true. This is because the size of the C1S submatrix relies on the position of the destructive column. 4. Over all, the number of consecutive blocks from the algorithm in Tables-(3, 4 and 5) show a good improvement, While the number of columns with C1P is not as expected.

Conclusions
We presented a metaheuristic method for solving the C1S problem. A basic GA is proposed with many new fitness functions and FF5 is chosen to solve the C1S problem. The minimum consecutive blocks or CBM in [3] is also solved using the GA. We applied our algorithm to a large number of randomly generated matrices and real-world instances. The results show that large submatrices with C1P can be found for matrices with small sizes. However, since the optimum solutions are not known, it is not possible to say how far the solutions returned by our approaches are actually resulted from them. Finally, we can say that the GA is a suitable algorithm for solving the C1S problem if we want to separate a small submatrix with C1S.