Introduction to Machine Learning: Class Notes 67577

2009·arXiv

Contents 1 Bayesian Decision Theory page 1 1.1 Independence Constraints 5 1.1.1 Example: Coin Toss 7 1.1.2 Example: Gaussian Density Estimation 7 1.2 Incremental Bayes Classiﬁer 9 1.3 Bayes Classiﬁer for 2-class Normal Distributions 10

Bayesian Decision Theory

During the next few lectures we will be looking at the inference from training data problem as a random process modeled by the joint probability distribution over input (measurements) and output (say class labels) variables. In general, estimating the underlying distribution is a daunting and unwieldy task, but there are a number of constraints or ”tricks of the trade” so to speak that under certain conditions make this task manageable and fairly effective.

To make things simple, we will assume a discrete world, i.e., that the values of our random variables take on a finite number of values. Consider for example two random variables X taking on k possible values and H taking on two values . The values of X could stand for a Body Mass Index (BMI) measurement of a person and H stands for the two possibilities standing for the ”person being over-weight” and as the possibility ”person of normal weight”. Given a BMI measurement we would like to estimate the probability of the person being over-weight.

The joint probability P(X, H) is a two dimensional array (2-way array) with 2k entries (cells). Each training example () falls into one of those cells, therefore ) holds the ratio between the number of hits into cell (i, j) and the total number of training examples (assuming the training data arrive i.i.d.). As a result

The projections of the array onto its vertical and horizontal axes by summing over columns or over rows is called marginalization and produces ) the sum over the j’th row is the probability i.e., the probability of a person being over-weight (or not) before we see any measurement — these are called priors. Likewise, is the probability ) which is the probability of receiving such a BMI measurement to begin with — this is often called evidence. Note

Fig. 1.1. Joint probability P(X, H) where X ranges over 5 discrete values and H over two values. Each entry contains the number of hits for the cell (joint probability ) is the number of hits divided by the total number of hits (22). See text for more details.

that, by definition, ) = 1. In Fig. 1.1 we have that 22 that is there is a higher prior probability of a person being over-weight than being of normal weight. Also is the highest meaning that we encounter BMI = with the highest probability.

tween the number of hits in cell (i, j) and the number of hits in the i’th column, i.e., the probability that the outcome is given the measurement . In Fig. 1.1 we have 7. Note that

Likewise, the conditional probability number of hits in cell (i, j) normalized by the number of hits in the j’th row and represents the probability of receiving BMI = given the class label (over-weight or not) of the person. In Fig. 1.1 we have 3/8 which is the probability of receiving BMI = given that the person is known to be of normal weight. Note that

The Bayes formula arises from:

from which we get:

The left hand side ) is called the posterior probability and is called the class conditional likelihood. The Bayes formula provides a way to estimate the posterior probability from the prior, evidence and class likelihood. It is useful in cases where it is natural to compute (or collect data of) the class likelihood, yet it is not quite simple to compute directly the posterior. For example, given a measurement ”12” we would like to estimate the probability that the measurement came from tossing a pair of dice or from spinning a roulette table. If x = 12 is our measurement, and stands for ”pair of dice” and for ”roulette” then it is natural to compute the class conditional: Computing the posterior directly is much more difficult. As another example, consider medical diagnosis. Once it is known that a patient suffers from some disease , it is natural to evaluate the probabilities ) of the emerging symptoms . As a result, in many inference problems it is natural to use the class conditionals as the basic building blocks and use the Bayes formula to invert those to obtain the posteriors.

The Bayes rule can often lead to unintuitive results — the one in particular is known as ”base rate fallacy” which shows how an nonuniform prior can influence the mapping from likelihoods to posteriors. On an intuitive basis, people tend to ignore priors and equate likelihoods to posteriors. The following example is typical: consider the ”Cancer test kit” problemwhich has the following features: given that the subject has Cancer ”C”, the probability of the test kit producing a positive decision ”+” is P(+ | C) = 0.98 (which means that 02) and the probability of the kit producing a negative decision ”-” given that the subject is healthy ”H” is (which means also that P(+ | H) = 0.03). The prior probability of Cancer in the population is P(C) = 0.01. These numbers appear at first glance as quite reasonable, i.e, there is a probability of 98% that the test kit will produce the correct indication given that the subject has Cancer. What we are actually interested in is the probability that the subject has Cancer given that the test kit generated a positive decision, i.e., P(C | +). Using Bayes rule:

which means that there is a 26.6% chance that the subject has Cancer given that the test kit produced a positive response — by all means a very poor performance.

If we draw the posteriors ) using the probability distribution array in Fig. 1.1 we will see that values of X smaller than a value which is in between . Therefore the decision which will minimize the probability of misclassification would

be to choose the class with the maximal posterior:

which is known as the Maximal A Posteriori (MAP) decision principle. Since P(x) is simply a normalization factor, the MAP principle is equivalent to:

In the case where information about the prior P(h) is not known or it is known that the prior is uniform, the we obtain the Maximum Likelihood (ML) principle:

The MAP principle is a particular case of a more general principle, known as ”proper Bayes”, where a loss is incorporated into the decision process. Let ) be the loss incurred by deciding on class when in fact the correct class. For example, the ”0/1” loss function is:

The least-squares loss function is: typically used when the outcomes are vectors in some high dimensional space rather than class labels. We define the expected risk:

The proper Bayes decision policy is to minimize the expected risk:

The MAP policy arises in the case ) is the 0/1 loss function:

Thus,

At this point we may pause and ask what have we obtained? well, not much. Clearly, the inference problem is captured by the joint probability distribution and we do not need all these formulas to see this. How do we obtain the necessary data to fill in the probability distribution array to begin with? Clearly without additional simplifying constraints the task is not practical as the size of these kind of arrays are exponential in the number of variables. There are three families of simplifying constraints used in the literature:

• statistical independence constraints, • parametric form of the class likelihood ) where the inference

becomes a density estimation problem, • structural assumptions — latent (hidden) variables, graphical models.

Today we will focus on the first of these simplifying constraints — statistical independence properties.

Consider two random variables X and Y . The variables are statistically independent ) meaning that information about the value of Y does not add anything about X. The independence condition is equivalent to the constraint: P(X, Y ) = P(X)P(Y ). This can be easily proven: if other hand, if P(X, Y ) = P(X)P(Y ) then

Let the values of X range over and the values of Y range over . The associated 2-way array, ) is represented by the outer product ) of two vectors P(X) = ()). In other words, the 2-way array viewed as a matrix is of rank 1 and is determined by k + l (minus 2 because the sum of each vector is 1) parameters rather than kl (minus 1) parameters.

Likewise, if statistically independent random variables where ranges over discrete and distinct values, then the n-way array ) is an outer-product of n vectors and is therefore determined by (minus n) parameters instead of (minus 1) parameters. Viewed as a tensor, the joint probability is a rank 1 tensor. The main point is that the statistical independence assumption reduced the representation of the multivariate joint distribution from exponential to linear size.

Since our variables are typically divided to measurement variables and an output/class variable H (or in general ), it is useful to introduce another, weaker form, of independence known as conditional independence. Variables X, Y are conditionally independent given H, denoted by ) meaning that given H, the value of Y does not add any information about X. This is equivalent to the condition P(X, Y | H) = P(X | H)P(Y | H). The proof goes as follows:

• If P(X | Y, H) = P(X | H), then

• If P(X, Y | H) = P(X | H)P(Y | H), then

Consider as an example, Joe and Mo live on opposite sides of the city. Joe goes to work by train and Mo by car. Let X be the event ”Joe is late to work” and Y be the event ”Mo is late for work”. Clearly X and Y are not independent because there could be other factors. For example, a train strike will cause Joe to be late, but because of the strike there would be extra traffic (people using their car instead of the train) thus causing Mo to be pate as well. Therefore, a third variable H standing for the event ”train strike” would decouple X and Y .

From a computational standpoint, the conditional independence assumption has a similar effect to the unconditional independence. Let X range over k distinct value, Y range over r distinct values and H range over s distinct values. Then P(X, Y, H) is a 3-way array of size that means that ), a 2-way ”slice” of the 3-way array along the H axis is represented by the outer-product of two vectors ). As a result the 3-way array is represented by 2) parameters instead of 1. Likewise, if the n-way array ) (which is a slice along the H axis of the (n + 1)-array )) is represented by an outer-product of n vectors, i.e., by parameters.

We will use the ML principle to estimate the bias of a coin. Let X be a random variable taking the value {0, 1} and H would be our hypothesis taking a real value in [0, 1] standing for the coin’s bias. If the coin’s bias is . We receive m i.i.d. examples . We wish to determine the value of q. Given that , the ML problem we must solve is:

Let 0 stand for the number of ’0’ instances, i.e., 1, ..., m}|. Therefore our ML problem becomes:

Taking the partial derivative with respect to q and setting it to zero:

produces the result:

So far we considered constraints induced by conditional independent state- ments among the random variables as a means to reduce the space and time complexity of the multivariate distribution array. Another approach would be to assume some form of parametric form governing the entries of the array — the most popular assumption is Gaussian distribution ) with mean vector and covariance matrix E. The parameters of the density function are denoted by ) and for every vector we have:

Assume we are given an i.i.d sample of and we would like to find the Bayes optimal

by maximizing the likelihood (here we are assuming that the the priors are equal, thus the maximum likelihood and the MAP would produce the same result). Because the sample was drawn i.i.d. we can assume that:

Let ) and since Log is monotonously increasing we have that

be recovered by taking derivatives with respect to = 0. We have:

We will start with a simple scenario where , i.e., all the covariances are zero and all the variances are equal to . After substitution (and removal of items which do not depend on

The partial derivative with respect to

from which we obtain:

The partial derivative with respect to

from which we obtain:

Note that the reason for dividing by n is due to the fact that

In the general case, E is a full rank symmetric matrix, then the derivative of eqn. (1.1) with respect to

and since is full rank we obtain . For the derivative with respect to E we note two auxiliary items:

Using the fact that ) we can transform for any vector z. Given that is symmetric, then:

Substituting we obtain:

from which we obtain:

Consider another application of conditional dependence which is the Bayes incremental rule. Suppose we have processed and computed somehow ). We are given a new measurement X and wish to compute (update) the posterior ). We will use the chain rule

P(X | Y, Z) = P(X, Y, Z)P(Y, Z ) = P(Z | X, Y )P(X | Y )P(Y )P(Z | Y )P(Y ) = P(Z | X, Y )P(X | Y )P(Z | Y )

to obtain:

from conditional independence, ). The term ) can expanded as follows:

After substitution we obtain:

The old posterior ) is now the prior for the updated formula. Consider the following example: We have a coin which could be either fair or biased towards Head at a probability of 0be the event that the coin is fair, and that the coin is biased. We start with prior probabilities 25 (we have a higher initial belief that the coin is fair). Suppose our first coin toss is a Head, i.e., Then,

and 286. Our posterior belief that the coin is fair has gone down after a Head toss. Assume we have another measurement then:

and 325, thus our belief that the coin is fair continues to go down after Head tosses.

For the last topic in this lecture consider the 2-class inference problem. We will encountered this problem in this course in the context of SVM and LDA. In the Bayes framework, if denotes the ”class member” variable with two possible outcomes, then the MAP decision policy calls for

making the decision based on data x:

or in other words the class would be chosen if The decision surface (as a function of x) is therefore described by:

The questions we ask here is what would the Bayes optimal decision surface be like if we assume that the two classes are normally distributed with different means and the same covariance matrix? What we will see is that under the condition of equal priors ) the decision surface is a hyperplane — and not only that, it is the same hyperplane produced by LDA.

, the the Bayes optimal decision surface is a hyperplane . In other words, the decision surface is described by:

Proof: The decision surface is described by which is equivalent to the statement that the ratio of the posteriors is 1, or equivalently that the log of the ratio is zero, and using Bayes formula we obtain:

In other words, the decision surface is described by

After expanding the two terms we obtain eqn. (1.2).

Maximum Likelihood/ Maximum Entropy Duality

In the previous lecture we defined the principle of Maximum Likelihood (ML): suppose we have random variables form a random sample from a discrete distribution whose joint probability distribution is where ) is a vector in the sample and is a parameter from some parameter space (which could be a discrete set of values — say class membership). When ) is considered as a function of it is called the likelihood function. The ML principle is to select the value of that maximizes the likelihood function over the observations (training set) If the observations are sampled i.i.d. (a common, not always valid, assumption), then the ML principle is to maximize:

which due to the product nature of the problem it becomes more convenient to maximize the log likelihood. We will take a closer look today at the ML principle by introducing a key element known as the relative entropy measure between distributions.

The ML principle states that the empirical distribution of an i.i.d. sequence of examples is the closest possible (in terms of relative entropy which would be defined later) to the true distribution. To make this statement clear let X be a set of symbols ) be the probability (belonging to a parametric family with parameter ) of drawing a symbol be a sequence of symbols drawn i.i.d. according to P. The occurrence frequency f(a) measures the number of draws of the symbol

and let the empirical distribution be defined by

The joint probability ) is equal to the product which according to the definitions above is equal to:

The ML principle is therefore equivalent to the optimization problem:

where denote the set of n-dimensional probability vectors (”probability simplex”). Let ). Since argmax) and given that ln the solution to this problem can be found by setting the partial derivative of the Lagrangian to zero:

where is the Lagrange multiplier associated with the equality constraint 0 are the Lagrange multipliers associated with the inequality constraints 0. We also have the complementary slackness condition that sets

After setting the partial derivative with respect to to zero we get:

Assume for now that . Then from complementary slackness we must have = 0 (because We are left therefore with the result . Following the constraint = 1 we obtain . As a result we obtain: = 0 we could use the convention 0 ln 0 = 0 and from continuity arrive to

We have arrived to the following theorem:

Theorem 1 The empirical distribution estimate ˆP is the unique Maximum Likelihood estimate of the probability model Q on the occurrence frequency f().

This seems like an obvious result but it actually runs deep because the result holds for a very particular (and non-intuitive at first glance) distance measure between non-negative vectors. Let dist(f, p) be some distance measure between the two vectors. The result above states that:

for some (family?) of distance measures dist(). It turns out that there is only onesuch distance measure, known as the relative-entropy, which satisfies the ML result stated above.

The relative-entropy (RE) measure D(x||y) between two non-negative vectors is defined as:

In the definition we use the convention that 0 ln = 0 and based on con- tinuity that 0 ln are also probability vectors, i.e., belong to is also known as the Kullback-Leibler divergence. The RE measure is not a distance metric as it is not symmetric, ), and does not satisfy the triangle inequality. Nevertheless, it has several interesting properties which make it a fundamental measure in statistical inference.

The relative entropy is always non-negative and is zero if and only if x = y. This comes about from the log-sum inequality:

Thus,

But ) which follows from the inequality ln(x + 1) > x/(x + 1) (which holds for We can state the following theorem:

be the occurrence frequency on a training sample. ˆis a ML estimate iff

Proof:

and

There are two (related) interesting points to make here. First, from the proof of Thm. 1 we observe that the non-negativity constraint not be enforced - as long as 0 (which holds by definition) the closest p to f under the constraint the fact that the closest point p to f comes out as a scaling of f (which is by definition the empirical distribution ˆP) arises because of the relative-entropy measure. For example, if we had used a least-squares distance measure the result would not be a scaling of f. In other words, we are looking for a projection of the vector f onto the probability simplex, i.e., the intersection of the hyperplane = 1 and the non-negative orthant 0. Under relative-entropy the projection is simply a scaling of f (and this is why we do not need to enforce non-negativity). Under least-sqaures, a projection onto the hyper-plane = 1 could take us out of the non-negative orthant (see Fig. 2.1 for illustration). So, relative-entropy is special in that regard — it not only provides the ML estimate, but also simplifies the optimization process(something which would be more noticeable when we handle a latent class model next lecture).

The relative-entropy measure is not symmetric thus we expect different out- comes of the optimization min) compared to min). The lat-

Fig. 2.1. Projection of a non-neagtaive vector f onto the hyperplane Under relative-entropy the projection ˆP is a scaling of f (and thus lives in the probability simplex). Under least-squares the projection lives outside of the probability simplex, i.e., could have negative coordinates.

ter of the two, i.e., minis some empirical evidence and Q is some model, provides the ML estimation. For example, in the next lecture we will consider Q the set of low-rank joint distributions (called latent class model) and see how the ML (via relative-entropy minimization) solution can be found.

Let denote the entropy function. With regard to min) we can state the following observation:

Claim 2

Proof:

which follows from the condition

In other words, the closest distribution to uniform is achieved by maximizing the entropy. To make this interesting we need to add constraints. Consider a linear constraint on . To be concrete, con- sider a die with six faces thrown many times and we wish to estimate the probabilities given only the . Say, the average is 3.5 which is what one would expect from an unbiased die. The ciple of insufficient reasoning calls for assuming uniformity unless there is additional information (a controversial assumption in some cases). In other words, if we have no information except that each 0 and that we should choose the uniform distribution since we have no reason to choose any other distribution. Thus, employing Laplace’s principle we would say that if the average is 3.5 then the most ”likely” distribution is the uniform. What if 2? This kind of problem can be stated as an optimization problem:

where 2. We have now two constraints and with the aid of Lagrange multipliers we can arrive to the result:

Note that because of the exponential 0 and again ”non-negativity comes for free”. Following the constraint = 1 we get exp1from which obtain:

where Z (a function of ) is a normalization factor and needs to be set by using (see later). There is nothing special about the uniform distribution, thus we could be seeking a probability vector p as close as possible to some prior probability under the constraints above:

with the result:

We could also consider adding more linear constraints on p of the form: . The result would be:

Probability distributions of this form are called Gibbs Distributions. In practical applications the linear constraints on p could arise from average information about the system such as temperature of a fluid (where the probabilities of the particles moving at various velocities), rainfall data or general environmental data (where represent the probability of finding animal colonies at discrete locations in a 3D map). A constraint of the form states that the expectation ] should be equal to the empirical distribution is either uniform or given as input. Let

and

We could therefore consider looking for the ML solution for the parameters of the Gibbs distribution:

where if ˆp is uniform then min ) can be replaced by max (because

As it turns out, the MaxEnt and ML are duals of each other and the intersection of the two sets contains only a single point which solves both problems.

Theorem 3 The following are equivalent:

In practice, the duality theorem is used to recover the parameters of the Gibbs distribution using the ML route (second line in the theorem above) — the algorithm for doing so is known as the iterative scaling algorithm (which we will not get into).

EM Algorithm: ML over Mixture of Distributions

In Lecture 2 we saw that the Maximum Likelihood (ML) principle over i.i.d. data is achieved by minimizing the relative entropy between a model Q and the occurrence-frequency of the training data. Specifically, let i.i.d. where each -tupple of symbols taken from an alphabet X having n different letters be the empirical joint distribution, i.e., an array with d dimensions where each axis has n entries, i.e., each entry ˆ, represents the (normalized) co-occurrence of the in the training set . We wish to find a joint distribution -array) which belongs to some model family of distributions Q closest as possible to ˆP in relative-entropy:

In this lecture we will focus on a model of distributions Q which represents mixtures of simple distributions H— known as latent class models. A latent class model arises when the joint probability ) we observe (i.e., from which ˆP is generated by observing samples ) is in fact a marginal of is a ”hidden” (or ”latent”) random variable which has k different discrete values

The idea is that given the value of the hidden variable H the problem of recovering the model ), which belongs to some family of joint distributions H, is a relatively simple problem. To make this idea clearer we consider the following example: Assume we have two coins. The first coin has a probability of heads (”0”) equal to p and the second coin has a probability of heads equal to q. At each trial we choose to toss coin 1 with probability and coin 2 with probability 1 . Once a coin has been chosen it is tossed 3 times, producing an observation given a set of such observations where each observation is a triplet of coin tosses (the same coin). Given D, we can construct the empirical distribution ˆP which is a 2 2 array defined as:

Let be a random variable associated with the observation that was generated by coin 1 and was generated by coin 2. If we knew the values of then our task would be simply to estimate two separate Bernoulli distributions by separating the triplets generated from coin 1 from those generated by coin 2. Since is not known, we have the marginal:

where (is a triplet coin toss and 0 number of heads (”0”) in the triplet of tosses. In other words, the likelihood P(x) of triplet of tosses ) is a linear combination (”mixture”) of two Bernoulli distributions. Let H stand for Bernoulli distributions:

where stands for the outer-product of with itself d times, i.e., an n- way array indexed by , and whose value there is equal to . The model family Q is a mixture of Bernoulli distributions:

where specifically for our coin-toss example becomes:

We see therefore that the eight entries of which minimizes over the set Q is determined by three parameters . For the coin-toss

example this looks like:

where . Trying to work out an algorithm for minimizing the unknown parameters would be somewhat ”unpleasant” (and even more so for other families of distributions H) because of the log-over-a-sum present in the optimization function — if we could somehow turn this into a sum-over-log our task would be much easier. We would then be able to turn the problem into a succession of problems over H rather than a single problem over Another point worth attention is the non-negativity of the output variables — simply minimizing the relative-entropy measure under the constraints of the class model Q would not guarantee a non-negative solution. As we shall see, breaking down the problem into a successions of problems over H would give us the ”non-negativity for free” feature.

The technique for turning the log-over-sum into a sum-over-log as part of finding the ML solution for a mixture model is known as the ExpectationMaximization (EM) algorithm introduced by Dempster, Laird and Rubin in 1977. It is based on two ideas: (i) introduce auxiliary variables, and (ii) use of Jensen’s inequality.

Let represent the training data where is taken from some instance space X which we leave unspecified. For now, we leave matters to be as general as possible and specifically we do not make independence assumptions on the data generation process.

The ML problem is to find a setting of parameters which maximizes the likelihood ), namely, we wish to maximize parameters , which is equivalent to maximizing the log-likelihood:

where y represents the hidden variables. We will denote Let ) be some (arbitrary) distribution of the hidden variables y conditioned on the parameters and the input sample 1. We define a ) as follows:

The inequality comes from Jensen’s inequality log when = 1. What we have obtained is an ”auxiliary” function satisfying

for all distributions ). The maximization of ) proceeds by interleaving the variables as we separately ascend on each set of variables. At the (t + 1) iteration we fix the current value of of the t’th iteration and maximize , and then maximize

The strategy of the EM algorithm is to maximize the lower bound with the hope that if we ascend on the lower bound function we will also ascend with respect to ). The claim below guarantees that an ascend on Q will also generate an ascend on L:

Claim 3 (Jordan-Bishop) The optimal at each step is

Proof: We will show that ) which proves the claim since , thus the best q-distribution we can

hope to find is one that makes the lower-bound meet

The proof provides also the validity for the approach of ascending along the lower bound ) because at the point the two functions coincide, i.e., the lower bound function at is equal to ) therefore if we continue and to ascend along as well— therefore, convergence is guaranteed. It can also be shown (but omitted here) that the point of convergence is a stationary point of shown originally by C.F. Jeff Wu in 1983 years after EM was introduced in 1977) under fairly general conditions. The second step of maximizing over then becomes:

This defines the EM algorithm. Often the ”Expectation” step is described as taking the expectation of:

followed by a Maximization step of finding that maximizes the expectation — hence the term EM for this algorithm.

Eqn. 3.8 describes a principle but not an algorithm because in general, without making assumptions on the statistical relationship between the data points and the hidden variable the problem presented in eqn. 3.8 is unwieldy. We will reduce eqn. 3.8 to something more manageable by making the i.i.d. assumption. This is detailed in the following section.

The EM optimization presented in eqn. 3.8 can be simplified if we assume the data points (and the hidden variable values) are i.i.d.

and

For any

this is because ) = 1. Substituting the simplifications above into eqn. 3.8 we obtain:

where y

We will apply the EM scheme to our running example of mixture of Bernoulli distributions. We wish to compute

and then maximize Q() with respect to

Q(θ, θ

where ). The values of known since ) are given from the previous iteration. The Bayes formula is used to compute

We wish to compute: max). The partial derivative with respect

from which we obtain the update formula of

The partial derivative with respect to p is:

from which we obtain the update formula:

Likewise the update rule for q is:

To conclude, we start with some initial ”guess” of the values of pute the values of and update iteratively the values of where at the end of each iteration the new values of are computed.

The Gaussian mixture model assumes that is a linear combination of Gaussian distributions

where

is Normally distributed with mean and covariance matrix be the i.i.d sample data and we wish to solve for the mean and covariances of the individual Gaussians (the ”factors”) and the mixing coefficients ). In order to make clear where the parameters are located we will write ) instead of are the mean and variance of the j’th factor. We denote by the collection of mixing coefficients be auxiliary variables per point and per factor y = j standing for:

The EM step (eqn. 3.9) is:

Note the constraint = 1. The update formula for is done through the use of Bayes formula:

where is a scaling factor so that The update formula for follow by taking partial derivatives of eqn. (3.10) and setting them to zero. Taking partial derivatives with respect

to we obtain the update rules:

In other words, the observations are weighted by before a Gaussian is fitted (k times, one for each factor).

The Gaussian mixture model is classically used for clustering applications. In a clustering application one receives a sample of points each point resides in . The task of the learner (in this case ”unsupervised” learning) is to group the m points into i = 1, ..., m stands for the required labeling. The clustering solution is an assignment of values to according to some clustering criteria.

In the Gaussian mixture model points are clustered together if they arise from the same Gaussian distribution. The EM algorithm provides a probabilistic assignment ) which we denoted above as

The multinomial mixture (the coins example we toyed with) is typically used for representing ”count” data, such as when representing text documents as high-dimensional vectors. A vector representation of a text document associates a word from a fixed vocabulary to a coordinate entry of the vector. The value of the entry represents the number of times that particular word appeared in the document. If we ignore the order in which the words appeared and count only their frequency, a set of documents set of words could be jointly represented by a co-occurence matrix contains the number of times word appeared in document . If we scale = 1 then we have a distribution P(w, d). This kind of representation of a set of documents is called ”bag of words”.

For purposes of search and filtering it is desired to reveal additional information about words and documents such as to which ”topic” a document belongs to or to which topics a word is associated with. This is similar to a clustering task where documents associated with the same topic are to be clustered together. This can be achieved by considering the topics as the value of a latent variable y:

where we made the assumption that (i.e., words and documents are conditionally independent given the topic). The conditional independent assumption gives rise to the multinomial mixture model. To be more specific, ley denote the k possible topics and let that = 1), then the latent class model becomes:

Note that P(w | y = j) is a vector which we denote as j) is also a vector we denote by stands for the outer-product of the two vectors, i.e., is a rank-1 matrix. The Maximum-Likelihood estimation problem is therefore to find vectors and scalars such that the empirical distribution represented by the unit scaled matrix G is as close as possible (in relative-entropy measure) to the low-rank matrix the constraints of non-negativity and are unit-scaled as well (

Let )) stand for the i’th example i = 1, ..., q where an example is a pair of word and document where is the index to the word alphabet and is the index to the document. The EM algorithm involves the following optimization step:

An update rule for ’th entry of ) is derived below: the derivative

of the Lagrangian is:

where N(r) stands for the frequency of the word in all the documents . Note that N(r) is the result of summing-up the that the vector N(1), ..., N(n) is the marginal the constraint = 1 we obtain the update rule:

Update rules for the remaining unknowns are similarly derived. Once EM has converged, then is the probability of the word to the j’th topic and is the probability that the s’th document comes from the

Support Vector Machines and Kernel Functions

In this lecture we begin the exploration of the 2-class hyperplane separation problem. We are given a training set of instances and class labels 1 (i.e., the training set is made up of “positive” and “negative” examples). We wish to find a hyperplane direction and an offset scalar 0 for positive examples and 0 for negative examples — which together means that the margins 0 are positive.

Assuming that such a hyperplane exists, clearly it is not unique. We therefore need to introduce another constraint so that we could find the most “sensible” solution among all (infinitley many) possible hyperplanes which separate the training data. Another issue is that the framework is very limited in the sense that for most real-world classification problems it is somewhat unlikely that there would exist a linear separating function to begin with. We therefore need to find a way to extend the framework to include non-linear decision boundaries at a reasonable cost. These two issues will be the focus of this lecture.

Regarding the first issue, since there is more than one separating hyper-plane (assuming the training data is linearly separable) then the question we need to ask ourselves is among all those solutions which of them has the best “generalization” properties? In other words, our goal in constructing a learning machine is not necessarily to do very well (or perfect) on the training data, because the training data is merely a sample of the instance space, and not necessarily a “representative” sample — it is simply a sample. Therefore, doing well on the sample (the training data) does not necessarily guarantee (or even imply) that we will do well on the entire instance space. The goal of constructing a learning machine is to maximize the performance on the test data (the instances we haven’t seen), which in turn means that we wish to generalize “good” classification performance on the training set onto the entire instance space.

A related issue to generalization is that the distribution used to generate the training data is unknown. Unlike the statistical inference material we had so far, this time we will not attempt to estimate the distribution. The reason one can derive optimal learning algorithms yet bypass the need for estimating distributions would be explained later in the course when PAClearning will be introduced. For now we will focus only on the algorithmic aspect of the learning problem.

The idea is to consider a subset of all hyperplanes which have a fixed margin where the margin is defined as the distance of the closest training point to the hyperplane:

The Support Vector Machine (SVM), first introduced by Vapnik and his colleagues in 1992, seeks a separating hyperplane which simultaneously minimizes the empirical error and maximizes the margin. The idea of maximizing the margin is intuitively appealing because a decision boundary which lies close to some of the training instances is less likely to generalize well because the learning machine will be susceptible to small perturbations of those instance vectors. A formal motivation for this approach is deferred to the PAC-learning material we will introduce later in the course.

4.1 Large Margin Classifier as a Quadratic Linear Programming

We would first like to set up the linear separating hyperplane as an optimization problem which is both consistent with the training data and maximizes the margin induce by the separating hyperplane over all possible consistent hyperplanes.

defined by

Since we are allowed to scale the parameters w, b at will (note that if 0) we can set the distance between the boundary points to the hyperplane to be 1by scaling w, b such the point(s) with smallest margin (closest to the hyperplane) will be normalized: = 1, therefore the margin is simply 2Fig. 5.1). Note that argmaxis equivalent to argmaxwhich in turn is equivalent to argminSince all positive points and negative points should be farther away from the boundary points we also have the separability constraints is a positive instance and is a negative instance. Both separability constraints can be combined: 1. Taken together, we have defined the following optimization problem:

This type of optimization problem has a quadratic criteria function and linear inequalities and is known in the literature as a Quadratic Linear Programming (QP) type of problem.

This particular QP, however, requires that the training data are linearly separable — a condition which may be unrealistic. We can relax this condition by introducing the concept of a “soft margin” in which the separability holds approximately with some error:

Where is some pre-defined weighting factor. The (non-negative) variables allow data points to be miss-classified thereby creating an approximate separation. Specifically, if is a positive instance (= 1) then the “soft” constraint becomes:

where if = 0 we are back to the original constraint where is either a boundary point or laying further away in the half space assigned to positive instances. When 0 the point can reside inside the margin or even in the half space assigned to negative instances. Likewise, if is a negative instance (1) then the soft constraint becomes:

Fig. 4.1. Separating hyperplane w, b with maximal margin. The boundary points are associated with non-vanishing Lagrange multipliers 0 and margin errors are associated with 0 where the criteria function encourages a small number of margin errors.

The criterion function penalizes (the -norm) for non-vanishing overall system will seek a solution with few as possible “margin errors” (see Fig. 5.1). Typically, when possible, an norm is preferable as the overly weighs high magnitude outliers which in some cases can dominate the energy function. Another note to make here is that strictly speaking the ”right thing” to do is to penalize the margin errors based on the norm , i.e., the number of non-zero entries, and drop the balancing parameter . This is because it does not matter how far away a point is from the hyperplane — all what matters is whether a point is classified correctly or not (see the definition of empirical error in Lecture 4). The problem with that is that the optimization problem would no longer be convex and non-convex problems are notoriously difficult to solve. Moreover, the class of convex optimization problems (as the one described in Eqn. 4.3) can be solved in polynomial time complexity.

So far we have described the problem formulation which when solved would provide a solution with “sensible” generalization properties. Although we can proceed using an off-the-shelf QLP solver, we will first pursue the ”dual” problem. The dual form will highlight some key properties of the approach and will enable us to extend the framework to handle non-linear decision surfaces at a very little cost. In the appendix we take a brief tour on the basic principles associated with constrained optimization, the Karush-Kuhn-Tucker (KKT) theorem and the dual form. Those are recommended to read before moving to the next section.

We return now to the primal problem (eqn. 6.3) representing the maximal margin separating hyperplane with margin errors:

We will now derive the Lagrangian Dual of this problem. By doing so a new key property will emerge facilitated by the fact that the criteria function ) (note there are no equality constraints thus there is no need for involves only inner-products of the training instance vectors . This property will form the key of mapping the original input space of dimension n to a higher dimensional space thereby allowing for non-linear decision surfaces for separating the training data.

Note that with this particular problem the strong duality conditions are satisfied because the criteria function and the inequality constraints form a convex set. The Lagrangian takes the following form:

Recall that

Since the minimum is obtained at the vanishing partial derivatives of the Lagrangian with respect to w, b, the next step would be to evaluate those

constraints and substitute them back into L() to obtain

From the first constraint (4.4) we obtain scribed by a linear combination of a subset of the training instances. The reason that not all instances participate in the linear superposition is due to the KKT conditions: 1, i.e., the instance is classified correctly and is not a boundary point, and conversely, when , i.e., when is a boundary point or when is a margin error (0) — note that for a margin error instance the value of would be the smallest possible required to reach an equality in the constraint because the criteria function penalizes large values of boundary points (and the margin errors) are called support vectors thus w is defined by the support vectors only. The third constraint (4.6) is equivalent to the constraint:

since 0. Also note that if 0, i.e., point is a margin-error point, then by KKT conditions we must have = 0. As a result . Therefore based on the values of alone we can make the following classifications:

Substituting these results/constraints back into the Lagrangian L() we obtain the dual problem:

The criterion function ) can be written in a more compact manner as follows: Let matrix whose entries are is the vector of (1is the vector (is the transpose (row vector). Note that M is positive 0 for all vectors = 0 — a property which will be important later.

The key feature of the dual problem is not so much that it is simpler than the primal (in fact it isn’t since the primal has no equality constraints) or that it has a more “elegant” feel, the key feature is that the problem is completely described by the inner products of the training instances i = 1, ..., m. This fact will be shown to be a crucial ingredient in the so called “kernel trick” for the computation of inner-products in high dimensional spaces using simple functions defined on pairs of training instances.

We ended with the dual formulation of the SVM problem and noticed that the input data vectors are represented by the Gram matrix M. In other words, only inner-products of the input vectors play a role in the dual formulation — there is no explicit use of or any other function of inner-products. This observation suggests the use of what is known as the ”kernel trick” to replace the inner-products by non-linear functions.

The common principle of kernel methods is to construct nonlinear variants of linear algorithms by substituting inner-products by nonlinear kernel functions. Under certain conditions this process can be interpreted as mapping of the original measurement vectors (so called ”input space”) onto some higher dimensional space (possibly infinitely high) commonly referred to as the ”feature space”. Mathematically, the kernel approach is defined as follows: let be vectors in the input space, say sider a mapping is an inner-product space. The kernel-trick is to calculate the inner-product in F using a kernel function ), while avoiding explicit mappings (evaluation of)

Common choices of kernel selection include the d’th order polynomial kernels and the Gaussian RBF kernels exp(). If an algorithm can be restated such that the input vectors appear in terms of inner-products only, one can substitute the inner-products by such a kernel function. The resulting kernel algorithm can be interpreted as running the original algorithm on the space F of mapped objects

We know that M of the dual form is positive semi-definite because M

can be written is ]. Therefore 0 for all choices of x (which means that the eigenvalues of M are non-negative). If the entries of M are to be replaced with the condition we must enforce on the function k() is that it is a positive definite kernel function. A positive definite function is defined such that for any set of vectors and for any values of q the matrix K whose entries are ) is positive semi-definite. Formally, the conditions for admissible kernels k() are known as Mercer’s conditions summarized below:

be symmetric and continuous. The following conditions are equivalent:

(iii) for all and for all q, the matrix is positive semi-definite.

Perhaps the non-obvious condition is No. 1 which allows for the feature map () to have infinitely many coordinates (a vector in Hilbert space). For example, as we shall see below, the kernel exp() is an inner- product of two vectors with infinitely many coordinates. We will consider next a number of popular kernels.

Let 0 is a natural number. Then, the corresponding feature map ) coordinates which take the value:

where!) is the multinomial coefficient (number of ways to distribute d balls into k bins where the j’th bin hold exactly balls):

The dimension of the vector space can be measured using the following combinatorial problem: how many arrangements of partitions to be placed among d items? the answer is). For example, k = d = 2 gives us :

where

The feature map ) contains all monomials whose power is lesser or equal to . This can be acheived by increasing the dimension to is used to fill the gap between Therefore the dimension of

(

Therefore, the entries of the vector ) take the values:

For example, k = d = 2 gives us :

where ). In this example, mapping from and hyperplanes to

Assume we would like to find a separating conic (Parabola, Hyperbola, Ellipse) function rather than a line in . The discussion so far suggests we construct the Gram matrix M in the dual form with the d = 2 polynomial kernel for some parameter of our choosing. The extra effort we will need to invest is negligible — simply replace every occurrence

The function (RBF) is a kernel function but with an infinite expansion. Without loss of generality let = 1, then we have:

From which we can see that the entries of the feature map

By adopting some kernel k() we are in fact mapping then proceed to solve for using some QLP solver. The QLP solution of the dual form will yield the solution for the Lagrange multipliers . We saw from eqn. (4.4) that we can express ) as a function of the (mapped) examples:

Rather than explicitly representing ) — a task which may be prohibitly expensive since in general the dimension of the feature space of a polynomial mapping is— we store all the support vectors (those input vectors with corresponding 0) and use them for the evaluation of test examples:

We see that the kernel trick enabled us to look for a non-linear separating surface by making an implicit mapping of the input space onto a higher dimensional feature space using the same dual form of the SVM formulation — the only change required was in the way the Gram matrix was constructed. The price paid for this convenience is to carry all the support vectors at the time of classification f(x).

A couple of notes may be worthwhile at this point. The constant b can be recovered from any of the support vectors. Say, is a positive support vector (but not a margin error, i.e., from which b can be recovered. The second note is that the number of support vectors is typically around 10% of the number of training examples (empirically). Thus the computational load during evaluation of f(x) may be relatively high. Approximations have been proposed in the literature by looking for a reduced number of support vectors (not necessarily aligned with the training set) — but this is beyond the scope of this course.

The kernel trick gained its popularity with the introduction of the SVM but since then has taken a life of its own and has been applied to principal component analysis (PCA), ridge regression, canonical correlation analysis (CCA), QR factorization and the list goes on. We will meet again with the kernel trick later on.

Spectral Analysis I: PCA, LDA, CCA

In this lecture (and the following one) we will focus on spectral methods for learning. Today we will focus on dimensionality reduction using Principle Component Analysis (PCA), multi-class learning using Linear Discriminant Analysis (LDA) and Canonical Correlation Analysis (CCA). In the next lecture we will focus on spectral clustering methods.

Dimensionality reduction appears when the dimension of the input vector is very large (imagine pixels in an image, for example) while the coordinate measurements are highly inter-dependent (again, imagine the redundancy present among neighboring pixels in an image). High dimensional data impose computational efficiency challenges and often translate to poor generalization abilities of the learning engine (see lectures on PAC). A dimensionality reduction can also be viewed as a feature extraction process where one takes as input a large feature set (the original measurements) and creates from them a much smaller number of new features which are then fed into the learning engine.

In this lecture we will focus on feature extraction from a very specific (and constrained) stanpoint. We would be looking for a mixing (linear combination) of the input coordinates such that we obtain a linear projection from . In doing so we wish to reduce the redundancy while preserving as much as possible the variance of the data. From a statistical standpoint this is achieved by transforming to a new set of variables, called principal components, which are uncorrelated so that the first few retain most of the variation present in all of the original coordinates. For example, in an image processing application the input images are highly redundant where neighboring pixel values are highly correlated. The purpose of feature extraction would be to transform the input image into a vector of output components with the least redundancy possible. Form a geometric standpoint, this is achieved by finding the ”closest” (in least squares sense) linear q-dimensional susbspace to the m sample points S. The new subspace is a lower dimensional ”best approximation” to the sample S. These two, equivalent, perspectives on data compression (dimensionality reduction) form the central idea of principal component analysis (PCA) which probably the oldest (going back to Pearson 1901) and best known of the techniques of multivariate analysis in statistics. The computation of PCA is very simple and the definition is straightforward, but has a wide variety of different applications, a number of different derivations, quite a number of different terminologies (especially outside the statistical literature) and is the basis for quite a number of variations on the basic technique.

We then extend the variance preserving approach for data representation for labeled data sets. We will describe the linear classifier approach (separating hyperplane) form the point of view of looking for a hyperplane such that when the data is projected onto it the separation is maximized (the distance between the class means is maximal) and the data within each class is compact (the variance/spread is minimized). The solution is also produced, just like PCA, by a spectral analysis of the data. This approach goes under the name of Fisher’s Linear Discriminant Analysis (LDA).

What is common between PCA and LDA is (i) the use of spectral matrix analysis — i.e., what can you do with eigenvalues and eigenvectors of matrices representing subspaces of the data? (ii) these techniques produce optimal results for normally distributed data and are very easy to implement. There is a large variety of uses of spectral analysis in statistical and learning literature including spectral clustering, Multi Dimensional Scaling (MDS) and data modeling in general. Another point to note is that this is the first time in the course where the type of data distribution plays a role in the analysis — the two techniques are defined for any distribution but are optimal only under the Gaussian distribution.

We will also describe a non-linear extension of PCA known as KernelPCA, but the focus would be mostly on PCA itself and its analysis from a couple of vantage points: (i) PCA as an optimal reconstruction after a dimension reduction, i.e., data compression, and (ii) PCA for redundancy reduction (decorrelation) of the output components.

Let be our sample data S of vectors in , arranged as columns of a matrix A. It will be convenient to assume that the data is centered, i.e., = 0. If the data is not centered we can always center it by computing the mean vector and replace the original data sample with the new sample . In a statistical sense, the coordinates of the vector are considered as random variables, thus a row in the matrix A is the sample of values of a particular random variable, drawn from some unknown probability distribution, associated with the row position. We wish to find vectors (arranged as columns of a matrix U), where ), such that the new feature measurements (who are the result of linear combinations of the original feature measurements x) have certain desirable properties.

The idea property to seek from the new coordinates y is statistical independence, i.e., ) which would mean that we have removed the redundancy of the original data x in the best possible manner. This goal, however, is too much to ask from a linear transformation and instead we would ask for a weaker property to hold: that the pairwise covariance ) = 0 vanishes, i.e., that the covariance matrix on the new coordinates is diagonal. A diagonal covariance insures some redundancy removal, but not as good as statistical independence. However, when the data is Normally distributed Σ) with mean and covariance Σ, then the transformation which diagonalizes the covariance matrix also guarantees statistical independence. Among all transformations that de-correlate the data we will seek the one that maximizes the spread (variance) of the sample data after being projected onto the new axes vectors.

The property we would like to maximize is that the projection of the sample data on the new axes is as spread as possible. To start this analysis, assume q = 1, i.e., the n components of the input vector x are reduced to a single output component . We are looking for a single vector whose direction maximizes the variance of the output component y.

Formally, we are looking for a unit vector u which maximizes (see Appendix A for basic statistical definitions and note that E[y] = 0 because ) = 0). In other words, the projected points onto the axis represented by the vector u are as spread as possible (in a least squares sense). In vector notation, the optimization problem takes the following form:

The Lagrangian of the problem is:

By taking the partial derivative = 0 we obtain the following necessary condition (see Appendix B):

which tells us that u is an eigenvector of the (symmetric and positive definite) matrix . There are n eigenvectors associated with we can easily convince ourselves that we are looking for the one associated with the maximal eigenvalue: substitute instead of in the criterion function and since the eigenvalues must be positive (since is positive definite), then the optimum is obtained for the maximal eigenvalue. The leading eigenvector is called the first principal axis of the data sample represented by the columns of the matrix is called the first principal component of the data sample.

For convenience, we denote as the leading eigenvector and eigenvalue of . Next, we look for with and which has maximum variance (and so on for Two random variables are uncorrelated if their covariance vanishes. By definition of covariance (see Appendix A) we obtain:

We can therefore use the condition = 0 to specify zero correlation between . The functional to be optimized becomes:

with the Lagrangian being:

By taking the partial derivative with respect to we obtain the necessary condition:

Multiply the equation by from the left:

and noting from above that = 0 we obtain result we obtain:

so once more we have that form an eigenvalue/eigenvector pair of As before, should be as large as possible from the remaining spectral decomposition. By induction, it can be shown that the remaining principal vectors are the decreasing order eigenvactors of and the variance of the i’th principal component

Taken together, the PCA is the solution of the following optimization problem:

It will be useful for later to write the optimization function in a more concise manner as follows. Let matrix whose columns are diagonal matrix and from above we have that . Using the fact that

trace(B) we can convert ) as follows:

Thus, PCA becomes the solution of the following optimization function:

The solution, as saw above, is that ] consists of the decreasing order eigenvectors of . At the optimum, ) is equal to trace(D) which is equal to the sum of eigenvalues

It is worthwhile noting that when , and the PCA transform is a change of basis in known as Karhunen-Loeve transform.

To conclude, the PCA transform looks for q orthogonal direction vectors (called the principal axes) such that the projection of input sample vectors onto the principal directions has the maximal spread, or equivalently that the variance of the output coordinates is maximal. The principal directions are the leading (with respect to descending eigenvalues) q eigenvectors of the matrix , the principal directions form a basis of with the property of de-correlating the data and maximizing the variance of the coordinates of the sample input vectors.

In the previous section we saw that PCA generates a new coordinate system where the coordinates in the new system are uncorrelated. This means that the covariance matrix over the principle components should be diagonal. In this section we will explore this perspective in more detail.

The covariance matrix of the sample data with zero mean is

therefore the matrix we derived above is a scaled version of the covariance of the sample data (see Appendix A). The scale factor 1/m was unimportant in the process above because the eigenvectors are of unit norm, thus any scale of would produce the same set of eigenvectors.

The off-diagonal entries of the covariance matrix represent the correlation (a measure of statistical dependence) between the i’th and j’th component vectors, i.e., the entries of the input vectors x. The existence of correlations among the components (features) of the input signal is a sign of redundancy, therefore from the point of view of transforming the input representation into one which is less redundant, we would like to find a transformation with an output representation y which is associated with a diagonal covariance matrix , i.e., the components of y are uncorrelated.

Formally, , therefore we wish to find an matrix for which is diagonal. If in addition, we would require that the variance of the output coordinates is maximized, i.e., ) is maximal (but then we need to constrain the length of the column vectors of U, i.e., set = 1) then we would get a unique solution for U where the columns are orthonormal and are defined as the first q eigenvectors of the covariance matrix . This is exactly the optimization problem defined by eqn. (5.1).

We see therefore that PCA “decorrelates” the input data. Decorrelation and statistical independence are not the same thing. If the coordinates are statistically independent then the covariance matrix is diagonal, but it does not follow that uncorrelated variables must be statistically independent — covariance is just one measure of dependence. In fact, the covariance is a measure of pairwise dependency only. However, it is a fact that uncorrelated variables are statistically independent if they have a multivariate normal distribution (a Gaussian). In other words, if the sample data x are drawn from a probability distribution p(x) which has Gaussian form, the PCA transforms the sample data into a statistically independent set of variables . The details are explained below.

Recall that a multivariate normal distribution of the random variables is defined as

Also recall that a linear combination of the variables produces also a normal distribution

therefore choose is a diagonal matrix ). We have in that case:

which can be written as a product of univariate normal distributions

which proves the assertion that decorrelated normally distributed variables are statistically independent.

A different, yet equivalent, perspective on the PCA transformation is as an optimal reconstruction (in a least squares sense) after a dimension reduction. We are given a sample data as before and we are looking for a small number of orthonormal principal vectors which define a q-dimensional linear subspace of approximate the original input vectors in a least squares sense. In other words, the projection ˆof the sample points onto the q-dimensional subspace should minimize over all possible q-dimensional subspaces of

Let U be the subspace spanned by the principal vectors (columns of U) and let projection matrix mapping a point projection ˆ. From the definition of projection, the vector be orthogonal to the subspace ) be the coordinates of ˆx with respect to the principal vectors, i.e., . Then, from orthogonality we have that (= 0 for all vectors . Since this is true for all = 0. Therefore, result the projection matrix P becomes:

satisfying . In the case the columns of U are orthonormal, we have . We are ready now to describe the optimization problem on U: we wish to find an orthonormal set of principal vectors, such that is minimized.

Note that is the square Frobenious norm of a matrix. The optimal reconstruction problem therefore becomes:

We will show now that:

which shows that the optimal reconstruction problem is solved by PCA (recall Eqn. 5.1). From the identity ), we have:

Expanding the right hand side gives us:

trace((A

The second and third term are equal (commutativity of trace) and is also equal to the 4th term due to commutativity of the trace and Taken together:

To conclude, we have proven that by taking the first q eigenvectors of we obtain a linear subspace which is as close as possible (in a least squares sense) to the original sample data. Hence, PCA can be viewed as a vehicle for optimal reconstruction after dimension reduction. The optimization problem whose solution is the leading q eigenvectors of is described in eqn. 5.1:

Consider the situation where n, the dimension of the input vectors, is relatively large compared to the number of sample vectors m. For example, consider input vectors representing 50 50 sized images of faces, i.e., n = 2500, where m = 100. In other words, we are looking for a small number of “face templates” (known as “eigenfaces”) which approximate well the original set of 100 face images. In this case, is very large, 25002500, yet the number of non-vanishing eigenvalues cannot be higher than 100. Given that the eigendecomposition process is ), the computational burden would be very high. However, it is possible to perform an eigendecomposition on 100 matrix) instead, as shown next.

Let the columns of Q be the first q < m eigenvectors of QD where D is diagonal containing the corresponding eigenvalues. After pre-multiplying both sides by A we obtain:

from which we conclude that AQ contains the first q eigenvectors (but unnormalized) of . We have therefore that

where we used the fact that . Note that eigenvalues of and are the same (because

We can take the case n >> m described in the previous section one step further and consider such large values of n which are practically uncomputable — a situation which results when mapping the original input vectors to a high dimensional space: For example, ) representing the d’th order monomials of the coordinates of which is exponential in d. The mappings of interest are those which are paired with a non-linear kernel function:

Performing PCA on )] is equivalent to finding the

non-linear surface in (the nature of the non-linearity depends on the choice of ()) which best approximates the original sample data The problem is that is not computable — however is computable because (

From the previous section, contains the first q eigenvectors of are computable). Since A itself is not computable we cannot represent U explicitly, but we can project a new vector ) onto the principal directions and obtain the principal components, i.e., the output vector ), as follows.

Given the principal components (entries of measure, for example, the ) and the projection ˆonto the linear subspace spanned by (without the need to explicitly compute the principal axes ), as follows.

We now extend the variance preserving approach for data representation for labeled data sets. We will focus on 2-class sets and look for a separating hyperplane:

such that x belongs to the first class if f(x) > 0 and x belongs to the second class if f(x) < 0. In the statistical literature this type of function is called a linear discriminant function. The decision boundary is given by the set of points satisfying f(x) = 0 which is a hyperplane. Fisher’s (1936) Linear Discriminant Analysis (LDA) is a variance preserving approach for finding a linear discriminant function.

Fig. 5.1. Linear discriminant analysis based on class centers alone is not sufficient. Seeking a projection which maximizes the distance between the projected centers will prefer the horizontal axis over the vertical, yet the two classes overlap on the horizontal axis. The projected distance along the vertical axis is smaller yet the classes are better separated. The conclusion is that the sample variance of the two classes must be taken into consideration as well.

We will then introduce another popular statistical technique called Canonical Correlation Analysis (CCA) for learning the mapping between input and output vectors using the notion ”angle” between subspaces.

What is common in the three techniques PCA, LDA and CCA is the use of spectral matrix analysis — i.e., what can you do with eigenvalues and eigenvectors of matrices representing subspaces of the data? These techniques produce optimal results for normally distributed data and are very easy to implement. There is a large variety of uses of spectral analysis in statistical and learning literature including spectral clustering, Multi Dimensional Scaling (MDS) and data modeling in general.

To appreciate the general idea behind Fisher’s LDA consider Fig. 5.1. Let the centers of classes one and two be denoted by respectively. A linear discriminant function is a projection onto a 1D subspace such that the classes would be separated the most in the 1D subspace. The obvious first step in this kind of analysis is to make sure that the projected centers ˆwould be separated as much as possible. We can easily see that the direction of the 1D subspace should be proportional to as follows:

The right-hand term is maximized when . As illustrated in Fig. 5.1, this type of consideration is not sufficient to capture separability in the projected subspace because the spread (variance) of the data points around their centers also play an important role. For example, the horizontal axis in the figure separates the centers better than the vertical axis but on the other hand does a worse job in separating the classes themselves because of the way the data points are spread around their centers. The argument in favor of separating the centers would work if the data points were living in a hyper-sphere around the centers, but will not be sufficient otherwise.

The basic idea behind Fisher’s LDA is to consider the sample covariance matrix of the individual classes as well as their centers, in the following way. The optimal 1D projection would that which maximizes the variance of the projected centers while minimizes the variance of the projected data points of each class separately. Mathematically, this idea can be implemented by maximizes the following ratio:

and likewise,

where ˆ

We will now formalize this approach and derive its solution. We will begin with a general description of a multiclass problem where the sample data points belong to q different classes, and later focus on the case of q = 2.

Let the sample data points S be members of number of points belonging to class is denoted by and the total number of the training set is denote the center of class denote the center of the complete training set S:

Let be the matrix associated with class whose columns consists of the mean shifted data points:

Then, is the covariance matrix associated with class (where ”w” stands for ”within”) be the sum of the class covariance matrices:

From the discussion in the previous section, it is wish to minimize. To see why this is so, note

Let B be the matrix holding the class centers:

and let (where ”b” stands for ”between”). From the discussion above it is which we wish to maximize. Taken together, we wish to maximize the ratio (called ”Rayleigh’s quotient”):

The necessary condition for optimality is:

From which we obtain the generalized eigensystem:

That is, w is the leading eigenvector of is invertible). The general case of finding q such axes involves finding the leading generalized eigenvectors of () — the derivation is out of scope of this lecture. Note that since is not symmetric there may be no real-value solution, which is a complication will not pursue further in this course. Instead we will focus now on the 2-class (q = 2) setting below.

The general derivation is simplified when there are only two classes. The covariance matrix becomes a rank-1 matrix:

As a result, is a vector in direction . Therefore, the solution for w from eqn. 5.2 is:

The decision boundary ) = 0 becomes:

This decision boundary will surface again in the course when we consider Bayseian inference. It will be shown that this decision boundary is the Maximum Likelihood solution in the case where the two classes are normally distributed with means and with the same covariance matrix

Both LDA and SVM search for a so called ”optimal” linear discriminant function, what is the difference? The heart of the matter lies in the definition of what constitutes a sufficient compact representation of the data. In LDA the assumption is that each class can be represented by its mean vector and its spread (i.e., covariance matrix). This is true for normally distributed data — but not true in general. This means that we should expect that LDA will produce the optimal discriminant linear function when each of the classes are normally distributed.

With SVM, on the other hand, there is no assumption on how the data is distributed. Instead, the emerging result is that the data is represented by the subset of data points which lie on the boundary between the two classes (the so called support vectors). Rather than making a parametric assumption on how the data can be captured (i.e., mean and covariance) the theory shows that the data can be captured by a special subset of points. The tools, as a result, are naturally more complex (quadratic linear programming versus spectral matrix analysis) — but the advantage is that optimality is guaranteed without making assumptions on the distribution of the data (i.e., distribution free). It can be shown that SVM and LDA would produce the same result if the class data is normally distributed.

CCA is a technique for learning a mapping using the notion of subspace similarity (an extension of the inner product between two vectors) from a training set of (a mapping, where y can be any point in as opposed to a discrete set of labels, is often referred to as a ”regression” (as opposed to ”classification”).

Like in PCA and LDA, the approach would be to look for projection axes such that the projection of the input and output vectors on those axes satisfy certain requirements — and like PCA and LDA the tools we would be using is matrix spectral analysis.

It will be convenient to stack our vectors as rows of an input matrix A and output matrix matrix whose rows are matrix whose rows are . Consider vectors and and project the input and output data onto them producing . The requirement we would like to place on the projection axes is that , or in other words that (is maximal. The requirement therefore is that the projection of the input points onto the u axis is similar to the projection of the output points onto the v axis. If we extend this notion to multiple axes necessarily orthogonal) and ) our requirement becomes that the new coordinates of the input points projected onto the subspace spanned by the u vectors are similar to the new coordinates of the output points projected onto the subspace spanned by the v vectors. In other words, we wish to find two q-dimensional subspaces one of the other of such that the two sets of projected points are as aligned as possible.

CCA goes a step further and makes the assumption that the input/output relationship is solely determined by the relation (angles) between the column spaces of A, B. In other words, the particular columns of A are not really important, what is important is the space spanned by the columns. Since g = Au is a point in (a linear combination of the columns of A) and h = Bv is a point in is the cosine angle, cos(between the two axes provided that we normalize the vectors g and h. If we continue this line of reasoning recursively, we obtain a set of angles 0 2), called ”principal angles”, between the two subspaces uniquely defined as:

subject to:

As a result, we obtain the following optimization function over axes u, v:

To solve this problem we first perform a ”QR” factorization of A and B. A ”QR” factorization of a matrix A is a Grahm-Schmidt process resulting in an orthonormal set of vectors arranged as the columns of a matrix column space is equal to the column space of A, and a matrix contains the coefficients of the linear combination of the columns of that . Since orthoganilzation is not unique, the Grahm-Schmidt process perfroms the orthogonalization such that is an upper-diagonal matrix. Likewise let . Because the column spaces of are the same, then for every u there exists a ˆoptimization problem now becomes:

The solution of this problem is when ˆare the leading singular vectors of . The singular value decomposition (SVD) of any matrix E is a decomposition where the columns of U are the leading eigenvectors of , the rows of are the leading eigenvectors of D is a diagonal matrix whose entries are the corresponding square eigenvalues (note that the eigenvalues of are the same). The SVD decomposition has the property that if we keep only the first q leading eigenvectors then is the closest (in least squares sense) rank q matrix to E.

Therefore, let ˆbe the SVD of using the first q eigenvectors. Then, our sought after axes ] is simply and likewise and the axes ] is equal to . The axes are called ”canon- ical vectors”, and the vectors (mutually orthogonal) are called ”variates”. The concept of principal angles is due to Jordan in 1875, where Hotelling in 1936 is the first to introduce the recursive definition above.

Given a new vector the resulting vector y can be found by solving the linear system (since our assumption is that in the new basis the coordinates of x and y are similar).

To conclude, the relationship between A and B is captured by creating similar variates, i.e., creating subspaces of dimension q such that the projections of the input vectors and the output vectors have similar coordinates. The process for obtaining the two q-dimensional subspaces is by performing a QR factorization of A and B followed by an SVD. Here again the spectral analysis of the input and output data matrices plays a pivoting role in the input/output association.

Spectral Analysis II: Clustering

In the previous lecture we ended up with the formulation:

and showed the solution G is the leading eigenvectors of the symmetric positive semi definite matrix (sample covariance matrix) with , those eigenvectors form a basis to a k-dimensional subspace of which is the closest (in norm sense) to the sample points . The axes (called principal axes) preserve the variance of the original data in the sense that the projection of the data points on the has maximum variance, projection on has the maximum variance over all vectors orthogonal to , etc. The spectral decomposition of the sample covariance matrix is a way to ”compress” the data by means of linear super-position of the original coordinates

We also ended with a ratio formulation:

where where scatter matrices defined such that is the variance of class centers (which we wish to maximize) and is the sum of within class variance (which we want to minimize). The solution w is the generalized eigenvector with maximal

In this lecture we will show additional applications where the search for leading eigenvectors plays a pivotal part of the solution. So far we have seen how spectral analysis relates to PCA and LDA and today we will focus on the classic Data Clustering problem of partitioning a set of points 2 classes, i.e., generating as output indicator variables . We will begin with ”K-means” algorithm for clustering and then move on to show how the optimization criteria relates to grapth-theoretic approaches (like Min-Cut, Ratio-Cut, Normalized Cuts) and spectral decomposition.

The K-means formulation (originally introduced by [4]) assumes that the clusters are defined by the distance of the points to their class centers only. In other words, the goal of clustering is to find those k mean vectors and provide the cluster assignment of each point in the set. The K-means algorithm is based on an interleaving approach where the cluster assignments are established given the centers and the centers are computed given the assignments. The optimization criterion is as follows:

Assume that are given from the previous iteration, then

and next assume that (cluster assignments) are given, then for any set we have that

In other words, given the estimated centers in the current round, the new assignments are computed by the closest center to each point given the updated assignments the new centers are estimated by taking the mean of each cluster. Since each step is guaranteed to reduce the optimization energy the process must converge — to some local optimum.

The drawback of the K-means algorithm is that the quality of the local optimum strongly depends on the initial guess (either the centers or the assignments). If we start with a wild guess for the centers it would be fairly unlikely that the process would converge to a good local minimum (i.e. one that is close to the global optimum). An alternative approach would be to define an approximate but simpler problem which has a closed form solution (such as obtained by computing eigenvectors of some matrix). The global optimum of the K-means is an NP-Complete problem (mentioned briefly in the next section).

Next, we will rewrite the K-means optimization criterion in matrix form and see that it relates to the spectral formulation (eqn. 6.1).

We rewrite eqn. 6.2 as follows [7]. Instead of carrying the class variables we define class sets and . The K-means optimization criterion seeks for the centers and the class sets:

Let and following the expansion of the squared norm and dropping we end up with an equivalent problem:

Next we substitute with its definition: (1and obtain a new equivalent formulation where the centers are eliminated form consideration:

which is more conveniently written as a maximization problem:

Since the resulting formulation involves only inner-products we could have replaced ) in eqn. 6.2 where the mapping ) is chosen such that ) can be replaced by some non-linear function known as the ”kernel trick” (discussed in previous lectures). Having the ability to map the input vectors onto some high-dimensional space before K-means is applied provides more flexibility and increases our chances of getting out a ”good” clustering from the global K-means solution (again, the local optimum depends on the initial conditions so it could be ”bad”). The RBF kernel is quite popular in this context with some pre-determined parameter. Note that can be interpreted loosely as the probability of to be clustered together.

Let symmetric positive-semi-definite matrix often referred to as the ”affinity” matrix. Let whose entries are for some class otherwise. In other words, if we sort the points according to cluster membership, then F is a block diagonal matrix with blocks block of 1’s scaled by 1. Then, Eqn. 6.3 can be written in terms of K as follows:

In order to form this as an optimization problem we need to represent the structure of F in terms of constraints. Let column-scaled indicator matrix: belongs to the j’th class) and = 0 otherwise. Let be the columns of G and it can be easily verified that 0) therefore Since trace(AB) = trace(BA) we can now write eqn. 6.4 in terms of G:

under conditions on G which we need to further spell out.

We will start with the necessary conditions. Clearly 0 (has non-negative entries). Because each point belongs to exactly one cluster we must have Furthermore we have that the rows and columns of 1, i.e., which means that F is doubly stochastic which translates to the constraint We have therefore three necessary conditions on , and (iii) The claim below asserts that these are also sufficient conditions:

Claim 4 The feasibility set of matrices G which satisfy the three conditions are of the form:

= 0 we have that has a single non-vanishing element in each row. It will be convenient to assume that the points are sorted according to the class membership, thus the columns of G have the non-vanishing entries in consecutive order and let be the number of non-vanishing entries in column vector of entries holding only the non-vanishing entries of the doubly stochastic constraint results that (j = 1, ..., k. Multiplying 1 from both sides yields (therefore

This completes the equivalence between the matrix formulation:

and the original K-means formulation of eqn. 6.2.

We have obtained the same optimization criteria as eqn. 6.1 with additional two constraints: G should be non-negative and should be doubly stochastic. The constraint comes from the requirement that each point is assigned to one class only. The doubly stochastic constraint comes from a ”class balancing” requirement which we will expand on below.

We will arrive to eqn. 6.5 from a graph-theoretic perspective. We start with representing the graph Min-Cut problem in matrix form, as follows. A convenient way to represent the data to be clustered is by an undirected graph with edge-weights where V = {1, ..., m} is the vertex set, is the edge set and is the positive weight function. Vertices of the graph correspond to data points , edges represent neighborhood relationships, and edge-weights represent the similarity (affinity) between pairs of linked vertices. The weight adjacency matrix K holds the weights where = 0 otherwise.

A cut in the graph is defined between two disjoint sets B = V , is the sum of edge-weights connecting the two sets: cut(A, B) = which is a measure of dissimilarity between the two sets. The Min-Cut problem is to find a minimal weight cut in the graph (can be solved in polynomial time through Max Network Flow solution). The following claim associates algebraic conditions on G with an indicator matrix:

Claim 5 The feasibility set of matrices G which satisfy the three conditions for some diagonal matrix D are of the form:

= 0 we have that = 0, i.e., G has a single non-vanishing element in each row. From G1 = 1 the single non-vanishing entry of each row must have the value of 1.

In the case of two classes (k = 2), the function ) is equal to . Therefore max) is equivalent to minimizing the cut: . As a result, the Min-Cut problem is equivalent to solving the optimization problem:

We seem to be close to eqn. 6.5 with the difference that G is orthogonal (instead of orthonormal) and the doubly-stochasitc constraint is replaced by G1 = 1. The difference can be bridged by considering a ”balancing” requirement. Min-Cut can produce an unbalanced partition where one set of vertices is very large and the other contains a spurious set of vertices having a small number of edges to the larger set. This is an undesirable outcome in the context of clustering. Consider a ”balancing” constraint which makes a strict requirement that all the k clusters have an equal number of points. We can relax the balancing constraint slightly by combining the balancing constraint with G1 = 1 into one single constraint is scaled doubly stochastic. Note that the two conditions propose the relaxed-balanced hard clustering scheme:

The scale m/k is a global scale that can be dropped without affecting the resulting solution, thus the Min-Cut with a relaxed balancing requirement becomes eqn. 6.5 which we saw is equivalent to K-means:

We saw above that the doubly-stochastic constraint has to do with a ”bal- ancing” desire. A further relaxation of the balancing desire is to perform the optimization in two steps: (i) replace the affinity matrix K with the closest (under some chosen error measure) doubly-stochastic matrix , (ii) find a solution to the problem:

because should come out close to and is doubly-stochastic, then should come out close to satisfying a doubly-stochastic constraint — this is the motivation behind the 2-step approach. Moreover, we drop the non-negativity constraint that the non-negativity constraint is crucial for the physical interpretation of G; nevertheless, for k = 2 clusters it is possible to make an interpretation, as we shall next. As a result we are left with a spectral decomposition problem of eqn. 6.1:

where the columns of G are the leading eigenvectors of . We will refer to the first step as a ”normalization” process and there are two popular normalizations in the literature — one leading to Ratio-Cuts and the other to Normalized-Cuts.

Let D = diag(K1) which is a diagonal matrix containing the row sums of K. The Ratio-Cuts normalization is to look for as the closest doubly-stochastic matrix to K by minimizing the norm — this turns out to be

Claim 6 (ratio-cut) Let K be a symmetric positive-semi-definite whose values are in the range [0, 1]. The closest doubly stochastic matrix the error norm is

for any matrix F, we must have:

Let

The Laplacian matrix of a graph is is an eigenvector of the Laplacian with eigenvalue is also an eigenvector of with eigenvalue 1 and since (= 0 then the smallest eigenvector v = 1 of the Laplacian is the largest of , and the second smallest eigenvector of the Laplacian (the ratio-cut result) corresponds to the second largest eigenvector of . Because the eigenvectors are orthogonal, the second eigenvector must have positive and negative entries (because the inner-product with 1 is zero) — thus the sign of the entries of the second eigenvector determines the class membership.

Ratio-Cuts, the second smallest eigenvector of the Laplacian an approximation due to Hall in the 70s [2] to the Min-Cut formulation. Let determine the class membership such that clustered together if have similar values. This leads to the following optimization problem:

The criterion function is equal to (1and the derivative of the Lagrangian (11) with respect to z gives rise to the necessary condition (and the Ratio-Cut scheme follows.

Normalized-Cuts looks for the closest doubly-stochastic matrix entropy error measure defined as:

We will encounter the relative entropy measure in more detail later in the course. We can show that must have the form ΛKΛ for some diagonal matrix Λ:

Claim 7 The closest doubly-stochastic matrix F under the relative-entropy error measure to a given non-negative symmetric matrix K, i.e., which minimizes:

has the form for some (unique) diagonal matrix Λ.

Proof: The Lagrangian of the problem is:

The derivative with respect to

from which we obtain:

Let D), then we have:

Since is symmetric we must have

Next, we can show that the diagonal matrix Λ can found by an iterative process where K is replaced by was defined above as diag(K1):

Claim 8 For any non-negative symmetric matrix , iterating the process converges to a doubly stochastic matrix.

The proof is based on showing that the permanent increases monotonically, i.e. ). Because the permanent is bounded the process must converge and if the permanent does not change (at the convergence point) the resulting matrix must be doubly stochastic. The resulting doubly stochastic matrix is the closest to K in relative-entropy.

Normalized-Cuts takes the result of the first iteration by replacing K with followed by the spectral decomposition (in case of k = 2 classes the partitioning information is found in the second leading eigenvector of — just like Ratio-Cuts but with a different in this manner is not the closest doubly-stochastic matrix to K but is fairly close (the first iteration is the dominant one in the process).

Normalized-Cuts, as the second leading eigenvector of is an approximation to a ”balanced” Min-Cut described first in [6]. Deriving it from first principles proceeds as follows:

Let be defined for any two subsets (not necessarily disjoint) of vertices. The normalized-cuts measures the cut cost as a fraction of the total edge connections to all the nodes in the graph:

A minimal Ncut partition will no longer favor small isolated points since the cut value would most likely be a large percentage of the total connections from that small set to all the other vertices. A related measure Nassoc(A, B) defined as:

reflects how tightly on average nodes within the group are connected to each other. Given that ) one can easily verify

that:

therefore the optimal bi-partition can be represented as maximizing A). The Nassoc naturally extends to k > 2 classes (partitions) as follows: Let be disjoint sets

We will now rewrite Nassoc in matrix form and establish equivalence to eqn. 6.7. Let ¯with the 1s indicating membership to the j’th class. Note that

therefore Note also that (1= 1, therefore ¯have that Taken together we have that maximizing Nassoc is equivalent to:

where . Note that this is exactly the K-means matrix setup of eqn. 6.5 where the doubly-stochastic constraint is relaxed into the replacement of . The constraint 0 is then dropped and the resulting solution for G is the k leading eigenvectors of

We have arrived via seemingly different paths to eqn. 6.8 which after we drop the constraint 0 we end up with a closed form solution consisting of the k leading eigenvectors of = 2 (two classes) one can easily verify that the partitioning information is fully contained in the second eigenvector. Let be the first leading eigenvectors of is an eigenvector with eigenvalue

In fact = 1 is the largest eigenvalue (left as an exercise) thus 0. Since is symmetric the contains positive and negative entries — those are interpreted as indicating class membership (positive to one class and negative to the other).

The case k > 2 is treated as an embedding (also known as Multi-Dimensional Scaling) by re-coordinating the points using the rows of G. In other

words, the i’th row of G is a representation of ditions where K is block diagonal (the distance between clusters is infinity) the rows associated with points clustered together are identical (i.e., the n original points are mapped to ) [5]. In practice, one performs the iterative K-means in the embedded space.

The Formal (PAC) Learning Model

We have see so far algorithms that explicitly estimate the underlying distribution of the data (Bayesian methods and EM) and algorithms that are in some sense optimal when the underlying distribution is Gaussian (PCA, LDA). We have also encountered an algorithm (SVM) that made no assumptions on the underlying distribution and instead tied the accuracy to the margin of the training data.

In this lecture and in the remainder of the course we will address the issue of ”accuracy” and ”generalization” in a more formal manner. Because the learner receives only a finite training sample, the learning function can do very well on the training set yet perform badly on new input instances. What we would like to establish are certain guarantees on the accuracy of the learner measured over all the instance space and not only on the training set. We will then use those guarantees to better understand what the largemargin principle of SVM is doing in the context of generalization.

In the remainder of this lecture we will refer to the following notations: the class of learning functions is denoted by C. A learning functions is often referred to as a ”concept” or ”hypothesis”. A target function function that has zero error on all input instances (such a function may not always exist).

In many learning situations of interest, we would like to assume that the learner receives m examples sampled by some fixed (yet unknown) distribution D and the learner must do its best with the training set in order to achieve the accuracy and confidence objectives. The Probably Approximate Correct (PAC) model, also known as the ”formal model”, first introduced by Valient in 1984, provides a probabilistic setting which formalizes the notions of accuracy and confidence.

The PAC model makes the following statistical assumption. We assume the learner receives a set which are sampled randomly and independently according to a distribution D over X. In other words, a random training set S of length m is distributed according to the product probability distribution . The distribution D is unknown, but we will see that one can obtain useful results by simply assuming that D is fixed — there is no need to attempt to recover D during the learning process. To recap, we make the following three assumptions: (i) D is unkown, (ii) D is fixed throughout the learning process, and (iii) the example instances are sampled independently of each other (are Identically and Independently Distributed — i.i.d.).

We distinguish between the ”realizable” case where a target concept is known to exist, and the unrealizable case, where there is no such guarantee. In the realizable case our training examples are i = 1, ..., m and D is defined over are given by the unrealizable case, is the distribution over (each element is a pair, one from X and the other from Y ).

We next define what is meant by the error induced by a concept function h(x). In the realizable case, given a function , the error of h is defined with respect to the distribution D:

where ind(F) is an indication function which returns ’1’ if the proposition F is true and ’0’ otherwise. The function err(h) is the probability that an instance x sampled according to D will be labeled incorrectly by h(x). Let 0 be a parameter given to the learner specifying the ”accuracy” of the learning process, i.e. we would like to achieve . Note that

In addition, we define a ”confidence” parameter 0, also given to the learner, which defines the probability that

or equivalently:

In other words, the learner is supposed to meet some accuracy criteria but is allowed to deviate from it by some small probability. Finally, the learning algorithm is supposed to be ”efficient” if the running time is polynomial in 1and the size of the concept target function () (measured by the number of bits necessary for describing it, for example).

We will say that an algorithm L learns a concept family C in the formal sense (PAC learnable) if for any and for every distribution D on the instance space X, the algorithm L generates efficiently a concept function such that the probability that is at least 1

The inclusion of the confidence value could seem at first unnatural. What we desire from the learner is to demonstrate a consistent performance regardless of the training sample Z. In other words, it is not enough that the learner produces a hypothesis h whose accuracy is above threshold, i.e., , for some training sample Z. We would like the accuracy performance to hold under all training samples (sampled from the distribution ) — since this requirement could be too difficult to satisfy, the formal model allows for some ”failures”, i.e, situations where training samples Z, as long as those failures are rare and the frequency of their occurrence is controlled (the parameter ) and can be as small as we like.

In the unrealizable case, there may be no function err(h) = 0, thus we need to define what we mean by the best a learning algorithm can achieve:

which is the best that can be done on the concept class C using functions that map between X and Y . Given the desired accuracy and confidence values the learner seeks a hypothesis such that:

We are ready now to formalize the discussion above and introduce the defi-nition of the formal learning model (Anthony & Bartlett [1], pp. 16):

Definition 1 (Formal Model) Let C be the concept class of functions that map from a set X to Y . A learning algorithm L is a function:

from the set of all training examples to C with the following property: given any there is an integer such that if any probability distribution is a training set of length m drawn randomly according to the product probability distribution with probability of at least 1 the hypothesis L is such that . We say that C is learnable (or PAC learnable) if there is a learning algorithm for C.

There are few points to emphasize. The sample size ) is a suf-ficient sample size for PAC learning C by L and is allowed to vary with . Decreasing the value of either makes the learning problem more difficult and in turn a larger sample size is required. Note however that on the distribution D! that is, a sufficient sample size can be given that will work for any distribution D — provided that D is fixed throughout the learning experience (both training and later for testing). This point is a crucial property of the formal model because if the sufficient sample size is allowed to vary with the distribution D then not only we would need to have some information about the distribution in order to set the sample complexity bounds, but also an adversary (supplying the training set) could control the rate of convergence of L to a solution (even if that solution can be proven to be optimal) and make it arbitrarily slow by suitable choice of D.

What makes the formal model work in a distribution-invariant manner is that it critically depends on the fact that in many interesting learning scenarios the concept class C is not too complex. For example, we will show later in the lecture that any finite concept class is learnable, and the sample complexity (in the realizable case) is

In the next lecture we will consider concept classes of infinite size and show that despite the fact that the class is infinite it still can be of low complexity!

Before we illustrate the concepts above with an example, there is another useful measure which is the empirical error (also known as the sample error) ˆerr(h) which is defined as the proportion of examples from Z on which h made a mistake:

(replace for the unrealizable case). The situation of bounding the true error err(h) by minimizing the sample error ˆerr(h) is very convenient — we will get to that later.

As an illustration of learnability we will consider the problem (introduced in Kearns & Vazirani [3]) of learning an axes-aligned rectangle from positive and negative examples. We will show that the problem is PAC-learnable and find out

In the rectangle learning game we are given a training set consisting of points in the 2D plane with a positive ’+’ or negative ’-’ label. The positive examples are sampled inside the target rectangle (parallel to the main axes) R and the negative examples are sampled outside of R. Given m examples sampled i.i.d according to some distribution D the learner is supposed to generate an approximate rectangle which is consistent with the training set (we are assuming that R exists) and which satisfies the accuracy and confidence constraints.

We first need to decide on a learning strategy. Since the solution is not uniquely defined given any training set Z, we need to add further constraints to guarantee a unique solution. We will choose as the axes-aligned concept which gives the tightest fit to the positive examples, i.e., the smallest area axes-aligned rectangle which contains the positive examples. If no positive examples are given then . We can also assume that Z contains at least three non-collinear positive examples in order to avoid complications associated with infinitesimal area rectangles. Note that we could have chosen other strategies, such as the middle ground between the tightest fit to the positive examples and the tightest fit (from below) to the negative examples, and so forth. Defining a strategy is necessary for the analysis below — the type of strategy is not critical though.

We next define the error ) on the concept generated by our learning strategy. We first note that with the strategy defined above we always have is the tightest fit solution which is consistent with the sample data (there could be a positive example outside of is not in the training set). We will define the ”weight” w(E) of a region E in the plane as

i.e., the probability that a random point sampled according to the distribution D will fall into the region. Therefore, the error associated with the concept

Fig. 7.1. Given the tightest-fit to positive examples strategy we have that The strip has weight 4 and the strip is defined as the upper strip covering the area between

and we wish to bound the error with probability of at least 1 after seeing m examples.

We will divide the region into four strips (see Fig.7.1) which overlap at the corners. We will estimate ) noting that the overlaps between the regions makes our estimates more pessimistic than they truly are (since we are counting the overlapping regions twice) thus making us lean towards the conservative side in our estimations.

Consider the upper strip ) then we are done. We are however interested in quantifying the probability that this is not the case. Assume and define a strip which starts from the upper axis of R and stretches to the extent such that have that . Furthermore:

the the label must be positive since the label is positive then given our learning strategy of fitting the tightest rectangle over the positive examples, then that

We have therefore that ) iff no point in appears in the sample (otherwise intersects with probability that a point sampled according to the distribution D will fall outside of . Given the independence assumption (examples are

drawn i.i.d.), we have:

Repeating the same analysis to regions and using the union bound ) we come to the conclusion that the probability that any of the four strips of has weight greater that 4 is at most 4(1 . In other words,

We can make the expression more convenient for manipulation by using the inequality (recall that 1 + (1from which it follows that (1 + and by taking the power of 0 we obtain (1 +

from which we obtain the bound:

To conclude, assuming that the learner adopts the tightest-fit to positive examples strategy and is given at least training examples in order to find the axes-aligned rectangle , we can assert that with probability 1 the error associated with (i.e., the probability that an (point will be classified incorrectly) is at most

We can see form the analysis above that indeed it applies to any distribution D where the only assumption we had to make is the independence of the draw. Also, the sample size m behaves well in the sense that if one desires a higher level of accuracy (smaller ) or a higher level of confidence (smaller ) then the sample size grows accordingly. The growth of m is linear in 1and linear in ln(1

In the previous section we illustrated the concept of learnability with a particular simple example. We will now focus on applying the learnability model to a more general family of learning examples. We will consider the family of all learning problems over finite concept classes example, the conjunction learning problem (over boolean formulas) with n literals contains only 3hypotheses because each variable can appear in the conjunction or not and if appears it could be negated or not. We have shown that n is the lower bound on the number of mistakes on the worst case analysis any on-line algorithm can achieve. With the definitions we have above on the formal model of learnability we can perform accuracy and sample complexity analysis that will apply to any learning problem over finite concept classes. This was first introduced by Valiant in 1984.

In the realizable case over , we will show that any algorithm L which returns a hypothesis which is consistent with the training set Z is a learning algorithm for C. In other words, any finite concept class is learnable and the learning algorithms simply need to generate consistent hypotheses. The sample complexity associated with the choice of can be shown as equal to:

the empirical error (the error obtained on Z) is a learning algorithm for C. The sample complexity can be shown as equal to:

Let be some consistent hypothesis with the training set Z (we know that such a hypothesis exists, in particular the target concept used for generating Z) and suppose that

Then, the probability (with respect to the product distribution agrees with on a random sample of length m is at most (1 the inequality we saw before

We wish to bound the error uniformly, i.e., that for all concepts . This requires the evaluation of:

There at most |C| such functions h, therefore using the Union-Bound the probability that some function in C has error larger than and is consistent

with on a random sample of length m is at most

For any positive , this probability is less than

This derivation can be summarized in the following theorem (Anthony & Bartlett [1], pp. 25):

Theorem 5 Let C be a finite set of functions from X to Y . Let L be an algorithm such that for any m and for any is a training sample , then the hypothesis h = L(Z) satisfies is a learning algorithm for C in the realizable case with sample complexity

In the realizable case an algorithm simply needs to generate a consistent hypothesize to be considered a learning algorithm in the formal sense. In the unrealizable situation (a target function might not exist) an algorithm which minimizes the empirical error, i.e., an algorithm L generates h = L(Z) having minimal sample error:

is a learning algorithm for C (assuming finite |C|). This is a particularly useful property given that the true errors of the functions in C are unknown. It seems natural to use the sample errors ˆerr(h) as estimates to the performance of L.

The fact that given a large enough sample (training set Z) then the sample error ˆerr(h) becomes close to the true error err(h) is somewhat of a restatement of the ”law of large numbers” of probability theory. For example, if we toss a coin many times then the relative frequency of ’heads’ approaches the true probability of ’head’ at a rate determined by the law of large numbers. We can bound the probability that the difference between the empirical error and the true error of some using Hoeffding’s inequality:

Claim 10 Let h be some function from X to Y = {0, 1}. Then

for any probability distribution and any positive integer m.

Proof: This is a straightforward application of Hoeffding’s inequality to Bernoulli variables. Hoeffding’s inequality says: Let X be a set, D a probability distribution on real-valued functions from X to an interval on the real line (

where

In our case (1

The Hoeffding bound almost does what we need, but not quite so. What we have is that for any given hypothesis , the empirical error is close to the true error with high probability. Recall that our goal is to minimize err(h) over all possible but we can access only ˆerr(h). If we can guarantee that the two are close to each other minimizing ˆwill approximately minimize err(h). Put formally, in order to ensure that L learns the class C, we must show that

In other words, we need to show that the empirical errors converge (at high probability) to the true errors . If that can be guaranteed, then with (high) probability 1 , for every

So, since the algorithm L running on training set Z returns h = L(Z) which minimizes the empirical error, we have:

which is what is needed in order that L learns C. Thus, what is left is to prove the following claim:

Claim 11

Proof: We will use the union bound. Finding the maximum over C is equivalent to taking the union of all the events:

using the union-bound and Claim 2, we have:

Finally, given that 2we obtain the sample complexity:

This discussion is summarized with the following theorem (Anthony & Bartlett [1], pp. 21):

Theorem 6 Let C be a finite set of functions from X to Y = {0, 1}. Let L be an algorithm such that for any m and for any training set i = 1, ..., m, then the hypothesis L(Z) satisfies:

Then L is a learning algorithm for C with sample complexity

Note that the main difference with the realizable case (Theorem 1) is the larger 1rather than 1. The realizable case requires a smaller training set since we are estimating a random quantity so the smaller the variance the less data we need.

The VC Dimension

The result of the PAC model (also known as the ”formal” learning model) is that if the concept class C is PAC-learnable then the learning strategy must simply consist of gathering a sufficiently large training sample S of size ), for given accuracy 0 and confidence 0 parameters, and finds a hypothesis which is consistent with S. The learning algorithm is then guaranteed to have a bounded error with probability 1 . The error measurement includes data not seen by the training phase.

This state of affair also holds (with some slight modifications on the sample complexity bounds) when there is no consistent hypothesis (the unrealizable case). In this case the learner simply needs to minimize the empirical error ˆerr(h) on the sample training data S, and if m is sufficiently large then the learner is guaranteed to have with probability 1 . The measure Opt(C) is defined as the minimal Note that in the realizable case Opt(C) = 0.

The property of bounding the true error err(h) by minimizing the sample error ˆerr(h) is very convenient. The fundamental question is under what conditions this type of generalization property applies? We saw in the previous lecture that a satisfactorily answer can be provided when the cardinality of the concept space is bounded, i.e. , which happens for Boolean concept space for example. In that lecture we have proven that:

is sufficient for guaranteeing a learning model in the formal sense, i.e., which has the generalization property described above.

In this lecture and the one that follows we have two goals in mind. First is to generalize the result of finite concept class cardinality to infinite cardinality — note that the bound above is not meaningful when Can we learn in the formal sense any non-trivial infinite concept class? (we already saw an example of a PAC-learnable infinite concept class which is the class of axes aligned rectangles). In order to answer this question we will need to a general measure of concept class complexity which will replace the cardinality term |C| in the sample complexity bound ). It is tempting to assume that the number of parameters which fully describe the concepts of C can serve as such a measure, but we will show that in fact one needs a more powerful measure called the Vapnik-Chervonenkis (VC) dimension. Our second goal is to pave the way and provide the theoretical foundation for the large margin principle algorithm (SVM) we derived in Lecture 4.

The basic principle behind the VC dimension measure is that although C may have infinite cardinality, the restriction of the application of concepts in C to a finite sample S has a finite outcome. This outcome is typically governed by an exponential growth with the size m of the sample S — but not always. The point at which the growth stops being exponential is when the ”complexity” of the concept class C has exhausted itself, in a manner of speaking.

We will assume C is a concept class over the instance space X — both of which can be infinite. We also assume that the concept class maps instances in X to {0, 1}, i.e., the input instances are mapped to ”positive” or ”negative” labels. A training sample S is drawn i.i.d according to some fixed but unknown distribution D and S consists of In our notations we will try to reserve to denote the target concept and concept. We begin with the following definition:

Definition 2

which is a set of vectors in

Π) is set whose members are m-dimensional Boolean vectors induced by functions of C. These members are often called dichotomies or behaviors on S induced or realized by C. If C makes a full realization then Πhave 2members. An equivalent description is a collection of subsets of S:

where each makes a partition of S into two sets — the positive and negative points. The set Π) contains therefore subsets of S (the positive points of S under h). A full realization will provide will use both descriptions of Π) as a collection of subsets of S and as a set of vectors interchangeably.

is considered shattered by C. In other words, S is shattered by C if C realizes all possible dichotomies of S.

Consider as an example a finite concept class applied to three instance vectors with the results:

Then,

With these definitions we are ready to describe the measure of concept class complexity.

Definition 4 (VC dimension) The VC dimension of C, noted as V Cdim(C), is the cardinality d of the largest set S shattered by C. If all sets S (arbitrarily large) can be shattered by

The VC dimension of a class of functions C is the point d at which all samples S with cardinality |S| > d are no longer shattered by C. As long as C shatters S it manifests its full ”richness” in the sense that one can obtain from S all possible results (dichotomies). Once that ceases to hold, i.e., when |S| > d, it means that C has ”exhausted” its richness (complexity). An infinite VC dimension means that C maintains full richness for all sample sizes. Therefore, the VC dimension is a combinatorial measure of a function class complexity.

Before we consider a number of examples of geometric concept classes and their VC dimension, it is important clarify the lower and upper bounds (existential and universal quantifiers) in the definition of VC dimension. The VC dimension is at least d if there exists some sample |S| = d which is shattered by C — this does not mean that all samples of size d are shattered by C. Conversely, in order to show that the VC dimension is at most d, one must show that no sample of size d + 1 is shattered. Naturally, proving an upper bound is more difficult than proving the lower bound on the VC dimension. The following examples are shown in a ”hand waiving” style and are not meant to form rigorous proofs of the stated bounds — they are shown for illustrative purposes only.

Intervals of the real line: The concept class C is governed by two parameters in the closed interval [0, 1]. A concept from this class will tag an input instance 0 < x < 1 as positive if and negative otherwise. The VC dimension is at least 2: select a sample of 2 points positioned in the open interval (0, 1). We need to show that there are values of which realize all the possible four dichotomies (+This is clearly possible as one can place the interval [] such the intersection with the interval [] is null, (thus producing ()), or to fully include [] (thus producing (+, +)) or to partially intersect [that are excluded (thus producing the remaining two dichotomies). To show that the VC dimension is at most 2, we need to show that any sample of three points on the line (0, 1) cannot be shattered. It is sufficient to show that one of the dichotomies is not realizable: the labeling (++) cannot be realizable by any interval [] — this is because if are labeled positive then by definition the interval [] must fully include the interval [] and since must be labeled positive as well. Thus V Cdim(C) = 2.

Axes-aligned rectangles in the plane: We have seen this concept class in the previous lecture — a point in the plane is labeled positive if it lies in an axes-aligned rectangle. The concept class C is thus governed by 4 parameters. The VC dimension is at least 4: consider a configuration of 4 input points arranged in a cross pattern (recall that we need only to show some sample S that can be shattered). We can place the rectangles (concepts of the class C) such that all 16 dichotomies can be realized (for example, placing the rectangle to include the vertical pair of points and exclude the horizontal pair of points would induce the labeling (+It is important to note that in this case, not all configurations of 4 points can be shattered — but to prove a lower bound it is sufficient to show the existence of a single shattered set of 4 points. To show that the VC dimension is at most 4, we need to prove that any set of 5 points cannot be shattered. For any set of 5 points there must be some point that is ”internal”, i.e., is neither the extreme left, right, top or bottom point of the five. If we label this internal point as negative and the remaining 4 points as positive then there is no axes-aligned rectangle (concept) which cold realize this labeling (because if the external 4 points are labeled positive then they must be fully within the concept rectangle, but then the internal point must also be included in the rectangle and thus labeled positive as well).

Separating hyperplanes: Consider first linear half spaces in the plane. The lower bound on the VC dimension is 3 since any three (non-collinear) points in can be shattered, i.e., all 8 possible labelings of the three points can be realized by placing a separating line appropriately. By having one of the points on one side of the line and the other two on the other side we can realize 3 dichotomies and by placing the line such that all three points are on the same side will realize the 4th. The remaining 4 dichotomies are realized by a sign flip of the four previous cases. To show that the upper bound is also 3, we need to show that no set of 4 points can be shattered. We consider two cases: (i) the four points form a convex region, i.e., lie on the convex hull defined by the 4 points, (ii) three of the 4 points define the convex hull and the 4th point is internal. In the first case, the labeling which is positive for one diagonal pair and negative to the other pair cannot be realized by a separating line. In the second case, a labeling which is positive for the three hull points and negative for the interior point cannot be realize. Thus, the VC dimension is 3 and in general the VC dimension for separating hyperplanes in

Union of a finite number of intervals on the line: This is an example of a concept class with an infinite VC dimension. For any sample of points on the line, one can place a sufficient number of intervals to realize any labeling.

The examples so far were simple enough that one might get the wrong impression that there is a correlation between the number of parameters required to describe concepts of the class and the VC dimension. As a counter example, consider the two parameter concept class:

which has an infinite VC dimension as one can show that for every set of m points on the line one can realize all possible labelings by choosing

a sufficiently large value of (which serves as the frequency of the sync function) and appropriate phase. We conclude this section with the following claim:

Theorem 7 The VC dimension of a finite concept class is bounded from above:

Proof: if V Cdim(C) = d then there exists at least 2functions in C because every function induces a labeling and there are at least 2Thus, from follows that

We saw that the VC dimension is a combinatorial measure of concept class complexity and we would like to have it replace the cardinality term in the sample complexity bound. The first result of interest is to show that if the VC dimension of the concept class is infinite then the class is not PAC learnable.

Theorem 8 Concept class is not learnable in the formal sense.

Proof: Assume the contrary that C is PAC learnable. Let L be the learning algorithm and m be the number of training examples required to learn the concept class with accuracy 9. That is, after seeing at least ) training examples, the learner generates a concept h which satisfies

Since the VC dimension is infinite there exist a sample set S with 2m instances which is shattered by C. Since the formal model (PAC) applies to any training sample we will use the set S as follows. We will define a probability distribution on the instance space X which is uniform on S (with probability ) and zero everywhere else.

Because S is shattered, then any target concept is possible so we will choose our target concept c in the following manner:

in other words, the labels ) are determined by a coin flip. The learner L selects an i.i.d. sample of m instances ¯S — which due to the structure of D means that the ¯and outputs a consistent hypothesis probability of error for each

The reason for that is because S is shattered by C, i.e., we can select any target concept for any labeling of S (the 2m examples) therefore we could select the labels of the m points not seen by the learner arbitrarily (by flipping a coin). Regardless of h, the probability of mistake is 0.5. The expectation on the error of h is:

This is because we have 2m points to sample (according to D as all other points have zero probability) from which the error on half of them is zero (as h is consistent on the training set ¯S) and the error on the remaining half is 0.5. Thus, the average error is 0.25. Note that E[err(h)] = 0.25 for any choice of as it is based on the sample size m. For any sample size m we can follow the construction above and generate the learning problem such that if the learner produces a consistent hypothesis the expectation of the error will be 0.25.

The result that E[err(h)] = 0.25 is not possible for the accuracy and confidence values we have set: with probability of at least 0.9 we have that 1 and with probability 0Taking the worst case of = 1 we come up with the average error:

We have therefore arrived to a contradiction that C is PAC learnable.

We next obtain a bound on the growth of when the sample size |S| = m is much larger than the VC dimension V Cdim(C) = d of the concept class. We will need few more definitions:

Definition 5 (Growth function)

The measure Π) is the maximum number of dichotomies induced by C for samples of size m. As long as . The question is what happens to the growth pattern of Πsee that the growth becomes polynomial — a fact which is crucial for the learnability of C.

Definition 6 For any natural numbers m, d we have the following definition:

By induction on m, d it is possible to prove the following:

Theorem 9

Proof: by induction on m, d. For details see [[3], pp. 56].

polynomial upper bound as follows.

Dividing both sides by

We need one more result before we are ready to present the main result of this lecture:

Theorem 10 (Sauer’s lemma) If V Cdim(C) = d, then for any m, Π

Proof: By induction on both d, m. For details see [[3], pp. 55–56].

Taken together, we have now a fairly interesting characterization on how the combinatorial measure of complexity of the concept class C scales up with the sample size m. When the VC dimension of C is infinite the growth is exponential, i.e., Πfor all values of m. On the other hand, when the concept class has a bounded VC dimension then the growth pattern undergoes a discontinuity from an exponential to a polynomial growth:

As a direct result of this observation, when m >> d is much larger than d the entropy becomes much smaller than m. Recall than from an information theoretic perspective, the entropy of a random variable Z with discrete values with probabilities is defined as:

where is a measure of ”information”, i.e., is large when is small (meaning that there is much information in the occurrence of an unlikely event) and vanishes when the event is certain entropy is therefore the expectation of information. Entropy is maximal for a uniform distribution . The entropy in information theory context can be viewed as the number of bits required for coding coding theory it can be shown that the entropy of a distribution provides the lower bound on the average length of any possible encoding of a uniquely decodable code fro which one symbol goes into one symbol. When the distribution is uniform we will need the maximal number of bits, i.e., one cannot compress the data. In the case of concept class C with VC dimension d, we see that one when all possible dichotomies are realized and thus one will need m bits (as there are 2dichotomies) for representing all the outcomes of the sample. However, when m >> d only a small fraction of the 2dichotomies can be realized, therefore the distribution of outcomes is highly non-uniform and thus one would need much less bits for coding the outcomes of the sample. The technical results which follow are therefore a formal way of expressing in a rigorous manner this simple truth — If it is possible to compress, then it is possible to learn. The crucial point is that learnability is a direct consequence of the ”phase transition” (from exponential to polynomial) in the growth of the number of dichotomies realized by the concept class.

In the next lecture we will continue to prove the ”double sampling” theorem which derives the sample size complexity as a function of the VC dimension.

The Double-Sampling Theorem

In this lecture will use the measure of VC dimension, which is a combinatorial measure of concept class complexity, to bound the sample size complexity.

9.1 A Polynomial Bound on the Sample Size m for PAC Learning

In this section we will follow the material presented in Kearns & Vazirani [3] pp. 57–61 and prove the following:

Theorem 11 (Double Sampling) Let C be any concept class of VC dimension d. Let L be any algorithm that when given a set S of m labeled examples , sampled i.i.d according to some fixed but unknown distribution D over the instance space X, of some concept as output a concept that is consistent with S. Then L is a learning algorithm in the formal sense provided that the sample size obeys:

for some constant

The idea behind the proof is to build an ”approximate” concept space which includes concepts arranged such that the distance between the approximate concepts h and the target concept c is at least — where distance is de-fined as the weight of the region in X which is in conflict with the target concept. To formalize this story we will need few more definitions. Unless specified otherwise, denotes the target concept and some concept.

Definition 7

is the region in instance space where both concepts do not agree — the error region. The probability that is equal to (by definition) err(h).

Definition 8

∆(c) is a set of error regions, one per concept over all concepts. The error regions are with respect to the target concept. The set ∆is the set of all error regions whose weight exceeds . Recall that weight is defined as the probability that a point sampled according to D will hit the region.

It will be important for later to evaluate the VC dimension of ∆(c). Unlike C, we are not looking for the VC dimension of a class of function but the VC dimension of a set of regions in space. Recall the definition of Πfrom the previous lecture: there were two equivalent definitions one based on a set of vectors each representing a labeling of the instances of S induced by some concept. The second, yet equivalent, definition is based on a set of subsets of S each induced by some concept (where the concept divides the sample points of S into positive and negative labeled points). So far it was convenient to work with the first definition, but for evaluating the VC dimension of ∆(c) it will be useful to consider the second definition:

that is, the collection of subsets of S induced by intersections with regions of ∆(c). An intersection between S and a region r is defined as the subset of points from S that fall into r. We can easily show that the VC dimensions of ) are equal:

Lemma 1

Proof: we have that the elements of Π) are susbsets of S, thus we need to show that for every S the cardinality of both sets is equal . To do that it is sufficient to show that for every element ) there is a unique corresponding element in Πbe the subset of S induced by the target concept c. The set s (a subset of S) is realized by some concept h (those points in S which were labeled positive by h). Therefore, the set ) is the subset of S containing the points that hit the region which is an element of Πthis is a one-to-one mapping we have that

, a sample set if every region in ∆is hit by at least one point of S:

In other words, if S hits all the error regions in ∆(c) whose weight exceeds -net. Consider as an example the concept class of intervals on the line [0, 1]. A concept is defined by an interval [] such that all points inside the interval are positive and all those outside are negative. Given is the target concept and is some concept, then the error region is the union of two intervals: consists of all points which are not in the interval of all points but which are not in h. Assume that the distribution D is uniform (just for the sake of this example) then, which is the length of the interval I. As a result, 2. The sample set

contains sample points from 0 to 1 with increments of 2. Therefore, every interval larger than must be hit by at least one point from S and by definition

It is important to note that if -net then we are guaranteed that be the consistent hypothesis with S (returned by the learning algorithm L). Becuase h is consistent, not been hit by S (recall that is the error region with respect to the target concept c, thus if h is consistent then it agrees with c over S and therefore S does not hit -net for ∆(c) we must have ) (recall that by definition S hits all error regions with weight larger than ). As a result, the error region must have a weight smaller than which means that

The conclusion is that if we can bound the probability that a random sample -net for ∆(c), then we have bounded the probability that a concept h consistent with . This is the goal of the proof of the double-sampling theorem which we are about to prove below:

dom sample of size m (sampled i.i.d. according to the unknown distribution D) and let A be the event that -net for ∆(c). From the preceding discussion our goal is to upper bound the probability for A to occur, i.e.,

If A occurs, i.e., -net, then by definition there must be some region ) which is not hit by . Note that ) for some concept h which is consistent with . At this point the space of possibilities is infinite, because the probability that we fail to hit random examples is at most (1. Thus the probability that we fail to hit ) is bounded from above by — which does not help us due to the fact that is infinite. The idea of the proof is to turn this into a finite space by using another sample, as follows.

Let be another random sample of size m. We will select m (for both ) to guarantee a high probability that many times. In fact we wish that with probability of at least 0.5:

We will use the Chernoff bound (lower tail) to obtain a bound on the right-hand side term. Recall that if we have m Bernoulli trials (coin tosses) with expectation and we consider the random variable with expectation (note that all 0 1 we have:

Considering the sampling of m examples that form as Bernoulli trials, we have that (since the probability that an example will hit r is at least 5. We obtain therefore:

which happens when To summarize what we have obtained so far, we have calculated the probability that times given that r was fixed using the previous sampling, i.e., given that does not form an -net. To formalize this, let B denote the combined event that does not form an 2 times. Then, we have shown that for

From this we can calculate prob(B):

which means that our original goal of bounding prob(A) is equivalent to finding a bound crucial point with the new goal is that to analyze the probability of the event B, we need only to consider a finite number of possibilities, namely to consider the regions of

This is because the occurrence of the event B is equivalent to saying that there is some ) such that 2 (i.e., the region r is hit at least 2 times) and . This is because Πcontains all the subsets of realized as intersections over all regions in ∆). Thus even though we have an infinite number of regions we still have a finite number of subsets. We wish therefore to analyze the following probability:

Let a random sample of 2m (note that since the sampling is i.i.d. it is equivalent to sampling separately) and r satisfying being fixed. Consider some random partitioning of consider then the problem of estimating the probability that This problem is equivalent to the following combinatorial question: we have 2m balls, each colored Red or Blue, with exaclty 2 Red balls. We divide the 2m balls into groups of equal size and we are interested in bounding the probability that all of the l balls fall in (that is, the probability that ). This in turn is equivalent to first dividing the 2m uncolored balls into groups and then randomly choose l of the balls to be colored Red and analyze the probability that all of the Red balls fall into . This probability is exactly

This probability was evaluated for a fixed S and r. Thus, the probability that this occurs for ) satisfying 2 (which is prob(B)) can be calculated by summing over all possible fixed r and applying the

union bound

from which it follows that:

Few comments are worthwhile at this point:

(i) It is possible to show that the upper bound on the sample complexity m is tight by showing that the lower bound on ) (see [[3], pp. 62]).

(ii) The treatment above holds also for the unrealizable case (target concept ) with slight modifications to the bound. In this context, the learning algorithm L must simply minimize the sample (empirical) error ˆerr(h) defined:

The generalization of the double-sampling theorem (Derroye’82) states that the empirical errors converge uniformly to the true errors:

Taken together, we have arrived to a fairly remarkable result. Despite the fact that the distribution D from which the training sample S is drawn from is unknown (but is known to be fixed), the learner simply needs to minimize the empirical error. If the sample size m is large enough the learner is guaranteed to have minimized the true errors for some accuracy and confidence parameters which define the sample size complexity. Equivalently,

Not only is the convergence is independent of D but also the rate of convergence is independent (namely, it does not matter where the optimal is located). The latter is very important because without it one could arbitrarily slow down the convergence rate by maliciously choosing D. The beauty of the results above is that D does not have an effect at all — one simply needs to choose the sample size to be large enough for the accuracy, confidence and VC dimension of the concept class to be learned over.

In Lecture 4 we discussed the large margin principle for finding an optimal separating hyperplane. It is natural to ask how does the PAC theory presented so far explains why a maximal margin hyperplane is optimal with regard to the formal sense of learning (i.e. to generalization from empirical errors to true errors)? We saw in the previous section that the sample complexity ) depends also on the VC dimension of the concept class — which is n + 1 for hyperplanes in . Thus, another natural question that may certainly arise is what is the gain in employing the ”kernel trick”? For a fixed m, mapping the input instance space X of dimension n to some higher (exponentially higher) feature space might simply mean that we are compromising the accuracy and confidence of the learner (since the VC dimension is equal to the instance space dimension plus 1).

Given a fixed sample size m, the best the learner can do is to minimize the empirical error and at the same time to try to minimize the VC dimension d of the concept class. The smaller d is, for a fixed m, the higher the accuracy and confidence of the learning algorithm. Likewise, the smaller d is, for a fixed accuracy and confidence values, the smaller sample size is required.

There are two possible ways to decrease d. First is to decrease the dimension n of the instance space X. This amounts to ”feature selection”, namely find a subset of coordinates that are the most ”relevant” to the learning task r perform a dimensionality reduction via PCA, for example. A second approach is to maximize the margin. Let the margin associated with the separating hyperplane h (i.e. consistent with the sample input vectors have a bounded norm, . It can be shown that the VC dimension of the concept class of hyperplanes with margin

Thus, if the margin is very small then the VC dimension remains n + 1. As the margin gets larger, there comes a point where and as a result the VC dimension decreases. Moreover, mapping the instance space X to some higher dimension feature space will not change the VC dimension as long as the margin remains the same. It is expected that the margin will not scale down or will not scale down as rapidly as the scaling up of dimension from image space to feature space.

To conclude, maximizing the margin (while minimizing the empirical error) is advantageous as it decreases the VC dimension of the concept class and causes the accuracy and confidence values of the learner to be largely immune to dimension scaling up while employing the kernel trick.

Appendix

98 Appendix

Let X, Y be two random variables and let f(x, y) be some function on ) be the probability of the event x and y occurring together. The expectation E[f(x, y)] is defined:

. The mean, variance and covariance are defined:

σ

In vector-matrix notation, let x represent the n random variables of i.e., is an instance vector and p(x) is the probability of the instance occurrence. Then the mean is a vector and the covariance matrix E are defined:

Note that the covariance matrix E is the linear superposition of rank-1 matrices (with coefficients p(x). The diagonal of E containes the variances of the variables For a uniform distribution and a sample data S consisting of m points, let the matrix whose columns consist of the points centered around the mean: . The (sample) covariance matrix is

A0.2 Derivatives of Matrix Operations: Scalar Functions of a Vector

The two most important examples of a scalar function of a vector x are the linear form and the quadratic form for some square matrix A.

where the derivative ) using the rule of products and noting that d(Ax) = Adx. Thus, is symmetric then

Consider first the general optimization with equality constraints which gives rise to the notion of Lagrange multipliers.

where is a vector function (each from . We want to derive a necessary and sufficient constraint for a point to be a local minimum subject to the k equality constraints h(x) = 0. Assume that point, meaning that the gradient vectors ) are linearly independent. Note that and the null space of this matrix:

defines the tangent plane at the point . We have the following fundamental theorem:

in other words, all vectors y spanning the tangent plane at the point also perpendicular to the gradient of

100 Appendix

curve on the surface h(x) = 0, i.e., h(x(t)) = 0. Let

(0) the tangent to the curve at . From the definition of tangency, the vector (0) is a local extremum of f(x), then

where the coefficients are called Lagrange Multipliers and the expression:

is called the Lagrangian of the optimization problem (0.1).

Consider next the general constrained optimization with inequality constraints (called “non-linear programming”):

where We will assume that the optimal solution regular point which has the following meaning: Let J be the set of indices ) = 0, then is a regular point if the gradient vectors are linearly independent. A basic result (we will not prove here) is the Karush-Kuhn-Tucker (KKT) theorem:

Let be a local minimum of the problem and suppose is a regular point. Then, there exist such that:

Fig. A0.1. Geometric interpreatation of Duality (see text).

Note that the condition ) = 0 is equivalent to the condition that ) = 0 (since 0 thus there sum cannot vanish unless each term vanishes) which in turn implies: The expression

is the Lagrangian of the problem (0.2) and the associated condition 0 is called the KKT condition.

The remaining concepts we need are the “duality” and the “Lagrangian Dual” problem.

102 Appendix

The optimization problem (0.2) is called the “Primal” problem. The Lagrangian Dual problem is defined as:

where