# l_1 and l_infinity regression using cvxpy
import numpy as np
import cvxpy as cvx
import matplotlib.pyplot as plt
# generate a synthetic dataset
# actual parameter values
theta1_act = 2
theta2_act = 5
# Number of points in dataset
N = 200
# Noise magnitude
mag = 30
np.random.seed(123)
# datapoints
x = np.arange(0,N)
y = theta1_act * x + theta2_act *np.ones([1,N]) + np.random.normal(0,mag,N)
plt.figure()
# Scatter plot of data
plt.scatter(x,y)
plt.show()
# generate a 200x2 matrix such that X^T*theta = x*theta_1 + theta_2
X = np.ones((N,2))
for i in range(N):
X[i,0] = x[i]
#define the parameter vector
theta = cvx.Variable(2)
#define the l1-norm slack variables
t = cvx.Variable(N)
one_vector = np.ones((N,1))
#l1 norm optimization problem
objective_1 = cvx.Minimize(one_vector.T@t)
constraints_1 = []
for i in range(N):
constraints_1.append(t[i] >= 0)
constraints_1.append(y[0,i] - X[i,]@theta <= t[i])
constraints_1.append(-y[0,i] + X[i,]@theta <= t[i])
prob1 = cvx.Problem(objective_1, constraints_1)
prob1.solve()
theta_1 = theta.value
opt_1 = prob1.value
print(prob1.status)
print(theta_1, opt_1)
#define the slack variable s
s = cvx.Variable(1)
#l_inf norm optimization problem
objective_inf = cvx.Minimize(s)
constraints_inf = [s >=0]
for i in range(N):
constraints_inf.append(y[0,i] - X[i,]@theta <= s)
constraints_inf.append(-y[0,i] + X[i,]@theta <= s)
prob_inf = cvx.Problem(objective_inf, constraints_inf)
prob_inf.solve()
print(prob_inf.status)
opt_inf = prob_inf.value
theta_inf = theta.value
print(theta_inf, opt_inf)
#define the equations of the lines
y1 = theta_1[0]*x + theta_1[1]
y_inf = theta_inf[0]*x + theta_inf[1]
plt.figure()
# Scatter plot of data
plt.scatter(x,y)
#plot each line of best fit
plt.plot(x, y1, '-r', label='l1-norm')
plt.plot(x, y_inf, '-b', label='$l_\infty$-norm')
plt.legend(loc='upper left')
plt.title("Lines of best fit for the synthetic dataset")
plt.show()
theta_act = np.array((theta1_act,theta2_act))
l1_err = np.linalg.norm(theta_act - theta_1)
l_inf_err = np.linalg.norm(theta1_act - theta_inf)
print("The l_1 optimization problem has an error of magnitude",
l1_err)
print("The l_\infty optimization problem has an error of magnitude",
l_inf_err, ".")