Skip to content

Instantly share code, notes, and snippets.

@chrisdembia
Last active August 29, 2015 14:01
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save chrisdembia/fb546144f64275bcdead to your computer and use it in GitHub Desktop.
Save chrisdembia/fb546144f64275bcdead to your computer and use it in GitHub Desktop.
ee364a_hw4_prob5p13.py
import numpy as np
import numpy.matlib as ml
# data for censored fitting problem.
np.random.seed(0)
n = 20 # dimension of x's
M = 25 # number of non-censored data points
K = 100 # total number of points
c_true = ml.randn(n, 1)
X_unsorted = ml.randn(n, K)
y_unsorted = X_unsorted.T * c_true + 0.1 * np.sqrt(n) * ml.randn(K, 1)
# Reorder measurements, then censor
sort_ind = np.argsort(y_unsorted, 0)
y = y_unsorted[sort_ind, 0]
X = X_unsorted[:, sort_ind]
D = 0.5 * (y[M - 1] + y[M])
y = y[0:M]
# In[34]:
from cvxpy import *
from cens_fit_data import *
#from numpy.linalg import norm
# In[35]:
c = Variable(n)
y_cens = Variable(K - M)
constraints = [y_cens >= D]
# In[36]:
uncensored_error = sum(square(y - X[:, 0:M].T * c))
censored_error = sum(square(y_cens - X[:, M::].T * c))
objective = Minimize(uncensored_error + censored_error)
problem_with_censored = Problem(objective, constraints)
problem_with_censored.solve()
#print('c Error when including censored data: {}'.format(np.linalg.norm(c_true - c.value) / np.linalg.norm(c_true)))
# In[28]:
error = norm2(vstack(y, y_cens) - X.T * c)
objective = Minimize(error)
problem_with_censored = Problem(objective, constraints)
problem_with_censored.solve(solver=CVXOPT)
# In[36]:
uncensored_error = (y - X[:, 0:M].T * c).T * (y - X[:, 0:M].T * c)
censored_error = (y_cens - X[:, M::].T * c).T * (y_cens - X[:, M::].T * c)
objective = Minimize(uncensored_error + censored_error)
problem_with_censored = Problem(objective, constraints)
problem_with_censored.solve()
# In[ ]:
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment