{'fraction_non_alphanumeric': 0.027992613803623518, 'fraction_numerical': 0.012470896780581798, 'mean_word_length': 4.870856245090338, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Optimal Decision Trees for the Algorithm Selection Problem: Integer Programming Based Approaches

Even though it is well known that for most relevant computational problems different algorithms may perform better on different classes of problem instances, most computational experiments still focus on determining a single best algorithm configuration based on aggregate results such as the average. In this paper, we propose Integer Programming based approaches to build decision trees for the Algorithm Selection Problem. These techniques allow to automatically: (i) find the most important problem features to determine problem classes; (ii) group the problems into classes and (iii) select the best algorithm configuration for each class. To evaluate this new approach, extensive computational experiments were executed using the linear programming algorithms implemented in the COIN-OR Branch&Cut solver in a comprehensive set of instances, including all MIPLIB benchmark instances. The results exceeded our initial expectations. While the single best parameter setting discovered decreased the total running time by 22%, our approach decreased the total running time by 40% in average in 10-fold cross validation experiments. These results indicate that our method generalizes quite well and does not overfit.


Introduction
From the observation that different algorithms may perform better on different classes of problems, Rice (1976) proposed a formal definition of the Algorithm Selection Problem (ASP). The main components of this problem are depicted in Fig. 1. Formally the ASP has the following input data: P : the problem space, a probably very large and diverse set of different problem instances; these instances have a number of characteristics, i.e. for linear programming, each possible constraint matrix defines a different problem instance; A : the algorithm space, the set of available algorithms to solve instances of problem P; since many algorithms have parameters that significantly change their behavior, differently from Rice (1976), we consider that each element in A is an algorithm with a specific parameter setting; thus, selecting the best algorithm also involves selecting the best parameter tuning for this algorithm; F : the feature space; ideally elements of F have a significantly lower dimension than elements of P, since not every problem instance influences the selection of the best algorithm; these features are also important to cluster problems having a common best algorithm, e.g., for linear programming some algorithms are known for performing well on problems with a dense constraint matrix; W : the criteria space, since algorithms can be evaluated with different criteria, such as processing time, memory consumption and simplicity, the evaluation of the execution results r = R n produced using an algorithm a to solve a problem instance p may be computed using a weight vector w ∈ W = [0, 1] n that describes the relative importance of each criterion.
The objective is to define a function S that, considering problem features, maps problem instances to the best performing algorithms. This function is a mapping function that always selects the best algorithm for every instance. Thus, if B is the ideal function, the objective is to define S minimizing: p∈P |w T r(B(f (p), w)) − w T r(S(f (p), w))| (1) The main motivation for solving the ASP is that usually there is no best algorithm in the general sense: even though some algorithms may perform better on average, usually some algorithms perform much better than others for some groups of instances. A "winner-take-all" approach will probably discard algorithms that perform poorly on average, even if they produce excellent results for a small, but still relevant, group of instances. This paper investigates the construction of S using decision trees. The use of decision trees to compute S was one of the suggestions included in the seminal paper of Rice (1976). To the best of our knowledge however, the use of decision trees for algorithm selection was mostly ignored in the literature. One recent exception is the work of Polyakovskiy et al. (2014) who evaluated many heuristics for the traveling thief problem and built a decision tree for algorithm recommendation. Polyakovskiy et al. (2014) did not report which algorithm was used to build this tree, but inform that the MatLab R Statistics Toolbox was used to produce an initial tree that was subsequently pruned to produce a compact tree. This is an important consideration: even though deep decision trees can achieve 100% of accuracy in the training dataset, they usually overfit, achieving low accuracy when predicting the class of new instances. The production of compact and accurate decision trees is an NP-Hard problem (Hyafil and Rivest, 1976). Thus, many greedy heuristics have been proposed, such as ID3 (Quinlan, 1986), C4.5 (Quinlan, 1993) and CART (Breiman et al., 1984). These heuristics recursively analyze each split in isolation and proceed recursively. Recently, (Bertsimas and Dunn, 2017) proposed Integer Programming for producing optimal decision trees for classification. Thus, the entire decision tree is evaluated to reach global optimality. Their results showed that much better classification trees were produced for an extensive test set. This result was somewhat unexpected since there is a popular belief that optimum decision trees could overfit at the expense of generalization. Trees are not only the organizational basis of many machine learning methods, but also an important structural information (Zhang et al., 2018). The main advantage of methods that produce a tree as result is the interpretability of the produced model, an important feature in some applications such as healthcare.
At this point it is important to clarify the relationship between decision trees for the ASP and decision trees for classification and regression, their most common applications. Even though the ASP can be seen as the classification problem of selecting the best algorithm for each instance, this modeling does not capture some important problem aspects: firstly, it is often the case that many algorithms may produce equivalent results for a given instance, this aspect can be fixed by working with multi-label classification algorithms (Tsoumakas and Katakis, 2007). Secondly, the evaluation of the individual decisions of a classification algorithm always returns zero or one for wrong and correct predictions, respectively. In the ASP each decision is evaluated according to a real number, which indicates how far the performance of the suggested algorithm is from the performance of the best algorithm, as stated in (1). Thus, the construction of optimal decision trees for the ASP can be seen as a generalization of the multilabel classification problem and is at least as hard. Another approach is to model the ASP as a regression problem, in which one tries to predict the performance of each algorithm for a given instance and the best one is selected (Xu et al., 2008;Leyton-Brown et al., 2003b;Battistutta et al., 2017). In this approach, the cost of determining the recommended algorithm for a new instance grows proportionally to the number of available algorithms and the performance of the classification algorithm. By contrast, a decision tree with limited depth can recommend an algorithm in constant time. Another shortcoming of regression based approaches is related to the loss of precision in solution evaluation: consider the case when the results produced by a regression algorithm are always the correct result plus some additional large constant value. Even though the ranking of the algorithms for a given instance would remain correct and the right algorithm would always be selected, this large constant would imply an (invalid) estimated error for the regression algorithm. The present paper therefore investigates the applicability of Integer Programming to build a mapping S with decision trees using the precise objective function (1) of the ASP.
It is also important to distinguish the ASP from the problem of discovering improved parameter settings. Popular software packages, such as irace (López-Ibáñez et al., 2016) and ParamILS (Hutter et al., 2014a) embed several heuristics to guide the search of improved parameters to a parameter setting that, ideally, performs well in a large set of instances. In this search, parameter settings with a poor performance in an initial sampling of instances may be discarded. In the ASP, even parameter settings with poor results for many instances may be worth investigating since they may be the best choice to a small group of instances with similar characteristics. Thus, exploring parameter settings for the ASP may be significantly more computationally expensive than finding the best on average parameter setting, requiring more exhaustive computational experiments. An intermediate approach was proposed by Kadioglu et al. (2010): initially the instance set is divided into clusters and then the parameter settings search begin. One shortcoming of this approach is the requirement of an a priori distance metric to cluster instances. It can be hard to decide which instance features may be more influential for the parameter setting phase before the results of an initial batch of experiments is available. Optimized decision trees for the ASP provide important information regarding which instance features are more influential to parameter settings since these parameters will appear in the first levels of the decision tree. Also, instances are automatically clustered in the leaves. It is important to observe that an iterative approach is possible: after instances are clustered using a decision tree for the ASP, a parallel search for better parameters to instance groups may be executed, generating a new input for the ASP and so on.
Another fundamental contextual information is that the ASP is a static tuning technique: no runtime information is considered to dynamically change some of the suggested algorithm/parameter settings, as in the so called reactive methods (Mascia et al., 2014). The static approach has the advantage that usually no considerable additional computational effort is required to get a recommended setting, but its success obviously depends on the construction of a sufficiently diverse, generally large, set of problem instances for the training dataset to cover all relevant cases. After the assembly of this dataset, a possibly large set of experiments must be performed to collect the results of many different algorithmic approaches for each instance. Finally, a clever recommendation algorithm must be trained to determine relevant features for recommending the best parameters for new problem instances. Misir and Sebag (2017) tackle the problem of recommending algorithms with incomplete data, i.e., if the experiments results matrix is sparse and only a few algorithms were executed for each problem instance. In this paper we consider the more computationally expensive case, where for the training dataset all problem instances were evaluated on all algorithms. This paper propose the construction of optimal decision trees for the ASP using Integer Programming techniques. To speed up the production of high quality feasible solutions, a variable neighborhood descent based mathematical programming heuristic was also developed. To validate our proposal, we selected the challenging task of improving the performance of the COIN-OR Linear Programming Solver -CLP, which is the Linear Programming (LP) solver employed within the COIN-OR Branch & Cut -CBC solver (Lougee-Heimer, 2003). CLP is currently considered the fastest open source LP solver (Mittelmann, 2018;Gearhart et al., 2013). The LP solver is the main component in Integer Programming solvers (Atamtürk and Savelsbergh, 2005) and it is executed at every node of the search tree. Mixed-Integer Programing is the most successful technique to optimally solve NP-Hard problems and has been applied to a large number of problems, from production planning (Pochet and Wolsey, 2006) to prediction of protein structures (Zhu, 2007).
To the best of our knowledge, this is the first time that mathematical programming based methods for solving the ASP are proposed and evaluated computationally in literature. As our results show, not only our algorithm produces more accurate predictions for the best algorithm on unknown instances, considering a 10-fold validation process (Section 4.3), but also has the important advantages of recommending algorithms in constant time and producing easily interpretable results. This paper is organized as follows: Section 2 presents the Integer Programming formulation for the construction of optimal decision trees for the ASP. Section 3 presents a variable neighborhood descent based mathematical programming heuristic. Our extensive computational experiments and their results are presented in Section 4 and, finally, Section 5 presents the results and 5 discusses the final remarks.

Optimal Decision Trees for the Algorithm Selection Problem
This section presents our integer programming model proposed for the construction of optimal decision trees for the ASP. The sets and parameters are described in Section 2.1. The corresponding decision variables are described in Section 2.2. The objective function associated with the problem and the constraints are described in Section 2.3. β penalty incorporated into the objective function when a leaf node contains a number of problem instances smaller than threshold τ ;

Input data
v p,f value of feature f for problem instance p; g l,n parameter that indicates which is the parent node of a given node n (considering a child node n where its parent is at its left); h l,n parameter that indicates which is the parent node of a given node n (considering a child node n where its parent is at its right); The maximum allowed tree depth is defined by d. To prevent overfitting, an additional cost (parameter β) is included into the objective function to penalize the occurrence of leaf nodes containing a number of problem instances smaller than threshold τ .
Parameters v p,f indicate the value of each feature f for each problem instance p. Parameters g l,n and h l,n indicate the parent node of a given node located at the left or right, respectively. Thus, if the parent of the node n is at left (n mod 2 = 0), then g l,n = (n + 1)/2 , otherwise g l,n = −1. Similarly, if the parent of the node n is at the right (n mod 2 = 1), then h l,n = (n + 1)/2 , otherwise, h l,n = −1.

Decision variables
The main decision variables x l,n,f,c are related to the feature and branch values at each branching node of the decision tree. From the choices defined by the model, the problem instances are grouped (variables y l,n,p ) according to the features and the cut-off points that were imposed on the branches. The model will determine the best algorithm for each group of problem instances placed on each leaf node (variables z n,a ). Variables u n are used to check if there are problem instances allocated to a given leaf node. This set will be linked to the set of variables m n -explained later in this section.
The next two sets of decision variables are used in the objective function. With the exception of set m n (m n ∈ Z + ), all other sets of variables are binary. To penalize leaf nodes with few instances, which could result in overfitting, variables m n are used to compute the number of problem instances that are missing for the leaf node n to reach a pre-established threshold of problem instances per leaf node, determined by the parameter τ . The set of decision variables w p,n,a is responsible for connecting the sets of decision variables y l,n,p and z n,a , i.e., to ensure that all problem instances allocated to a particular leaf node have the same recommended algorithm and that this algorithm is exactly the one corresponding to z n,a . In addition, the connection between the set w p,n,a and y l,n,p ensures that the problem instances allocated to leaf nodes respect branching decisions on parent nodes. m n = number of problem instances missing from the leaf node n to reach a pre-established threshold of problem instances per leaf node.

Objective function and constraints
The objective of our model is to construct a tree of determined maximum depth that minimizes the distance of the performance obtained using the recommended algorithm from the ideal performance for each problem p. Here we consider that this non-negative value is already computed in r. There is an additional cost involved in the objective function to penalize the occurrence of leaf nodes with only a few supporting instances. Follows the objective function (2) and the set of constraints (3-19) of our model: Equations 3 ensure that each internal node of the tree must have exactly one feature and branching value selected. Each problem instance must be allocated to exactly one leaf node and one algorithm (Equations 4) and each leaf node must have exactly one associated algorithm (Equations 5). Inequalities 6 guarantee that the recommended algorithm for a leaf node is the same as the algorithm of the problem instances allocated to this node. Inequalities 7 guarantee that variables y d,n,p and w p,n,a have values associated with the problem instances belonging to each leaf node.
Constraint set 8 ensures that variables u n are 1 if and only if there is at least one problem instance associated with leaf node n. Constraints 9 ensure that variable m n is set to the number of problem instances missing from the leaf node n to reach the threshold τ . If m n = 0, then the leaf node n contains at least τ problem instances.
Constraints 10 ensure that any problem instance allocated in a particular node must belong to the associated parent node. Finally, constraints 11 and 12 ensure that problem instances allocated in a particular node respect the feature and branching values selected at the parent node. Constraints 11 are generated when g l,n = −1 and v p,f ≤ c and ensure that problem instance p cannot be allocated at node n of level l (y l,n,p = 0), when feature f and branching value c are chosen for its parent node (x (l−1),g l,n ,f,c = 1). Similarly, Constraints 12 are generated when h l,n = −1 and v p,f > c and ensure that problem instance p cannot be allocated at node n of level l (y l,n,p = 0), when feature f and branching value c are chosen for its parent node (x (l−1),h l,n ,f,c = 1). Constraints 13-19 are related to the domain of the decision variables defined in the model.

VND to accelerate the discovery of better solutions
The model proposed in Section 2 can be optimized by standalone Mixed-Integer Programming (MIP) solvers and in finite time the optimal solution for the ASP will be produced. Despite the continuous evolution of MIP solvers (Johnson et al., 2000;Gamrath et al., 2015), the optimization of large MIP models in restricted times, in the general case, is still challenging. Thus, we performed some scalability tests (Section 4.2) to check how practical it is the use of the our complete model to create optimal decision trees for the ASP in datasets of different sizes in limited times. Since our objective is to produce a method that can tackle large datasets of experiment results, we also propose a mathematical programming heuristic (Fischetti and Fischetti, 2016) based on Variable Neighborhood Descent (VND) (Mladenović and Hansen, 1997) to speed the production of feasible solutions. VND is a local search method that consists of exploring the solution space through systematic change of neighborhood structures. Its success is based on the fact that different neighborhood structures do not usually have the same local minimum. Algorithm 1 shows the pseudo-code for our approach called VND-ASP. The overall method operates as follows: a Greedy Randomized algorithm is employed to generate initial feasible solutions (GRC-ASP). Multiple runs of this constructive algorithm are used to construct an Elite Set of solutions. The best solution of this set is used in our VND local search (MIPSearch), where parts of this solution are fixed and the remaining parts are optimized with the MIP model presented previously. Also, a subset Q including at most q algorithms is built. It includes all algorithms that appear in the elite set E and additional algorithms selected with a MIP model as follows: a covering like model to select q algorithms is solved where each instance should be covered by at least q < q algorithms, minimizing the cost of covering each problem instance with the selected algorithm.
Algorithm 1 VND-ASP (r, h, l, D, d, α, m, n, Q, A, P) Input. matrix r: algorithm performance matrix; set D: all different branching values for all features (F ×C f ); set A: set of algorithms; set P: set of problems. Parameters. h: matheuristic execution timeout; l: MIP search execution timeout; d: maximum depth; α: represents a continuous value between 0.1 and 1.0 that controls the greedy or random construction of the solution; m: maximum number of trees; n: maximum number of iterations without improvement; q: defines the number of algorithms of the set A, q : minimum number of algorithms to cover each problem instance. Our matheuristic has the following parameters: h indicates the time limit for running the entire algorithm; l indicates the time limit for a single exploration of a sub-problem in the MIP based neighborhood search; d is the maximum tree depth; parameter α ∈ [0, 1] controls the randomization of the greedy algorithm; parameter m controls the maximum size of the set of trees E; executions of the GRC-ASP algorithm are restricted to at most n iterations without updated in the elite set. Set D represents the different features and cut-off points of the problem instances.
The list J contains all sub-problems that can be obtained by fixing solution components in all neighborhoods. We shuffle these sub-problems in a list so that there is no priority for searching first in one neighborhood relative to another. This strategy is inspired by (Souza et al., 2010), where several neighborhood orders were in tested in a VND algorithm and the randomized order obtained better results. Whenever the incumbent solution is updated, the list J is shuffled again and starts to explore it from the beginning.
In the following sections we describe in more details the algorithm used to generate initial feasible solutions (Section 3.1) and the MIP based neighborhoods employed in our algorithm (Section 3.2).

Constructive Algorithm
The initial solution s is obtained from the best solution of the set of trees E. This set is obtained using a hybrid approach inspired by Quinlan (1993)'s C4.5 algorithm to generate a decision tree and the Greedy Randomized Constructive (GRC) (Resende and Ribeiro, 2014) search. GRC search to a certain randomness in the greedy criterion adopted by the C4.5 algorithm. Algorithm 2 shows the hybrid approach. Lines 26-34 of Algorithm 2 (GRC-ASP) show the adaptation made in the C4.5 algorithm for use of the restricted candidate list of the GRC search.
Another adaptation considers the metric to split the nodes. Algorithm C4.5 uses information gain metric. This metric aims to choose the attribute that minimizes the impurity of the data. In a data set, it is a measure of the lack of homogeneity of the input data in relation to its classification. In our case, we used the performance degradation metric to split the nodes. This metric searches for the attribute that minimizes the degradation of performance obtained using the recommended algorithm from the ideal performance for each problem p.

Neighborhoods
Since optimizing the complete MIP model may be too expensive, five neighborhoods have been designed to be employed in a fix-and-optimize context. Each neighborhood defines a set of subproblems to be optimized. These neighborhoods are explained in what follows, together with examples shown in Figs. 2-6. Examples consider a decision tree with 4 levels: variables in gray are fixed and variables highlighted in black will be optimized. These neighborhoods are explored in a Variable Neighborhood Descent using the MIP solver in a fixand-optimize strategy. Not only the neighborhoods, but the sub-problems of all neighboorhoods, are explored in a random order.
We will explain how the decision variables x l,n,f,c , y l,n,p and z n,a are optimized in each of the neighborhoods. The following neighborhoods were developed: • Neighborhood N 1 : optimizes the selection of the feature and the cutoff point of an internal node and consequently optimizes the allocation of problems in child nodes as well as optimizes the choice of the recommended algorithm in the respective leaf nodes.
• Neighborhood N 2 : optimizes the selection of the feature and the cutoff point of an internal node (this node cannot be the root node) and optimizes the choice of the feature and the cut-off point of the associated parent node, it consequently optimizes the allocation of problems in child nodes of associated parent node as well as the choice of the recommended algorithm in the respective leaf nodes.
• Neighborhood N 3 : optimizes the selection of the feature and the cutoff point of all nodes at one level (the level of the root node cannot be chosen) of the decision tree. Consequently optimizes the allocation of the problems in the nodes of the subsequent levels to the chosen level, it in addition to the choice of the recommended algorithm in the respective leaf nodes.
• Neighborhood N 4 : optimizes the selection of the feature and the cutoff point of the root node and the choice of the feature and the cut-off point of an internal node, so that this node is at least at the third level of the tree (l = 2). Consequently it optimizes the allocation of problems in all nodes of the tree as well as the choice of the recommended algorithm in the respective leaf nodes.
• Neighborhood N 5 : optimizes the selection of the feature and the cutoff point of a particular path from the root node to one of the tree's leaf nodes. Consequently it both optimizes the allocation of problems in all nodes of the tree and the choice of the recommended algorithm in the respective leaf nodes.

Experiments
The computational experiments are divided into two groups: the first group (Section 4.2) is dedicated to check the scalability of our integer programming model using a standalone MIP solver in datasets with different numbers of problem instances and algorithms. The second group of experiments (Section 4.4) uses cross-validation and divides the complete base into 10 partitions, where training and test partitions are used. Section 4.1 presents the initial configurations for these two experiments.

Computational experiments dataset
This section describes how the experiments database was created. Details of the sets of problem instances and algorithms will be presented in Sections 4.1.1 and 4.1.2.

Problem instances
Computational experiments were performed for a diverse set of 1004 problem instances including the MIPLIB 3, 2003, 2010 and 2017 (Koch et al., 2011) benchmark sets. Additional instances from Nurse Rostering (Santos et al., 2016), School Timetabling (Fonseca et al., 2017) and Graph Drawing (Silva and Santos, 2017) were also included. We extracted 37 features (|F| = 37) associated to variables, constraints and coefficients in the constraint matrix for problem instances. These features are similar to the ones used in (Hutter et al., 2014b) with the notable exception that features that are computationally expensive to extract were discarded to ensure that our approach would incur no overhead when incorporated into an algorithm. When building the problem instances dataset, special care was taken to ensure that no application was overrepresented. Table 1 shows the minimum (min), maximum (max), average (avg) and standard deviation (sd) of each feature over the complete set of problem instances. The density feature was computed as ( nz rows×cols ) × 100. Fig. 7 summarizes the 37 features for ASP, grouped by features related to variables, constraints and coefficients in the constraint matrix.

Number of variables: cols. 2 -4.
Number of variables of type: bin, int and cont.
Number of non-zeros and density: nz and density.

Available algorithms
The definition of the solution method for the LP solver in CBC involves selecting the algorithm, such as dual simplex or the barrier and defining several parameters, such as the perturbation value and the pre-solve effort. Overall 532 different algorithm configurations were evaluated for each one of the 1004 problem instances . A timeout T = 4000 was set for each execution. The computational results matrix r was filled with the execution time for regular executions, i.e. executions that finished before the time limit and provided correct results. Executions for a given problem instance p and algorithm a that crashed, exceeded the time limit or produced wrong results were penalized by setting r pa = 8000. This large batch of experiments was executed in computers with 32 Gb of RAM and 10 Intel R i9-7900X processing cores. Tasks where scheduled in parallel (7 threads simultaneously) with the GNU Parallel package (Tange, 2011).

Experiments to evaluate scalability of the integer programming model
To evaluate the performance and the scalability of the proposed formulation in a standalone MIP solver, models for generating trees with different depths (d = {1, . . . , 3}) with datasets of different sizes built by randomly selecting subsets of results of the complete experimental results of the COIN-OR CBC solver were solved with the state-of-the-art CPLEX 12.9 MIP solver on a computer with 32GB of RAM and 6 Intel R i7-4960X cores. In this experiment, we measured the final gap reported by the solver between the best lower and upper bound at the end of execution with one hour time limit. These experiments considered generating trees minimum number of 10 instances per leaf node (τ = 10) with a penalty of (β = 50) for leaf nodes violating this constraint. Fig. 8 shows the performance of our integer programming model, considering bases with 50 algorithms and problem instances ranging from 50 to 500. 20  As it can be seen, optimal decision trees were generated for models with up to 200 problem instances. For larger datasets the execution terminated with increasingly larger gaps for the produced bounds, at the point that for models with more than 300 instances a feasible solution was not found in the time limit. For the model with 500 problem instances not even the LP relaxation of the MIP model was computed in the time limit and no lower bound was available. Thus, for the complete dataset, experiments in the next subsections were performed only with the proposed VND-ASP heuristic.

Experiments with the complete dataset
To create optimized decision trees considering the entire experimental results dataset is quite challenging: there are more than half a million observations 1 , far beyond the limits indicated in the previous section. Thus, only experiments with our mathematical programming heuristics were conducted for this dataset. Fig. 9 presents the decision tree constructed with VNS-ASP using the following parameters: maximum tree depth d = 3, total time limit h = 72000, MIP search timeout l = 4000, elite set size m = 20, initial algorithms sub-setsize q = 100, q = 20, minimum number of instances per leaf node τ = 50 and penalty cost β = 500 . The estimated performance improvement with this decision tree is 61%, a remarkable improvement. Please note, however, that this improvement does not reflects the expected performance improvement of this tree in unknown instances, which is the really important estimate. The estimated results of the decision trees produced with our method on unknown instances is computed in the next section in 10-fold cross validation experiments.
An inspection in the contents of our decision tree shows that the range of the coefficients in the constraint matrix plays an important role for determining the best algorithm. The feature selected for the root node aRatioLSA is computed as the ratio between the largest and the smallest absolute non-zero values in the constraint matrix. Each leaf node has a set of instances allocated to it, depending on the the decision on all parent nodes and a recommended algorithm, which is the algorithm with better results on these instances. As an example, for the left-most branch of the tree, the best algorithm configuration select used the Primal simplex algorithm setting the "idiot" parameter to value 7 considering 127 LP problems allocated to this node.

Experiment using cross-validation on the complete base of problem instances
To evaluate the predictive power of our method, i.e. the expected performance on unknown instances, a 10-fold cross validation experiment was performed: a randomly shuffled complete dataset was divided into 10 partitions and at each iteration 9 of the partitions were used to create the decision tree (training dataset) and the remaining partition used for evaluating the decision tree (test dataset). Each partition had 480928 examples (904 problem instances × 532 available algorithms), with the exception of the last four partitions that contained 480396 examples (903 problem instances × 532 available algorithms).
The results of the cross-validation are given in Fig. 10. This figure shows the average performance degradation considering the ideal performance to solve the LP relaxation of all problem instances (the lower bound). Results of VND-ASP with maximum tree depth 4 (VND-ASP(D=5)) and 5 (VND-ASP(D=5)) are included. The remaining parameters of VND-ASP are the same described in the previous subsection. We also compare our results with the results produced with different configurations of the Random Forest (RF) algorithm implemented in Weka (Hall et al., 2009) (RF(T=1,. . .,T=200)). Default CBC settings (Default) and results obtained selecting a single best algorithm (SBA) are also included. As it can be seen, our results indicated performance gains of 40%, comparing with the default settings and are noticeably better than the results that can be obtained by selecting a single best algorithm were gains of 22% are estimated. VND-ASP results are also mostly superior to the ones obtained with RF, with the exception of RF(T=50), were equivalent results were obtained. We believe that this is a very positive result, since the result produced by our algorithm (a single tree) is more easily interpretable than the result produced by RF. More importantly, algorithms can be recommended much faster (in constant time), with our approach since the processing cost does not depends on the number of available algorithms as in RF, where a series of regression problems need to be solved to recommend an algorithm for each new instance.

Discussion and closing remarks
This paper presented a new mathematical programming formulation to solve the Algorithm Selection Problem. This formulation produces globally optimal decision trees with limited depth. The main advantage of this approach is that, even though the construction of the tree itself can be computationally expensive, once the tree is built algorithm recommendations can be done in constant time. A dataset with experimental results of many linear programming solver configurations of the COIN-OR Branch-&-Cut linear programming solver (CLP) was built solving a comprehensive set of instances from various applications. This initial batch of experiments itself already revealed improved parameter settings for the LP solver with the discovery of a new algorithm configuration that performed 22% faster than default CLP settings.
For the production of optimal decision trees, scalability tests were performed in subsets of the dataset to evaluate the feasibility of solving the proposed formulation in standalone MIP solvers considering different problem sizes. Since at some point the resulting MIP model becomes too hard to be optimized exactly, a mathematical programming based VND local search heuristic was also proposed to handle larger datasets.
To evaluate the predictive power of our method, a 10-fold cross validation experiment was conducted. The results were very promising: executions with the recommended parameter settings were 40% faster than CLP default settings, almost doubling the improvement that could be obtained using a single best parameter setting. Our results are comparable with the results obtained after tuning the Random Forest algorithm, with the advantage that the predictive model produced by our method, a single tree, is easily interpretable and, more importantly, the algorithm recommendation cost does not depends on the number of available algorithms.
Future directions include evaluating stronger alternative integer programming formulations for this problem: as the scalability test showed, there is still a significant gap between the lower and upper bounds produced for the larger datasets. The positive results for the ASP are also a good indication that the application of a our methodology to classification and regression problems.