Page MenuHomec4science

project1.py
No OneTemporary

File Metadata

Created
Sun, Feb 23, 23:29

project1.py

import numpy as np
from proj1_helpers import *
from functions import *
def sigmoid(t):
"""apply sigmoid function on t."""
#print(np.e**(-t))
return 1/(1+np.e**(-t))
def calculate_log_loss(y, tx, w):
"""compute the cost by negative log likelihood."""
return np.sum(np.log(1 + np.e**(tx@w)) - y*(tx@w))
def calculate_log_gradient(y, tx, w):
"""compute the gradient of loss."""
return tx.T@(sigmoid(tx@w)-y)
def calculate_hessian(y, tx, w):
"""return the hessian of the loss function."""
S = sigmoid(tx@w)*(1-sigmoid(tx@w))
S = np.identity(len(S))*S
#a_shape_before = S.shape
#a_shape_after = S[np.logical_not(np.isnan(S))].shape
#print(a_shape_before,a_shape_after)
return (tx.T@S@tx)
def penalized_logistic_regression(y, tx, w, lambda_):
"""return the loss, gradient, and hessian."""
# ***************************************************
# INSERT YOUR CODE HERE
# return loss, gradient, and hessian: TODO
# ***************************************************
loss = np.sum(np.log(1 + np.exp(tx@w)) - y*tx@w)
#print(calculate_log_loss(y,tx,w))
gradient = tx.T@(sigmoid(tx@w) - y) + lambda_*w
hessian = calculate_hessian(y,tx,w) + lambda_*np.identity(len(w))
return loss,gradient,hessian
def learning_by_penalized_gradient(y, tx, w, gamma, lambda_):
"""
Do one step of gradient descent, using the penalized logistic regression.
Return the loss and updated w.
"""
# ***************************************************
# INSERT YOUR CODE HERE
# return loss, gradient: TODO
# ***************************************************
loss,gradient,hessian = penalized_logistic_regression(y,tx,w,lambda_)
# ***************************************************
# INSERT YOUR CODE HERE
# update w: TODO
# ***************************************************
w_new = w-gamma*np.linalg.pinv(hessian)@gradient
#w_new = w-gamma*gradient
return loss, w_new
def logistic_regression(y, tx, w):
"""return the loss, gradient, and hessian."""
loss = calculate_log_loss(y,tx,w)
gradient = calculate_log_gradient(y,tx,w)
hessian = calculate_hessian(y,tx,w)
return loss, gradient, hessian
def learning_by_newton_method(y, tx, w):
"""
Do one step on Newton's method.
return the loss and updated w.
"""
loss,gradient,hessian = logistic_regression(y,tx,w)
# ***************************************************
#print("hessian",np.min(hessian),np.max(hessian))
w_new = w - np.linalg.pinv(hessian)@gradient
return loss, w_new
def normalize(x):
minimum = np.min(x)
maximum = np.max(x)
range_x = maximum - minimum
x -= minimum
x = x/range_x
return x
def normalize2(x_te, x_tr):
minimum = min(np.min(x_te),np.min(x_tr))
maximum = max(np.max(x_te),np.max(x_tr))
range_x = maximum - minimum
x_te -= minimum
x_te = x_te/range_x
return x_te
y_tr, x_tr_noisy, ids_tr = load_csv_data("train.csv")
y_te,x_te,ids_te = load_csv_data("test.csv")
y_tr = np.expand_dims(y_tr,axis=1)
y_tr = (y_tr+1)/2
print("1")
x_tr1 = real_mean(x_tr_noisy) #format_data(x_tr_noisy)
x_te1 = real_mean2(x_te,x_tr_noisy)
x_tr = normalize2(x_tr1,x_te1)
degree = 5
tx_tr = build_poly_matrix(x_tr,degree)
#tx_tr = np.c_[np.ones((y_tr.shape[0], 1)),x_tr]
#print(x_tr)
#x_tr = x_tr_noisy
#x_tr = x_tr.T[:10].T
#print("2")
"""
N = x_tr_noisy.shape[0] - len(np.unique(np.where(x_tr_noisy==-999)[0]))
D = 30 #x_tr_noisy.shape[1] - len(np.delete(x_tr,np.unique(np.where(x_tr_noisy==-999)[1])))
x_tr1 = real_mean(x_tr_noisy) #format_data(x_tr_noisy)
x_te1 = real_mean2(x_te,x_tr)
x_tr = normalize2(x_tr1,x_te1)
x_tr2 = np.delete(x_tr,np.unique(np.where(x_tr_noisy==-999)[0]),axis=0)#.reshape(N,D)
y_tr2 = np.delete(y_tr,np.unique(np.where(x_tr_noisy==-999)[0]),axis=0)#.reshape(N,1)
tx_tr2 = np.c_[np.ones((y_tr2.shape[0], 1)),x_tr2]
"""
#print("TX : ",tx_tr)
#print(x_tr_noisy)
w = np.zeros((tx_tr.shape[1], 1))
#print("ok")
def one_iteration(y,tx,w,start,stop):
max_iter = 100
threshold = 1e-8
gamma= 0.3
lambda_ = 0.0000004
losses = []
# start the logistic regression
for iter in range(max_iter):
# get loss and update w.
#print(tx_tr)
#print("x: ",tx_tr.shape)
#print("y: ",y_tr.shape)
#print("w: ",w.shape)
"""
print("AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA")
print(tx_tr[:2])
print("BBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBB")
print(y_tr[:2])
print("CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC")"""
loss, w = learning_by_penalized_gradient(y[start:stop], tx[start:stop], w, gamma, lambda_)
#print("AAAH BON",y_tr[:1],tx_tr[:1])
# log info
if iter % 1 == 0:
print("Current iteration={i}, the loss={l}".format(i=iter, l=loss))
#converge criterion
losses.append(loss)
if len(losses) > 1 and np.abs(losses[-1] - losses[-2]) < threshold:
break
return w.T
step = 20000
shuffled_indices = np.arange(len(y_tr))
np.random.shuffle(shuffled_indices)
y_to_use = y_tr[shuffled_indices]
tx_to_use = tx_tr[shuffled_indices]
w_i= one_iteration(y_to_use,tx_to_use,w,0,step)
w_list = np.array(w_i)
print(1,"/",int(250000/step), calculate_loss(y_to_use, tx_to_use, w_i.T))
for i in range(step,len(y_to_use),step):
w_i = one_iteration(y_to_use,tx_to_use,w,i,i+step)
#print(i," : ",w_i)
#print(w_list)
print(int(i/step)+1,"/",int(250000/step), calculate_loss(y_to_use, tx_to_use, w_i.T))
w_list = np.vstack((w_list,w_i))
#print(len(w_list))
#print(w_list)
w_final = np.mean(w_list,axis=0)
#print(w_final)
print("loss={l}".format(l=calculate_loss(y_tr, tx_tr, w_final)))

Event Timeline