def simple_predict(x, theta):
    """Computes the vector of prediction y_hat from two non-empty numpy.ndarray.
     Args:
      x: has to be an numpy.ndarray, a vector of dimension m * 1.
     theta: has to be an numpy.ndarray, a vector of dimension 2 * 1.
      Returns:
      y_hat as a numpy.ndarray, a vector of dimension m * 1.
      None if x or theta are empty numpy.ndarray. None if x or theta dimensions are not appropriate.
     Raises: This function should not raise any Exception.
    """
    if len(x) == 0:
        return None
    answer = []
    for i in range(len(x)):
        temp = theta[0]
        for j in range(len(x[0])):
            temp += (x[i][j] * theta[j + 1])
        answer.append(temp)
    return (answer)

def simple_predict(x, theta):
    """Computes the vector of prediction y_hat from two non-empty numpy.ndarray.
     Args:
      x: has to be an numpy.ndarray, a vector of dimension m * 1.
     theta: has to be an numpy.ndarray, a vector of dimension 2 * 1.
      Returns:
      y_hat as a numpy.ndarray, a vector of dimension m * 1.
      None if x or theta are empty numpy.ndarray. None if x or theta dimensions are not appropriate.
     Raises: This function should not raise any Exception.
    """
    def add_intercept(x):
        if len(x) == 0 or x.ndim >= 3:
            return None
        if x.ndim == 1:
            return np.vstack((np.ones(len(x)), x)).T
        else:
            return np.insert(x, 0, 1, axis=1)
    if len(x) == 0:
        return None
    return (add_intercept(x) @ theta)

import numpy as np

def fit_(x, y, theta, alpha, n_cycle):
    """ Description: Fits the model to the training dataset contained in x and y.
    Args:
    x:
    has to be a numpy.ndarray, a vector of dimension m * n:  (number of training examples, number of features).
    y:
    has to be a numpy.ndarray, a vector of dimension m * 1: (number of training examples, 1).
    theta: has to be a numpy.ndarray, a vector of dimension (n + 1) * 1: (number of features + 1, 1).
    alpha: has to be a float, the learning rate
    max_iter: has to be an int, the number of iterations done during the gradient descent
    Returns:
    new_theta: numpy.ndarray, a vector of dimension (number of features + 1, 1).
    None if there is a matching dimension problem.
    Raises: This function should not raise any Exception. """
    def gradient(x, y, theta):
        def add_intercept(x):
            if len(x) == 0 or x.ndim >= 3:
                return None
            if x.ndim == 1:
                return np.vstack((np.ones(len(x)), x)).T
            else:
                return np.insert(x, 0, 1, axis=1)
        def predict_(x, theta):
            if x.ndim == 1:
                x = x[:,np.newaxis]
            if len(x) == 0:
                return None
            return add_intercept(x) @ theta
        if len(x) == 0 or len(y) == 0 or len(theta) == 0:
            return None
        return add_intercept(x).T @ (predict_(x, theta) - y) / len(x)


    if len(x) == 0 or len(y) == 0 or len(theta) == 0:
        return None
    theta_r = np.array(theta, dtype=float)
    while n_cycle:
        for i, v in enumerate(gradient(x, y, theta_r)):
            theta_r[i] -= (alpha * v)
        n_cycle -= 1
    return theta_r

X = np.array([[1., 1., 2., 3.], [5., 8., 13., 21.], [34., 55., 89., 144.]])
Y = np.array([[23.], [48.], [218.]])
mylr = MyLinearRegression([[1.], [1.], [1.], [1.], [1]])
print(mylr.cost_(X,Y))
mylr.fit_(X, Y)
print(mylr.thetas)
print(mylr.predict_(X))
print(mylr.cost_elem_(X,Y))
print(mylr.cost_(X,Y))

csv_data = pd.read_csv("../resources/spacecraft_data.csv")
print(csv_data["Age"])
x = np.array(csv_data["Age"]).reshape(-1,1)
y = np.array(csv_data["Sell_price"]).reshape(-1,1)
plt.scatter(x, y,color="navy", s=10)
plt.xlabel('x1: age (in years)')
plt.ylabel('y : Sell_price(in kiloeuros)')
plt.title('Plot of the selling prices of spacecrafts with respect to their age')
mylr = MyLinearRegression([[1.], [1.]])
mylr.fit_(x, y)
print(mylr.thetas)
print(mylr.cost_(x, y))
plt.scatter(x, mylr.predict_(x), color="cornflowerblue", s=7)
plt.legend(['Sell_price', 'Predicted sell_price'])
plt.grid()
plt.show()

mylr = MyLinearRegression([[338.94317973],
 [-22.67763021],
 [  5.84252624],
 [ -2.59281776]])
mylr.fit_(x, y)
plt.scatter(x_dist, y,color="purple", s=10)
plt.xlabel('x1: Terameters (in Tmeters)')
plt.ylabel('y : Sell_price(in kiloeuros)')
plt.title('Plot of the selling prices of spacecrafts with respect to their dist')
plt.scatter(x_dist, mylr.predict_(x), color="violet", s=7)
plt.legend(['Sell_price', 'Predicted sell_price'])
plt.grid()
plt.show()

def add_polynomial_features(x, power):
    """Add polynomial features to vector x by raising its values up to the power given in argument.
    Args: x: has to be an numpy.ndarray, a vector of dimension m * 1.
    power: has to be an int, the power up to which the components of vector x are going to be raised., →
    Returns: The matrix of polynomial features as a numpy.ndarray,
        of dimension m * n, containg he polynomial feature values for all training examples.
    None if x is an empty numpy.ndarray.
    Raises: This function should not raise any Exception. """
    temp = x.copy()
    for i in range(2, power + 1):
        temp = np.append(temp, np.power(x, i), axis=1 )
    return temp

csv_data = pd.read_csv("../resources/are_blue_pills_magics.csv")
y_n = []
x = np.array(csv_data["Micrograms"]).reshape(-1,1)
x = zscore(x)

y = np.array(csv_data["Score"]).reshape(-1,1)
y = zscore(y)

x2 = add_polynomial_features(x, 2)
x3 = add_polynomial_features(x, 3)
x4 = add_polynomial_features(x, 4)
x5 = add_polynomial_features(x, 5)
x6 = add_polynomial_features(x, 6)
x7 = add_polynomial_features(x, 7)
x8 = add_polynomial_features(x, 8)
x9 = add_polynomial_features(x, 9)

mylr2 = MyLinearRegression([[88.85],[-9.0 ], [0]])
mylr3 = MyLinearRegression([[88.85],[-9.0 ], [0], [0]])
mylr4 = MyLinearRegression([[88.85],[-9.0 ], [0], [0], [0]])

mylr5 = MyLinearRegression([[88.85],[-9.0 ], [0], [0], [0.0], [0.0]])
mylr6 = MyLinearRegression([[88.85],[-9.0 ], [0], [0.0], [0.0], [0.0], [0.0]])

mylr7 = MyLinearRegression([[88.85],[-9.0 ], [0.0], [0.0], [0.0], [0.0], [0.0], [0.0]])
mylr8 = MyLinearRegression([[88.85],[-9.0 ], [0.0], [0.0], [0.0], [0.0], [0.0], [0.0], [0.0]])

mylr9 = MyLinearRegression([[88.85],[-9.0 ], [0.0], [0.0], [0.0], [0.0], [0.0], [0.0], [0.0], [0.0]])
mylr2.fit_(x2, y)
y_n.append(mylr2.cost_(x2,y))

mylr3.fit_(x3, y)
y_n.append(mylr3.cost_(x3,y))

mylr4.fit_(x4, y)
y_n.append(mylr4.cost_(x4,y))

mylr5.fit_(x5, y)
y_n.append(mylr5.cost_(x5,y))

mylr6.fit_(x6, y)
y_n.append(mylr6.cost_(x6,y))

mylr7.fit_(x7, y)
y_n.append(mylr7.cost_(x7,y))

mylr8.fit_(x8, y)
y_n.append(mylr8.cost_(x8,y))

mylr9.fit_(x9, y)
y_n.append(mylr9.cost_(x9,y))

plt.bar(range(2, 10), y_n)

plt.xlabel('n exponent')
plt.ylabel('cost')
plt.title('train nine separate Linear Regression models with polynomial hypotheses with degrees ranging from 2 to 10.')
plt.show()

x = np.arange(1,11).reshape(-1,1)
y = np.array([[ 1.39270298], [ 3.88237651], [ 4.37726357], [ 4.63389049], [ 7.79814439], [ 6.41717461], [ 8.63429886], [ 8.19939795], [10.37567392], [10.68238222]])
plt.scatter(x,y)
plt.show()

from polynomial_model import add_polynomial_features
from mylinearregression import MyLinearRegression as MyLR
# Build the model:
x_ = add_polynomial_features(x, 3)
my_lr = MyLR(np.ones(4).reshape(-1,1))
my_lr.fit_(x_, y)
## To get a smooth curve, we need a lot of data points
continuous_x = np.arange(1,10.01, 0.01).reshape(-1,1)
x_ = add_polynomial_features(continuous_x, 3)
y_hat = my_lr.predict_(x_)
plt.scatter(x,y)
# print(my_lr.thetas)
plt.plot(continuous_x, y_hat, color='orange')
plt.show()

고찰
우선 데이터를 정규화 했을 때,
x값이 0.1이라 하자
x^2, x^3 .. x^9
0.01 0.001 0.000000001
이 된다. 상대적으로 x^2보다 x^9이 aplha(learning rate)에 민감하게 반응할텐데 왜 pdf에서는 alpha를 동일하게 사용하는지 모르겠다 각각마다 learning rate를 다르게 주면 더 정확하게 측정할 수 있을 것 같다.
alpha가 너무 크면 overshoot가 일어나 발산이되고, alpha가 너무 작으면 수렴까지의 속도가 너무 오래걸린다.

따라서 작은 learning rate부터 알고리즘 연산을 해서 점점 큰 learning rate를 선택해야한다.
 ex)
0.001 -> 0.003 -> 0.01 -> 0.03 -> 0.1 -> 0.3 -> 1

이 문제에서 출제자의 의도는 다항식의 degree 가 높아질 수록 cost는 확실히 낮아지지만 이는 우리가 가지고 있는 데이터에 fit된것 뿐이고 전체 데이터와는 오차가 큰
overfit이 일어나는 상황을 의도해서
다항식의 차수가 높다고 해서 좋은 모델이 아니라는것을 알려주려고 하는 것 같다.
하지만 overfit인 상황을 만들기 위해선 많은 시간이 필요하다 ...

train을 빠르게 하기위해 반복문을 최적화 시켰다.
  while n_cycle:
            self.thetas -= self.alpha * ((self.gradient(x, y)))
            # for i, v in enumerate(self.gradient(x, y)):
            #     self.thetas[i] -= (self.alpha * v)
            n_cycle -= 1
        return self.thetas

x = np.arrange(10)
x[[8,9]]
>> [8,9]

x = np.arange(10)
y = np.arange(10)

s = np.arange(x.shape[0])
np.random.shuffle(s)
x = x[s]
y = y[s]

def data_spliter(x, y, proportion):
    """Shuffles and splits the dataset (given by x and y) into a training and a test set,
    while respecting the given proportion of examples to be kept in the traning set., →
    Args:
        x: has to be an numpy.ndarray, a matrix of dimension m * n.
    y: has to be an numpy.ndarray, a vector of dimension m * 1. proportion: has to be a float,
    the proportion of the dataset that will be assigned to the training set., →
    Returns: (x_train, x_test, y_train, y_test) as a tuple of numpy.ndarray None if x or y is an empty numpy.ndarray.
    None if x and y do not share compatible dimensions.
    Raises: This function should not raise any Exception. """
    n = int(float(x.shape[0] * proportion))
    idxs = np.arange(x.shape[0])
    np.random.shuffle(idxs)
    x_train = x[0:n]
    x_test = x[n:]
    y_train = y[0:n]
    y_test = y[n:]
    return (x_train, x_test, y_train, y_test)

csv_data = pd.read_csv("../resources/are_blue_pills_magics.csv")
y_n = []
x = np.array(csv_data["Micrograms"]).reshape(-1,1)
x = zscore(x)

y = np.array(csv_data["Score"]).reshape(-1,1)
y = zscore(y)

temp = data_spliter(x, y, 0.5)
x_train = temp[0]
x_test = temp[1]
y_train = temp[2]
y_test = temp[3]

x2 = add_polynomial_features(x_train, 2)
x3 = add_polynomial_features(x_train, 3)
x4 = add_polynomial_features(x_train, 4)
x5 = add_polynomial_features(x_train, 5)
x6 = add_polynomial_features(x_train, 6)
x7 = add_polynomial_features(x_train, 7)
x8 = add_polynomial_features(x_train, 8)
x9 = add_polynomial_features(x_train, 9)

x2_test = add_polynomial_features(x_test, 2)
x3_test = add_polynomial_features(x_test, 3)
x4_test = add_polynomial_features(x_test, 4)
x5_test = add_polynomial_features(x_test, 5)
x6_test = add_polynomial_features(x_test, 6)
x7_test = add_polynomial_features(x_test, 7)
x8_test = add_polynomial_features(x_test, 8)
x9_test = add_polynomial_features(x_test, 9)

mylr2 = MyLinearRegression([[88.85],[-9.0 ], [1]])
mylr3 = MyLinearRegression([[88.85],[-9.0 ], [1], [1]])
mylr4 = MyLinearRegression([[88.85],[-9.0 ], [1], [1], [1]])

mylr5 = MyLinearRegression([[88.85],[-9.0 ], [1], [1], [1], [1]])
mylr6 = MyLinearRegression([[88.85],[-9.0 ], [1], [1], [1], [1], [1]])

mylr7 = MyLinearRegression([[88.85],[-9.0 ], [1], [1], [1], [1], [1], [1]])
mylr8 = MyLinearRegression([[88.85],[-9.0 ], [1], [1], [1], [1], [1], [1], [1]])

mylr9 = MyLinearRegression([[88.85],[-9.0 ], [1], [1], [1], [1], [1], [1], [1], [1]])

mylr2.fit_(x2, y_train)
mylr2.alpha = 0.00001
mylr2.fit_(x2, y_train)
mylr2.alpha = 0.00003
mylr2.fit_(x2, y_train)
mylr2.alpha = 0.0001
mylr2.fit_(x2, y_train)
mylr2.alpha = 0.0003
mylr2.fit_(x2, y_train)
y_n.append(mylr2.cost_(x2_test,y_test))

mylr3.fit_(x3, y_train)
mylr3.alpha = 0.00001
mylr3.fit_(x3, y_train)
mylr3.alpha = 0.00003
mylr3.fit_(x3, y_train)
mylr3.alpha = 0.0001
mylr3.fit_(x3, y_train)
mylr3.alpha = 0.0003
mylr3.fit_(x3, y_train)
mylr3.alpha = 0.001
mylr3.fit_(x3, y_train)
y_n.append(mylr3.cost_(x3_test,y_test))

mylr4.fit_(x4, y_train)
mylr4.alpha = 0.00001
mylr4.fit_(x4, y_train)
mylr4.alpha = 0.00003
mylr4.fit_(x4, y_train)
mylr4.alpha = 0.0001
mylr4.fit_(x4, y_train)
mylr4.alpha = 0.0003
mylr4.fit_(x4, y_train)
mylr4.alpha = 0.001
mylr4.fit_(x4, y_train)
y_n.append(mylr4.cost_(x4_test,y_test))

mylr5.fit_(x5, y_train)
mylr5.alpha = 0.00001
mylr5.fit_(x5, y_train)
mylr5.alpha = 0.00003
mylr5.fit_(x5, y_train)
mylr5.alpha = 0.0001
mylr5.fit_(x5, y_train)
y_n.append(mylr5.cost_(x5_test,y_test))

mylr6.fit_(x6, y_train)
mylr6.alpha = 0.00001
mylr6.fit_(x6, y_train)
mylr6.alpha = 0.00003
mylr6.fit_(x6, y_train)
mylr6.alpha = 0.0001
mylr6.fit_(x6, y_train)
mylr6.alpha = 0.0003
mylr6.fit_(x6, y_train)
mylr6.alpha = 0.001
mylr6.fit_(x6, y_train)
y_n.append(mylr6.cost_(x6_test,y_test))

mylr7.fit_(x7, y_train)
mylr7.alpha = 0.00001
mylr7.fit_(x7, y_train)
mylr7.alpha = 0.00003
mylr7.fit_(x7, y_train)
mylr7.alpha = 0.0001
mylr7.fit_(x7, y_train)
mylr7.alpha = 0.0003
mylr7.fit_(x7, y_train)
mylr7.alpha = 0.001
mylr7.fit_(x7, y_train)
y_n.append(mylr7.cost_(x7_test,y_test))

mylr8.fit_(x8, y_train)
mylr8.alpha = 0.000003
mylr8.fit_(x8, y_train)
mylr8.alpha = 0.00001
mylr8.fit_(x8, y_train)
mylr8.alpha = 0.00003
mylr8.fit_(x8, y_train)
mylr8.alpha = 0.0001
mylr8.fit_(x8, y_train)
mylr8.alpha = 0.0003
mylr8.fit_(x8, y_train)
mylr8.alpha = 0.001
mylr8.fit_(x8, y_train)
y_n.append(mylr8.cost_(x8_test,y_test))

mylr9.alpha = 0.0000001
mylr9.fit_(x9, y_train)
mylr9.alpha = 0.0000003
mylr9.fit_(x9, y_train)
mylr9.alpha = 0.0000008
mylr9.fit_(x9, y_train)
mylr9.alpha = 0.000001
mylr9.fit_(x9, y_train)
mylr9.alpha = 0.000003
mylr9.fit_(x9, y_train)
mylr9.alpha = 0.000008
mylr9.fit_(x9, y_train)
mylr9.alpha = 0.00001
mylr9.fit_(x9, y_train)
mylr9.alpha = 0.00003
mylr9.fit_(x9, y_train)
mylr9.alpha = 0.00008
mylr9.fit_(x9, y_train)
mylr9.alpha = 0.0001
mylr9.fit_(x9, y_train)
mylr9.alpha = 0.0003
mylr9.fit_(x9, y_train)
mylr9.alpha = 0.0008
mylr9.fit_(x9, y_train)
mylr9.alpha = 0.001
mylr9.fit_(x9, y_train)
mylr9.alpha = 0.003
mylr9.fit_(x9, y_train)
y_n.append(mylr9.cost_(x9_test,y_test))
print(mylr9.thetas)
print(mylr9.cost_(x9_test,y_test))
plt.bar(range(2, 10), y_n)

plt.xlabel('n exponent')
plt.ylabel('cost')
plt.title('train nine separate Linear Regression models with polynomial hypotheses with degrees ranging from 2 to 10.')
plt.show()