The 1990’s saw two major advances in machine learning: the support vector machine (SVM) and the AdaBoost algorithm. Two fundamental ideas behind these algorithms are especially far-reaching. The first one is that we can transform many classical linear algorithms into highly flexible nonlinear algorithms by using kernel functions. The second one is that we can make accurate predictions by building an ensemble of models without much fine-tuning for each, rather than carefully fine-tuning a single model.
In this expository article, I first present the main ideas behind kernel methods (Section 2) and ensemble methods (Section 3) by reviewing four prototypical algorithms: the support vector machine (SVM, e.g., Cristianini and Shawe-Taylor 2000)
, kernel principal component analysis(kPCA, Schölkopf et al. 1998), AdaBoost (Freund and Schapire 1996), and random forest (Breiman 2001). I then illustrate the influence of these ideas on my own research (Section 4) by highlighting two recent algorithms that I have invented with my collaborators: LAGO (Zhu et al. 2006), a fast kernel machine for rare target detection; and Darwinian evolution in parallel universes (Zhu and Chipman 2006), an ensemble method for variable selection.
To better focus on the main ideas and not be distracted by the technicalities, I shall limit myself mostly to the two-class classification problem, although the SVM, AdaBoost and random forest can all deal with multi-class classification and regression problems as well. Technical details that do not affect the understanding of the main ideas are also omitted.
I begin with kernel methods. Even though the idea of kernels is fairly old, it is the support vector machine (SVM) that ignited a new wave of research in this area over the past 10 to 15 years.
In a two-class classification problem, we have predictor vectors and class labels ,
. SVM seeks an optimal hyperplane to separate the two classes.
A hyperplane in consists of all that satisfy the linear equation:
Given and , a hyperplane is called a separating hyperplane if there exists such that
Clearly, a hyperplane can be reparameterized by scaling, e.g.,
for any scalar . In particular, we can scale the hyperplane so that (1) becomes
that is, scaled so that . A separating hyperplane satisfying condition (2) is called a canonical separating hyperplane (CSHP).
If two classes are perfectly separable, then there exist an infinite number of separating hyperplanes. Figure 1 shows two competing hyperplanes in such a situation. The SVM is based on the notion that the “best” canonical separating hyperplane to separate two classes is the one that is the farthest away from the training points. This notion is formalized mathematically by the margin of a hyperplane — hyperplanes with larger margins are better. In particular, the margin of a hyperplane is equal to
where is the signed distance between observation and the hyperplane; see Figure 1 for an illustration. Figure 1 also shows to a certain extent why large margins are good on an intuitive level; there is also an elaborate set of theories to justify this (see, e.g., Vapnik 1995).
It can be shown (e.g., Hastie et al. 2001, Section 4.5) that is equal to
To find the “best” CSHP with the largest margin, we are interested in solving the following optimization problem:
The extra variables are introduced to relax the separability condition (2) because, in general, we can’t assume the two classes are always perfectly separable. The term acts as a penalty to control the degree of such relaxation, and is a tuning parameter.
The main message from the brief introduction above is this: SVM tries to find the best CSHP; it is therefore a linear classifier. The usual immediate response to this message is: So what? How does this make the SVM much different from and superior to classical logistic regression?
Equivalently, the constrained optimization problem above can be written as (e.g., Hastie et al. 2001, Exercise 12.1)
For statisticians, the objective function in (6
) has the familiar form of a loss function plus a penalty term. For the SVM, the loss function is, and it is indeed very similar to the binomial log-likelihood used by logistic regression (e.g., Hastie et al. 2001, Figure 12.4). But the usual logistic regression model does not include the penalty term . This is the familiar ridge penalty and often stabilizes the solution, especially in high-dimensional problems. Indeed, this gives the SVM an advantage.
However, one can’t possibly expect a linear classifier to succeed in general situations, no matter how optimal the hyperplane is. So, why is the SVM such a sensational success?
2.2 The “kernel trick”
Cristianini and Shawe-Taylor (2000, Chapters 5 and 6) provided detailed derivations to show that the optimal looks like this:
where “SV” denotes the set of “support vectors” with strictly positive; the coefficients , are solutions to the (dual) problem:
This means the resulting hyperplane can be written as
The key point here is the following: In order to obtain , one solves (7)-(8), a problem that depends on the predictors only through their inner-products ; once the ’s are obtained, the ultimate decision function (9) is also just a function of inner-products in the predictor space.
Therefore, one can make SVM a lot more general simply by defining a “different kind of inner-product,” say, , in place of . The function is a called a kernel function, where is a hyper-parameter, which is often determined empirically by cross-validation. Then, (7) becomes
and the decision function (9) becomes
The boundary is linear in the space of where is such that
but generally it is nonlinear in the original predictor space (unless one picks a linear kernel function). Mercer’s theorem (Mercer 1909) guarantees the existence of such as long as is a non-negative definite kernel function. The beauty here is that we don’t even need to define the mapping explicitly; all we have to do is to pick a kernel function . This makes the SVM very general.
2.3 Kernelization of linear algorithms
That we can apply a linear method in a different space is, of course, not a new idea to statisticians at all. For example, we all know how to fit a high-order polynomial using linear regression — simply add the termsto the regression equation!
The idea that we don’t need to explicitly create these high-order terms is perhaps somewhat less familiar. Actually, it is not really a new idea, either; it is less familiar only in the sense that students usually don’t learn about it in “Regression Analysis 101.”
However, the SVM does deserve some credit in this regard. Even though the basic idea of kernels is fairly old, it is the SVM that has revived it and brought it back into the spotlight for applied statisticians. The basic idea is as follows.
A typical data matrix we encounter in statistics, , is , stacking observations as -dimensional row vectors. That is,
It is easy to see that
is an matrix of pairwise inner-products. Therefore, if a linear algorithm can be shown to depend on the data matrix only through
then it can be easily “kernelized” — we simply replace each inner-product entry of with , where is a desired kernel function.
2.4 Kernel PCA
Kernel principal component analysis (kPCA; Schölkopf et al. 1998) is a successful example of “kernelizing” a well-known classic linear algorithm. To focus on the main idea, let us assume that the data matrix is already centered so that each column has mean zero. Let
Then, the (ordered) eigenvectors of, say , are the principal components. Being eigenvectors, they satisfy the equations
This shows that can be represented in the form of — by letting , to be specific. We will plug into (15
) and reparameterize the eigenvalue problem in terms of.
For , this leads to
If we left-multiply both sides by , we get
which shows that can be obtained by solving a problem that depends on the data matrix only through the inner-product matrix .
Schölkopf et al. (1998) explained why, in the context of kPCA, it is sufficient to reduce (17) to ; I do not go into this detail here. Once we obtain the ’s, suppose we’d like to project new data onto a few leading principal components, e.g., . We immediately find that
and it is easily seen that is just a matrix of pairwise inner products between each new and old observations.
Hence, it becomes clear that both finding and projecting onto principal components depend on just the inner-products and, according to Section 2.3, PCA can be “kernelized” easily. Figure 2 shows a toy example. There are some spherical data in
. The data being spherical, all directions have equal variance and there are no meaningful principal components in the traditional sense. But by using a Gaussian kernel — equation (18) below with — in place of all the inner-products, the first kernel principal direction obtained gives a meaningful order of how far each observation is away from the origin. In this case, kernel PCA has successfully discovered the (only) underlying pattern in the data, one that is impossible to detect with classical PCA.
2.5 Discussion: Kernel methods are like professional cameras
Any acute reader must have noticed that, so far, I have never really discussed the kernel function explicitly. This is not an accident. It is often claimed that one important advantage of these kernel methods lies in their modularity: to solve a different problem, just use a different kernel function. Any discussion about kernel functions, therefore, is best carried out in the context of a specific problem.
Of course, to be effective in practice, we must use the right kernel function. What’s more, we must choose the right hyper-parameter as well, and the performance of the method can be quite sensitive to these choices in practice. These are no trivial tasks and often require a considerable amount of data analytic experience as well as knowledge of the specific application area.
In this regard, these kernel-based algorithms are very much like professional cameras. They are capable of producing great pictures even under very difficult conditions, but you need to give them to a professional photographer. If you give them to an amateur or novice, you can’t expect great pictures. The photographer must know how to select the right lens, set the right shutter speed, and use the right aperture for any given condition. If any of these parameters is not set appropriately, the result could be a disaster. But that does not mean the camera itself is a poor piece of equipment; it simply means one must be adequately trained to operate it. Much of the power of these professional cameras lies precisely in the fact that they allow a knowledgeable and experienced user to control exactly how each single picture should be taken.
2.5.1 Example: Spam data
As a very simple illustration, let us try to see how well the SVM can
predict on the spam data set, available at
are a total of observations, each with a binary response and
predictors. For more details about this data set, refer to the
aforementioned web site. I use an
R package called
fit SVMs and use the kernel function
A random sample of observations are used as training data and the remaining observations are used as test data. Using different values of and , a series of SVMs are fitted on the training data and then applied to the test data. The total number of misclassification errors on the test data are recorded and plotted for each pair of ; see Figure 3(a). Here, is the penalty parameter in equation (4).
Figure 3(a) shows that the performance of SVM using this particular kernel function is very sensitive to the parameter but not as sensitive to the parameter . Given , the prediction performance of SVM is often quite stable for a wide range of ’s, but bad choices of can lead to serious deteriorations in the prediction performance. Therefore, if one uses the SVM without carefully tuning the parameter , the result can be disastrous.
I now turn to ensemble methods. Again, I shall mainly focus on the two-class classification problem with predictor vectors and class labels , .
AdaBoost constructs a collection of classifiers rather than one single classifier. The entire collection makes up an ensemble, and it is the ensemble — not any single classifier alone — that makes the final classification.
Table 1 contains an exact description of the AdaBoost algorithm. Here is a description of the algorithm in plain English: Start by assigning equal weights to all observations in the training data. Sequentially build a series of classifiers. At each step, fit a classifier, say , to the training data using the current weights. Calculate the (properly weighted) right-to-wrong ratio of this classifier; call it . For those observations incorrectly classified by , inflate their weights by a factor of . With the new weights, build the next classifier. In the end, each classifier in the ensemble will cast a vote; its vote is to be weighted by the logarithm of its right-to-wrong ratio, .
For people hearing about this algorithm for the very first time, AdaBoost certainly has a very strong mystical flavor to it. Intuitively, we can perhaps appreciate to some extent that the right-to-wrong ratio must be important for any classifier, but it is not at all clear why we should reweight incorrectly classified observations by this ratio each time, nor is it immediately clear why the final vote from each individual member of the ensemble should be weighted by the logarithm of this ratio.
This is no easy mystery to untangle. Friedman et al. (2000) gave a very nice argument and revealed that the AdaBoost algorithm actually minimizes an exponential loss function using a forward stagewise approach. In particular, AdaBoost chooses the best and one step at a time to minimize
which they showed to be very similar to maximizing the binomial log-likelihood. This particular interpretation has not only untangled the AdaBoost mystery (at least to some extent), but also led to many new (and sometimes better) versions of boosting algorithms.
3.2 Random forest
Professor Leo Breiman came up with the same basic idea of using a collection or an ensemble of models to make predictions, except he constructed his ensemble in a slightly different manner. Breiman called his ensembles random forests; details are given in Table 2.
The history behind Breiman’s random forest is very interesting. In 1996, he first proposed an ensemble algorithm called Bagging (Breiman 1996), which is essentially the random forest algorithm with just the bootstrap step (Table 2, step 1a). In 2001, he added the random subset step (Table 2, step 1b) and created random forest (Breiman 2001).
Why did he add the extra random subset step?
3.3 Breiman’s theorem
Breiman (2001) proved a remarkable theoretical result. First, he gave a formal definition of random forests: The set
is called a random forest.
This definition requires some explanation. Here, is a classifier completely parameterized by . For example, if is a classification tree, then the parameter
specifies all the splits and the estimates in the terminal nodes. Next, the statement “” means that each is generated independently and identically from some underlying random mechanism, .
To be specific, in Breiman’s implementation, iid sampling from the random mechanism consists of: (i) iid sampling from the empirical distribution (the bootstrap step), and (ii) iid sampling from the set (the random subset step).
Breiman then proved that the prediction error of a random forest, , satisfies the inequality
where is the mean correlation between any two members of the forest (ensemble) and , the mean strength of a typical member of the forest (ensemble). This result — including the exact definitions of and — is fairly technical; details can be found in Breiman (2001). Moreover, the actual bound itself is often useless. For example, if and , then one gets
but of course the error rate is less than 100%.
So, why is this result significant?
3.4 The secret of ensembles
The fundamental idea of using an ensemble classifier rather than a single classifier is nothing short of being revolutionary. It also is remarkable that building these ensembles is often relatively mindless. Take Breiman’s random forest, for example. There is no need to prune the individual trees.
Clearly, there are many different ways to build an ensemble, AdaBoost and Breiman’s random forest being two primary examples. What’s the most effective way?
Recall the formal definition of random forests. The random mechanism that generates the individual members of the forest is unspecified. You are free to pick any mechanism you want. Surely some mechanisms are bound to be more effective than others. What’s the most effective mechanism?
Breiman’s result is significant because it tells us what makes a good random forest. Breiman’s theorem (19) tells us that a good random forest should have a small and a large . That is, we should try to reduce the correlation between individual classifiers within the ensemble and make each individual classifier as accurate as possible.
This explains why Breiman added the random subset step into his original Bagging algorithm: extra randomness is needed to reduce the correlation between individual trees; the bootstrap step alone is not enough!
From this, we can write
Therefore, we can see that
which means the misclassification error of under the new weights is exactly 50% — the worst possible error.
The next classifier, , is built using these new weights, so it is set up to work with a (weighted) dataset that the current classifier, , cannot classify. This is sometimes referred to as “decoupling” in the boosting literature — the classifier is decoupled from .
In Breiman’s language, we can say that the adaptive and hitherto mysterious reweighting mechanism in AdaBoost is actually aiming to reduce the correlation between consecutive members of the ensemble.
3.5 Discussion: Ensemble methods are like foolproof cameras
Compared with kernel methods, ensemble methods are very much like foolproof cameras. They are relatively easy for the less experienced users to operate. This does not mean they don’t have any tuning parameters; they do. Even when using a foolproof camera, one must still make a few decisions, e.g., whether or not to turn on the flash, and so on. But relatively speaking, the number of decisions one has to make is much more limited and these decisions are also relatively easy to make.
For example, in Breiman’s random forest, the size of the subset, (Table 2, step 1b), is an important tuning parameter. If is too large, it will cause to be too large. In the extreme case of , all the trees in the forest will be searching over the entire set of variables in order to make splits, and they will be identical — since the tree-growing algorithm is deterministic conditional on the data. On the other hand, if is too small, it will cause to be too small. In the extreme case of , all the trees will essentially be making random splits, and they will not be very good classifiers. There is plenty of empirical evidence to suggest, however, that the parameter is still relatively easy to choose in practice. Moreover, the parameter is not as sensitive as the complexity parameter of a kernel function (also see Section 3.5.1 below). Translation: Even if you are a bit off, the consequences will not be quite so disastrous.
I have had many occasions working with graduate students trying to make predictions using the SVM and Breiman’s random forest. They almost always produce much better predictions with the random forest, even on problems that are well-suited for the SVM! Sometimes, their SVMs actually perform worse than linear logistic regression. Certainly, there are many cases in practice where one would not expect the SVM to be much superior to linear logistic regression, e.g., when the true decision boundary is in fact linear. But if used correctly, the SVM should at least be comparable with linear logistic regression; there is no reason why it ever would be much worse. These experiences remind me over and over again just how difficult it can be for a novice to use the SVM.
But, as I stated in Section 2.5, you can’t blame the professional camera if you don’t know how to use it properly. There is always a tradeoff. With limited flexibility, even a fully-experienced professional photographer won’t be able to produce images of the highest professional quality with just a foolproof camera, especially under nonstandard and difficult conditions. That’s why professional cameras are still on the market. But we have to admit: most consumers are amateur photographers and, more often than not, they are taking pictures under fairly standard conditions. That’s why the demand for foolproof cameras far exceeds that for professional cameras. I think the demand for statistical tools follows a similar pattern.
3.5.1 Example: Spam data (continued)
As a simple illustration, let us take a look at how well the random forest
can predict on the spam data set. I use exactly the same set-up as in
Section 2.5.1 and the
randomForest package in
Using different values of and , a series of random forests are
fitted on the training data and then applied to the test data. The total
number of misclassification errors on the test data are recorded and
plotted; see Figure 3(b). Here, we can see that the
performance of random forests is more sensitive to the parameter than
to the parameter . Given , the prediction performance of random
forests is fairly stable as long as is sufficiently large, e.g.,
in this case. But it is important to use an that is neither
too small nor too big, e.g., in this case.
However, if we compare panels (a) and (b) in Figure 3, we can see that choosing the right for SVM is much more critical than choosing the right for random forest; performance deterioration is much more serious for bad choices of than for bad choices of .
It is also clear from Figure 3 that, for this particular data set, an SVM with kernel function (18) is not competitive against a random forest, even if well tuned. In order to be competitive, it is necessary to use a different kernel function. I do not pursue this possibility here because getting the SVM to work for this data set is far from the main point of our discussion, but this example does demonstrate that choosing the right kernel function
and picking the right hyperparameterare very important, and that an ensemble method such as the random forest can be somewhat easier to use in this regard.
I now share a few personal perspectives on statistical learning research. Here, I am working with a particular definition of the word “perspective” from the American Heritage Dictionary: a subjective evaluation of relative significance [emphasis added].
4.1 Statistical learning research
My discussions in Sections 2.5 and 3.5 have led me to ask the following question: If I were the president of a big camera manufacturing company, how would I run such a business? Other than standard business divisions such as accounting and human resources, I see three main lines of operation:
(Consulting and Consumer Outreach) Advise and teach photographers how to use various products and how to use the right equipment to produce great pictures under various difficult conditions. This is my consulting and consumer outreach division.
(High-end R&D) Understand the need of professional photographers and manufacture new, specialized equipment still lacking on the market. This is my R&D division for my high-end consumers.
(Mass R&D) Build the next-generation foolproof camera. This is my R&D division for my mass consumers.
I see a great deal of parallelism in statistical learning research. For statistical learning research, the consulting and consumer outreach division applies different learning methods to solve various difficult real-world problems; the high-end R&D division develops new, specialized algorithms for analyzing new types of data or data with special characteristics; and the mass R&D division develops better off-the-shelf learning algorithms.
With this particular point of view in mind, I end this article by briefly describing two personal learning products: a new kernel method from my high-end R&D division, and a new ensemble method from my mass R&D division.
4.2 A high-end R&D product: LAGO
Consider a two-class problem in which the class of interest () is very rare; most observations belong to a majority, background class (). Given a set of unlabelled observations, the goal is to rank those belonging to ahead of the rest.
Of course, one can use any classifier to do this as long as the classifier is capable of producing not only a class label but also an estimated posterior probabilityor a classification score. For example, the SVM does not estimate posterior probabilities, but the final decision function (11) is a classification score which can be used (at least operationally) to rank unlabelled observations — whether this is effective or not is a separate issue.
The final decision function produced by SVM (11) is of the form
where is a kernel function. For example, we can take to be diagonal and let be the Gaussian kernel
where is the determinant of .
The function (20
) is sometimes called a (single-layer) radial basis function network (RBFnet). Generally speaking, to construct an RBFnet one must compute and specify three ingredients:
the location parameter of each kernel function — together, they make up the set ;
the shape parameter of each kernel function; and
the coefficient in front of each kernel function.
Typically, one first specifies and and then estimates the ’s by least-squares or maximum likelihood. Often, one sets and treats the parameter as a global tuning parameter — this is what SVM does. Determining the ’s or the best set
from training data, however, is an NP-hard combinatorial optimization problem in general.
The SVM can be viewed as an algorithm for determining the set and the ’s simultaneously (Schölkopf et al. 1997); the set is simply the set of all support vectors. In order to do so, SVM solves a quadratic programming instead of a combinatorial optimization problem.
The product from my R&D division is an algorithm called LAGO (Zhu et al. 2006). The decision function constructed by LAGO for ranking unlabelled observations is as follows:
The parameter is a global tuning parameter. In the simplest case, we take
where is the average distance between the kernel center, , and its -nearest neighbors from , i.e.,
The notation “” denotes the -nearest neighbors of from ; and is a distance function, e.g., .
By comparing (22)-(23) with (20), we can easily see that LAGO can also be viewed as an algorithm for constructing an RBFnet, just like the SVM. In particular, the three ingredients of the RBFnet are specified as follows:
Every is a training observation from the rare class, .
Each kernel function is spherical with radius proportional to the average distance between its center and its the -nearest neighbors from .
Simply set and .
Here we see that the only computation needed is the calculation of ; all other ingredients are completely determined a priori. The calculation of , of course, is considerably simpler than quadratic programming, making LAGO many times faster and simpler than the SVM. Instead of solving an optimization problem to find support vectors, LAGO fully exploits the special nature of these rare-class detection problems and simply uses all training observations from the rare class as its “support vectors,” a significant shortcut. Our empirical experiences show that the shortcut is highly worthwhile. We find that LAGO almost always performs as well as and sometimes even better than the SVM for these rare-class classification and detection problems.
) can be viewed as a kernel density estimate ofadjusted locally by a factor that is approximately inversely proportional to , i.e., . The resulting ranking function is thus approximately a monotonic transformation of the posterior probability that item belongs to the rare class.
The only nontrivial calculation performed by the algorithm, equation (24), is somewhat special and nonstandard. The original idea came from a Chinese board game called GO. Consider the two black stones labelled A and B in Figure 4. A GO player will tell you that B controls more territories on the board than A. Why? Because, when compared with B, A is closer to more enemy (white) stones. Therefore, imagine two classes fighting for control over a common space. Given an observation from , if we want to use a kernel function to describe its effective control over the entire space, we should use a large kernel radius if its nearby neighbors from are a long distance away and a small kernel radius if its nearby neighbors from are a short distance away. Equation (24) captures this basic principle.
4.2.3 eLAGO versus sLAGO
That is, the kernel function was chosen to be elliptical rather than spherical. To distinguish the two, we call (25) eLAGO and (23), sLAGO. For many real-world rare-class problems, the dataset often contains a limited amount of information because is very rare. As such, the extra flexibility afforded by eLAGO is seldom needed in practice.
4.2.4 Discussion: LAGO is a specialized kernel method
LAGO is a kernel method, much like the SVM. There are two tuning parameters, and . Experiments similar to those described in Section 2.5.1 and Figure 3 have shown that the performance of LAGO is not very sensitive to and much more sensitive to . In practice, it often suffices to fix .
LAGO is not a general-purpose method; it is a specialized algorithm for a special learning problem, namely rare-class classification and detection. Its main advantages are its speed and simplicity. Discussions in Section 2.5 have made it clear that these kernel methods must be carefully tuned, e.g., using empirical procedures such as cross-validation. This means that, in practice, one almost always has to run these algorithms repeatedly for many times. One may be tempted to think that, if one algorithm takes 10 minutes to run and another takes 1 minute, the difference is still “negligible” for all practical purposes, but such ten-fold differences are often significantly magnified if one has to run these two algorithms repeatedly for many times.
Apart from these practical matters such as time savings, the more important lesson from this research lies in the basic principles behind the construction of LAGO (22). Here, we see that it always pays to exploit the special nature of an underlying problem. For these rare-class problems, there is only limited amount of useful information in the training data. LAGO fully exploits this fact by immediately zooming into the useful information (i.e., ) and making a few quick local adjustments based on — equation (24).
4.3 A mass R&D product: Darwinian evolution in parallel universes
Let us now consider a different problem, the variable selection problem. Given potential predictors, which combination is the best for predicting ? Let be the space of all possible subsets of . The typical approach is as follows: First, define a proper evaluation criterion,
Preferably should be a fair measure of . Common examples of include the Akaike information criterion (AIC, Akaike 1973), the Bayesian information criterion (BIC, Schwarz 1978), and generalized cross-validation (GCV), to name a few. Then, use a search algorithm to find the best which optimizes .
4.3.1 Two challenges: computation and criterion
There are two main challenges. The first one is computation. With potential predictors, the size of is . This gets large very quickly. For example, take and suppose we can evaluate a billion () subsets per second. How long will it take us to evaluate all of them? The answer is about 40,000 billion years:
This may seem serious, but it actually is not the problem we shall be concerned about here. Everyone must face this problem; there is no way out — just yet. For moderately large
, exhaustive search is impossible; stepwise or heuristic search algorithms must be used.
The second challenge is more substantial, especially for statisticians, and that’s the question of what makes a good evaluation criterion, . It is well-known that both the AIC and the BIC are problematic in practice. Roughly speaking, with finite data, the AIC tends to favor subsets that are too large, while the BIC tends to favor ones that are too small. For classic linear models, both the AIC and the BIC have the form:
where is the size of , or the number of variables included. The AIC uses whereas the BIC uses , being the sample size. Therefore, it appears that is too small and is too big. But if this is the case, surely there must be a magic somewhere in between? So why not find out what it is? While this logic is certainly quite natural, it by no means implies that the task is easy.
4.3.2 Darwinian evolution in parallel universes
The product from my R&D division is a very simple yet surprisingly effective method for variable selection by using Darwinian evolution in parallel universes (Zhu and Chipman 2006).
Here is how the algorithm works in a nutshell. Create a number of parallel universes. In each universe, run an evolutionary algorithm using the (apparently incorrect) AIC as the objective function forjust a few generations — the evolutionary algorithm is a heuristic stochastic search algorithm that mimics Darwin’s “natural selection” to optimize any given objective function (see, e.g., Goldberg 1989). Whatever it is, there will be a current best solution in each universe when we stop. For example, the current best subset in universe 1 may be ; in universe 2, it may be ; in universe 3, perhaps ; and so on. These form an ensemble. Now take a majority vote and select those variables that show up in significantly more universes than the rest. In the example here, this would be — and that’s the answer.
4.3.3 Explanation with a toy example
Why does this simple strategy work? A small toy example is enough to illustrate the gist of the idea. Generate
In other words, there are 10 potential predictors but the true model contains only 3 of them: , and . With just 10 variables, there are altogether subsets, and we can still afford to exhaustively compute the AIC for each one of them. Figure 5 plots the AIC versus the size for all possible subsets. A number of characteristic observations can be made:
The subset that has the smallest AIC is wrong; it includes a few variables too many.
On the AIC scale, many subsets are very close to each other and it is hard to tell them apart.
Let’s separate the subsets into two groups. Group I consists of those that include all the true variables — they are labelled with circles () in the plot. Group II consists of those that miss out on at least one of the true variables — they are labelled with crosses (), pluses (), and triangles (). Then, on the AIC scale, a significant gap exists between these two groups.
Having made these observations, we are now ready to explain why parallel evolution works. The large gap between group I and group II (observation 3) means that members from group I are significantly superior and hence easily favored by evolution. Therefore, after evolving for just a few generations, the current best subset in each universe is likely a member from group I. They are the ones that make up our ensemble. What do they have in common? They all include the 3 true variables.
But in order for majority vote to be effective in selecting the right variables, it is necessary that the true variables are the only thing that these ensemble members have in common. That’s why we can’t run the evolution for too long in each universe. With a short evolution, since members of group I are hard to distinguish from each other on the AIC scale (observation 2), the random nature of evolution will cause each universe to settle on different members from this group. If, on the other hand, we run the evolution for too long, the current best subsets from different universes will start to develop something else in common — they will all start to converge to the minimum AIC solution, which includes spurious variables (observation 1).
Figure 6 illustrates how parallel evolution works on this toy example. After running the evolutionary algorithm for just generations in each universe, we measure the importance of a variable by how often it shows up across the parallel universes. The correct solution for this example is . When a single universe is used (), we get the wrong solution — a spurious variable, namely variable , also shows up. But as more and more parallel universes are used, only the truly important variables, i.e., variables 2, 5 and 8 in this case, can “survive” the majority vote. We can see from Figure 6 that when as few as universes are used, the correct solution is already easily discernible: out of the 10 universes, variables 2, 5, and 8 each showed up at least 9 times; variable 6 showed up 4 times; and all other variables showed up at most twice.
4.3.4 Discussion: Parallel evolution is an easy-to-use ensemble method
Parallel evolution for variable selection is a successful example of using ensembles in a very different context. By using an ensemble, we can significantly “boost up” the performance of an apparently wrong variable selection criterion such as the AIC. The procedure is very easy to use. Most importantly, it is trivial to adapt this principle to general variable selection problems regardless of whether the underlying model is a classic linear model, a generalized linear model, a generalized additive model, a Cox proportional hazard model, or any other model for which the question of variable selection is meaningful. As such, it is not unfair to call parallel evolution a first-generation, foolproof, off-the-shelf variable selector.
A number of smart statisticians have questioned whether it is necessary to use the evolutionary algorithm. For example, one can apply Breiman’s Bagging principle and create an ensemble as follows: Draw a bootstrap sample of the data. Using the bootstrap sample, run a stepwise algorithm to optimize the AIC and choose a subset. Do this many times, and we get an ensemble of subsets. Take majority vote. Clearly, this would also work. I have experimented with this idea and found that it is not as effective; the probability of selecting the right subset of variables decreases significantly in simulation experiments. Why? Breiman’s theorem (Section 3.3) points us to an answer. Because bootstrapping alone does not create enough diversity within the ensemble. These subsets share too many things in common with the minimal AIC solution.
4.4 Section Summary
In this section, I have discussed a new kernel-based algorithm for rare target detection, LAGO, and a new ensemble method for variable selection based on parallel evolution. In doing so, a more general formulation of LAGO is presented (Section 4.2) using much better mathematical notation, e.g., equation (22). A simpler version, sLAGO, is given for the first time. Better explanations are also given for why parallel evolution (Section 4.3) works, e.g., Figure 5. Many people have the incorrect understanding that parallel evolution is merely a better search algorithm for variable selection. This is simply not true. In Section 4.3, it is emphasized that, instead of a better search algorithm, parallel evolution is actually an ensemble method that boosts up the performance of an apparently incorrect search criterion such as the AIC.
So, what have we learned? First of all, we learned that, by using kernel functions, we can use many linear algorithms such as separating hyperplanes and principal component analysis to find nonlinear patterns (Section 2). This easily can be done as long as the underlying linear algorithm can be shown to depend on the data only through pairwise inner-products, i.e., . Then, we simply can replace the inner-product with a kernel function . However, even though such a framework is straightforward, we also learned that it is important in practice to use the right kernel function and to carefully select the hyperparameter (Section 2.5). We saw that this is not necessarily an easy task (Section 2.5.1).
We then learned about ensemble methods (Section 3). The fundamental idea there is to use a collection of perhaps not-so-well-tuned models rather than a single model that often requires careful fine-tuning. This usually makes ensemble methods easier to use for non-experts. I then emphasized that, even for ensembles, it is necessary to perform some fine-tuning (Section 3.5) — this typically involves creating the right amount of diversity in the ensemble (Section 3.4 and 3.5). However, we saw that fine-tuning an ensemble algorithm is often easier than fine-tuning a kernel-based algorithm (Section 3.5.1).
I then argued (Section 4.1) that kernel methods and ensemble methods need to co-exist in practice. In particular, non-experts may tend to prefer ensemble methods because they are easier to use, whereas experts may tend to prefer kernel methods because they provide more flexibility for solving nonstandard and difficult problems (Sections 2.5 and 3.5). Hence, it is important for researchers in statistical machine learning to advance both types of methodology. I then presented some of my own research on both fronts: LAGO, a fast kernel machine for rare target detection; and Darwinian evolution in parallel universes, an ensemble method for variable selection.
This expository article is based on two workshops (September 2006, Sydney, Australia; October 2006, Adelaide, Australia) and a graduate course (January-April 2007, University of Waterloo, Canada) that bear the same title. The workshops were made possible by an AusCan scholarship from the Statistical Societies of Australia and Canada. I am indebted to Dr. Peter Ricci of the Australian Taxation Office for initiating these workshops. I have received many valuable comments and feedbacks from the participants of the workshops as well as from the students and faculty members who painstakingly sat through my graduate course. I thank Dr. Stanley S. Young for providing me with temporary office space at the National Institute of Statistical Sciences (NISS), where this article was first written. I also thank Dr. Dianna Xu and Dr. Dylan Small for providing me with temporary office spaces at Bryn Mawr College and the University of Pennsylvania, respectively, so that I could work on the revision of this article with little distraction. My research is partially supported by the Natural Science and Engineering Research Council (NSERC) of Canada, Canada’s National Program on Complex Data Structures (NPCDS) and the Mathematics of Information Technology And Complex Systems (MITACS) network. I’d like to acknowledge the work by Alexandra Laflamme-Sanders and Dandi Qiao to make LAGO and parallel evolution into R libraries; these libraries will be publically available from http://cran.r-project.org/. Finally, I’d like to thank the anonymous referees for helping me improve this manuscript, and especially referee number one for giving me the strongest support and encouragement possible.
- Akaike (1973) Akaike, H. (1973). Information theory and an extension of the maximum likelihood principle. Second International Symposium on Information Theory, pages 267–281.
- Breiman (1996) Breiman, L. (1996). Bagging predictors. Machine Learning, 24(2), 123–140.
- Breiman (2001) Breiman, L. (2001). Random forests. Machine Learning, 45(1), 5–32.
- Cristianini and Shawe-Taylor (2000) Cristianini, N. and Shawe-Taylor, J. (2000). An Introduction to Support Vector Machines and Other Kernel-based Learning Methods. Cambridge University Press.
- Freund and Schapire (1996) Freund, Y. and Schapire, R. (1996). Experiments with a new boosting algorithm. In Machine Learning: Proceedings of the Thirteenth International Conference, pages 148–156, San Francisco. Morgan Kauffman.
- Friedman et al. (2000) Friedman, J. H., Hastie, T. J., and Tibshirani, R. J. (2000). Additive logistic regression: A statistical view of boosting (with discussion). The Annals of Statistics, 28(2), 337–407.
- Goldberg (1989) Goldberg, D. E. (1989). Genetic Algorithms in Search, Optimization and Machine Learning. Addison-Wesley.
- Hastie et al. (2001) Hastie, T. J., Tibshirani, R. J., and Friedman, J. H. (2001). The Elements of Statistical Learning: Data-mining, Inference and Prediction. Springer-Verlag.
- Mercer (1909) Mercer, J. (1909). Functions of positive and negative type and their connection with the theory of integral equations. Philosophical Transactions of the Royal Society of London, Series A, 209, 415–446.
- Schölkopf et al. (1997) Schölkopf, B., Sung, K. K., Burges, C. J. C., Girosi, F., Niyogi, P., Poggio, T., and Vapnik, V. (1997). Comparing support vector machines with Gaussian kernels to radial basis function classifiers. IEEE Transactions on Signal Processing, 45(11), 2758–2765.
- Schölkopf et al. (1998) Schölkopf, B., Smola, A., and Müller, K. B. (1998). Nonlinear component analysis as a kernel eigenvalue problem. Neural Computation, 10, 1299–1319.
- Schwarz (1978) Schwarz, G. (1978). Estimating the dimension of a model. The Annals of Statistics, 6, 461–464.
Vapnik, V. N. (1995).
The Nature of Statistical Learning Theory. Springer-Verlag.
- Zhu and Chipman (2006) Zhu, M. and Chipman, H. A. (2006). Darwinian evolution in parallel universes: A parallel genetic algorithm for variable selection. Technometrics, 48, 491–502.
- Zhu et al. (2006) Zhu, M., Su, W., and Chipman, H. A. (2006). LAGO: A computationally efficient approach for statistical detection. Technometrics, 48, 193–205.