Exercises

1. Compute an SVM with Quadratic Programming Consider the following two-dimensional binary classification dataset:

   \(x_1\)	\(x_2\)	Class \(y\)
   2	2	+1
   2	0	+1
   0	0	-1
   0	2	-1

   Compute the SVM for this dataset with a quadratic solver for example cvxpy. This Python package allows for a implementing an optimization objective over the objective function and the constraints. The typical workflow looks like this:

   import cvxpy as cp
   a = np.array([2.0, -1.0]) # constant vector
   
   # Define optimization variables
   x = cp.Variable(2)      # vector with 2 components
   c = cp.Variable()       # scalar
   
   # Define the objective function
   objective = cp.Minimize(...)
   
   # Define the constraints
   constraints = [ x>=0, 
       cp.sum_squares(x) +c <= 10, # x_1^2 + x_2^2 + c <= 10
       a @ x >=3, # 2*x_1 -x_2 >=3
       ...
   ]
   
   # 4. Construct the optimization problem and solve it
   problem = cp.Problem(objective, constraints)
   problem.solve()
   
   # 6. Get the result
   print(problem.status)
   print(problem.value)
   print(x.value)
   print(c.value)
   

   Tasks

   1. Write down the primal and the dual optimization problem explicitly for the four observations.

   2. Solve the quadratic optimization problem using a quadratic solver.

   3. Report the optimal values of \(w_1\), \(w_2\), and \(b\).

   1. For the hard-margin SVM, the primal problem is defined as follows. We get one constraint for each data point, stating that the data point needs to be on the correct side of the decision boundary. The primal objectie is stated in Eq. (39) and for the given dataset, it translates to the following objective:

      \[\begin{align*} \min_{w_1,w_2,b} &\frac12(w_1^2+w_2^2)\\ \text{s.t. } &2w_1+2w_2+b\geq 1\\ & 2w_1+b\geq 1\\ & -b\geq 1\\ & -2w_2-b\geq 1 \end{align*}\]

      To write the dual problem, we introduce one Lagrange multiplier \(\lambda_i\geq0\) for each observation. According to Eq. (42) the dual objective is defined as follows:

      \[\begin{align*} \min_{\lambda_1,\ldots,\lambda_4}\quad &\frac14 \left( 8\lambda_1^2 +4\lambda_1\lambda_2 -4\lambda_1\lambda_4 +4\lambda_2^2 +4\lambda_4^2 \right) - \sum_{i=1}^4\lambda_i\\ \text{s.t. }\quad& \lambda_1+\lambda_2-\lambda_3-\lambda_4=0,\\ &\lambda_i\geq0 \text{ for all }1\leq i\leq 4. \end{align*}\]

      Alternatively, we can define the matrix \(Q\) with entries \(Q_{ij}=y_i y_j x_i^\top x_j\)

      \[\begin{split} Q= \begin{pmatrix} 8 & 4 & 0 & -4\\ 4 & 4 & 0 & 0\\ 0 & 0 & 0 & 0\\ -4 & 0 & 0 & 4 \end{pmatrix}. \end{split}\]

      and use Eq. (43) to define our dual objective.

   2. We solve the dual quadratic objective with the Python package:

      import numpy as np
      import cvxpy as cp
      
      X = np.array([
          [2., 2.],
          [2., 0.],
          [0., 0.],
          [0., 2.]
      ])
      
      y = np.array([1., 1., -1., -1.])
      
      # Construct Q_ij = y_i y_j x_i^T x_j
      Q = np.outer(y, y) * (X @ X.T)
      
      # One dual variable for each observation
      lambdas = cp.Variable(4)
      
      # Dual objective
      objective = cp.Minimize(
          0.25 * cp.quad_form(lambdas, Q)
          - cp.sum(lambdas)
      )
      
      # Dual constraints
      constraints = [
          y @ lambdas == 0,
          lambdas >= 0       
      ]
      
      problem = cp.Problem(objective, constraints)
      problem.solve(solver=cp.OSQP)
      
      print("Status:", problem.status)
      print("Optimal objective:", problem.value)
      print("lambda =", lambdas.value)
      

   3. As solution for the Lagrange multipliers we get \(\lambda = (\frac13,\frac23,\frac23,\frac13)\). According to Eq. (44) we can recover the normal vector as

      \[\begin{split} \begin{align*} \mathbf{w}&=\frac12\sum_{i=1}^4\lambda_i y_i \mathbf{x}_i\\ &= \frac12\left(\frac13 \begin{pmatrix} 2\\2\end{pmatrix} + \frac23 \begin{pmatrix} 2\\0\end{pmatrix} - \frac13 \begin{pmatrix} 0\\2\end{pmatrix}\right)\\ &= \begin{pmatrix} 1\\0\end{pmatrix} \end{align*} \end{split}\]

      Using Python, we can compute the normal vector \(\mathbf{w}\) as (lambdas.value * y).T @ X. For computing the bias term, we use again Eq. (44), stating that it can be recovered for any data point having a nonzero Lagrange multiplier. Since all data points have a nonzero Lagrange multiplier, we can choose. It’s most simple to use the third data point, which is equal to the zero vector. This way, we obtain

      \[ b= y_3 - \mathbf{w}^\top \mathbf{x}_3 = -1. \]