1. # --------------------------------------------------------------------------
2. # ------------  Metody Systemowe i Decyzyjne w Informatyce  ----------------
3. # --------------------------------------------------------------------------
4. #  Zadanie 1: Regresja liniowa
5. #  autorzy: A. Gonczarek, J. Kaczmar, S. Zareba
6. #  2017
7. # --------------------------------------------------------------------------
8. from _operator import inv
9. from numpy.linalg import inv
10. import numpy as np
11. from utils import polynomial
12.  
13.  
14. def mean_squared_error(x, y, w):
15.     '''
16.     :param x: ciag wejsciowy Nx1
17.     :param y: ciag wyjsciowy Nx1
18.     :param w: parametry modelu (M+1)x1
19.     :return: blad sredniokwadratowy pomiedzy wyjsciami y
20.     oraz wyjsciami uzyskanymi z wielowamiu o parametrach w dla wejsc x
21.     '''
22.     error = 0.0
23.     for i in range(x.size):
24.         error += (y[i]-polynomial(x[i], w))**2
25.     return float(error/x.size)
26.  
27. def design_matrix(x_train, M):
28.     '''
29.     :param x_train: ciag treningowy Nx1
30.     :param M: stopien wielomianu 0,1,2,...
31.     :return: funkcja wylicza Design Matrix Nx(M+1) dla wielomianu rzedu M
32.     '''
33.     matrix = np.empty((x_train.size, (M+1)))
34.     for i in range(x_train.size):
35.         for j in range(M + 1):
36.             matrix[i, j] = x_train[i]**j
37.     return matrix
38.  
39. def least_squares(x_train, y_train, M):
40.     '''
41.     :param x_train: ciag treningowy wejscia Nx1
42.     :param y_train: ciag treningowy wyjscia Nx1
43.     :param M: rzad wielomianu
44.     :return: funkcja zwraca krotke (w,err), gdzie w sa parametrami dopasowanego wielomianu, a err blad sredniokwadratowy
45.     dopasowania
46.     '''
47.     fi = design_matrix(x_train, M)
48.     w = inv((fi.transpose() @ fi)) @ fi.transpose() @ y_train
49.     return w, mean_squared_error(x_train, y_train, w)
50.  
51. def regularized_least_squares(x_train, y_train, M, regularization_lambda):
52.     '''
53.     :param x_train: ciag treningowy wejscia Nx1
54.     :param y_train: ciag treningowy wyjscia Nx1
55.     :param M: rzad wielomianu
56.     :param regularization_lambda: parametr regularyzacji
57.     :return: funkcja zwraca krotke (w,err), gdzie w sa parametrami dopasowanego wielomianu zgodnie z kryterium z regularyzacja l2,
58.     a err blad sredniokwadratowy dopasowania
59.     '''
60.     fi = design_matrix(x_train, M)
61.     fi_matrix = (fi.transpose() @ fi)
62.     lambda_matrix = np.eye(fi_matrix.shape[0]) * regularization_lambda
63.     w = inv(fi_matrix + lambda_matrix) @ fi.transpose() @ y_train
64.     return w, mean_squared_error(x_train, y_train, w)
65.  
66. def model_selection(x_train, y_train, x_val, y_val, M_values):
67.     '''
68.     :param x_train: ciag treningowy wejscia Nx1
69.     :param y_train: ciag treningowy wyjscia Nx1
70.     :param x_val: ciag walidacyjny wejscia Nx1
71.     :param y_val: ciag walidacyjny wyjscia Nx1
72.     :param M_values: tablica stopni wielomianu, ktore maja byc sprawdzone
73.     :return: funkcja zwraca krotke (w,train_err,val_err), gdzie w sa parametrami modelu, ktory najlepiej generalizuje dane,
74.     tj. daje najmniejszy blad na ciagu walidacyjnym, train_err i val_err to bledy na sredniokwadratowe na ciagach treningowym
75.     i walidacyjnym
76.     '''
77.  
78.     w_score = []
79.     w_error_score = []
80.     for i in M_values:
81.         w, train_err = least_squares(x_train, y_train, i)
82.         val_err = mean_squared_error(x_val, y_val, w)
83.         w_score.append((w, train_err, val_err))
84.         w_error_score.append(val_err)
85.  
86.     min_index = w_error_score.index(min(w_error_score))
87.     return w_score[min_index]
88.  
89.  
90. def regularized_model_selection(x_train, y_train, x_val, y_val, M, lambda_values):
91.     '''
92.     :param x_train: ciag treningowy wejscia Nx1
93.     :param y_train: ciag treningowy wyjscia Nx1
94.     :param x_val: ciag walidacyjny wejscia Nx1
95.     :param y_val: ciag walidacyjny wyjscia Nx1
96.     :param M: stopien wielomianu
97.     :param lambda_values: lista ze wartosciami roznych parametrow regularyzacji
98.     :return: funkcja zwraca krotke (w,train_err,val_err,regularization_lambda), gdzie w sa parametrami modelu, ktory najlepiej generalizuje dane,
99.     tj. daje najmniejszy blad na ciagu walidacyjnym. Wielomian dopasowany jest wg kryterium z regularyzacja. train_err i val_err to
100.     bledy na sredniokwadratowe na ciagach treningowym i walidacyjnym. regularization_lambda to najlepsza wartosc parametru regularyzacji
101.     '''
102.     w_score = []
103.     w_error_score = []
104.     for i in range(len(lambda_values)):
105.         w, train_err = regularized_least_squares(x_train, y_train, M, lambda_values[i])
106.         val_err = mean_squared_error(x_val, y_val, w)
107.         w_score.append((w, train_err, val_err, lambda_values[i]))
108.         w_error_score.append(val_err)
109.  
110.     min_index = w_error_score.index(min(w_error_score))
111.     return w_score[min_index]
112.