Large scale semi-supervised linear support vector machines

[0001] The invention is a method for large scale semi-supervised linear support vector machines. Specifically, the invention provides a system and method that uses a modified finite Newton method for fast training of linear SVMs using a mix of labeled and unlabeled examples.

[0002] Many problems in information processing involve the selection or classification of items in a large data set. For example, web-based companies such as Yahoo! have to frequently classify web pages as belonging to one group or the other, e.g., as commercial or non-commercial.

[0003] Currently, large amounts of data can be cheaply and automatically collected. However, labeling of the data typically involves expensive and fallible human participation. For example, a single web-crawl by a search engine, such as Yahoo or Google, indexes billions of webpages. Only a very small fraction of these web-pages can be hand-labeled by human editorial teams and assembled into topic directories. The remaining web-pages form a massive collection of unlabeled documents.

[0004] The modified finite Newton algorithm, described in co-pending application Ser. No. 10/949,821, entitled “A Method And Apparatus For Efficient Training Of Support Vector Machines,” filed Sep. 24, 2004, the entirety of which is incorporated herein by reference, describes a method for training linear support vector machines (SVMs) on sparse datasets with a potentially very large number of examples and features. Such datasets are generated frequently in domains like document classification. However, the system and method described in that application incorporates only labeled data in a finite Newton algorithm (abbreviated L2-SVM-MFN). Large scale learning is often realistic only in a semi-supervised setting where a small set of labeled examples is available together with a large collection of unlabeled data.

[0005] A system and method that provides for extension of linear SVMs for semi-supervised classification problems involving large, and possibly very high-dimensional but sparse, partially labeled datasets is desirable. In many information retrieval and data mining applications, linear methods are strongly preferred because of their ease of implementation, interpretability and empirical performance. The preferred embodiments of the system and method described herein clearly address this and other needs.

[0006] In one preferred embodiment, a transductive support vector machine (TSVM) system and method for linear semi-supervised classification on large, sparse datasets is provided. The system and method exploits data sparsity and linearity to provide superior scalability. According to another preferred embodiment, a multiple switching heuristic further improves TSVM training significantly for large scale applications.

[0007] In another preferred embodiment, a system and method for semi-supervised SVMs is disclosed that uses global optimization using mean field methods. The method generates a family of objective functions with non-convexity that is controlled by an annealing parameter. The global minimizer is tracked with respect to this parameter. This method alleviates the problem of local minima in the TSVM optimization procedure that results in significantly better operation in many applications, while being computationally judicious.

[0008] According to another preferred embodiment, a computerized system and method for semi-supervised learning is provided. A training set example elements is received. Elements of the training set that are determined to fall within a classification group are labeled. The training set thereby has labeled elements and unlabelled elements. The system uses selected labeled elements and unlabeled elements as examples in a large scale semi-supervised support vector machine to select a linear classifier. The system receives unclassified data items. Using the selected linear classifier, the received unclassified data elements are classified.

[0009] FIG. 1 is a block diagram illustrating components of a search engine in which one embodiment operates;

[0010] FIG. 2 is an example of a news web page that can be categorized using one embodiment;

[0011] FIG. 3 is a flow diagram illustrating steps performed by the system according to one embodiment associated with search engine relevance ranking;

[0012] FIG. 4 is a table listing conjugate gradient scheme iterations according to one embodiment;

[0013] FIG. 5 is a table outlining a full line search routine method performed by one embodiment;

[0014] FIG. 6 is a table having an abridged pseudo-code listing for a method performed by one embodiment;

[0015] FIG. 7 is a graph illustrating an L₂ loss function over unlabelled examples for transductive SVM;

[0016] FIG. 8 is a table outlining a transductive SVM algorithm performed by one embodiment;

[0017] FIG. 9 is a graph illustrating an L₂ loss function over unlabelled examples for transductive SVM;

[0018] FIG. 10 is a table illustrating steps for optimizing a weight vector w according to one embodiment;

[0019] FIG. 11 is a graph of a root function used by an embodiment;

[0020] FIG. 12 is table illustrating steps for optimizing a belief probability p according to one embodiment; and

[0021] FIG. 13 is table illustrating steps for a method of mean field annealing according to one embodiment.

[0022] A preferred embodiment of a large scale semi-supervised linear support vector machine, constructed in accordance with the claimed invention, is directed towards a modified finite Newton method for fast training of linear SVMs. According to one embodiment, the system uses two efficient and scalable methods for training semi-supervised SVMs. The first method provides a fast implementation of a variant of linear a transductive SVM (TSVM). The second method is based on a mean field annealing (MFA) approach, designed to alleviate the problem of local minimum in the TSVM optimization procedure.

[0023] In one embodiment, as an example, and not by way of limitation, an improvement in Internet search engine labeling of web pages is provided. The World Wide Web is a distributed database comprising billions of data records accessible through the Internet. Search engines are commonly used to search the information available on computer networks, such as the World Wide Web, to enable users to locate data records of interest. A search engine system 100 is shown in FIG. 1. Web pages, hypertext documents, and other data records from a source 101, accessible via the Internet or other network, are collected by a crawler 102. The crawler 102 collects data records from the source 101. For example, in one embodiment, the crawler 102 follows hyperlinks in a collected hypertext document to collect other data records. The data records retrieved by crawler 102 are stored in a database 108. Thereafter, these data records are indexed by an indexer 104. Indexer 104 builds a searchable index of the documents in database 108. Common prior art methods for indexing may include inverted files, vector spaces, suffix structures, and hybrids thereof. For example, each web page may be broken down into words and respective locations of each word on the page. The pages are then indexed by the words and their respective locations. A primary index of the whole database 108 is then broken down into a plurality of sub-indices and each sub-index is sent to a search node in a search node cluster 106.

[0024] To use search engine 100, a user 112 typically enters one or more search terms or keywords, which are sent to a dispatcher 110. Dispatcher 110 compiles a list of search nodes in cluster 106 to execute the query and forwards the query to those selected search nodes. The search nodes in search node cluster 106 search respective parts of the primary index produced by indexer 104 and return sorted search results along with a document identifier and a score to dispatcher 110. Dispatcher 110 merges the received results to produce a final result set displayed to user 112 sorted by relevance scores.

[0025] As a part of the indexing process, or for other reasons, most search engine companies have a frequent need to classify web pages as belonging to one “group” or another. For example, a search engine company may find it useful to determine if a web page is of a commercial nature (selling products or services), or not. As another example, it may be helpful to determine if a web page contains a news article about finance or another subject, or whether a web page is spam related or not. Such web page classification problems are binary classification problems (x versus not x). To develop a classifier that can do such distinguishing usually takes a large sample of web pages labeled by editorial teams and for use as a training set for designing a classifier.

[0026] Referring to FIG. 2, there is shown an example of a web page that has been categorized. In this example, the web page is categorized as a “Business” related web page, as indicated by the topic indicator 225 at the top of the page. Other category indicators 225 are shown. Thus, if a user had searched for business categorized web pages, then the web page of FIG. 2 would be listed, having been categorized as such.

[0027] With reference to FIG. 3, a flow diagram illustrates the steps performed according to one embodiment. In step 300, a training set of documents is provided. In one embodiment, this training set includes randomly selected unlabelled web pages. In step 302, a sub-set of the training set is labeled as meeting the criteria for inclusion in or exclusion from a group (classification group). For example, the group may include sports web pages in the case of the Internet web-page example discussed above. This step may involve, by way of example, and not by way of limitation, human review of sample web pages for labeling of web pages that are determined to fall within the group (for example, sports web pages) and also possibly labeling some pages as outside the group (for example, non-sports pages), in the labeled training set. Unlike previous systems and methods, it is not strictly necessary to label examples of data that does not fall within the criteria (for example, non-sports pages). In one embodiment, non-labeled examples (e.g. unlabelled sports and non-sports pages) can be used by the system, but it is not necessary. Once the sub-set is labeled, the whole training set of labeled and unlabelled data is available to the system.

[0028] In step 304, a modified finite Newton method for fast training of large scale semi-supervised linear SVMs is used to find a classifier. Even though the labeled training set can be any set of examples from within and outside the group, it is preferable for the method to retain knowledge of the fraction of the group in the entire distribution. For example, in the case of the group consisting of sports web pages, if this group forms 1% of the entire collection of web pages, then this number is preferably conveyed as an input to the method as compared to unlabeled examples. In the description of the method below, this number, written as a fraction (for example, 0.01) is called r.

[0029] In Step 308, new data, such as a web-page of a document to be classified is input into the system. In Step 310, the semi-supervised SVM obtained in Step 304 is used to classify new data.

[0030] A modified version of L₂-SVM-MFN algorithm is now discussed to suit the description of two semi-supervised extensions discussed thereafter. The first extension provides an efficient implementation of a variant of Linear Transductive SVM (TSVM). The second extension is based on a mean field annealing (MFA) algorithm designed to track the global optimum of the TSVM objective function.

[0031] Modified Finite Newton Linear L₂-SVM

[0032] The modified finite Newton L₂-SVM method (L₂-SVM-MFN) for Linear SVMs is ideally suited to sparse datasets with large number of examples and possibly large number of features. In a typical application, such as document classification, many training documents are collected and processed into a format that is convenient for mathematical manipulations. For example, each document may be represented as a collection of d features associated with a vocabulary of d words. These features may simply indicate the presence or absence of a word (binary features), or measure the frequency of a word suitably normalized by its importance (term frequency inverse document frequency (TFIDF) features). Even though the vocabulary might be large, only a very small number of words appear in any document relative to the vocabulary size. Thus, each document is sparsely represented as a bag of words. A label is then manually assigned to each document identifying a particular category to which it belongs (e.g. “commercial” or not). The task of a classification algorithm (e.g. SVM) is to produce a classifier that can reliably identify the category of new documents based on training documents.

[0033] Given a binary classification problem with l labeled examples {x_i,y_i}_i=1^l where the input patterns x_iεR^d (e.g. a document) and the labels y_iε{+1,−1}, L₂-SVM-MFN provides an efficient primal solution to the following SVM optimization problem: $\begin{array}{cc}{w}^{*}=\underset{w\in {ℝ}^d}{\arg \min }\frac{1}{2}\sum _{i=1}^t{\max \left[0,1-y_i\left(w^t{x}_i\right)\right]}^2+\frac{\lambda }{2}{w}^2& \left(1\right)\end{array}$
Here, λ is a real-valued regularization parameter and sign(w*^Tx) is the final classifier.

[0034] This objective function differs from the standard SVM problem in some respects. First, instead of using the hinge loss as the data fitting term, the square of the hinge loss (or the so-called quadratic soft margin loss function) is used. This makes the objective function continuously differentiable, allowing easier applicability of gradient techniques. Secondly, the bias term (“b”) is also regularized. In the problem formulation of equation 1, it is assumed that an additional component in the weight vector and a constant feature in the example vectors have been added to indirectly incorporate the bias. This formulation combines the simplicity of a least squares aspect with algorithmic advantages associated with SVMs. It should also be noted that the methods discussed herein can be applied to other loss functions.

[0035] In one embodiment, a version of L2-SVM-MFN where a weighted quadratic soft margin loss function is used. $\begin{array}{cc}{w}^{*}=\underset{w\in {ℝ}^d}{\arg \min }f\left(w\right)=\underset{w\in {ℝ}^d}{\arg \min }\frac{1}{2}\sum _{i\in j\left(w\right)}{c}_i{d}_i^2\left(w\right)+\frac{\lambda }{2}{w}^2& \left(2\right)\end{array}$

[0036] In equation 2, equation 1 is re-written in terms of a partial summation of d_i(w)=w^Tx_i−y_i over an index set j(w)={i:y_i(w^Tx_i)<1}. Additionally, the loss associated with the i^th example has a cost c_i. f(w) refers to the objective function being minimized, evaluated at a candidate solution w. It is noted that if the index set j(w) were independent of w and ran over all data points, this would simply be the objective function for weighted linear regularized least squares (RLS).

[0037] In one embodiment, f is a strictly convex, piecewise quadratic, continuously differentiable function having a unique minimizer. The gradient of f at w is given by: $\nabla f\left(w\right)=\lambda w+\sum _{i\in j\left(w\right)}{c}_i{d}_i\left(w\right)x_i=\lambda w+X_{j\left(w\right)}^T{C}_{j\left(w\right)}\left[X_{j\left(w\right)}w-Y_{j\left(w\right)}\right]$
where X_j(w) is a matrix whose rows are the feature vectors of training points corresponding to the index set j(w), Y_j(w) is a column vector containing labels for these points, and C_j(w) is a diagonal matrix that contains the costs c_i for these points along its diagonal.

[0038] In one embodiment, L₂-SVM-MFN is a primal algorithm that uses the Newton's method for unconstrained minimization of a convex function. The classical Newton's method is based on a second order approximation of the objective function, and involves updates of the following kind:
w^(k+1)=w^(k)+δ^(k)n^(k)  (3)
where the step size δ^kεR, and the Newton direction n^kεR^d is given by:
n^(k)=−[∇²f(w^(k))]⁻¹∇f(w^(k))

[0039] ∇f(w^(k)) is the gradient vector and ∇²f(w^(k)) is the Hessian matrix of f at w^(k). However, the Hessian does not exist everywhere, since f is not twice differentiable at those weight vectors w where w^Tx_i=y_i for some index i. For this reason, a finite Newton method works around this issue through a generalized definition of the Hessian matrix. The modified finite Newton procedure proceeds as follows. The step w^(k)=w^(k)+n^(k) in the Newton direction can be seen to be given by solving the following linear system associated with a weighted linear regularized least squares problem over the data subset defined by the indices j(w^(k)):
[λI+X_j(w(k))C_j(w(k))X_j(w(k))] w^(k)=X_j(w(k))^TC_j(w(k))Y_j(w(k)))  (4)
where I is the identity matrix. Once w^(k) is obtained, w^(k+1) is obtained from equation 3 by setting w^(k+1)=w^(k)+δ^k( w^(k)−w^(k)) after performing an exact line search for δ^k, i.e. by exactly solving a one-dimensional minimization problem: ${\delta }^{\left(k\right)}=\underset{\delta \ge 0}{\mathrm{argmin}}f\left(w^{\left(k\right)}+\delta \left({\stackrel{\_}{w}}^{\left(k\right)}-w^{\left(k\right)}\right)\right)$

[0040] The modified finite Newton procedure has the property of finite convergence to the optimal solution. The key features that bring scalability and numerical robustness to L₂-SVM-MFN are: (a) Solving the regularized least squares system of equation 4 by a numerically well-behaved conjugate gradient scheme referred to as conjugate gradient for least squares method (CGLS), which is designed for large, sparse data matrices X. (b) Due to the one-sided nature of margin loss functions, these systems are solved over restricted index sets j(w), which can be much smaller than the whole dataset. This also allows additional heuristics to be developed such as terminating CGLS early when working with a crude starting guess, such as 0, and allowing the line search step to yield a point where the index set j(w) is small. Subsequent optimization steps then work on smaller subsets of the data.

[0041] CGLS

[0042] The CGLS procedure solves large, sparse, weighted regularized least squares problems of the following form:
[λI+X^TCX]β=X^TCY  (5)
The key computational issue in equation 5 is to avoid the construction of the large and dense matrix X^TCX, and work only with the sparse matrix X and the diagonal cost matrix (stored as a vector) C.

[0043] Starting with a guess solution, β₀, Conjugate Gradient performs iterations of the form:
β^(j+1)=β^(j)+γ^(j)p^(j)
where p^(j) is a search direction and γ^(j)εR gives the step in that direction. The residual vector (the difference vector between LHS and RHS of equation 5 for a candidate β, which is also the gradient of the associated quadratic form evaluated at β) is therefore updated as:
r^(j+1)=X^TCY−[λI+X^TCX]β^(j+1)=X^Tz^(j+1)−λβ^(j+1)
The following intermediate vectors are introduced: $\begin{array}{c}{z}^{\left(j+1\right)}=C\left(Y-X{\beta }^{\left(j+1\right)}\right)\\ =C\left(Y-\left[X{\beta }^{\left(j\right)}+{\gamma }^{\left(j\right)}{\mathrm{Xp}}^{\left(j\right)}\right]\right)\\ =z^{\left(j\right)}-{\gamma }^{\left(j\right)}{\mathrm{Cq}}^{\left(j\right)}\\ \mathrm{where}{q}^{\left(j\right)}={\mathrm{Xp}}^{\left(j\right)}\end{array}$
The optimal choice of γ^(j) is given by: ${\gamma }^{\left(j\right)}=\frac{{r^{\left(j\right)}}^2}{p^{\left(j\right)T}\left(\gamma I+X^T\mathrm{CX}\right)p^{\left(j\right)}}=\frac{{r^{\left(j\right)}}^2}{\lambda {p^{\left(j\right)}}^2+q^{\left(j\right)T}{\mathrm{Cq}}^{\left(j\right)}}$
Finally, the search directions are updated as: $p^{\left(j+1\right)}=r^{\left(j+1\right)}+w^{\left(j\right)}{p}^{\left(j\right)}$ $\mathrm{where}{w}^{\left(j\right)}=\frac{{r^{\left(j+1\right)}}^2}{{r^{\left(j\right)}}^2}$

[0044] The CGLS iterations are terminated when the norm of the gradient r^(j+1) becomes small enough relative to the norm of the iterate z^(j+1) or if the number of iterations exceed a certain maximum allowable number.

[0045] The CGLS iterations are listed in the table of FIG. 4. The data matrix X is only involved in the computations through matrix vector multiplication for computing the iterates q^(j) and r^(j). This forms the dominant expense in each iteration (the product with C simply scales each element of a vector). If there are no non-zero elements in the data matrix, this has O(n₀) cost, where n₀ is the number of non-zero elements in the data matrix. As a subroutine of L₂-SVM-MFN, CGLS is typically called on a small subset of the full data set. The total cost of CGLS is O(t_cglsn₀) where t_cgls is the number of iterations, which depends on the practical rank of X and is typically found to be very small relative to the dimensions of X (number of examples and features). The memory requirements are also minimal: only five vectors need to be maintained, including the outputs over the currently active set of data points.

[0046] Finally, another feature of CGLS is worth noting. Suppose the solution β of a regularizer least squares problem is available, i.e. the linear system in equation 5 has been solved using CGLS. If there is a need to solve a perturbed linear system, it is greatly advantageous in many settings to start the CG iterations for the new system with β as the initial guess. This is often called seeding. If the starting residual is small, CGLS can converge much faster than with a guess of 0 vector. The utility of this feature depends on the nature and degree of perturbation. In L2-SVM-MFN, the candidate solution w^(k) obtained after line search in iteration k is seeded for the CGLS computation of w^k in the next iteration. Also, in tuning λ over a range of values, it is computationally valuable to seed the solution for a particular λ onto the next value. For the transductive SVM implementation with L₂-SVM-MFN, solutions are seeded across linear systems with slightly perturbed label vectors and data matrices.

[0047] Line Search

[0048] Given the vectors w, w in some iteration of L₂-SVM-MFN, the line search step includes solving: ${\delta }^{*}=\underset{\delta \ge 0}{\arg \min }\varphi \left(\delta \right)=f\left(w_{\delta }\right)$
where w_δ=w+δ( w−w).

[0049] The one-dimensional function φ(δ) is the restriction of the objective function ƒ on the ray from w onto w. Hence, like ƒ, φ(δ) is also a continuously differentiable, strictly convex, piecewise quadratic function with a unique minimizer δ* given by φ′(δ8)=0. Thus, one needs to find the root of the piecewise linear function $\begin{array}{cc}{\varphi }^{\prime }\left(\delta \right)=\lambda w_{\delta }^T\left(\stackrel{\_}{w}-w\right)+\sum _{i\in j\left(w_{\delta }\right)}{c}_i{d}_i\left(w_{\delta }\right)\left({\stackrel{\_}{o}}_i-o_i\right)& \left(6\right)\end{array}$
where o=Xw, o=X w.

[0050] The linear pieces of φ′ are defined over those intervals where j(w_δ) remains constant. Thus, the break points occur at a certain set of values δ_i where $w_{{\delta }_j}^T{x}_i=y_i.$ for some data point indexed by, i.e. ${\delta }_i=\frac{y_i-o_i}{{\stackrel{\_}{o}}_i-o_i}=\frac{1-y_i{o}_i}{y_i\left({\stackrel{\_}{o}}_j-o_i\right)}.$

[0051] Among these values, only those indices i where δ_i≧0 i.e. if iεj(w) (then y_io_i<1), so y_i( o_i−o_i)>0 or if i∉j(w) (then y_io_i>1) are considered, so y_i( o_i−o_i)<0. When δ is increased past a δ_i, in the former case the index i leaves j(w) and in the latter case it enters j(w). Reordering the indices so that δ_i, are sorted in a non-decreasing order as δ_j1,δ_j2 . . . , the root is then easily checked in each interval (δ_jk,δ_jk+1), k=1,2 . . . by keeping track of the slope of the linear piece in that interval. The slope is constant for each interval and non-decreasing as the search progresses through these ordered intervals. The interval in which the slope becomes non-negative for the first time brackets the root. Defining the extension of the linear piece in the interval $\left({\delta }_{j_k},{\delta }_{j_{k+1}}\right)\mathrm{as}{\varphi }_k^{\prime }\left(\delta \right)=\lambda w_{\delta }^T\left(\stackrel{\_}{w}-w\right)+\sum _{i\in j\left(w{\delta }_{j_k}\right)}{c}_i{d}_i\left(w_{\delta }\right)\left({\stackrel{\_}{o}}_i-o_i\right),$
the slope and the root computations are conveniently done by keeping track of $\begin{array}{c}L={\varphi }_k^{\prime }\left(0\right)\\ =\lambda w^T\left(\stackrel{\_}{w}-w\right)+\\ \sum _{i\in j\left(w_{{\delta }_{j_k}}\right)}{c}_i\left(o_i-y_i\right)\left({\stackrel{\_}{o}}_i-o_i\right)\mathrm{and}\\ R={\varphi }_k^{\prime }\left(1\right)\\ =\lambda {\stackrel{\_}{w}}^T\left(\stackrel{\_}{w}-w\right)+\\ \sum _{i\in j\left(w_{{\delta }_{j_k}}\right)}{c}_i\left({\stackrel{\_}{o}}_i-y_i\right)\left({\stackrel{\_}{o}}_i-o_i\right).\end{array}$
The full line search routine is outlined in the table of FIG. 5.

[0052] The table of FIG. 6 provides an abridged pseudo-code for L₂-SVM-MFN. Its computational complexity therefore is O(t_mfn t_cglsn₀) where t_mfn is the number of outer iterations of CGLS calls and line search, and t_cgls is the average number of CGLS iterations. These depend on the data set and the tolerance desired in the stopping criterion, but are typically very small. Therefore the complexity is found to be linear in the number of entries in the data matrix.

[0053] Semi-Supervised Linear SVMs

[0054] It is now assumed that there are l labeled examples {x_i,y_i}_i=1^l and u unlabeled examples {x′_j}_j=1^u with x_i,x′_jεR^d and y_iε{−1,+1}. The goal is to construct a linear classifier sign(w^Tx) that utilizes unlabeled data, typically in situations where l<<u. Semi-supervised algorithms provide L₂-SVM-MFN the capability of dealing with unlabeled data.

[0055] Transductive SVM

[0056] Transductive SVM (TSVM) appends an additional term in the SVM objective function whose role is to drive the classification hyperplane towards low data density regions. The following optimization problem is setup for standard TSVM: $w^{*}=\underset{w\in {ℝ}^d,{\left\{y_j^{\prime }\in \left\{-1,+1\right\}\right)}_{j=1}^u}{\arg \min }\frac{\lambda }{2}{w}^2+\frac{1}{2l}\sum _{i=1}^l\max \left[0,1-y_i\left(w^T{x}_i\right)\right]+\frac{{\lambda }^{\prime }}{2u}\sum _{j=1}^u\max \left[0,1-y_j^{\prime }\left(w^T{x}_j^{\prime }\right)\right]$ $\mathrm{subject}\mathrm{to}\text{:}\frac{1}{u}\sum _{j=1}^u\max \left[0,\mathrm{sign}\left(w^T{x}_j^{\prime }\right)\right]=r$
The labels on the unlabeled data, y′₁ . . . y′_u are {+1,−1}-valued variables in the optimization problem. In other words, TSVM seeks a hyperplane w and a labeling of the unlabeled examples, so that the SVM objective function is minimized, subject to the constraint that a fraction r of the unlabeled data be classified positive. SVM margin maximization in the presence of unlabeled examples can be interpreted as an implementation of the cluster assumption. In the optimization problem above, λ′ is a user-provided parameter that provides control over the influence of unlabeled data. If there is enough labeled data, λ,λ′,r can be tuned by cross-validation. An initial estimate of r can be made from the fraction of labeled examples that belong to the positive class and subsequent fine tuning can be done based on performance on a validation set.

[0057] In one method, this optimization is implemented by first using an inductive SVM to label the unlabeled data and then iteratively switching labels and retraining SVMs to improve the objective function. The TSVM algorithm is a wrapper around an SVM training procedure. In one existing software implemented system, the training of SVMs in the inner loops of TSVM uses dual decomposition techniques. In sparse, linear settings significant speed improvements can be obtained with L2-SVM-MFN over this software system. Thus, by implementing TSVM with L₂-SVM-MFN, improvements are made for semi-supervised learning on large, sparse datasets. The L₂-SVM-MFN retraining steps in the inner loop of TSVM are typically executed extremely fast by using seeding techniques. Additionally, in one embodiment of TSVM, more than one pair of labels may be switched in each iteration.

[0058] To develop the TSVM implementation with L₂-SVM-MFN, the objective function is considered corresponding to equation 7 below, but with the L₂ loss function: $w^{*}=\underset{w\in {ℝ}^d,{\left\{y_j^{\prime }\in \left\{-1,+1\right\}\right)}_{j=1}^u}{\arg \min }\frac{\lambda }{2}{w}^2+\frac{1}{2l}\sum _{i=1}^l{\max \left[0,1-y_i\left(w^T{x}_i\right)\right]}^2+\frac{{\lambda }^{\prime }}{2u}\sum _{j=1}^u{\max \left[0,1-y_j^{\prime }\left(w^T{x}_j^{\prime }\right)\right]}^2$ $\mathrm{subject}\mathrm{to}\text{:}\frac{1}{u}\sum _{j=1}^u\max \left[0,\mathrm{sign}\left(w^T{x}_j^{\prime }\right)\right]=r$

[0059] Note that this objective function above can also be equivalently written in terms of the following loss over each unlabeled example x:
min(max[0,1−(w^Tx)]², max[0,1+(w^Tx)]²)=max[0,1−|w^Tx|]²

[0060] A value of the label variable y is selected that minimizes the loss on the unlabeled example x, and rewritten in terms of the absolute value of the output of the classifier on x. This loss function is shown in FIG. 7. It should be noted that the L₁ and L₂ loss terms over unlabeled examples are very similar on the interval [−1,+1]. The non-convexity of this loss function implies that the TSVM training procedure is susceptible to local optima issues. Next, a mean field annealing procedure is outlined that can overcome this problem.

[0061] The TSVM algorithm with L₂-SVM-MFN is outlined in the table of FIG. 8. A classifier is obtained by first running L₂-SVM-MFN on just the labeled examples. Temporary labels are assigned to the unlabeled data by thresholding the soft outputs of this classifier so that the fraction of the total number of unlabeled examples that are temporarily labeled positive equals the parameter r.

[0062] Then starting from a small value of λ′, the unlabeled data is gradually brought in by increasing λ′ by a factor of R where R is set to 2 in the outer loop. This gradual increase of the influence of the unlabeled data is a way to protect TSVM from being immediately trapped in a local minimum. An inner loop identifies up to S pairs of unlabeled examples with positive and negative temporary labels such that switching these labels would decrease the objective function. L₂-SVM-MFN is then retrained with the switched labels.

[0063] Transductive L₂-SVM-MFN with multiple-pair switching converges in a finite number of steps. Previous transductive support vector machine implementations used single switching (S=1) of labels. However, in the implementation of the present embodiment, larger values for S (S>1, i.e., multiple switching) are used which leads to a significant improvement in efficiency without any loss in performance.

[0064] L₂-SVM-MFN on large sparse datasets combined with the efficiency gained from seeding w in the re-training steps (after switching labels or after increasing λ′) is very effective. Consider an iteration in Loop 2 of TSVM where a new pair of labels has been switched, and the solution w from the last retraining of L₂-SVM-MFN (marked as Re-training 2 in FIG. 8) is available for seeding. When the last L₂-SVM-MFN converged, its solution w is given by the linear systems:
[λI+X_I(w)^TC_I(w)X_I(w)]w=X_I(w)^TC_I(w)Y
where Y is the current label vector. When labels Y_i,Y_j are switched, back at the top of loop 2, the label vector is updated as:
Y=Y+2e_ij
where e_ij is a vector whose elements zero everywhere except in the i^th and the j^th position which are +1 and −1 or −1 and +1 respectively. Note also that if i,jεj(w) the re-training of L₂-SVM-MFN with w as the starting guess immediately encounters a call CGLS to solve the following perturbed system:
[λI+X_j(w)^TC_j(w)X_j(w)] w=X_j(w)^TC_j(w)[Y+2e_ij]
The starting residual vector r⁰ is given by: $\begin{array}{c}{r}^0=X_{j\left(w\right)}^T{C}_{j\left(w\right)}\left[Y+2e_{\mathrm{ij}}\right]-\left[\lambda I+X_{j\left(w\right)}^T{C}_{j\left(w\right)}{X}_{j\left(w\right)}\right]w\\ =r\left(w\right)+2X_{j\left(w\right)}^T{C}_{j\left(w\right)}{e}_{\mathrm{ij}}\\ \le \epsilon +2{\lambda }^{\prime }x_i-x_j\end{array}$
where r(w) in the second step is the final residual of w which fell below ε at the convergence of the last re-training. In applications, such as text categorization, TFIDF feature vectors are often length normalized and have positive entries. Therefore, ∥x_i−x_j∥≧√{square root over (2)}. This gives the following bound on the starting residual:
r⁰≧ε+2√{square root over (2)}λ′
which is much smaller than a bound of n√{square root over (n)}λ′ with a zero starting vector. Seeding is quite effective for Loop 1 as well, where λ′ is changed. With the two additional loops, the complexity of Transductive L₂-TSVM-MFN becomes O(n_switches t_mfn t_cglsn₀), where n_switches is the number of label switches. The outer loop executes a fixed number of times; the inner loop calls L₂-TSVM-MFN n_switches times. Typically, n_switches is expected to strongly depend on the data set and also on the number of labeled examples. Since it is difficult to apriori estimate the number of switches, this is an issue that is best understood from empirical observations.

[0065] Mean Field Annealing

[0066] The transductive SVM loss function over the unlabeled examples can be seen from FIG. 7 to be non-convex. This makes the TSVM optimization procedure susceptible to local minimum issues causing a loss in its performance in many situations. A new algorithm that is based on mean field annealing can be used to overcome this problem while also being computationally very attractive for large scale applications.

[0067] Mean field annealing (MFA) is an established tool for combinatorial optimization that approaches the problem from information theoretic principles. The discrete variables in the optimization problem are relaxed to continuous probability variables and a non-negative temperature parameter T is used to track the global optimum.

[0068] First, the TSVM objective function is re-written as follows: $w^{*}=\underset{w\in {ℝ}^d,{\left\{{\mu }_j\in \left\{0,1\right\}\right\}}_{j=1}^u}{\arg \min }\frac{\lambda }{2}{w}^2+\frac{1}{2l}\sum _{i=1}^l{\max \left[0,1-y_i\left(w^T{x}_i\right)\right]}^2+\frac{{\lambda }^{\prime }}{2u}\sum _{j=1}^u\max \left[{\mu }_j{\max \left[0,1-\left(w^T{x}_j^{\prime }\right)\right]}^2+\left(1-{\mu }_j\right){\max \left[0,1+\left(w^T{x}_j^{\prime }\right)\right]}^2\right)$

[0069] Binary valued variables μ_j=(1+y_j)/2 are introduced. Let p_jε[0,1] denote the probability that the unlabeled example belongs to the positive class. The Ising model 10 of Mean field annealing motivates the following objective function, where the binary variables μ_j are relaxed to probability variables p_j, and include entropy terms for the distributions defined by p_j: $\begin{array}{cc}{w}_T^{*}=\underset{w\in {ℝ}^d,{\left\{p_j\in \left[01\right]\right\}}_{j=1}^u}{\arg \min }\frac{\lambda }{2}{w}^2+\frac{1}{2l}\sum _{i=1}^l{\max \left[0,1-y_i\left(w^T{x}_i\right)\right]}^2+\frac{{\lambda }^{\prime }}{2u}\sum _{j=1}^u\left(p_j{\max \left[0,1-\left(w^T{x}_j^{\prime }\right)\right]}^2+\left(1-p_j\right){\max \left[0,1+\left(w^T{x}_j^{\prime }\right)\right]}^2\right)+\frac{T}{2u}\sum _{j=1}^u\left(p_j\log p_j+\left(1-p_j\right)\log \left(1-p_j\right)\right)& \left(8\right)\end{array}$

[0070] The “temperature” T parameterizes a family of objective functions. The objective function for a fixed T is minimized under the following class balancing constraints: $\begin{array}{cc}\frac{1}{u}\sum _{j=1}^u{p}_j=r& \left(9\right)\end{array}$

[0071] where r is the fraction of the number of unlabeled examples belonging to the positive class. As in TSVM, r is treated as a user-provided parameter. It may also be estimated from the labeled examples.

[0072] The solution to the optimization problem above is tracked as the temperature parameter T is lowered to 0. The final solution is given as: $\begin{array}{cc}{w}^{*}=\underset{T\to 0}{\lim }{w}_T^{*}& \left(10\right)\end{array}$
In practice, the system (indexer 104) monitors the value of the objective function in the optimization path and returns the solution corresponding to the minimum value achieved.

[0073] To develop an intuition for this method, the loss term is considered in the objective function associated with an unlabeled example as a function of the output of the classifier. This loss term is based on calculations to be described below. FIG. 9 plots this loss term for various values of T. As the temperature is decreased, the loss function deforms from a squared-loss shape where a global optimum is easier to achieve, to the TSVM loss function in FIG. 7. At high temperatures a global optimum is easier to obtain. The global minimizer is then slowly tracked as the temperature is lowered towards zero.

[0074] The optimization is done in stages, starting with high values of T and then gradually decreasing T towards 0. For each T, the problem in equations 8,9 is optimized by alternating the minimization over w and p=[p₁ . . . p_u] respectively. Fixing p, the optimization over w is done by L₂-SVM-MFN. Fixing w, the optimization over p can also be done easily as described below. Both these problems involve convex optimization and can be done exactly and efficiently. Details of these optimization steps follow.

[0075] Optimizing w

[0076] Described are the steps to efficiently implement the L₂-SVM-MFN loop for optimizing w keeping p fixed. The call to L₂-SVM-MFN is made on the data {circumflex over (X)}=[X^TX^′TX^′T]^T whose first l rows are formed by the labeled examples, and the next 2u rows are formed by the unlabeled examples appearing as two repeated blocks. The associated label vector and cost matrix are given by $\begin{array}{cc}C=\mathrm{diag}\left[\stackrel{\stackrel{l}{︷}}{\frac{1}{l}\dots \frac{1}{l}},\stackrel{\hat{Y}=\left[y_1,y_2\dots y_l\stackrel{\stackrel{u}{︷}}{1,1,\dots 1,}\stackrel{\stackrel{u}{︷}}{-1,-1,\dots -1}\right]}{\stackrel{\stackrel{u}{︷}}{\frac{{\lambda }^{\prime }{p}_1}{u}\dots \frac{{\lambda }^{\prime }{p}_u}{u}}\stackrel{\stackrel{u}{︷}}{\frac{{\lambda }^{\prime }\left(1-p_1\right)}{u}\dots \frac{{\lambda }^{\prime }\left(1-p_u\right)}{u}}}\right]& \left(11\right)\end{array}$
Even though each unlabeled data contributes two terms to the objective function, effectively only one term contributes to the complexity. This is because matrix-vector products, which form the dominant expense in L₂-SVM-MFN, are performed only on unique rows of a matrix. The output may be duplicated for duplicate rows.

[0077] In fact, the CGLS calls in L₂-SVM-MFN can be re-written so that the unlabeled examples appear only once in the data matrix. The CGLS call at some iteration where the active index set is j=j(w) for some current candidate weight vector w is as follows:
[λI+{circumflex over (X)}_j^TC_j{circumflex over (X)}_j] w={circumflex over (X)}^TC_jŶ_j  (12)

[0078] Note that if |w^Tx′_j|≧1, the unlabeled example x′_j appears as one row in the data matrix {circumflex over (X)}_j with label given by −sign(w^Tx′_j). If |w^Tx′_j|<1, the unlabeled example x′_j appears as two identical rows x′_j with both labels. Let j_lε1 . . . l be the indices of the labeled examples in the active set, j′₁ε1 . . . u be the indices of unlabeled examples with |w^Tx′_j|≧1 and j′₂ε1 . . . u be the indices of unlabeled examples with |w^Tx′_j|<1. Note that the index of every unlabeled example appears in one of these sets i.e., j′₁∪j′₂=1 . . . u. Equation 12 may be re-written as: $\left[\lambda I+\frac{1}{l}\sum _{i\in j_l}{x}_i^T{x}_i+\frac{{\lambda }^{\prime }}{u}\sum _{j\in j_1^{\prime }}{c}_j{x}_j^{\prime T}{x}_j+\frac{{\lambda }^{\prime }}{u}\sum _{j\in j_2^{\prime }}{x}_j^{\prime T}{x}_j\right]\stackrel{\_}{w}=\frac{1}{l}\sum _{i\in j_l}{y}_i{x}_i-\frac{{\lambda }^{\prime }}{u}\sum _{j\in j_1^{\prime }}{c}_j\mathrm{sign}\left(w^T{x}_j\right)x_j+\frac{{\lambda }^{\prime }}{u}\sum _{j\in j_2^{\prime }}\left(2p-1\right)x_j$
where c_j=p_j if sign(w^Tx′_j)=−1 and c_j=1−p_j if sign(w^Tx′_j)=1. Rewriting in matrix notation, an equivalent linear system is obtained that can be solved by CGLS:
[λI+{tilde over (X)}^T{tilde over (C)}{tilde over (X)}] w={tilde over (X)}^T{tilde over (C)}{tilde over (Y)}  (13)
where X=[X_ji^TX^′], C is a diagonal matrix and Y is the vector of effectively active labels. Each of these data objects have l+u rows. These are given by: $\begin{array}{cc}{\stackrel{\_}{C}}_{\mathrm{jj}}=\frac{1}{l},{\stackrel{\_}{Y}}_j=y_{i}j\in 1\dots j_l{\stackrel{\_}{C}}_{\left(j+j_l\right)\left(j+j_l\right)}=\frac{{\lambda }^{\prime }{p}_j}{u},Y_{j+j_l}=1\mathrm{if}j\in 1\dots u,j\in j_1^{\prime },\mathrm{sign}\left(w^T{x}_j^{\prime }\right)=-1{\stackrel{\_}{C}}_{\left(j+j_l\right)\left(j+j_l\right)}=\frac{{\lambda }^{\prime }\left(1-p_j\right)}{u},{\stackrel{\_}{Y}}_{j+j_l}=-1\mathrm{if}j\in 1\dots uj\in j_1^{\prime },\mathrm{sign}\left(w^T{x}_j^{\prime }\right)=1{\stackrel{\_}{C}}_{\left(j+j_l\right)\left(j+j_l\right)}=\frac{{\lambda }^{\prime }}{u},{\stackrel{\_}{Y}}_{j+j_l}=\left(2p_j-1\right)\mathrm{if}j\in 1\dots uj\in j_2^{\prime }& \left(14\right)\end{array}$
Thus, CGLS needs to only operate on data matrices with one instance of each unlabeled example using a suitably modified cost matrix and label vector.

[0079] After the CGLS step, the optimality conditions are checked. The optimality conditions can be re-written as:
∀iεj_l y_i o_i≦1+ε
∀iεj_l^c y_i o_i≧1−ε
∀jεj′₁ | o′_j|≧1−ε
∀jεj′₂ | o′_j|≦1+ε
For the subsequent line search step, appropriate output and label vectors are reassembled to call the routine in the table of FIG. 5. The steps for optimizing w are outlined in the table of FIG. 10.

[0080] Optimizing p

[0081] For the latter problem of optimizing p for a fixed w, the Lagrangian is constructed: $ℒ=\frac{{\lambda }^{\prime }}{2u}\sum _{j=1}^u\left(p_j{\max \left[0,1-\left(w^T{x}_j^{\prime }\right)\right]}^2+\left(1+p_j\right){\max \left[0,1\left(w^T{x}_j^{\prime }\right)\right]}^2+\frac{T}{2u}\underset{j=1}{\overset{u}{\sum (}}{p}_j\log p_j+\left(1-p_j\right)\log \left(1-p_j\right)\right)-v\left[\frac{1}{u}\sum _{j=1}^u{p}_j-r\right]$
Differentiating the Lagrangian with respect to p_j, the following is obtained: $\frac{\partial ℒ}{\partial p_j}=\frac{{\lambda }^{\prime }}{2u}\left({\max \left[0,1-\left(w^T{x}_j^{\prime }\right)\right]}^2-{\max \left[0,1+\left(w^T{x}_j^{\prime }\right)\right]}^2\right)+\frac{T}{2u}\log \frac{p_j}{1-p_j}-\frac{v}{u}=0$
Defining
g_j=λ′(max[0,1−(w^Tx′_j)]²−max[0,1+(w^Tx′_j)]²),
the expression for p_j is given by: $\begin{array}{cc}{p}_j=\frac{1}{1+e^{\frac{g_j-2v}{T}}}& \left(15\right)\end{array}$
Substituting this expression in the balance constraint in equation 9, a one-dimensional non-linear equation in v is as follows: $\frac{1}{u}\sum _{j=1}^u\frac{1}{1+e^{\frac{g_i-2v}{T}}}=r$
The value of 2 v is computed exactly by using a hybrid combination of Newton-Raphson iterations and the bisection method to find the root of the function $B\left(v\right)=\frac{1}{u}\sum _{i=l+1}^{l+u}\frac{1}{1+e^{\frac{g_i-v}{T}}}-r$
This method is rapid due to the quadratic convergence properties of Newton-Raphson iterations and fail-safe due to bisection steps. Note that the root exists and is unique, since one can see that B(v)→1−r>0 when v→∞, B(v)→−r<0 when v→−∞, and B(v) is a continuous, monotonically increasing function, as shown in FIG. 11. The root finding begins by bracketing the root in the interval [v₋,v₊] so that B(v₋)<0,B(v₊)>0 where v₋,v₊ are given by: $v_{-}=\min \left(g_1\dots g_i\right)-T\log \frac{1-r}{r}$ $v_{+}=\max \left(g_1\dots g_u\right)-T\log \frac{1-r}{r}$

[0082] The hybrid root finding algorithm performs Newton-Raphson iterations with a starting guess of $\frac{v_{-}+v_{+}}{2}$
and invokes the bisection method whenever an iterate goes outside the brackets. The steps for optimizing p are outlined in table if FIG. 12.

[0083] Stopping Criteria

[0084] For a fixed T, this alternate minimization proceeds until some stopping criterion is satisfied. A natural criterion is the known mean Kullback-Liebler divergence (relative entropy) KL(p,q) between current values of p_i and the values, say q_i, at the end of last iteration. Thus the stopping criterion for fixed T is: $\begin{array}{cc}\mathrm{KL}\left(p,q\right)=\sum _{j=1}^u{p}_j\log \frac{p_j}{q_j}+\left(1-p_j\right)\log \frac{1-p_j}{1-q_j}<u\epsilon & \left(16\right)\end{array}$
A good value for ε is 10⁻⁶. The temperature may be decreased in the outer loop until the total entropy falls below a threshold T, which is also taken to be 10⁻⁶: $\begin{array}{cc}H\left(p\right)=-\sum _{j=l}^u\left(p_j\log p_j+\left(1-p_j\right)\log \left(1-p_j\right)\right)<u\epsilon & \left(17\right)\end{array}$
The TSVM objective function is monitored as the optimization proceeds. $\begin{array}{cc}J\left(w\right)=\frac{\lambda }{2}{w}^2+\frac{1}{2l}\sum _{i=1}^l{\max \left[0,1-y_i\left(w^T{x}_i\right)\right]}^2+\frac{{\lambda }^{\prime }}{2u}\sum _{j=1}^u{\max \left[0,1-w^T{x}_j^{\prime }\right]}^2& \left(18\right)\end{array}$

[0085] The weight vector corresponding to the minimum transductive cost is returned as the solution. The steps of mean field annealing with L₂-SVM-MFN are outlined in the table of FIG. 13.

[0086] Application to One-Class Problems

[0087] The transductive SVM and mean field annealing algorithms described above can be deployed in a variety of settings involving different amounts of labeled and unlabeled data. Many settings present the task of identifying members of a certain class as opposed to distinguishing between well-specified classes. For example, in order to identify web documents concerning sports, it is much easier to label sports documents than to label the diverse and ill-characterized set of non-sports documents. In such problems, labeled examples come from a single class, very often with large amounts of unlabeled data containing some instances of the class of interest and many instances of the “others” class.

[0088] Being a special case of semi-supervised learning, the problem of one-class learning with unlabeled data can be addressed by the algorithms developed in this paper. These algorithms implement the cluster assumption under constraints on class ratios. For one-class problems, unlabeled data is expected to be helpful in biasing the classification hyperplane to pass through a low density region keeping clusters of the relevant class on one side and instances of the “others” class on the other side. The fraction r can first be crudely estimated using a small random sample and then finely adjusted based on validation performance of the trained classifier.

[0089] The various embodiments described above are provided by way of illustration only and should not be construed to limit the invention. Those skilled in the art will readily recognize various modifications and changes that may be made to the claimed invention without following the example embodiments and applications illustrated and described herein, and without departing from the true spirit and scope of the claimed invention, which is set forth in the following claims.