import numpy as np
import matplotlib.pyplot as plt
import cvxpy as cp

n = 10
onev = np.ones(n)
y = 10*onev
y[-1] = 11

plt.figure(figsize=(14,8), dpi=100)
#plt.stem(y, use_line_collection=True)
plt.stem(y)
plt.xlabel('index')
plt.ylabel(r'$y_i$')
plt.grid()
plt.show()

yhat = cp.Variable()

loss_2norm = cp.norm(yhat-y, 2)
loss_1norm = cp.norm(yhat-y, 1)
loss_infty = cp.norm(yhat-y, 'inf')

obj = cp.Minimize(loss_2norm)
prob = cp.Problem(obj)
prob.solve()
yhat_2norm = yhat.value

obj = cp.Minimize(loss_1norm)
prob = cp.Problem(obj)
prob.solve()
yhat_1norm = yhat.value

obj = cp.Minimize(loss_infty)
prob = cp.Problem(obj)
prob.solve()
yhat_infty = yhat.value

plt.figure(figsize=(14,8), dpi=100)
plt.stem(y)
plt.axhline(yhat_2norm, color='r', alpha=0.5, linewidth=4, \
            label=r'${\arg\min}_\hat{y} ||{\bf 1}\hat{y}-y||_2$')
plt.axhline(yhat_1norm, color='g', alpha=0.5, linewidth=4, \
            label=r'${\arg\min}_\hat{y} ||{\bf 1}\hat{y}-y||_1$')
plt.axhline(yhat_infty, color='c', alpha=0.5, linewidth=4, \
            label=r'$${\arg\min}_\hat{y} ||{\bf 1}\hat{y}-y||_\infty$')
plt.grid(), plt.legend()
plt.xlabel('index')
plt.ylabel(r'$y_i$')
plt.ylim(9,12) 
plt.show()


import numpy as np
import matplotlib.pyplot as plt
import cvxpy as cp

m = 100
n = 1000

np.random.seed(7030)
A = np.random.randn(m, n)
x = np.random.rand(n)
x *= np.random.rand(n)>0.99
y_clean = np.dot(A, x)
y = y_clean + 0.1*np.random.randn(m)

plt.figure(figsize=(14,8), dpi=100)
plt.subplot(211)
plt.plot(x)
plt.grid()
plt.xlabel('index'), plt.ylabel(r'$x$')
plt.subplot(212)
plt.plot(y_clean, label='clean')
plt.plot(y, alpha=0.8, label='corrupted')
plt.grid(), plt.legend()
plt.xlabel('index')
plt.ylabel(r'$y$')
plt.show()

xhat = cp.Variable(n)

loss_2norm = cp.norm(xhat, 2)
obj = cp.Minimize(loss_2norm)
constr = [ A@xhat==y ]
prob = cp.Problem(obj, constr)
prob.solve()
xhat_2norm = xhat.value

plt.figure(figsize=(14,8), dpi=100)
plt.plot(x, label='true')
plt.plot(xhat_2norm, alpha=0.8, label=r'$\ell_2$ reconstructed')
plt.grid(), plt.legend()
plt.xlabel('index'), plt.ylabel(r'$x$')
plt.show()

xhat = cp.Variable(n)

loss_2norm = cp.sum_squares(A@xhat-y) 
regu_2norm = cp.norm(xhat,2)
obj = cp.Minimize(loss_2norm + 10*regu_2norm)
prob = cp.Problem(obj)
prob.solve()
xhat_2norm = xhat.value

plt.figure(figsize=(14,8), dpi=100)
plt.plot(x, label='true')
plt.plot(xhat_2norm, alpha=0.8, label=r'$\ell_2$ reconstructed')
plt.grid(), plt.legend()
plt.xlabel('index'), plt.ylabel(r'$x$')
plt.show()

xhat = cp.Variable(n)

loss_2norm = cp.sum_squares(A@xhat-y) 
regu_1norm = cp.norm(xhat,1)
obj = cp.Minimize(loss_2norm + 10*regu_1norm)
prob = cp.Problem(obj)
prob.solve()
xhat_1norm = xhat.value

plt.figure(figsize=(14,8), dpi=100)
plt.plot(x, label='true')
plt.plot(xhat_1norm, alpha=0.8, label=r'$\ell_1$ reconstructed')
plt.grid(), plt.legend()
plt.xlabel('index'), plt.ylabel(r'$x$')
plt.show()