-
| Grade Trend | Hours slept | Hours Studied | Grade |
| Above | Low | High | Above |
| Below | High | Low | Below |
| Above | Low | High | Above |
| Above | High | High | Above |
| Below | Low | High | Below |
| Above | Low | Low | Below |
| Below | High | High | Below |
| Below | Low | High | Below |
| Above | Low | Low | Below |
| Above | High | High | Above |
A Support Vector Machine (SVM) is a very powerful and versatile +Machine Learning method, capable of performing linear or nonlinear +classification, regression, and even outlier detection. It is one of +the most popular models in Machine Learning, and anyone interested in +Machine Learning should have it in their toolbox. SVMs are +particularly well suited for classification of complex but small-sized or +medium-sized datasets. +
+ +The case with two well-separated classes only can be understood in an +intuitive way in terms of lines in a two-dimensional space separating +the two classes (see figure below). +
+ +The basic mathematics behind the SVM is however less familiar to most of us. +It relies on the definition of hyperplanes and the +definition of a margin which separates classes (in case of +classification problems) of variables. It is also used for regression +problems. +
+ +With SVMs we distinguish between hard margin and soft margins. The +latter introduces a so-called softening parameter to be discussed +below. We distinguish also between linear and non-linear +approaches. The latter are the most frequent ones since it is rather +unlikely that we can separate classes easily by say straight lines. +
@@ -440,7 +278,7 @@ MathJax.Hub.Config({
-
In computations we will translate all classes into numbers. Being -these binary classes, they can easily be split into ones and zeros. +
The theory behind support vector machines (SVM hereafter) is based on +the mathematical description of so-called hyperplanes. Let us start +with a two-dimensional case. This will also allow us to introduce our +first SVM examples. These will be tailored to the case of two specific +classes, as displayed in the figure here based on the usage of the petal data.
-See handwritten notes for Thursday November 11
+We assume here that our data set can be well separated into two +domains, where a straight line does the job in the separating the two +classes. Here the two classes are represented by either squares or +circles. +
+ + +from sklearn import datasets
+from sklearn.svm import SVC, LinearSVC
+from sklearn.linear_model import SGDClassifier
+from sklearn.preprocessing import StandardScaler
+import matplotlib
+import matplotlib.pyplot as plt
+plt.rcParams['axes.labelsize'] = 14
+plt.rcParams['xtick.labelsize'] = 12
+plt.rcParams['ytick.labelsize'] = 12
+
+
+iris = datasets.load_iris()
+X = iris["data"][:, (2, 3)] # petal length, petal width
+y = iris["target"]
+
+setosa_or_versicolor = (y == 0) | (y == 1)
+X = X[setosa_or_versicolor]
+y = y[setosa_or_versicolor]
+
+
+
+C = 5
+alpha = 1 / (C * len(X))
+
+lin_clf = LinearSVC(loss="hinge", C=C, random_state=42)
+svm_clf = SVC(kernel="linear", C=C)
+sgd_clf = SGDClassifier(loss="hinge", learning_rate="constant", eta0=0.001, alpha=alpha,
+ max_iter=100000, random_state=42)
+
+scaler = StandardScaler()
+X_scaled = scaler.fit_transform(X)
+
+lin_clf.fit(X_scaled, y)
+svm_clf.fit(X_scaled, y)
+sgd_clf.fit(X_scaled, y)
+
+print("LinearSVC: ", lin_clf.intercept_, lin_clf.coef_)
+print("SVC: ", svm_clf.intercept_, svm_clf.coef_)
+print("SGDClassifier(alpha={:.5f}):".format(sgd_clf.alpha), sgd_clf.intercept_, sgd_clf.coef_)
+
+# Compute the slope and bias of each decision boundary
+w1 = -lin_clf.coef_[0, 0]/lin_clf.coef_[0, 1]
+b1 = -lin_clf.intercept_[0]/lin_clf.coef_[0, 1]
+w2 = -svm_clf.coef_[0, 0]/svm_clf.coef_[0, 1]
+b2 = -svm_clf.intercept_[0]/svm_clf.coef_[0, 1]
+w3 = -sgd_clf.coef_[0, 0]/sgd_clf.coef_[0, 1]
+b3 = -sgd_clf.intercept_[0]/sgd_clf.coef_[0, 1]
+
+# Transform the decision boundary lines back to the original scale
+line1 = scaler.inverse_transform([[-10, -10 * w1 + b1], [10, 10 * w1 + b1]])
+line2 = scaler.inverse_transform([[-10, -10 * w2 + b2], [10, 10 * w2 + b2]])
+line3 = scaler.inverse_transform([[-10, -10 * w3 + b3], [10, 10 * w3 + b3]])
+
+# Plot all three decision boundaries
+plt.figure(figsize=(11, 4))
+plt.plot(line1[:, 0], line1[:, 1], "k:", label="LinearSVC")
+plt.plot(line2[:, 0], line2[:, 1], "b--", linewidth=2, label="SVC")
+plt.plot(line3[:, 0], line3[:, 1], "r-", label="SGDClassifier")
+plt.plot(X[:, 0][y==1], X[:, 1][y==1], "bs") # label="Iris-Versicolor"
+plt.plot(X[:, 0][y==0], X[:, 1][y==0], "yo") # label="Iris-Setosa"
+plt.xlabel("Petal length", fontsize=14)
+plt.ylabel("Petal width", fontsize=14)
+plt.legend(loc="upper center", fontsize=14)
+plt.axis([0, 5.5, 0, 2])
+
+plt.show()
+
-
See handwritten notes for Thursday November 11
-The aim of the SVM algorithm is to find a hyperplane in a +\( p \)-dimensional space, where \( p \) is the number of features that +distinctly classifies the data points. +
+ +In a \( p \)-dimensional space, a hyperplane is what we call an affine subspace of dimension of \( p-1 \). +As an example, in two dimension, a hyperplane is simply as straight line while in three dimensions it is +a two-dimensional subspace, or stated simply, a plane. +
+ +In two dimensions, with the variables \( x_1 \) and \( x_2 \), the hyperplane is defined as
+$$ +b+w_1x_1+w_2x_2=0, +$$ + +where \( b \) is the intercept and \( w_1 \) and \( w_2 \) define the elements of a vector orthogonal to the line +\( b+w_1x_1+w_2x_2=0 \). +In two dimensions we define the vectors \( \boldsymbol{x} =[x1,x2] \) and \( \boldsymbol{w}=[w1,w2] \). +We can then rewrite the above equation as +
+ +$$ +\boldsymbol{x}^T\boldsymbol{w}+b=0. +$$@@ -426,7 +276,7 @@ MathJax.Hub.Config({
-
See handwritten notes for Thursday November 11
-We limit ourselves to two classes of outputs \( y_i \) and assign these classes the values \( y_i = \pm 1 \). +In a \( p \)-dimensional space of say \( p \) features we have a hyperplane defines as +
+$$ +b+wx_1+w_2x_2+\dots +w_px_p=0. +$$ +If we define a +matrix \( \boldsymbol{X}=\left[\boldsymbol{x}_1,\boldsymbol{x}_2,\dots, \boldsymbol{x}_p\right] \) +of dimension \( n\times p \), where \( n \) represents the observations for each feature and each vector \( x_i \) is a column vector of the matrix \( \boldsymbol{X} \), +
+$$ +\boldsymbol{x}_i = \begin{bmatrix} x_{i1} \\ x_{i2} \\ \dots \\ \dots \\ x_{ip} \end{bmatrix}. +$$ -For final tree, see the above handwritten notes
+If the above condition is not met for a given vector \( \boldsymbol{x}_i \) we have
+$$ +b+w_1x_{i1}+w_2x_{i2}+\dots +w_px_{ip} >0, +$$ + +if our output \( y_i=1 \). +In this case we say that \( \boldsymbol{x}_i \) lies on one of the sides of the hyperplane and if +
+$$ +b+w_1x_{i1}+w_2x_{i2}+\dots +w_px_{ip} < 0, +$$ + +for the class of observations \( y_i=-1 \), +then \( \boldsymbol{x}_i \) lies on the other side. +
+ +Equivalently, for the two classes of observations we have
+$$ +y_i\left(b+w_1x_{i1}+w_2x_{i2}+\dots +w_px_{ip}\right) > 0. +$$ + +When we try to separate hyperplanes, if it exists, we can use it to construct a natural classifier: a test observation is assigned a given class depending on which side of the hyperplane it is located.
@@ -428,7 +288,7 @@ MathJax.Hub.Config({
- -
Let us try to develop our intuition about SVMs by limiting ourselves to a two-dimensional +plane. To separate the two classes of data points, there are many +possible lines (hyperplanes if you prefer a more strict naming) +that could be chosen. Our objective is to find a +plane that has the maximum margin, i.e the maximum distance between +data points of both classes. Maximizing the margin distance provides +some reinforcement so that future data points can be classified with +more confidence. +
- -# Common imports
-import numpy as np
-import pandas as pd
-import matplotlib.pyplot as plt
-from sklearn.tree import DecisionTreeClassifier
-from sklearn.model_selection import train_test_split
-from sklearn.tree import export_graphviz
-from sklearn.preprocessing import StandardScaler, OneHotEncoder
-from sklearn.compose import ColumnTransformer
-from IPython.display import Image
-from pydot import graph_from_dot_data
-import os
-
-# Where to save the figures and data files
-PROJECT_ROOT_DIR = "Results"
-FIGURE_ID = "Results/FigureFiles"
-DATA_ID = "DataFiles/"
-
-if not os.path.exists(PROJECT_ROOT_DIR):
- os.mkdir(PROJECT_ROOT_DIR)
-
-if not os.path.exists(FIGURE_ID):
- os.makedirs(FIGURE_ID)
-
-if not os.path.exists(DATA_ID):
- os.makedirs(DATA_ID)
-
-def image_path(fig_id):
- return os.path.join(FIGURE_ID, fig_id)
-
-def data_path(dat_id):
- return os.path.join(DATA_ID, dat_id)
-
-def save_fig(fig_id):
- plt.savefig(image_path(fig_id) + ".png", format='png')
-
-infile = open(data_path("grades.csv"),'r')
-
-# Read the experimental data with Pandas
-from IPython.display import display
-grades = pd.read_csv(infile,names = ('Trend','Sleep','Studied','Grade'))
-grades = pd.DataFrame(grades)
-
-# Features and targets
-X = grades.loc[:, grades.columns != 'Grade'].values
-y = grades.loc[:, grades.columns == 'Grade'].values
-
-# Create the encoder.
-encoder = OneHotEncoder(handle_unknown="ignore")
-# Assume for simplicity all features are categorical.
-encoder.fit(X)
-# Apply the encoder.
-X = encoder.transform(X)
-print(X)
-# Then do a Classification tree
-tree_clf = DecisionTreeClassifier(max_depth=2)
-tree_clf.fit(X, y)
-print("Train set accuracy with Decision Tree: {:.2f}".format(tree_clf.score(X,y)))
-#transfer to a decision tree graph
-export_graphviz(
- tree_clf,
- out_file="DataFiles/grade.dot",
- rounded=True,
- filled=True
-)
-cmd = 'dot -Tpng DataFiles/grade.dot -o DataFiles/grades.png'
-os.system(cmd)
-
-What a linear classifier attempts to accomplish is to split the +feature space into two half spaces by placing a hyperplane between the +data points. This hyperplane will be our decision boundary. All +points on one side of the plane will belong to class one and all points +on the other side of the plane will belong to the second class two. +
+Unfortunately there are many ways in which we can place a hyperplane +to divide the data. Below is an example of two candidate hyperplanes +for our data sample. +
@@ -508,7 +272,7 @@ os.system(cmd)
-
import os
-from sklearn.datasets import load_breast_cancer
-from sklearn.tree import DecisionTreeClassifier
-from sklearn.model_selection import train_test_split
-from sklearn.metrics import confusion_matrix
-from sklearn.tree import export_graphviz
+Let us define the function
+$$
+f(x) = \boldsymbol{w}^T\boldsymbol{x}+b = 0,
+$$
-from IPython.display import Image
-from pydot import graph_from_dot_data
-import pandas as pd
-import numpy as np
+as the function that determines the line \( L \) that separates two classes (our two features), see the figure here.
+Any point defined by \( \boldsymbol{x}_i \) and \( \boldsymbol{x}_2 \) on the line \( L \) will satisfy \( \boldsymbol{w}^T(\boldsymbol{x}_1-\boldsymbol{x}_2)=0 \).
-cancer = load_breast_cancer()
-X = pd.DataFrame(cancer.data, columns=cancer.feature_names)
-print(X)
-y = pd.Categorical.from_codes(cancer.target, cancer.target_names)
-y = pd.get_dummies(y)
-print(y)
-X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=1)
-tree_clf = DecisionTreeClassifier(max_depth=5)
-tree_clf.fit(X_train, y_train)
-
-export_graphviz(
- tree_clf,
- out_file="DataFiles/cancer.dot",
- feature_names=cancer.feature_names,
- class_names=cancer.target_names,
- rounded=True,
- filled=True
-)
-cmd = 'dot -Tpng DataFiles/cancer.dot -o DataFiles/cancer.png'
-os.system(cmd)
-
-The signed distance \( \delta \) from any point defined by a vector \( \boldsymbol{x} \) and a point \( \boldsymbol{x}_0 \) on the line \( L \) is then
+$$ +\delta = \frac{1}{\vert\vert \boldsymbol{w}\vert\vert}(\boldsymbol{w}^T\boldsymbol{x}+b). +$$@@ -473,7 +265,7 @@ os.system(cmd)
-
# Common imports
-import numpy as np
-from sklearn.model_selection import train_test_split
-from sklearn.tree import DecisionTreeClassifier
-from sklearn.datasets import make_moons
-from sklearn.tree import export_graphviz
-from pydot import graph_from_dot_data
-import pandas as pd
-import os
+How do we find the parameter \( b \) and the vector \( \boldsymbol{w} \)? What we could
+do is to define a cost function which now contains the set of all
+misclassified points \( M \) and attempt to minimize this function
+
-np.random.seed(42)
-X, y = make_moons(n_samples=100, noise=0.25, random_state=53)
-X_train, X_test, y_train, y_test = train_test_split(X,y,random_state=0)
-tree_clf = DecisionTreeClassifier(max_depth=5)
-tree_clf.fit(X_train, y_train)
+$$
+C(\boldsymbol{w},b) = -\sum_{i\in M} y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b).
+$$
-export_graphviz(
- tree_clf,
- out_file="DataFiles/moons.dot",
- rounded=True,
- filled=True
-)
-cmd = 'dot -Tpng DataFiles/moons.dot -o DataFiles/moons.png'
-os.system(cmd)
-
-We could now for example define all values \( y_i =1 \) as misclassified in case we have \( \boldsymbol{w}^T\boldsymbol{x}_i+b < 0 \) and the opposite if we have \( y_i=-1 \). Taking the derivatives gives us
+$$ +\frac{\partial C}{\partial b} = -\sum_{i\in M} y_i, +$$ + +and
+$$ +\frac{\partial C}{\partial \boldsymbol{w}} = -\sum_{i\in M} y_ix_i. +$$@@ -464,7 +270,7 @@ os.system(cmd)
-
Scikit-Learn has also another way to visualize the trees which is very useful, here with the Iris data.
+We can now use the Newton-Raphson method or different variants of the gradient descent family (from plain gradient descent to various stochastic gradient descent approaches) to solve the equations
+$$ +b \leftarrow b +\eta \frac{\partial C}{\partial b}, +$$ +and
+$$ +\boldsymbol{w} \leftarrow \boldsymbol{w} +\eta \frac{\partial C}{\partial \boldsymbol{w}}, +$$ - -from sklearn.datasets import load_iris
-from sklearn import tree
-X, y = load_iris(return_X_y=True)
-tree_clf = tree.DecisionTreeClassifier()
-tree_clf = tree_clf.fit(X, y)
-# and then plot the tree
-tree.plot_tree(tree_clf)
-
-where \( \eta \) is our by now well-known learning rate.
@@ -450,7 +262,7 @@ tree.plot_tree(tree_clf)
-
Alternatively, the tree can also be exported in textual format with the function exporttext. -This method doesn’t require the installation of external libraries and is more compact: +
The equations we discussed above can be coded rather easily (the +framework is similar to what we developed for logistic +regression). We are going to set up a simple case with two classes only and we want to find a line which separates them the best possible way.
-from sklearn.datasets import load_iris
-from sklearn.tree import DecisionTreeClassifier
-from sklearn.tree import export_text
-iris = load_iris()
-decision_tree = DecisionTreeClassifier(random_state=0, max_depth=2)
-decision_tree = decision_tree.fit(iris.data, iris.target)
-r = export_text(decision_tree, feature_names=iris['feature_names'])
-print(r)
+
-
Two algorithms stand out in the set up of decision trees:
-We discuss both algorithms with applications here. The popular library -Scikit-Learn uses the CART algorithm. For classification problems -you can use either the gini index or the entropy to split a tree -in two branches. +
There are however problems with this approach, although it looks +pretty straightforward to implement. When running the above code, we see that we can easily end up with many diffeent lines which separate the two classes. +
+ +For small +gaps between the entries, we may also end up needing many iterations +before the solutions converge and if the data cannot be separated +properly into two distinct classes, we may not experience a converge +at all.
@@ -427,9 +259,6 @@ in two branches.
-
For classification, the CART algorithm splits the data set in two subsets using a single feature \( k \) and a threshold \( t_k \). -This could be for example a threshold set by a number below a certain circumference of a malign tumor. +
A better approach is rather to try to define a large margin between +the two classes (if they are well separated from the beginning).
-How do we find these two quantities? -We search for the pair \( (k,t_k) \) that produces the purest subset using for example the gini factor \( G \). -The cost function it tries to minimize is then +
Thus, we wish to find a margin \( M \) with \( \boldsymbol{w} \) normalized to +\( \vert\vert \boldsymbol{w}\vert\vert =1 \) subject to the condition
+ $$ -C(k,t_k) = \frac{m_{\mathrm{left}}}{m}G_{\mathrm{left}}+ \frac{m_{\mathrm{right}}}{m}G_{\mathrm{right}}, +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq M \hspace{0.1cm}\forall i=1,2,\dots, p. $$ -where \( G_{\mathrm{left/right}} \) measures the impurity of the left/right subset and \( m_{\mathrm{left/right}} \) - is the number of instances in the left/right subset -
+All points are thus at a signed distance from the decision boundary defined by the line \( L \). The parameters \( b \) and \( w_1 \) and \( w_2 \) define this line.
-Once it has successfully split the training set in two, it splits the subsets using the same logic, then the subsubsets -and so on, recursively. It stops recursing once it reaches the maximum depth (defined by the -\( max\_depth \) hyperparameter), or if it cannot find a split that will reduce impurity. A few other -hyperparameters control additional stopping conditions such as the \( min\_samples\_split \), -\( min\_samples\_leaf \), \( min\_weight\_fraction\_leaf \), and \( max\_leaf\_nodes \). +
We seek thus the largest value \( M \) defined by
+$$ +\frac{1}{\vert \vert \boldsymbol{w}\vert\vert}y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq M \hspace{0.1cm}\forall i=1,2,\dots, n, +$$ + +or just
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq M\vert \vert \boldsymbol{w}\vert\vert \hspace{0.1cm}\forall i. +$$ + +If we scale the equation so that \( \vert \vert \boldsymbol{w}\vert\vert = 1/M \), we have to find the minimum of +\( \boldsymbol{w}^T\boldsymbol{w}=\vert \vert \boldsymbol{w}\vert\vert \) (the norm) subject to the condition +
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq 1 \hspace{0.1cm}\forall i. +$$ + +We have thus defined our margin as the invers of the norm of +\( \boldsymbol{w} \). We want to minimize the norm in order to have a as large as +possible margin \( M \). Before we proceed, we need to remind ourselves +about Lagrangian multipliers.
@@ -438,10 +284,6 @@ hyperparameters control additional stopping conditions such as the \( min\_sampl
-
The CART algorithm for regression works is similar to the one for classification except that instead of trying to split the -training set in a way that minimizes say the gini or entropy impurity, it now tries to split the training set in a way that minimizes our well-known mean-squared error (MSE). The cost function is now +
Consider a function of three independent variables \( f(x,y,z) \) . For the function \( f \) to be an +extreme we have
$$ -C(k,t_k) = \frac{m_{\mathrm{left}}}{m}\mathrm{MSE}_{\mathrm{left}}+ \frac{m_{\mathrm{right}}}{m}\mathrm{MSE}_{\mathrm{right}}. +df=0. $$ -Here the MSE for a specific node is defined as
+A necessary and sufficient condition is
$$ -\mathrm{MSE}_{\mathrm{node}}=\frac{1}{m_\mathrm{node}}\sum_{i\in \mathrm{node}}(\overline{y}_{\mathrm{node}}-y_i)^2, +\frac{\partial f}{\partial x} =\frac{\partial f}{\partial y}=\frac{\partial f}{\partial z}=0, $$ -with
+due to
$$ -\overline{y}_{\mathrm{node}}=\frac{1}{m_\mathrm{node}}\sum_{i\in \mathrm{node}}y_i, +df = \frac{\partial f}{\partial x}dx+\frac{\partial f}{\partial y}dy+\frac{\partial f}{\partial z}dz. $$ -the mean value of all observations in a specific node.
+In many problems the variables \( x,y,z \) are often subject to constraints (such as those above for the margin) +so that they are no longer all independent. It is possible at least in principle to use each +constraint to eliminate one variable +and to proceed with a new and smaller set of independent varables. +
-Without any regularization, the regression task for decision trees, -just like for classification tasks, is prone to overfitting. +
The use of so-called Lagrangian multipliers is an alternative technique when the elimination +of variables is incovenient or undesirable. Assume that we have an equation of constraint on +the variables \( x,y,z \) +
+$$ +\phi(x,y,z) = 0, +$$ + +resulting in
+$$ +d\phi = \frac{\partial \phi}{\partial x}dx+\frac{\partial \phi}{\partial y}dy+\frac{\partial \phi}{\partial z}dz =0. +$$ + +Now we cannot set anymore
+$$ +\frac{\partial f}{\partial x} =\frac{\partial f}{\partial y}=\frac{\partial f}{\partial z}=0, +$$ + +if \( df=0 \) is wanted +because there are now only two independent variables! Assume \( x \) and \( y \) are the independent +variables. +Then \( dz \) is no longer arbitrary.
@@ -437,11 +293,6 @@ just like for classification tasks, is prone to overfitting.
-
The example we will look at is a classical one in many Machine -Learning applications. Based on various meteorological features, we -have several so-called attributes which decide whether we at the end -will do some outdoor activity like skiing, going for a bike ride etc -etc. The table here contains the feautures outlook, temperature, -humidity and wind. The target or output is whether we ride -(True=1) or whether we do something else that day (False=0). The -attributes for each feature are then sunny, overcast and rain for the -outlook, hot, cold and mild for temperature, high and normal for -humidity and weak and strong for wind. +
However, we can add to
+$$ +df = \frac{\partial f}{\partial x}dx+\frac{\partial f}{\partial y}dy+\frac{\partial f}{\partial z}dz, +$$ + +a multiplum of \( d\phi \), viz. \( \lambda d\phi \), resulting in
+$$ +df+\lambda d\phi = (\frac{\partial f}{\partial z}+\lambda +\frac{\partial \phi}{\partial x})dx+(\frac{\partial f}{\partial y}+\lambda\frac{\partial \phi}{\partial y})dy+ +(\frac{\partial f}{\partial z}+\lambda\frac{\partial \phi}{\partial z})dz =0. +$$ + +Our multiplier is chosen so that
+$$ +\frac{\partial f}{\partial z}+\lambda\frac{\partial \phi}{\partial z} =0. +$$ + +We need to remember that we took \( dx \) and \( dy \) to be arbitrary and thus we must have
+$$ +\frac{\partial f}{\partial x}+\lambda\frac{\partial \phi}{\partial x} =0, +$$ + +and
+$$ +\frac{\partial f}{\partial y}+\lambda\frac{\partial \phi}{\partial y} =0. +$$ + +When all these equations are satisfied, \( df=0 \). We have four unknowns, \( x,y,z \) and +\( \lambda \). Actually we want only \( x,y,z \), \( \lambda \) needs not to be determined, +it is therefore often called +Lagrange's undetermined multiplier. +If we have a set of constraints \( \phi_k \) we have the equations
+$$ +\frac{\partial f}{\partial x_i}+\sum_k\lambda_k\frac{\partial \phi_k}{\partial x_i} =0. +$$ -The table here summarizes the various attributes and
-| Day | Outlook | Temperature | Humidity | Wind | Ride |
| 1 | Sunny | Hot | High | Weak | 0 |
| 2 | Sunny | Hot | High | Strong | 1 |
| 3 | Overcast | Hot | High | Weak | 1 |
| 4 | Rain | Mild | High | Weak | 1 |
| 5 | Rain | Cool | Normal | Weak | 1 |
| 6 | Rain | Cool | Normal | Strong | 0 |
| 7 | Overcast | Cool | Normal | Strong | 1 |
| 8 | Sunny | Mild | High | Weak | 0 |
| 9 | Sunny | Cool | Normal | Weak | 1 |
| 10 | Rain | Mild | Normal | Weak | 1 |
| 11 | Sunny | Mild | Normal | Strong | 1 |
| 12 | Overcast | Mild | High | Strong | 1 |
| 13 | Overcast | Hot | Normal | Weak | 1 |
| 14 | Rain | Mild | High | Strong | 0 |
@@ -452,12 +283,6 @@ humidity and weak and strong for wind.
-
In order to solve the above problem, we define the following Lagrangian function to be minimized
+$$ +{\cal L}(\lambda,b,\boldsymbol{w})=\frac{1}{2}\boldsymbol{w}^T\boldsymbol{w}-\sum_{i=1}^n\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)-1\right], +$$ +where \( \lambda_i \) is a so-called Lagrange multiplier subject to the condition \( \lambda_i \geq 0 \).
- -# Common imports
-import numpy as np
-import pandas as pd
-import matplotlib.pyplot as plt
-from sklearn.tree import DecisionTreeClassifier
-from sklearn.model_selection import train_test_split
-from sklearn.tree import export_graphviz
-from sklearn.preprocessing import StandardScaler, OneHotEncoder
-from sklearn.compose import ColumnTransformer
-from IPython.display import Image
-from pydot import graph_from_dot_data
-import os
+Taking the derivatives with respect to \( b \) and \( \boldsymbol{w} \) we obtain
+$$
+\frac{\partial {\cal L}}{\partial b} = -\sum_{i} \lambda_iy_i=0,
+$$
-# Where to save the figures and data files
-PROJECT_ROOT_DIR = "Results"
-FIGURE_ID = "Results/FigureFiles"
-DATA_ID = "DataFiles/"
+and
+$$
+\frac{\partial {\cal L}}{\partial \boldsymbol{w}} = 0 = \boldsymbol{w}-\sum_{i} \lambda_iy_i\boldsymbol{x}_i.
+$$
-if not os.path.exists(PROJECT_ROOT_DIR):
- os.mkdir(PROJECT_ROOT_DIR)
+Inserting these constraints into the equation for \( {\cal L} \) we obtain
+$$
+{\cal L}=\sum_i\lambda_i-\frac{1}{2}\sum_{ij}^n\lambda_i\lambda_jy_iy_j\boldsymbol{x}_i^T\boldsymbol{x}_j,
+$$
-if not os.path.exists(FIGURE_ID):
- os.makedirs(FIGURE_ID)
-
-if not os.path.exists(DATA_ID):
- os.makedirs(DATA_ID)
-
-def image_path(fig_id):
- return os.path.join(FIGURE_ID, fig_id)
-
-def data_path(dat_id):
- return os.path.join(DATA_ID, dat_id)
-
-def save_fig(fig_id):
- plt.savefig(image_path(fig_id) + ".png", format='png')
-
-infile = open(data_path("rideclass.csv"),'r')
-
-# Read the experimental data with Pandas
-from IPython.display import display
-ridedata = pd.read_csv(infile,names = ('Outlook','Temperature','Humidity','Wind','Ride'))
-ridedata = pd.DataFrame(ridedata)
-
-# Features and targets
-X = ridedata.loc[:, ridedata.columns != 'Ride'].values
-y = ridedata.loc[:, ridedata.columns == 'Ride'].values
-
-# Create the encoder.
-encoder = OneHotEncoder(handle_unknown="ignore")
-# Assume for simplicity all features are categorical.
-encoder.fit(X)
-# Apply the encoder.
-X = encoder.transform(X)
-print(X)
-# Then do a Classification tree
-tree_clf = DecisionTreeClassifier(max_depth=2)
-tree_clf.fit(X, y)
-print("Train set accuracy with Decision Tree: {:.2f}".format(tree_clf.score(X,y)))
-#transfer to a decision tree graph
-export_graphviz(
- tree_clf,
- out_file="DataFiles/ride.dot",
- rounded=True,
- filled=True
-)
-cmd = 'dot -Tpng DataFiles/cancer.dot -o DataFiles/cancer.png'
-os.system(cmd)
-
-subject to the constraints \( \lambda_i\geq 0 \) and \( \sum_i\lambda_iy_i=0 \). +We must in addition satisfy the Karush-Kuhn-Tucker (KKT) condition +
+$$ +\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) -1\right] \hspace{0.1cm}\forall i. +$$ +When \( \lambda_i > 0 \), the vectors \( \boldsymbol{x}_i \) are called support vectors. They are the vectors closest to the line (or hyperplane) and define the margin \( M \).
@@ -502,13 +278,6 @@ os.system(cmd)
-
The above functions (gini, entropy and misclassification error) are -important components of the so-called CART algorithm. We will discuss -this algorithm below after we have discussed the information gain -algorithm ID3. +
We can rewrite
+$$ +{\cal L}=\sum_i\lambda_i-\frac{1}{2}\sum_{ij}^n\lambda_i\lambda_jy_iy_j\boldsymbol{x}_i^T\boldsymbol{x}_j, +$$ + +and its constraints in terms of a matrix-vector problem where we minimize w.r.t. \( \lambda \) the following problem
+$$ +\frac{1}{2} \boldsymbol{\lambda}^T\begin{bmatrix} y_1y_1\boldsymbol{x}_1^T\boldsymbol{x}_1 & y_1y_2\boldsymbol{x}_1^T\boldsymbol{x}_2 & \dots & \dots & y_1y_n\boldsymbol{x}_1^T\boldsymbol{x}_n \\ +y_2y_1\boldsymbol{x}_2^T\boldsymbol{x}_1 & y_2y_2\boldsymbol{x}_2^T\boldsymbol{x}_2 & \dots & \dots & y_1y_n\boldsymbol{x}_2^T\boldsymbol{x}_n \\ +\dots & \dots & \dots & \dots & \dots \\ +\dots & \dots & \dots & \dots & \dots \\ +y_ny_1\boldsymbol{x}_n^T\boldsymbol{x}_1 & y_ny_2\boldsymbol{x}_n^T\boldsymbol{x}_2 & \dots & \dots & y_ny_n\boldsymbol{x}_n^T\boldsymbol{x}_n \\ +\end{bmatrix}\boldsymbol{\lambda}-\mathbb{1}\boldsymbol{\lambda}, +$$ + +subject to \( \boldsymbol{y}^T\boldsymbol{\lambda}=0 \). Here we defined the vectors \( \boldsymbol{\lambda} =[\lambda_1,\lambda_2,\dots,\lambda_n] \) and +\( \boldsymbol{y}=[y_1,y_2,\dots,y_n] \).
-In the example here we have converted all our attributes into numerical values \( 0,1,2 \) etc.
- - - -# Split a dataset based on an attribute and an attribute value
-def test_split(index, value, dataset):
- left, right = list(), list()
- for row in dataset:
- if row[index] < value:
- left.append(row)
- else:
- right.append(row)
- return left, right
-
-# Calculate the Gini index for a split dataset
-def gini_index(groups, classes):
- # count all samples at split point
- n_instances = float(sum([len(group) for group in groups]))
- # sum weighted Gini index for each group
- gini = 0.0
- for group in groups:
- size = float(len(group))
- # avoid divide by zero
- if size == 0:
- continue
- score = 0.0
- # score the group based on the score for each class
- for class_val in classes:
- p = [row[-1] for row in group].count(class_val) / size
- score += p * p
- # weight the group score by its relative size
- gini += (1.0 - score) * (size / n_instances)
- return gini
-
-# Select the best split point for a dataset
-def get_split(dataset):
- class_values = list(set(row[-1] for row in dataset))
- b_index, b_value, b_score, b_groups = 999, 999, 999, None
- for index in range(len(dataset[0])-1):
- for row in dataset:
- groups = test_split(index, row[index], dataset)
- gini = gini_index(groups, class_values)
- print('X%d < %.3f Gini=%.3f' % ((index+1), row[index], gini))
- if gini < b_score:
- b_index, b_value, b_score, b_groups = index, row[index], gini, groups
- return {'index':b_index, 'value':b_value, 'groups':b_groups}
-
-dataset = [[0,0,0,0,0],
- [0,0,0,1,1],
- [1,0,0,0,1],
- [2,1,0,0,1],
- [2,2,1,0,1],
- [2,2,1,1,0],
- [1,2,1,1,1],
- [0,1,0,0,0],
- [0,2,1,0,1],
- [2,1,1,0,1],
- [0,1,1,1,1],
- [1,1,0,1,1],
- [1,0,1,0,1],
- [2,1,0,1,0]]
-
-split = get_split(dataset)
-print('Split: [X%d < %.3f]' % ((split['index']+1), split['value']))
-
-diff --git a/doc/pub/week45/html/._week45-bs037.html b/doc/pub/week45/html/._week45-bs037.html index 8e9de3ab8..86a82c657 100644 --- a/doc/pub/week45/html/._week45-bs037.html +++ b/doc/pub/week45/html/._week45-bs037.html @@ -37,175 +37,10 @@ doconce format html week45.do.txt --html_style=bootstrap --pygments_html_style=d @@ -303,83 +174,44 @@ MathJax.Hub.Config({ Contents @@ -391,37 +223,36 @@ MathJax.Hub.Config({
-
The ID3 algorithm learns decision trees by constructing -them in a top down way, beginning with the question which attribute should be tested at the root of the tree? +
Solving the above problem, yields the values of \( \lambda_i \). +To find the coefficients of your hyperplane we need simply to compute
+$$ +\boldsymbol{w}=\sum_{i} \lambda_iy_i\boldsymbol{x}_i. +$$ -The ID3 algorithm selects which attribute to test at each node in the -tree. -
+With our vector \( \boldsymbol{w} \) we can in turn find the value of the intercept \( b \) (here in two dimensions) via
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1, +$$ -We would like to select the attribute that is most useful for classifying -examples. -
+resulting in
+$$ +b = \frac{1}{y_i}-\boldsymbol{w}^T\boldsymbol{x}_i, +$$ -What is a good quantitative measure of the worth of an attribute?
+or if we write it out in terms of the support vectors only, with \( N_s \) being their number, we have
+$$ +b = \frac{1}{N_s}\sum_{j\in N_s}\left(y_j-\sum_{i=1}^n\lambda_iy_i\boldsymbol{x}_i^T\boldsymbol{x}_j\right). +$$ -Information gain measures how well a given attribute separates the -training examples according to their target classification. -
+With our hyperplane coefficients we can use our classifier to assign any observation by simply using
+$$ +y_i = \mathrm{sign}(\boldsymbol{w}^T\boldsymbol{x}_i+b). +$$ -The ID3 algorithm uses this information gain measure to select among the candidate -attributes at each step while growing the tree. -
+Below we discuss how to find the optimal values of \( \lambda_i \). Before we proceed however, we discuss now the so-called soft classifier.
@@ -440,15 +271,6 @@ attributes at each step while growing the tree.
-
import matplotlib.pyplot as plt
-import numpy as np
-from sklearn.model_selection import train_test_split
-from sklearn.datasets import load_breast_cancer
-from sklearn.svm import SVC
-from sklearn.linear_model import LogisticRegression
-from sklearn.tree import DecisionTreeClassifier
+Till now, the margin is strictly defined by the support vectors. This defines what is called a hard classifier, that is the margins are well defined.
-# Load the data
-cancer = load_breast_cancer()
+Suppose now that classes overlap in feature space, as shown in the
+figure here. One way to deal with this problem before we define the
+so-called kernel approach, is to allow a kind of slack in the sense
+that we allow some points to be on the wrong side of the margin.
+
-X_train, X_test, y_train, y_test = train_test_split(cancer.data,cancer.target,random_state=0)
-print(X_train.shape)
-print(X_test.shape)
-# Logistic Regression
-logreg = LogisticRegression(solver='lbfgs')
-logreg.fit(X_train, y_train)
-print("Test set accuracy with Logistic Regression: {:.2f}".format(logreg.score(X_test,y_test)))
-# Support vector machine
-svm = SVC(gamma='auto', C=100)
-svm.fit(X_train, y_train)
-print("Test set accuracy with SVM: {:.2f}".format(svm.score(X_test,y_test)))
-# Decision Trees
-deep_tree_clf = DecisionTreeClassifier(max_depth=None)
-deep_tree_clf.fit(X_train, y_train)
-print("Test set accuracy with Decision Trees: {:.2f}".format(deep_tree_clf.score(X_test,y_test)))
-#now scale the data
-from sklearn.preprocessing import StandardScaler
-scaler = StandardScaler()
-scaler.fit(X_train)
-X_train_scaled = scaler.transform(X_train)
-X_test_scaled = scaler.transform(X_test)
-# Logistic Regression
-logreg.fit(X_train_scaled, y_train)
-print("Test set accuracy Logistic Regression with scaled data: {:.2f}".format(logreg.score(X_test_scaled,y_test)))
-# Support Vector Machine
-svm.fit(X_train_scaled, y_train)
-print("Test set accuracy SVM with scaled data: {:.2f}".format(logreg.score(X_test_scaled,y_test)))
-# Decision Trees
-deep_tree_clf.fit(X_train_scaled, y_train)
-print("Test set accuracy with Decision Trees and scaled data: {:.2f}".format(deep_tree_clf.score(X_test_scaled,y_test)))
-
-We introduce thus the so-called slack variables \( \boldsymbol{\xi} =[\xi_1,x_2,\dots,x_n] \) and +modify our previous equation +
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1, +$$ +to
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1-\xi_i, +$$ + +with the requirement \( \xi_i\geq 0 \). The total violation is now \( \sum_i\xi \). +The value \( \xi_i \) in the constraint the last constraint corresponds to the amount by which the prediction +\( y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1 \) is on the wrong side of its margin. Hence by bounding the sum \( \sum_i \xi_i \), +we bound the total amount by which predictions fall on the wrong side of their margins. +
+ +Misclassifications occur when \( \xi_i > 1 \). Thus bounding the total sum by some value \( C \) bounds in turn the total number of +misclassifications. +
@@ -472,16 +271,6 @@ deep_tree_clf.fit(X_train_scaled, y_train)
-
from __future__ import division, print_function, unicode_literals
+This has in turn the consequences that we change our optmization problem to finding the minimum of
+$$
+{\cal L}=\frac{1}{2}\boldsymbol{w}^T\boldsymbol{w}-\sum_{i=1}^n\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)-(1-\xi_)\right]+C\sum_{i=1}^n\xi_i-\sum_{i=1}^n\gamma_i\xi_i,
+$$
-# Common imports
-import numpy as np
-import os
+subject to
+$$
+y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1-\xi_i \hspace{0.1cm}\forall i,
+$$
-# to make this notebook's output stable across runs
-np.random.seed(42)
+with the requirement \( \xi_i\geq 0 \).
-# To plot pretty figures
-import matplotlib
-import matplotlib.pyplot as plt
-from matplotlib.colors import ListedColormap
-plt.rcParams['axes.labelsize'] = 14
-plt.rcParams['xtick.labelsize'] = 12
-plt.rcParams['ytick.labelsize'] = 12
+Taking the derivatives with respect to \( b \) and \( \boldsymbol{w} \) we obtain
+$$
+\frac{\partial {\cal L}}{\partial b} = -\sum_{i} \lambda_iy_i=0,
+$$
+and
+$$
+\frac{\partial {\cal L}}{\partial \boldsymbol{w}} = 0 = \boldsymbol{w}-\sum_{i} \lambda_iy_i\boldsymbol{x}_i,
+$$
-from sklearn.svm import SVC
-from sklearn import datasets
-from sklearn.tree import DecisionTreeClassifier
-from sklearn.datasets import make_moons
-from sklearn.tree import export_graphviz
+and
+$$
+\lambda_i = C-\gamma_i \hspace{0.1cm}\forall i.
+$$
-Xm, ym = make_moons(n_samples=100, noise=0.25, random_state=53)
+Inserting these constraints into the equation for \( {\cal L} \) we obtain the same equation as before
+$$
+{\cal L}=\sum_i\lambda_i-\frac{1}{2}\sum_{ij}^n\lambda_i\lambda_jy_iy_j\boldsymbol{x}_i^T\boldsymbol{x}_j,
+$$
-deep_tree_clf1 = DecisionTreeClassifier(random_state=42)
-deep_tree_clf2 = DecisionTreeClassifier(min_samples_leaf=4, random_state=42)
-deep_tree_clf1.fit(Xm, ym)
-deep_tree_clf2.fit(Xm, ym)
+but now subject to the constraints \( \lambda_i\geq 0 \), \( \sum_i\lambda_iy_i=0 \) and \( 0\leq\lambda_i \leq C \).
+We must in addition satisfy the Karush-Kuhn-Tucker condition which now reads
+
+$$
+\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) -(1-\xi_)\right]=0 \hspace{0.1cm}\forall i,
+$$
+$$
+\gamma_i\xi_i = 0,
+$$
-def plot_decision_boundary(clf, X, y, axes=[0, 7.5, 0, 3], iris=True, legend=False, plot_training=True):
- x1s = np.linspace(axes[0], axes[1], 100)
- x2s = np.linspace(axes[2], axes[3], 100)
- x1, x2 = np.meshgrid(x1s, x2s)
- X_new = np.c_[x1.ravel(), x2.ravel()]
- y_pred = clf.predict(X_new).reshape(x1.shape)
- custom_cmap = ListedColormap(['#fafab0','#9898ff','#a0faa0'])
- plt.contourf(x1, x2, y_pred, alpha=0.3, cmap=custom_cmap)
- if not iris:
- custom_cmap2 = ListedColormap(['#7d7d58','#4c4c7f','#507d50'])
- plt.contour(x1, x2, y_pred, cmap=custom_cmap2, alpha=0.8)
- if plot_training:
- plt.plot(X[:, 0][y==0], X[:, 1][y==0], "yo", label="Iris-Setosa")
- plt.plot(X[:, 0][y==1], X[:, 1][y==1], "bs", label="Iris-Versicolor")
- plt.plot(X[:, 0][y==2], X[:, 1][y==2], "g^", label="Iris-Virginica")
- plt.axis(axes)
- if iris:
- plt.xlabel("Petal length", fontsize=14)
- plt.ylabel("Petal width", fontsize=14)
- else:
- plt.xlabel(r"$x_1$", fontsize=18)
- plt.ylabel(r"$x_2$", fontsize=18, rotation=0)
- if legend:
- plt.legend(loc="lower right", fontsize=14)
-plt.figure(figsize=(11, 4))
-plt.subplot(121)
-plot_decision_boundary(deep_tree_clf1, Xm, ym, axes=[-1.5, 2.5, -1, 1.5], iris=False)
-plt.title("No restrictions", fontsize=16)
-plt.subplot(122)
-plot_decision_boundary(deep_tree_clf2, Xm, ym, axes=[-1.5, 2.5, -1, 1.5], iris=False)
-plt.title("min_samples_leaf = {}".format(deep_tree_clf2.min_samples_leaf), fontsize=14)
-plt.show()
-
-and
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) -(1-\xi_) \geq 0 \hspace{0.1cm}\forall i. +$$@@ -494,18 +288,6 @@ plt.show()
A Support Vector Machine (SVM) is a very powerful and versatile +Machine Learning method, capable of performing linear or nonlinear +classification, regression, and even outlier detection. It is one of +the most popular models in Machine Learning, and anyone interested in +Machine Learning should have it in their toolbox. SVMs are +particularly well suited for classification of complex but small-sized or +medium-sized datasets. +
+ +The case with two well-separated classes only can be understood in an +intuitive way in terms of lines in a two-dimensional space separating +the two classes (see figure below). +
+ +The basic mathematics behind the SVM is however less familiar to most of us. +It relies on the definition of hyperplanes and the +definition of a margin which separates classes (in case of +classification problems) of variables. It is also used for regression +problems. +
+ +With SVMs we distinguish between hard margin and soft margins. The +latter introduces a so-called softening parameter to be discussed +below. We distinguish also between linear and non-linear +approaches. The latter are the most frequent ones since it is rather +unlikely that we can separate classes easily by say straight lines. +
+The theory behind support vector machines (SVM hereafter) is based on +the mathematical description of so-called hyperplanes. Let us start +with a two-dimensional case. This will also allow us to introduce our +first SVM examples. These will be tailored to the case of two specific +classes, as displayed in the figure here based on the usage of the petal data. +
+ +We assume here that our data set can be well separated into two +domains, where a straight line does the job in the separating the two +classes. Here the two classes are represented by either squares or +circles. +
+ + +from sklearn import datasets
+from sklearn.svm import SVC, LinearSVC
+from sklearn.linear_model import SGDClassifier
+from sklearn.preprocessing import StandardScaler
+import matplotlib
+import matplotlib.pyplot as plt
+plt.rcParams['axes.labelsize'] = 14
+plt.rcParams['xtick.labelsize'] = 12
+plt.rcParams['ytick.labelsize'] = 12
+
+
+iris = datasets.load_iris()
+X = iris["data"][:, (2, 3)] # petal length, petal width
+y = iris["target"]
+
+setosa_or_versicolor = (y == 0) | (y == 1)
+X = X[setosa_or_versicolor]
+y = y[setosa_or_versicolor]
+
+
+
+C = 5
+alpha = 1 / (C * len(X))
+
+lin_clf = LinearSVC(loss="hinge", C=C, random_state=42)
+svm_clf = SVC(kernel="linear", C=C)
+sgd_clf = SGDClassifier(loss="hinge", learning_rate="constant", eta0=0.001, alpha=alpha,
+ max_iter=100000, random_state=42)
+
+scaler = StandardScaler()
+X_scaled = scaler.fit_transform(X)
+
+lin_clf.fit(X_scaled, y)
+svm_clf.fit(X_scaled, y)
+sgd_clf.fit(X_scaled, y)
+
+print("LinearSVC: ", lin_clf.intercept_, lin_clf.coef_)
+print("SVC: ", svm_clf.intercept_, svm_clf.coef_)
+print("SGDClassifier(alpha={:.5f}):".format(sgd_clf.alpha), sgd_clf.intercept_, sgd_clf.coef_)
+
+# Compute the slope and bias of each decision boundary
+w1 = -lin_clf.coef_[0, 0]/lin_clf.coef_[0, 1]
+b1 = -lin_clf.intercept_[0]/lin_clf.coef_[0, 1]
+w2 = -svm_clf.coef_[0, 0]/svm_clf.coef_[0, 1]
+b2 = -svm_clf.intercept_[0]/svm_clf.coef_[0, 1]
+w3 = -sgd_clf.coef_[0, 0]/sgd_clf.coef_[0, 1]
+b3 = -sgd_clf.intercept_[0]/sgd_clf.coef_[0, 1]
+
+# Transform the decision boundary lines back to the original scale
+line1 = scaler.inverse_transform([[-10, -10 * w1 + b1], [10, 10 * w1 + b1]])
+line2 = scaler.inverse_transform([[-10, -10 * w2 + b2], [10, 10 * w2 + b2]])
+line3 = scaler.inverse_transform([[-10, -10 * w3 + b3], [10, 10 * w3 + b3]])
+
+# Plot all three decision boundaries
+plt.figure(figsize=(11, 4))
+plt.plot(line1[:, 0], line1[:, 1], "k:", label="LinearSVC")
+plt.plot(line2[:, 0], line2[:, 1], "b--", linewidth=2, label="SVC")
+plt.plot(line3[:, 0], line3[:, 1], "r-", label="SGDClassifier")
+plt.plot(X[:, 0][y==1], X[:, 1][y==1], "bs") # label="Iris-Versicolor"
+plt.plot(X[:, 0][y==0], X[:, 1][y==0], "yo") # label="Iris-Setosa"
+plt.xlabel("Petal length", fontsize=14)
+plt.ylabel("Petal width", fontsize=14)
+plt.legend(loc="upper center", fontsize=14)
+plt.axis([0, 5.5, 0, 2])
+
+plt.show()
+
+The aim of the SVM algorithm is to find a hyperplane in a +\( p \)-dimensional space, where \( p \) is the number of features that +distinctly classifies the data points. +
+ +In a \( p \)-dimensional space, a hyperplane is what we call an affine subspace of dimension of \( p-1 \). +As an example, in two dimension, a hyperplane is simply as straight line while in three dimensions it is +a two-dimensional subspace, or stated simply, a plane. +
+ +In two dimensions, with the variables \( x_1 \) and \( x_2 \), the hyperplane is defined as
+
+$$
+b+w_1x_1+w_2x_2=0,
+$$
+
+
+
where \( b \) is the intercept and \( w_1 \) and \( w_2 \) define the elements of a vector orthogonal to the line +\( b+w_1x_1+w_2x_2=0 \). +In two dimensions we define the vectors \( \boldsymbol{x} =[x1,x2] \) and \( \boldsymbol{w}=[w1,w2] \). +We can then rewrite the above equation as +
+ +
+$$
+\boldsymbol{x}^T\boldsymbol{w}+b=0.
+$$
+
+
We limit ourselves to two classes of outputs \( y_i \) and assign these classes the values \( y_i = \pm 1 \). +In a \( p \)-dimensional space of say \( p \) features we have a hyperplane defines as +
+
+$$
+b+wx_1+w_2x_2+\dots +w_px_p=0.
+$$
+
+
+
If we define a +matrix \( \boldsymbol{X}=\left[\boldsymbol{x}_1,\boldsymbol{x}_2,\dots, \boldsymbol{x}_p\right] \) +of dimension \( n\times p \), where \( n \) represents the observations for each feature and each vector \( x_i \) is a column vector of the matrix \( \boldsymbol{X} \), +
+
+$$
+\boldsymbol{x}_i = \begin{bmatrix} x_{i1} \\ x_{i2} \\ \dots \\ \dots \\ x_{ip} \end{bmatrix}.
+$$
+
+
+
If the above condition is not met for a given vector \( \boldsymbol{x}_i \) we have
+
+$$
+b+w_1x_{i1}+w_2x_{i2}+\dots +w_px_{ip} >0,
+$$
+
+
+
if our output \( y_i=1 \). +In this case we say that \( \boldsymbol{x}_i \) lies on one of the sides of the hyperplane and if +
+
+$$
+b+w_1x_{i1}+w_2x_{i2}+\dots +w_px_{ip} < 0,
+$$
+
+
+
for the class of observations \( y_i=-1 \), +then \( \boldsymbol{x}_i \) lies on the other side. +
+ +Equivalently, for the two classes of observations we have
+
+$$
+y_i\left(b+w_1x_{i1}+w_2x_{i2}+\dots +w_px_{ip}\right) > 0.
+$$
+
+
+
When we try to separate hyperplanes, if it exists, we can use it to construct a natural classifier: a test observation is assigned a given class depending on which side of the hyperplane it is located.
+Let us try to develop our intuition about SVMs by limiting ourselves to a two-dimensional +plane. To separate the two classes of data points, there are many +possible lines (hyperplanes if you prefer a more strict naming) +that could be chosen. Our objective is to find a +plane that has the maximum margin, i.e the maximum distance between +data points of both classes. Maximizing the margin distance provides +some reinforcement so that future data points can be classified with +more confidence. +
+ +What a linear classifier attempts to accomplish is to split the +feature space into two half spaces by placing a hyperplane between the +data points. This hyperplane will be our decision boundary. All +points on one side of the plane will belong to class one and all points +on the other side of the plane will belong to the second class two. +
+ +Unfortunately there are many ways in which we can place a hyperplane +to divide the data. Below is an example of two candidate hyperplanes +for our data sample. +
+Let us define the function
+
+$$
+f(x) = \boldsymbol{w}^T\boldsymbol{x}+b = 0,
+$$
+
+
+
as the function that determines the line \( L \) that separates two classes (our two features), see the figure here.
+ +Any point defined by \( \boldsymbol{x}_i \) and \( \boldsymbol{x}_2 \) on the line \( L \) will satisfy \( \boldsymbol{w}^T(\boldsymbol{x}_1-\boldsymbol{x}_2)=0 \).
+ +The signed distance \( \delta \) from any point defined by a vector \( \boldsymbol{x} \) and a point \( \boldsymbol{x}_0 \) on the line \( L \) is then
+
+$$
+\delta = \frac{1}{\vert\vert \boldsymbol{w}\vert\vert}(\boldsymbol{w}^T\boldsymbol{x}+b).
+$$
+
+
How do we find the parameter \( b \) and the vector \( \boldsymbol{w} \)? What we could +do is to define a cost function which now contains the set of all +misclassified points \( M \) and attempt to minimize this function +
+ +
+$$
+C(\boldsymbol{w},b) = -\sum_{i\in M} y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b).
+$$
+
+
+
We could now for example define all values \( y_i =1 \) as misclassified in case we have \( \boldsymbol{w}^T\boldsymbol{x}_i+b < 0 \) and the opposite if we have \( y_i=-1 \). Taking the derivatives gives us
+
+$$
+\frac{\partial C}{\partial b} = -\sum_{i\in M} y_i,
+$$
+
+
+
and
+
+$$
+\frac{\partial C}{\partial \boldsymbol{w}} = -\sum_{i\in M} y_ix_i.
+$$
+
+
We can now use the Newton-Raphson method or different variants of the gradient descent family (from plain gradient descent to various stochastic gradient descent approaches) to solve the equations
+
+$$
+b \leftarrow b +\eta \frac{\partial C}{\partial b},
+$$
+
+
+
and
+
+$$
+\boldsymbol{w} \leftarrow \boldsymbol{w} +\eta \frac{\partial C}{\partial \boldsymbol{w}},
+$$
+
+
+
where \( \eta \) is our by now well-known learning rate.
+The equations we discussed above can be coded rather easily (the +framework is similar to what we developed for logistic +regression). We are going to set up a simple case with two classes only and we want to find a line which separates them the best possible way. +
+ + +
+
+There are however problems with this approach, although it looks +pretty straightforward to implement. When running the above code, we see that we can easily end up with many diffeent lines which separate the two classes. +
+ +For small +gaps between the entries, we may also end up needing many iterations +before the solutions converge and if the data cannot be separated +properly into two distinct classes, we may not experience a converge +at all. +
+A better approach is rather to try to define a large margin between +the two classes (if they are well separated from the beginning). +
+ +Thus, we wish to find a margin \( M \) with \( \boldsymbol{w} \) normalized to +\( \vert\vert \boldsymbol{w}\vert\vert =1 \) subject to the condition +
+ +
+$$
+y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq M \hspace{0.1cm}\forall i=1,2,\dots, p.
+$$
+
+
+
All points are thus at a signed distance from the decision boundary defined by the line \( L \). The parameters \( b \) and \( w_1 \) and \( w_2 \) define this line.
+ +We seek thus the largest value \( M \) defined by
+
+$$
+\frac{1}{\vert \vert \boldsymbol{w}\vert\vert}y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq M \hspace{0.1cm}\forall i=1,2,\dots, n,
+$$
+
+
+
or just
+
+$$
+y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq M\vert \vert \boldsymbol{w}\vert\vert \hspace{0.1cm}\forall i.
+$$
+
+
+
If we scale the equation so that \( \vert \vert \boldsymbol{w}\vert\vert = 1/M \), we have to find the minimum of +\( \boldsymbol{w}^T\boldsymbol{w}=\vert \vert \boldsymbol{w}\vert\vert \) (the norm) subject to the condition +
+
+$$
+y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq 1 \hspace{0.1cm}\forall i.
+$$
+
+
+
We have thus defined our margin as the invers of the norm of +\( \boldsymbol{w} \). We want to minimize the norm in order to have a as large as +possible margin \( M \). Before we proceed, we need to remind ourselves +about Lagrangian multipliers. +
+Consider a function of three independent variables \( f(x,y,z) \) . For the function \( f \) to be an +extreme we have +
+
+$$
+df=0.
+$$
+
+
+
A necessary and sufficient condition is
+
+$$
+\frac{\partial f}{\partial x} =\frac{\partial f}{\partial y}=\frac{\partial f}{\partial z}=0,
+$$
+
+
+
due to
+
+$$
+df = \frac{\partial f}{\partial x}dx+\frac{\partial f}{\partial y}dy+\frac{\partial f}{\partial z}dz.
+$$
+
+
+
In many problems the variables \( x,y,z \) are often subject to constraints (such as those above for the margin) +so that they are no longer all independent. It is possible at least in principle to use each +constraint to eliminate one variable +and to proceed with a new and smaller set of independent varables. +
+ +The use of so-called Lagrangian multipliers is an alternative technique when the elimination +of variables is incovenient or undesirable. Assume that we have an equation of constraint on +the variables \( x,y,z \) +
+
+$$
+\phi(x,y,z) = 0,
+$$
+
+
+
resulting in
+
+$$
+d\phi = \frac{\partial \phi}{\partial x}dx+\frac{\partial \phi}{\partial y}dy+\frac{\partial \phi}{\partial z}dz =0.
+$$
+
+
+
Now we cannot set anymore
+
+$$
+\frac{\partial f}{\partial x} =\frac{\partial f}{\partial y}=\frac{\partial f}{\partial z}=0,
+$$
+
+
+
if \( df=0 \) is wanted +because there are now only two independent variables! Assume \( x \) and \( y \) are the independent +variables. +Then \( dz \) is no longer arbitrary. +
+However, we can add to
+
+$$
+df = \frac{\partial f}{\partial x}dx+\frac{\partial f}{\partial y}dy+\frac{\partial f}{\partial z}dz,
+$$
+
+
+
a multiplum of \( d\phi \), viz. \( \lambda d\phi \), resulting in
+
+$$
+df+\lambda d\phi = (\frac{\partial f}{\partial z}+\lambda
+\frac{\partial \phi}{\partial x})dx+(\frac{\partial f}{\partial y}+\lambda\frac{\partial \phi}{\partial y})dy+
+(\frac{\partial f}{\partial z}+\lambda\frac{\partial \phi}{\partial z})dz =0.
+$$
+
+
+
Our multiplier is chosen so that
+
+$$
+\frac{\partial f}{\partial z}+\lambda\frac{\partial \phi}{\partial z} =0.
+$$
+
+
+
We need to remember that we took \( dx \) and \( dy \) to be arbitrary and thus we must have
+
+$$
+\frac{\partial f}{\partial x}+\lambda\frac{\partial \phi}{\partial x} =0,
+$$
+
+
+
and
+
+$$
+\frac{\partial f}{\partial y}+\lambda\frac{\partial \phi}{\partial y} =0.
+$$
+
+
+
When all these equations are satisfied, \( df=0 \). We have four unknowns, \( x,y,z \) and +\( \lambda \). Actually we want only \( x,y,z \), \( \lambda \) needs not to be determined, +it is therefore often called +Lagrange's undetermined multiplier. +If we have a set of constraints \( \phi_k \) we have the equations +
+
+$$
+\frac{\partial f}{\partial x_i}+\sum_k\lambda_k\frac{\partial \phi_k}{\partial x_i} =0.
+$$
+
+
In order to solve the above problem, we define the following Lagrangian function to be minimized
+
+$$
+{\cal L}(\lambda,b,\boldsymbol{w})=\frac{1}{2}\boldsymbol{w}^T\boldsymbol{w}-\sum_{i=1}^n\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)-1\right],
+$$
+
+
+
where \( \lambda_i \) is a so-called Lagrange multiplier subject to the condition \( \lambda_i \geq 0 \).
+ +Taking the derivatives with respect to \( b \) and \( \boldsymbol{w} \) we obtain
+
+$$
+\frac{\partial {\cal L}}{\partial b} = -\sum_{i} \lambda_iy_i=0,
+$$
+
+
+
and
+
+$$
+\frac{\partial {\cal L}}{\partial \boldsymbol{w}} = 0 = \boldsymbol{w}-\sum_{i} \lambda_iy_i\boldsymbol{x}_i.
+$$
+
+
+
Inserting these constraints into the equation for \( {\cal L} \) we obtain
+
+$$
+{\cal L}=\sum_i\lambda_i-\frac{1}{2}\sum_{ij}^n\lambda_i\lambda_jy_iy_j\boldsymbol{x}_i^T\boldsymbol{x}_j,
+$$
+
+
+
subject to the constraints \( \lambda_i\geq 0 \) and \( \sum_i\lambda_iy_i=0 \). +We must in addition satisfy the Karush-Kuhn-Tucker (KKT) condition +
+
+$$
+\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) -1\right] \hspace{0.1cm}\forall i.
+$$
+
+
+
+
When \( \lambda_i > 0 \), the vectors \( \boldsymbol{x}_i \) are called support vectors. They are the vectors closest to the line (or hyperplane) and define the margin \( M \).
+We can rewrite
+
+$$
+{\cal L}=\sum_i\lambda_i-\frac{1}{2}\sum_{ij}^n\lambda_i\lambda_jy_iy_j\boldsymbol{x}_i^T\boldsymbol{x}_j,
+$$
+
+
+
and its constraints in terms of a matrix-vector problem where we minimize w.r.t. \( \lambda \) the following problem
+
+$$
+\frac{1}{2} \boldsymbol{\lambda}^T\begin{bmatrix} y_1y_1\boldsymbol{x}_1^T\boldsymbol{x}_1 & y_1y_2\boldsymbol{x}_1^T\boldsymbol{x}_2 & \dots & \dots & y_1y_n\boldsymbol{x}_1^T\boldsymbol{x}_n \\
+y_2y_1\boldsymbol{x}_2^T\boldsymbol{x}_1 & y_2y_2\boldsymbol{x}_2^T\boldsymbol{x}_2 & \dots & \dots & y_1y_n\boldsymbol{x}_2^T\boldsymbol{x}_n \\
+\dots & \dots & \dots & \dots & \dots \\
+\dots & \dots & \dots & \dots & \dots \\
+y_ny_1\boldsymbol{x}_n^T\boldsymbol{x}_1 & y_ny_2\boldsymbol{x}_n^T\boldsymbol{x}_2 & \dots & \dots & y_ny_n\boldsymbol{x}_n^T\boldsymbol{x}_n \\
+\end{bmatrix}\boldsymbol{\lambda}-\mathbb{1}\boldsymbol{\lambda},
+$$
+
+
+
subject to \( \boldsymbol{y}^T\boldsymbol{\lambda}=0 \). Here we defined the vectors \( \boldsymbol{\lambda} =[\lambda_1,\lambda_2,\dots,\lambda_n] \) and +\( \boldsymbol{y}=[y_1,y_2,\dots,y_n] \). +
+Solving the above problem, yields the values of \( \lambda_i \). +To find the coefficients of your hyperplane we need simply to compute +
+
+$$
+\boldsymbol{w}=\sum_{i} \lambda_iy_i\boldsymbol{x}_i.
+$$
+
+
+
With our vector \( \boldsymbol{w} \) we can in turn find the value of the intercept \( b \) (here in two dimensions) via
+
+$$
+y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1,
+$$
+
+
+
resulting in
+
+$$
+b = \frac{1}{y_i}-\boldsymbol{w}^T\boldsymbol{x}_i,
+$$
+
+
+
or if we write it out in terms of the support vectors only, with \( N_s \) being their number, we have
+
+$$
+b = \frac{1}{N_s}\sum_{j\in N_s}\left(y_j-\sum_{i=1}^n\lambda_iy_i\boldsymbol{x}_i^T\boldsymbol{x}_j\right).
+$$
+
+
+
With our hyperplane coefficients we can use our classifier to assign any observation by simply using
+
+$$
+y_i = \mathrm{sign}(\boldsymbol{w}^T\boldsymbol{x}_i+b).
+$$
+
+
+
Below we discuss how to find the optimal values of \( \lambda_i \). Before we proceed however, we discuss now the so-called soft classifier.
+Till now, the margin is strictly defined by the support vectors. This defines what is called a hard classifier, that is the margins are well defined.
+ +Suppose now that classes overlap in feature space, as shown in the +figure here. One way to deal with this problem before we define the +so-called kernel approach, is to allow a kind of slack in the sense +that we allow some points to be on the wrong side of the margin. +
+ +We introduce thus the so-called slack variables \( \boldsymbol{\xi} =[\xi_1,x_2,\dots,x_n] \) and +modify our previous equation +
+
+$$
+y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1,
+$$
+
+
+
to
+
+$$
+y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1-\xi_i,
+$$
+
+
+
with the requirement \( \xi_i\geq 0 \). The total violation is now \( \sum_i\xi \). +The value \( \xi_i \) in the constraint the last constraint corresponds to the amount by which the prediction +\( y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1 \) is on the wrong side of its margin. Hence by bounding the sum \( \sum_i \xi_i \), +we bound the total amount by which predictions fall on the wrong side of their margins. +
+ +Misclassifications occur when \( \xi_i > 1 \). Thus bounding the total sum by some value \( C \) bounds in turn the total number of +misclassifications. +
+This has in turn the consequences that we change our optmization problem to finding the minimum of
+
+$$
+{\cal L}=\frac{1}{2}\boldsymbol{w}^T\boldsymbol{w}-\sum_{i=1}^n\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)-(1-\xi_)\right]+C\sum_{i=1}^n\xi_i-\sum_{i=1}^n\gamma_i\xi_i,
+$$
+
+
+
subject to
+
+$$
+y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1-\xi_i \hspace{0.1cm}\forall i,
+$$
+
+
+
with the requirement \( \xi_i\geq 0 \).
+ +Taking the derivatives with respect to \( b \) and \( \boldsymbol{w} \) we obtain
+
+$$
+\frac{\partial {\cal L}}{\partial b} = -\sum_{i} \lambda_iy_i=0,
+$$
+
+
+
and
+
+$$
+\frac{\partial {\cal L}}{\partial \boldsymbol{w}} = 0 = \boldsymbol{w}-\sum_{i} \lambda_iy_i\boldsymbol{x}_i,
+$$
+
+
+
and
+
+$$
+\lambda_i = C-\gamma_i \hspace{0.1cm}\forall i.
+$$
+
+
+
Inserting these constraints into the equation for \( {\cal L} \) we obtain the same equation as before
+
+$$
+{\cal L}=\sum_i\lambda_i-\frac{1}{2}\sum_{ij}^n\lambda_i\lambda_jy_iy_j\boldsymbol{x}_i^T\boldsymbol{x}_j,
+$$
+
+
+
but now subject to the constraints \( \lambda_i\geq 0 \), \( \sum_i\lambda_iy_i=0 \) and \( 0\leq\lambda_i \leq C \). +We must in addition satisfy the Karush-Kuhn-Tucker condition which now reads +
+
+$$
+\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) -(1-\xi_)\right]=0 \hspace{0.1cm}\forall i,
+$$
+
+
+
+$$
+\gamma_i\xi_i = 0,
+$$
+
+
+
and
+
+$$
+y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) -(1-\xi_) \geq 0 \hspace{0.1cm}\forall i.
+$$
+
+
A Support Vector Machine (SVM) is a very powerful and versatile +Machine Learning method, capable of performing linear or nonlinear +classification, regression, and even outlier detection. It is one of +the most popular models in Machine Learning, and anyone interested in +Machine Learning should have it in their toolbox. SVMs are +particularly well suited for classification of complex but small-sized or +medium-sized datasets. +
+ +The case with two well-separated classes only can be understood in an +intuitive way in terms of lines in a two-dimensional space separating +the two classes (see figure below). +
+ +The basic mathematics behind the SVM is however less familiar to most of us. +It relies on the definition of hyperplanes and the +definition of a margin which separates classes (in case of +classification problems) of variables. It is also used for regression +problems. +
+ +With SVMs we distinguish between hard margin and soft margins. The +latter introduces a so-called softening parameter to be discussed +below. We distinguish also between linear and non-linear +approaches. The latter are the most frequent ones since it is rather +unlikely that we can separate classes easily by say straight lines. +
+ +The theory behind support vector machines (SVM hereafter) is based on +the mathematical description of so-called hyperplanes. Let us start +with a two-dimensional case. This will also allow us to introduce our +first SVM examples. These will be tailored to the case of two specific +classes, as displayed in the figure here based on the usage of the petal data. +
+ +We assume here that our data set can be well separated into two +domains, where a straight line does the job in the separating the two +classes. Here the two classes are represented by either squares or +circles. +
+ + +from sklearn import datasets
+from sklearn.svm import SVC, LinearSVC
+from sklearn.linear_model import SGDClassifier
+from sklearn.preprocessing import StandardScaler
+import matplotlib
+import matplotlib.pyplot as plt
+plt.rcParams['axes.labelsize'] = 14
+plt.rcParams['xtick.labelsize'] = 12
+plt.rcParams['ytick.labelsize'] = 12
+
+
+iris = datasets.load_iris()
+X = iris["data"][:, (2, 3)] # petal length, petal width
+y = iris["target"]
+
+setosa_or_versicolor = (y == 0) | (y == 1)
+X = X[setosa_or_versicolor]
+y = y[setosa_or_versicolor]
+
+
+
+C = 5
+alpha = 1 / (C * len(X))
+
+lin_clf = LinearSVC(loss="hinge", C=C, random_state=42)
+svm_clf = SVC(kernel="linear", C=C)
+sgd_clf = SGDClassifier(loss="hinge", learning_rate="constant", eta0=0.001, alpha=alpha,
+ max_iter=100000, random_state=42)
+
+scaler = StandardScaler()
+X_scaled = scaler.fit_transform(X)
+
+lin_clf.fit(X_scaled, y)
+svm_clf.fit(X_scaled, y)
+sgd_clf.fit(X_scaled, y)
+
+print("LinearSVC: ", lin_clf.intercept_, lin_clf.coef_)
+print("SVC: ", svm_clf.intercept_, svm_clf.coef_)
+print("SGDClassifier(alpha={:.5f}):".format(sgd_clf.alpha), sgd_clf.intercept_, sgd_clf.coef_)
+
+# Compute the slope and bias of each decision boundary
+w1 = -lin_clf.coef_[0, 0]/lin_clf.coef_[0, 1]
+b1 = -lin_clf.intercept_[0]/lin_clf.coef_[0, 1]
+w2 = -svm_clf.coef_[0, 0]/svm_clf.coef_[0, 1]
+b2 = -svm_clf.intercept_[0]/svm_clf.coef_[0, 1]
+w3 = -sgd_clf.coef_[0, 0]/sgd_clf.coef_[0, 1]
+b3 = -sgd_clf.intercept_[0]/sgd_clf.coef_[0, 1]
+
+# Transform the decision boundary lines back to the original scale
+line1 = scaler.inverse_transform([[-10, -10 * w1 + b1], [10, 10 * w1 + b1]])
+line2 = scaler.inverse_transform([[-10, -10 * w2 + b2], [10, 10 * w2 + b2]])
+line3 = scaler.inverse_transform([[-10, -10 * w3 + b3], [10, 10 * w3 + b3]])
+
+# Plot all three decision boundaries
+plt.figure(figsize=(11, 4))
+plt.plot(line1[:, 0], line1[:, 1], "k:", label="LinearSVC")
+plt.plot(line2[:, 0], line2[:, 1], "b--", linewidth=2, label="SVC")
+plt.plot(line3[:, 0], line3[:, 1], "r-", label="SGDClassifier")
+plt.plot(X[:, 0][y==1], X[:, 1][y==1], "bs") # label="Iris-Versicolor"
+plt.plot(X[:, 0][y==0], X[:, 1][y==0], "yo") # label="Iris-Setosa"
+plt.xlabel("Petal length", fontsize=14)
+plt.ylabel("Petal width", fontsize=14)
+plt.legend(loc="upper center", fontsize=14)
+plt.axis([0, 5.5, 0, 2])
+
+plt.show()
+
+The aim of the SVM algorithm is to find a hyperplane in a +\( p \)-dimensional space, where \( p \) is the number of features that +distinctly classifies the data points. +
+ +In a \( p \)-dimensional space, a hyperplane is what we call an affine subspace of dimension of \( p-1 \). +As an example, in two dimension, a hyperplane is simply as straight line while in three dimensions it is +a two-dimensional subspace, or stated simply, a plane. +
+ +In two dimensions, with the variables \( x_1 \) and \( x_2 \), the hyperplane is defined as
+$$ +b+w_1x_1+w_2x_2=0, +$$ + +where \( b \) is the intercept and \( w_1 \) and \( w_2 \) define the elements of a vector orthogonal to the line +\( b+w_1x_1+w_2x_2=0 \). +In two dimensions we define the vectors \( \boldsymbol{x} =[x1,x2] \) and \( \boldsymbol{w}=[w1,w2] \). +We can then rewrite the above equation as +
+ +$$ +\boldsymbol{x}^T\boldsymbol{w}+b=0. +$$ + + +We limit ourselves to two classes of outputs \( y_i \) and assign these classes the values \( y_i = \pm 1 \). +In a \( p \)-dimensional space of say \( p \) features we have a hyperplane defines as +
+$$ +b+wx_1+w_2x_2+\dots +w_px_p=0. +$$ + +If we define a +matrix \( \boldsymbol{X}=\left[\boldsymbol{x}_1,\boldsymbol{x}_2,\dots, \boldsymbol{x}_p\right] \) +of dimension \( n\times p \), where \( n \) represents the observations for each feature and each vector \( x_i \) is a column vector of the matrix \( \boldsymbol{X} \), +
+$$ +\boldsymbol{x}_i = \begin{bmatrix} x_{i1} \\ x_{i2} \\ \dots \\ \dots \\ x_{ip} \end{bmatrix}. +$$ + +If the above condition is not met for a given vector \( \boldsymbol{x}_i \) we have
+$$ +b+w_1x_{i1}+w_2x_{i2}+\dots +w_px_{ip} >0, +$$ + +if our output \( y_i=1 \). +In this case we say that \( \boldsymbol{x}_i \) lies on one of the sides of the hyperplane and if +
+$$ +b+w_1x_{i1}+w_2x_{i2}+\dots +w_px_{ip} < 0, +$$ + +for the class of observations \( y_i=-1 \), +then \( \boldsymbol{x}_i \) lies on the other side. +
+ +Equivalently, for the two classes of observations we have
+$$ +y_i\left(b+w_1x_{i1}+w_2x_{i2}+\dots +w_px_{ip}\right) > 0. +$$ + +When we try to separate hyperplanes, if it exists, we can use it to construct a natural classifier: a test observation is assigned a given class depending on which side of the hyperplane it is located.
+ + +Let us try to develop our intuition about SVMs by limiting ourselves to a two-dimensional +plane. To separate the two classes of data points, there are many +possible lines (hyperplanes if you prefer a more strict naming) +that could be chosen. Our objective is to find a +plane that has the maximum margin, i.e the maximum distance between +data points of both classes. Maximizing the margin distance provides +some reinforcement so that future data points can be classified with +more confidence. +
+ +What a linear classifier attempts to accomplish is to split the +feature space into two half spaces by placing a hyperplane between the +data points. This hyperplane will be our decision boundary. All +points on one side of the plane will belong to class one and all points +on the other side of the plane will belong to the second class two. +
+ +Unfortunately there are many ways in which we can place a hyperplane +to divide the data. Below is an example of two candidate hyperplanes +for our data sample. +
+ +Let us define the function
+$$ +f(x) = \boldsymbol{w}^T\boldsymbol{x}+b = 0, +$$ + +as the function that determines the line \( L \) that separates two classes (our two features), see the figure here.
+ +Any point defined by \( \boldsymbol{x}_i \) and \( \boldsymbol{x}_2 \) on the line \( L \) will satisfy \( \boldsymbol{w}^T(\boldsymbol{x}_1-\boldsymbol{x}_2)=0 \).
+ +The signed distance \( \delta \) from any point defined by a vector \( \boldsymbol{x} \) and a point \( \boldsymbol{x}_0 \) on the line \( L \) is then
+$$ +\delta = \frac{1}{\vert\vert \boldsymbol{w}\vert\vert}(\boldsymbol{w}^T\boldsymbol{x}+b). +$$ + + +How do we find the parameter \( b \) and the vector \( \boldsymbol{w} \)? What we could +do is to define a cost function which now contains the set of all +misclassified points \( M \) and attempt to minimize this function +
+ +$$ +C(\boldsymbol{w},b) = -\sum_{i\in M} y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b). +$$ + +We could now for example define all values \( y_i =1 \) as misclassified in case we have \( \boldsymbol{w}^T\boldsymbol{x}_i+b < 0 \) and the opposite if we have \( y_i=-1 \). Taking the derivatives gives us
+$$ +\frac{\partial C}{\partial b} = -\sum_{i\in M} y_i, +$$ + +and
+$$ +\frac{\partial C}{\partial \boldsymbol{w}} = -\sum_{i\in M} y_ix_i. +$$ + + +We can now use the Newton-Raphson method or different variants of the gradient descent family (from plain gradient descent to various stochastic gradient descent approaches) to solve the equations
+$$ +b \leftarrow b +\eta \frac{\partial C}{\partial b}, +$$ + +and
+$$ +\boldsymbol{w} \leftarrow \boldsymbol{w} +\eta \frac{\partial C}{\partial \boldsymbol{w}}, +$$ + +where \( \eta \) is our by now well-known learning rate.
+ +The equations we discussed above can be coded rather easily (the +framework is similar to what we developed for logistic +regression). We are going to set up a simple case with two classes only and we want to find a line which separates them the best possible way. +
+ + +
+
+There are however problems with this approach, although it looks +pretty straightforward to implement. When running the above code, we see that we can easily end up with many diffeent lines which separate the two classes. +
+ +For small +gaps between the entries, we may also end up needing many iterations +before the solutions converge and if the data cannot be separated +properly into two distinct classes, we may not experience a converge +at all. +
+ +A better approach is rather to try to define a large margin between +the two classes (if they are well separated from the beginning). +
+ +Thus, we wish to find a margin \( M \) with \( \boldsymbol{w} \) normalized to +\( \vert\vert \boldsymbol{w}\vert\vert =1 \) subject to the condition +
+ +$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq M \hspace{0.1cm}\forall i=1,2,\dots, p. +$$ + +All points are thus at a signed distance from the decision boundary defined by the line \( L \). The parameters \( b \) and \( w_1 \) and \( w_2 \) define this line.
+ +We seek thus the largest value \( M \) defined by
+$$ +\frac{1}{\vert \vert \boldsymbol{w}\vert\vert}y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq M \hspace{0.1cm}\forall i=1,2,\dots, n, +$$ + +or just
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq M\vert \vert \boldsymbol{w}\vert\vert \hspace{0.1cm}\forall i. +$$ + +If we scale the equation so that \( \vert \vert \boldsymbol{w}\vert\vert = 1/M \), we have to find the minimum of +\( \boldsymbol{w}^T\boldsymbol{w}=\vert \vert \boldsymbol{w}\vert\vert \) (the norm) subject to the condition +
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq 1 \hspace{0.1cm}\forall i. +$$ + +We have thus defined our margin as the invers of the norm of +\( \boldsymbol{w} \). We want to minimize the norm in order to have a as large as +possible margin \( M \). Before we proceed, we need to remind ourselves +about Lagrangian multipliers. +
+ +Consider a function of three independent variables \( f(x,y,z) \) . For the function \( f \) to be an +extreme we have +
+$$ +df=0. +$$ + +A necessary and sufficient condition is
+$$ +\frac{\partial f}{\partial x} =\frac{\partial f}{\partial y}=\frac{\partial f}{\partial z}=0, +$$ + +due to
+$$ +df = \frac{\partial f}{\partial x}dx+\frac{\partial f}{\partial y}dy+\frac{\partial f}{\partial z}dz. +$$ + +In many problems the variables \( x,y,z \) are often subject to constraints (such as those above for the margin) +so that they are no longer all independent. It is possible at least in principle to use each +constraint to eliminate one variable +and to proceed with a new and smaller set of independent varables. +
+ +The use of so-called Lagrangian multipliers is an alternative technique when the elimination +of variables is incovenient or undesirable. Assume that we have an equation of constraint on +the variables \( x,y,z \) +
+$$ +\phi(x,y,z) = 0, +$$ + +resulting in
+$$ +d\phi = \frac{\partial \phi}{\partial x}dx+\frac{\partial \phi}{\partial y}dy+\frac{\partial \phi}{\partial z}dz =0. +$$ + +Now we cannot set anymore
+$$ +\frac{\partial f}{\partial x} =\frac{\partial f}{\partial y}=\frac{\partial f}{\partial z}=0, +$$ + +if \( df=0 \) is wanted +because there are now only two independent variables! Assume \( x \) and \( y \) are the independent +variables. +Then \( dz \) is no longer arbitrary. +
+ +However, we can add to
+$$ +df = \frac{\partial f}{\partial x}dx+\frac{\partial f}{\partial y}dy+\frac{\partial f}{\partial z}dz, +$$ + +a multiplum of \( d\phi \), viz. \( \lambda d\phi \), resulting in
+$$ +df+\lambda d\phi = (\frac{\partial f}{\partial z}+\lambda +\frac{\partial \phi}{\partial x})dx+(\frac{\partial f}{\partial y}+\lambda\frac{\partial \phi}{\partial y})dy+ +(\frac{\partial f}{\partial z}+\lambda\frac{\partial \phi}{\partial z})dz =0. +$$ + +Our multiplier is chosen so that
+$$ +\frac{\partial f}{\partial z}+\lambda\frac{\partial \phi}{\partial z} =0. +$$ + +We need to remember that we took \( dx \) and \( dy \) to be arbitrary and thus we must have
+$$ +\frac{\partial f}{\partial x}+\lambda\frac{\partial \phi}{\partial x} =0, +$$ + +and
+$$ +\frac{\partial f}{\partial y}+\lambda\frac{\partial \phi}{\partial y} =0. +$$ + +When all these equations are satisfied, \( df=0 \). We have four unknowns, \( x,y,z \) and +\( \lambda \). Actually we want only \( x,y,z \), \( \lambda \) needs not to be determined, +it is therefore often called +Lagrange's undetermined multiplier. +If we have a set of constraints \( \phi_k \) we have the equations +
+$$ +\frac{\partial f}{\partial x_i}+\sum_k\lambda_k\frac{\partial \phi_k}{\partial x_i} =0. +$$ + + +In order to solve the above problem, we define the following Lagrangian function to be minimized
+$$ +{\cal L}(\lambda,b,\boldsymbol{w})=\frac{1}{2}\boldsymbol{w}^T\boldsymbol{w}-\sum_{i=1}^n\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)-1\right], +$$ + +where \( \lambda_i \) is a so-called Lagrange multiplier subject to the condition \( \lambda_i \geq 0 \).
+ +Taking the derivatives with respect to \( b \) and \( \boldsymbol{w} \) we obtain
+$$ +\frac{\partial {\cal L}}{\partial b} = -\sum_{i} \lambda_iy_i=0, +$$ + +and
+$$ +\frac{\partial {\cal L}}{\partial \boldsymbol{w}} = 0 = \boldsymbol{w}-\sum_{i} \lambda_iy_i\boldsymbol{x}_i. +$$ + +Inserting these constraints into the equation for \( {\cal L} \) we obtain
+$$ +{\cal L}=\sum_i\lambda_i-\frac{1}{2}\sum_{ij}^n\lambda_i\lambda_jy_iy_j\boldsymbol{x}_i^T\boldsymbol{x}_j, +$$ + +subject to the constraints \( \lambda_i\geq 0 \) and \( \sum_i\lambda_iy_i=0 \). +We must in addition satisfy the Karush-Kuhn-Tucker (KKT) condition +
+$$ +\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) -1\right] \hspace{0.1cm}\forall i. +$$ + +When \( \lambda_i > 0 \), the vectors \( \boldsymbol{x}_i \) are called support vectors. They are the vectors closest to the line (or hyperplane) and define the margin \( M \).
+ +We can rewrite
+$$ +{\cal L}=\sum_i\lambda_i-\frac{1}{2}\sum_{ij}^n\lambda_i\lambda_jy_iy_j\boldsymbol{x}_i^T\boldsymbol{x}_j, +$$ + +and its constraints in terms of a matrix-vector problem where we minimize w.r.t. \( \lambda \) the following problem
+$$ +\frac{1}{2} \boldsymbol{\lambda}^T\begin{bmatrix} y_1y_1\boldsymbol{x}_1^T\boldsymbol{x}_1 & y_1y_2\boldsymbol{x}_1^T\boldsymbol{x}_2 & \dots & \dots & y_1y_n\boldsymbol{x}_1^T\boldsymbol{x}_n \\ +y_2y_1\boldsymbol{x}_2^T\boldsymbol{x}_1 & y_2y_2\boldsymbol{x}_2^T\boldsymbol{x}_2 & \dots & \dots & y_1y_n\boldsymbol{x}_2^T\boldsymbol{x}_n \\ +\dots & \dots & \dots & \dots & \dots \\ +\dots & \dots & \dots & \dots & \dots \\ +y_ny_1\boldsymbol{x}_n^T\boldsymbol{x}_1 & y_ny_2\boldsymbol{x}_n^T\boldsymbol{x}_2 & \dots & \dots & y_ny_n\boldsymbol{x}_n^T\boldsymbol{x}_n \\ +\end{bmatrix}\boldsymbol{\lambda}-\mathbb{1}\boldsymbol{\lambda}, +$$ + +subject to \( \boldsymbol{y}^T\boldsymbol{\lambda}=0 \). Here we defined the vectors \( \boldsymbol{\lambda} =[\lambda_1,\lambda_2,\dots,\lambda_n] \) and +\( \boldsymbol{y}=[y_1,y_2,\dots,y_n] \). +
+ +Solving the above problem, yields the values of \( \lambda_i \). +To find the coefficients of your hyperplane we need simply to compute +
+$$ +\boldsymbol{w}=\sum_{i} \lambda_iy_i\boldsymbol{x}_i. +$$ + +With our vector \( \boldsymbol{w} \) we can in turn find the value of the intercept \( b \) (here in two dimensions) via
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1, +$$ + +resulting in
+$$ +b = \frac{1}{y_i}-\boldsymbol{w}^T\boldsymbol{x}_i, +$$ + +or if we write it out in terms of the support vectors only, with \( N_s \) being their number, we have
+$$ +b = \frac{1}{N_s}\sum_{j\in N_s}\left(y_j-\sum_{i=1}^n\lambda_iy_i\boldsymbol{x}_i^T\boldsymbol{x}_j\right). +$$ + +With our hyperplane coefficients we can use our classifier to assign any observation by simply using
+$$ +y_i = \mathrm{sign}(\boldsymbol{w}^T\boldsymbol{x}_i+b). +$$ + +Below we discuss how to find the optimal values of \( \lambda_i \). Before we proceed however, we discuss now the so-called soft classifier.
+ +Till now, the margin is strictly defined by the support vectors. This defines what is called a hard classifier, that is the margins are well defined.
+ +Suppose now that classes overlap in feature space, as shown in the +figure here. One way to deal with this problem before we define the +so-called kernel approach, is to allow a kind of slack in the sense +that we allow some points to be on the wrong side of the margin. +
+ +We introduce thus the so-called slack variables \( \boldsymbol{\xi} =[\xi_1,x_2,\dots,x_n] \) and +modify our previous equation +
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1, +$$ + +to
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1-\xi_i, +$$ + +with the requirement \( \xi_i\geq 0 \). The total violation is now \( \sum_i\xi \). +The value \( \xi_i \) in the constraint the last constraint corresponds to the amount by which the prediction +\( y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1 \) is on the wrong side of its margin. Hence by bounding the sum \( \sum_i \xi_i \), +we bound the total amount by which predictions fall on the wrong side of their margins. +
+ +Misclassifications occur when \( \xi_i > 1 \). Thus bounding the total sum by some value \( C \) bounds in turn the total number of +misclassifications. +
+ +This has in turn the consequences that we change our optmization problem to finding the minimum of
+$$ +{\cal L}=\frac{1}{2}\boldsymbol{w}^T\boldsymbol{w}-\sum_{i=1}^n\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)-(1-\xi_)\right]+C\sum_{i=1}^n\xi_i-\sum_{i=1}^n\gamma_i\xi_i, +$$ + +subject to
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1-\xi_i \hspace{0.1cm}\forall i, +$$ + +with the requirement \( \xi_i\geq 0 \).
+ +Taking the derivatives with respect to \( b \) and \( \boldsymbol{w} \) we obtain
+$$ +\frac{\partial {\cal L}}{\partial b} = -\sum_{i} \lambda_iy_i=0, +$$ + +and
+$$ +\frac{\partial {\cal L}}{\partial \boldsymbol{w}} = 0 = \boldsymbol{w}-\sum_{i} \lambda_iy_i\boldsymbol{x}_i, +$$ + +and
+$$ +\lambda_i = C-\gamma_i \hspace{0.1cm}\forall i. +$$ + +Inserting these constraints into the equation for \( {\cal L} \) we obtain the same equation as before
+$$ +{\cal L}=\sum_i\lambda_i-\frac{1}{2}\sum_{ij}^n\lambda_i\lambda_jy_iy_j\boldsymbol{x}_i^T\boldsymbol{x}_j, +$$ + +but now subject to the constraints \( \lambda_i\geq 0 \), \( \sum_i\lambda_iy_i=0 \) and \( 0\leq\lambda_i \leq C \). +We must in addition satisfy the Karush-Kuhn-Tucker condition which now reads +
+$$ +\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) -(1-\xi_)\right]=0 \hspace{0.1cm}\forall i, +$$ + +$$ +\gamma_i\xi_i = 0, +$$ + +and
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) -(1-\xi_) \geq 0 \hspace{0.1cm}\forall i. +$$ +A Support Vector Machine (SVM) is a very powerful and versatile +Machine Learning method, capable of performing linear or nonlinear +classification, regression, and even outlier detection. It is one of +the most popular models in Machine Learning, and anyone interested in +Machine Learning should have it in their toolbox. SVMs are +particularly well suited for classification of complex but small-sized or +medium-sized datasets. +
+ +The case with two well-separated classes only can be understood in an +intuitive way in terms of lines in a two-dimensional space separating +the two classes (see figure below). +
+ +The basic mathematics behind the SVM is however less familiar to most of us. +It relies on the definition of hyperplanes and the +definition of a margin which separates classes (in case of +classification problems) of variables. It is also used for regression +problems. +
+ +With SVMs we distinguish between hard margin and soft margins. The +latter introduces a so-called softening parameter to be discussed +below. We distinguish also between linear and non-linear +approaches. The latter are the most frequent ones since it is rather +unlikely that we can separate classes easily by say straight lines. +
+ +The theory behind support vector machines (SVM hereafter) is based on +the mathematical description of so-called hyperplanes. Let us start +with a two-dimensional case. This will also allow us to introduce our +first SVM examples. These will be tailored to the case of two specific +classes, as displayed in the figure here based on the usage of the petal data. +
+ +We assume here that our data set can be well separated into two +domains, where a straight line does the job in the separating the two +classes. Here the two classes are represented by either squares or +circles. +
+ + +from sklearn import datasets
+from sklearn.svm import SVC, LinearSVC
+from sklearn.linear_model import SGDClassifier
+from sklearn.preprocessing import StandardScaler
+import matplotlib
+import matplotlib.pyplot as plt
+plt.rcParams['axes.labelsize'] = 14
+plt.rcParams['xtick.labelsize'] = 12
+plt.rcParams['ytick.labelsize'] = 12
+
+
+iris = datasets.load_iris()
+X = iris["data"][:, (2, 3)] # petal length, petal width
+y = iris["target"]
+
+setosa_or_versicolor = (y == 0) | (y == 1)
+X = X[setosa_or_versicolor]
+y = y[setosa_or_versicolor]
+
+
+
+C = 5
+alpha = 1 / (C * len(X))
+
+lin_clf = LinearSVC(loss="hinge", C=C, random_state=42)
+svm_clf = SVC(kernel="linear", C=C)
+sgd_clf = SGDClassifier(loss="hinge", learning_rate="constant", eta0=0.001, alpha=alpha,
+ max_iter=100000, random_state=42)
+
+scaler = StandardScaler()
+X_scaled = scaler.fit_transform(X)
+
+lin_clf.fit(X_scaled, y)
+svm_clf.fit(X_scaled, y)
+sgd_clf.fit(X_scaled, y)
+
+print("LinearSVC: ", lin_clf.intercept_, lin_clf.coef_)
+print("SVC: ", svm_clf.intercept_, svm_clf.coef_)
+print("SGDClassifier(alpha={:.5f}):".format(sgd_clf.alpha), sgd_clf.intercept_, sgd_clf.coef_)
+
+# Compute the slope and bias of each decision boundary
+w1 = -lin_clf.coef_[0, 0]/lin_clf.coef_[0, 1]
+b1 = -lin_clf.intercept_[0]/lin_clf.coef_[0, 1]
+w2 = -svm_clf.coef_[0, 0]/svm_clf.coef_[0, 1]
+b2 = -svm_clf.intercept_[0]/svm_clf.coef_[0, 1]
+w3 = -sgd_clf.coef_[0, 0]/sgd_clf.coef_[0, 1]
+b3 = -sgd_clf.intercept_[0]/sgd_clf.coef_[0, 1]
+
+# Transform the decision boundary lines back to the original scale
+line1 = scaler.inverse_transform([[-10, -10 * w1 + b1], [10, 10 * w1 + b1]])
+line2 = scaler.inverse_transform([[-10, -10 * w2 + b2], [10, 10 * w2 + b2]])
+line3 = scaler.inverse_transform([[-10, -10 * w3 + b3], [10, 10 * w3 + b3]])
+
+# Plot all three decision boundaries
+plt.figure(figsize=(11, 4))
+plt.plot(line1[:, 0], line1[:, 1], "k:", label="LinearSVC")
+plt.plot(line2[:, 0], line2[:, 1], "b--", linewidth=2, label="SVC")
+plt.plot(line3[:, 0], line3[:, 1], "r-", label="SGDClassifier")
+plt.plot(X[:, 0][y==1], X[:, 1][y==1], "bs") # label="Iris-Versicolor"
+plt.plot(X[:, 0][y==0], X[:, 1][y==0], "yo") # label="Iris-Setosa"
+plt.xlabel("Petal length", fontsize=14)
+plt.ylabel("Petal width", fontsize=14)
+plt.legend(loc="upper center", fontsize=14)
+plt.axis([0, 5.5, 0, 2])
+
+plt.show()
+
+The aim of the SVM algorithm is to find a hyperplane in a +\( p \)-dimensional space, where \( p \) is the number of features that +distinctly classifies the data points. +
+ +In a \( p \)-dimensional space, a hyperplane is what we call an affine subspace of dimension of \( p-1 \). +As an example, in two dimension, a hyperplane is simply as straight line while in three dimensions it is +a two-dimensional subspace, or stated simply, a plane. +
+ +In two dimensions, with the variables \( x_1 \) and \( x_2 \), the hyperplane is defined as
+$$ +b+w_1x_1+w_2x_2=0, +$$ + +where \( b \) is the intercept and \( w_1 \) and \( w_2 \) define the elements of a vector orthogonal to the line +\( b+w_1x_1+w_2x_2=0 \). +In two dimensions we define the vectors \( \boldsymbol{x} =[x1,x2] \) and \( \boldsymbol{w}=[w1,w2] \). +We can then rewrite the above equation as +
+ +$$ +\boldsymbol{x}^T\boldsymbol{w}+b=0. +$$ + + +We limit ourselves to two classes of outputs \( y_i \) and assign these classes the values \( y_i = \pm 1 \). +In a \( p \)-dimensional space of say \( p \) features we have a hyperplane defines as +
+$$ +b+wx_1+w_2x_2+\dots +w_px_p=0. +$$ + +If we define a +matrix \( \boldsymbol{X}=\left[\boldsymbol{x}_1,\boldsymbol{x}_2,\dots, \boldsymbol{x}_p\right] \) +of dimension \( n\times p \), where \( n \) represents the observations for each feature and each vector \( x_i \) is a column vector of the matrix \( \boldsymbol{X} \), +
+$$ +\boldsymbol{x}_i = \begin{bmatrix} x_{i1} \\ x_{i2} \\ \dots \\ \dots \\ x_{ip} \end{bmatrix}. +$$ + +If the above condition is not met for a given vector \( \boldsymbol{x}_i \) we have
+$$ +b+w_1x_{i1}+w_2x_{i2}+\dots +w_px_{ip} >0, +$$ + +if our output \( y_i=1 \). +In this case we say that \( \boldsymbol{x}_i \) lies on one of the sides of the hyperplane and if +
+$$ +b+w_1x_{i1}+w_2x_{i2}+\dots +w_px_{ip} < 0, +$$ + +for the class of observations \( y_i=-1 \), +then \( \boldsymbol{x}_i \) lies on the other side. +
+ +Equivalently, for the two classes of observations we have
+$$ +y_i\left(b+w_1x_{i1}+w_2x_{i2}+\dots +w_px_{ip}\right) > 0. +$$ + +When we try to separate hyperplanes, if it exists, we can use it to construct a natural classifier: a test observation is assigned a given class depending on which side of the hyperplane it is located.
+ + +Let us try to develop our intuition about SVMs by limiting ourselves to a two-dimensional +plane. To separate the two classes of data points, there are many +possible lines (hyperplanes if you prefer a more strict naming) +that could be chosen. Our objective is to find a +plane that has the maximum margin, i.e the maximum distance between +data points of both classes. Maximizing the margin distance provides +some reinforcement so that future data points can be classified with +more confidence. +
+ +What a linear classifier attempts to accomplish is to split the +feature space into two half spaces by placing a hyperplane between the +data points. This hyperplane will be our decision boundary. All +points on one side of the plane will belong to class one and all points +on the other side of the plane will belong to the second class two. +
+ +Unfortunately there are many ways in which we can place a hyperplane +to divide the data. Below is an example of two candidate hyperplanes +for our data sample. +
+ +Let us define the function
+$$ +f(x) = \boldsymbol{w}^T\boldsymbol{x}+b = 0, +$$ + +as the function that determines the line \( L \) that separates two classes (our two features), see the figure here.
+ +Any point defined by \( \boldsymbol{x}_i \) and \( \boldsymbol{x}_2 \) on the line \( L \) will satisfy \( \boldsymbol{w}^T(\boldsymbol{x}_1-\boldsymbol{x}_2)=0 \).
+ +The signed distance \( \delta \) from any point defined by a vector \( \boldsymbol{x} \) and a point \( \boldsymbol{x}_0 \) on the line \( L \) is then
+$$ +\delta = \frac{1}{\vert\vert \boldsymbol{w}\vert\vert}(\boldsymbol{w}^T\boldsymbol{x}+b). +$$ + + +How do we find the parameter \( b \) and the vector \( \boldsymbol{w} \)? What we could +do is to define a cost function which now contains the set of all +misclassified points \( M \) and attempt to minimize this function +
+ +$$ +C(\boldsymbol{w},b) = -\sum_{i\in M} y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b). +$$ + +We could now for example define all values \( y_i =1 \) as misclassified in case we have \( \boldsymbol{w}^T\boldsymbol{x}_i+b < 0 \) and the opposite if we have \( y_i=-1 \). Taking the derivatives gives us
+$$ +\frac{\partial C}{\partial b} = -\sum_{i\in M} y_i, +$$ + +and
+$$ +\frac{\partial C}{\partial \boldsymbol{w}} = -\sum_{i\in M} y_ix_i. +$$ + + +We can now use the Newton-Raphson method or different variants of the gradient descent family (from plain gradient descent to various stochastic gradient descent approaches) to solve the equations
+$$ +b \leftarrow b +\eta \frac{\partial C}{\partial b}, +$$ + +and
+$$ +\boldsymbol{w} \leftarrow \boldsymbol{w} +\eta \frac{\partial C}{\partial \boldsymbol{w}}, +$$ + +where \( \eta \) is our by now well-known learning rate.
+ +The equations we discussed above can be coded rather easily (the +framework is similar to what we developed for logistic +regression). We are going to set up a simple case with two classes only and we want to find a line which separates them the best possible way. +
+ + +
+
+There are however problems with this approach, although it looks +pretty straightforward to implement. When running the above code, we see that we can easily end up with many diffeent lines which separate the two classes. +
+ +For small +gaps between the entries, we may also end up needing many iterations +before the solutions converge and if the data cannot be separated +properly into two distinct classes, we may not experience a converge +at all. +
+ +A better approach is rather to try to define a large margin between +the two classes (if they are well separated from the beginning). +
+ +Thus, we wish to find a margin \( M \) with \( \boldsymbol{w} \) normalized to +\( \vert\vert \boldsymbol{w}\vert\vert =1 \) subject to the condition +
+ +$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq M \hspace{0.1cm}\forall i=1,2,\dots, p. +$$ + +All points are thus at a signed distance from the decision boundary defined by the line \( L \). The parameters \( b \) and \( w_1 \) and \( w_2 \) define this line.
+ +We seek thus the largest value \( M \) defined by
+$$ +\frac{1}{\vert \vert \boldsymbol{w}\vert\vert}y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq M \hspace{0.1cm}\forall i=1,2,\dots, n, +$$ + +or just
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq M\vert \vert \boldsymbol{w}\vert\vert \hspace{0.1cm}\forall i. +$$ + +If we scale the equation so that \( \vert \vert \boldsymbol{w}\vert\vert = 1/M \), we have to find the minimum of +\( \boldsymbol{w}^T\boldsymbol{w}=\vert \vert \boldsymbol{w}\vert\vert \) (the norm) subject to the condition +
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) \geq 1 \hspace{0.1cm}\forall i. +$$ + +We have thus defined our margin as the invers of the norm of +\( \boldsymbol{w} \). We want to minimize the norm in order to have a as large as +possible margin \( M \). Before we proceed, we need to remind ourselves +about Lagrangian multipliers. +
+ +Consider a function of three independent variables \( f(x,y,z) \) . For the function \( f \) to be an +extreme we have +
+$$ +df=0. +$$ + +A necessary and sufficient condition is
+$$ +\frac{\partial f}{\partial x} =\frac{\partial f}{\partial y}=\frac{\partial f}{\partial z}=0, +$$ + +due to
+$$ +df = \frac{\partial f}{\partial x}dx+\frac{\partial f}{\partial y}dy+\frac{\partial f}{\partial z}dz. +$$ + +In many problems the variables \( x,y,z \) are often subject to constraints (such as those above for the margin) +so that they are no longer all independent. It is possible at least in principle to use each +constraint to eliminate one variable +and to proceed with a new and smaller set of independent varables. +
+ +The use of so-called Lagrangian multipliers is an alternative technique when the elimination +of variables is incovenient or undesirable. Assume that we have an equation of constraint on +the variables \( x,y,z \) +
+$$ +\phi(x,y,z) = 0, +$$ + +resulting in
+$$ +d\phi = \frac{\partial \phi}{\partial x}dx+\frac{\partial \phi}{\partial y}dy+\frac{\partial \phi}{\partial z}dz =0. +$$ + +Now we cannot set anymore
+$$ +\frac{\partial f}{\partial x} =\frac{\partial f}{\partial y}=\frac{\partial f}{\partial z}=0, +$$ + +if \( df=0 \) is wanted +because there are now only two independent variables! Assume \( x \) and \( y \) are the independent +variables. +Then \( dz \) is no longer arbitrary. +
+ +However, we can add to
+$$ +df = \frac{\partial f}{\partial x}dx+\frac{\partial f}{\partial y}dy+\frac{\partial f}{\partial z}dz, +$$ + +a multiplum of \( d\phi \), viz. \( \lambda d\phi \), resulting in
+$$ +df+\lambda d\phi = (\frac{\partial f}{\partial z}+\lambda +\frac{\partial \phi}{\partial x})dx+(\frac{\partial f}{\partial y}+\lambda\frac{\partial \phi}{\partial y})dy+ +(\frac{\partial f}{\partial z}+\lambda\frac{\partial \phi}{\partial z})dz =0. +$$ + +Our multiplier is chosen so that
+$$ +\frac{\partial f}{\partial z}+\lambda\frac{\partial \phi}{\partial z} =0. +$$ + +We need to remember that we took \( dx \) and \( dy \) to be arbitrary and thus we must have
+$$ +\frac{\partial f}{\partial x}+\lambda\frac{\partial \phi}{\partial x} =0, +$$ + +and
+$$ +\frac{\partial f}{\partial y}+\lambda\frac{\partial \phi}{\partial y} =0. +$$ + +When all these equations are satisfied, \( df=0 \). We have four unknowns, \( x,y,z \) and +\( \lambda \). Actually we want only \( x,y,z \), \( \lambda \) needs not to be determined, +it is therefore often called +Lagrange's undetermined multiplier. +If we have a set of constraints \( \phi_k \) we have the equations +
+$$ +\frac{\partial f}{\partial x_i}+\sum_k\lambda_k\frac{\partial \phi_k}{\partial x_i} =0. +$$ + + +In order to solve the above problem, we define the following Lagrangian function to be minimized
+$$ +{\cal L}(\lambda,b,\boldsymbol{w})=\frac{1}{2}\boldsymbol{w}^T\boldsymbol{w}-\sum_{i=1}^n\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)-1\right], +$$ + +where \( \lambda_i \) is a so-called Lagrange multiplier subject to the condition \( \lambda_i \geq 0 \).
+ +Taking the derivatives with respect to \( b \) and \( \boldsymbol{w} \) we obtain
+$$ +\frac{\partial {\cal L}}{\partial b} = -\sum_{i} \lambda_iy_i=0, +$$ + +and
+$$ +\frac{\partial {\cal L}}{\partial \boldsymbol{w}} = 0 = \boldsymbol{w}-\sum_{i} \lambda_iy_i\boldsymbol{x}_i. +$$ + +Inserting these constraints into the equation for \( {\cal L} \) we obtain
+$$ +{\cal L}=\sum_i\lambda_i-\frac{1}{2}\sum_{ij}^n\lambda_i\lambda_jy_iy_j\boldsymbol{x}_i^T\boldsymbol{x}_j, +$$ + +subject to the constraints \( \lambda_i\geq 0 \) and \( \sum_i\lambda_iy_i=0 \). +We must in addition satisfy the Karush-Kuhn-Tucker (KKT) condition +
+$$ +\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) -1\right] \hspace{0.1cm}\forall i. +$$ + +When \( \lambda_i > 0 \), the vectors \( \boldsymbol{x}_i \) are called support vectors. They are the vectors closest to the line (or hyperplane) and define the margin \( M \).
+ +We can rewrite
+$$ +{\cal L}=\sum_i\lambda_i-\frac{1}{2}\sum_{ij}^n\lambda_i\lambda_jy_iy_j\boldsymbol{x}_i^T\boldsymbol{x}_j, +$$ + +and its constraints in terms of a matrix-vector problem where we minimize w.r.t. \( \lambda \) the following problem
+$$ +\frac{1}{2} \boldsymbol{\lambda}^T\begin{bmatrix} y_1y_1\boldsymbol{x}_1^T\boldsymbol{x}_1 & y_1y_2\boldsymbol{x}_1^T\boldsymbol{x}_2 & \dots & \dots & y_1y_n\boldsymbol{x}_1^T\boldsymbol{x}_n \\ +y_2y_1\boldsymbol{x}_2^T\boldsymbol{x}_1 & y_2y_2\boldsymbol{x}_2^T\boldsymbol{x}_2 & \dots & \dots & y_1y_n\boldsymbol{x}_2^T\boldsymbol{x}_n \\ +\dots & \dots & \dots & \dots & \dots \\ +\dots & \dots & \dots & \dots & \dots \\ +y_ny_1\boldsymbol{x}_n^T\boldsymbol{x}_1 & y_ny_2\boldsymbol{x}_n^T\boldsymbol{x}_2 & \dots & \dots & y_ny_n\boldsymbol{x}_n^T\boldsymbol{x}_n \\ +\end{bmatrix}\boldsymbol{\lambda}-\mathbb{1}\boldsymbol{\lambda}, +$$ + +subject to \( \boldsymbol{y}^T\boldsymbol{\lambda}=0 \). Here we defined the vectors \( \boldsymbol{\lambda} =[\lambda_1,\lambda_2,\dots,\lambda_n] \) and +\( \boldsymbol{y}=[y_1,y_2,\dots,y_n] \). +
+ +Solving the above problem, yields the values of \( \lambda_i \). +To find the coefficients of your hyperplane we need simply to compute +
+$$ +\boldsymbol{w}=\sum_{i} \lambda_iy_i\boldsymbol{x}_i. +$$ + +With our vector \( \boldsymbol{w} \) we can in turn find the value of the intercept \( b \) (here in two dimensions) via
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1, +$$ + +resulting in
+$$ +b = \frac{1}{y_i}-\boldsymbol{w}^T\boldsymbol{x}_i, +$$ + +or if we write it out in terms of the support vectors only, with \( N_s \) being their number, we have
+$$ +b = \frac{1}{N_s}\sum_{j\in N_s}\left(y_j-\sum_{i=1}^n\lambda_iy_i\boldsymbol{x}_i^T\boldsymbol{x}_j\right). +$$ + +With our hyperplane coefficients we can use our classifier to assign any observation by simply using
+$$ +y_i = \mathrm{sign}(\boldsymbol{w}^T\boldsymbol{x}_i+b). +$$ + +Below we discuss how to find the optimal values of \( \lambda_i \). Before we proceed however, we discuss now the so-called soft classifier.
+ +Till now, the margin is strictly defined by the support vectors. This defines what is called a hard classifier, that is the margins are well defined.
+ +Suppose now that classes overlap in feature space, as shown in the +figure here. One way to deal with this problem before we define the +so-called kernel approach, is to allow a kind of slack in the sense +that we allow some points to be on the wrong side of the margin. +
+ +We introduce thus the so-called slack variables \( \boldsymbol{\xi} =[\xi_1,x_2,\dots,x_n] \) and +modify our previous equation +
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1, +$$ + +to
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1-\xi_i, +$$ + +with the requirement \( \xi_i\geq 0 \). The total violation is now \( \sum_i\xi \). +The value \( \xi_i \) in the constraint the last constraint corresponds to the amount by which the prediction +\( y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1 \) is on the wrong side of its margin. Hence by bounding the sum \( \sum_i \xi_i \), +we bound the total amount by which predictions fall on the wrong side of their margins. +
+ +Misclassifications occur when \( \xi_i > 1 \). Thus bounding the total sum by some value \( C \) bounds in turn the total number of +misclassifications. +
+ +This has in turn the consequences that we change our optmization problem to finding the minimum of
+$$ +{\cal L}=\frac{1}{2}\boldsymbol{w}^T\boldsymbol{w}-\sum_{i=1}^n\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)-(1-\xi_)\right]+C\sum_{i=1}^n\xi_i-\sum_{i=1}^n\gamma_i\xi_i, +$$ + +subject to
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b)=1-\xi_i \hspace{0.1cm}\forall i, +$$ + +with the requirement \( \xi_i\geq 0 \).
+ +Taking the derivatives with respect to \( b \) and \( \boldsymbol{w} \) we obtain
+$$ +\frac{\partial {\cal L}}{\partial b} = -\sum_{i} \lambda_iy_i=0, +$$ + +and
+$$ +\frac{\partial {\cal L}}{\partial \boldsymbol{w}} = 0 = \boldsymbol{w}-\sum_{i} \lambda_iy_i\boldsymbol{x}_i, +$$ + +and
+$$ +\lambda_i = C-\gamma_i \hspace{0.1cm}\forall i. +$$ + +Inserting these constraints into the equation for \( {\cal L} \) we obtain the same equation as before
+$$ +{\cal L}=\sum_i\lambda_i-\frac{1}{2}\sum_{ij}^n\lambda_i\lambda_jy_iy_j\boldsymbol{x}_i^T\boldsymbol{x}_j, +$$ + +but now subject to the constraints \( \lambda_i\geq 0 \), \( \sum_i\lambda_iy_i=0 \) and \( 0\leq\lambda_i \leq C \). +We must in addition satisfy the Karush-Kuhn-Tucker condition which now reads +
+$$ +\lambda_i\left[y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) -(1-\xi_)\right]=0 \hspace{0.1cm}\forall i, +$$ + +$$ +\gamma_i\xi_i = 0, +$$ + +and
+$$ +y_i(\boldsymbol{w}^T\boldsymbol{x}_i+b) -(1-\xi_) \geq 0 \hspace{0.1cm}\forall i. +$$ +