Matma_AI_cyber/Projekt_2/fixed/algorytm.ipynb

Algorytm najszybszego spadku dla regresji wielomianowej

Algorytm przyjmuje zbiór danych - x oraz y i próbuje wyzaczyć funkcję wilomianową, która najlepiej przewiduje wartości y na podstawie x. Wynikiem jest wyznaczenie współczynników wielomianu.

Importy

from matplotlib import pyplot as plt
import numpy as np
import random

Początkowe współczynniki

Tworzymy losowe początkowe współczynniki wielomianu - od nich algorytm rozpocznie dopasowanie. Oraz oryginalne współczynniki na podstawie których zostanie wyznaczony zbiór danych. Aby zmienić generowany zbiór danych należy zmienić tablicę coeffs.

coeffs = [2, -5, 4]

k = int(input("Podaj stopień wielomianu: "))
rand_coeffs = [ random.randrange(-10, 10) for i in range(k+1) ]

Wyznaczenie wartości wielomianu

Funkcja na podstawie współczynników oraz x wyznacza wartość y wielomianu.

def eval_polynomial(coeffs, x, k):
    y = 0

    for i in range(k+1):
        a = coeffs[i]*x**(k-i)
        y += a
    return y

Wartości wielomianu z szumem

Funkcja jest analogiczna do poprzedniej - wyznacza wartość wielomianu na podstawie wpółczynników oraz x, ale dodatkowo dodaje szum do wyjściowych wartości - funkcja zostanie użyta przy generowaniu danych.

def eval_polynomial_jitter(coeffs, x, j, k):
    y = eval_polynomial(coeffs, x, k)
    
    interval_min = y-j
    interval_max = y+j
    
    return random.uniform(interval_min, interval_max)

Wygenerowanie danych

Kod generuje zbiór danych. Na podstawie wartości x od -10 do 10 i losowych współczynników wielomianu generuje wartości y z szumem. Parametr j określa jak moco dane będą zaszumione

hundred_xs=np.random.uniform(-10,10,100)
print(hundred_xs)

j=50
x_y_pairs = []
for x in hundred_xs:
    y  = eval_polynomial_jitter(coeffs, x, j, len(coeffs)-1)
    x_y_pairs.append((x,y))
    
xs = []
ys = []
for a,b in x_y_pairs:
    xs.append(a)
    ys.append(b)
    
plt.figure(figsize=(20,10))
plt.plot(xs, ys, 'g+')
plt.title('Original data')
plt.show()

[ 9.52487237  4.31728957 -2.86004387 -4.06680766 -7.10993755 -6.3659089
  3.06102533  6.27143777  9.1380008   2.19180977  4.94214734  8.06342966
  4.58156633 -8.14890179  5.20520688 -8.59541259  3.6293768  -8.99512129
 -9.31440353 -5.42600609  7.82342339 -1.93271916 -0.81554086 -9.6994592
  5.13719864  9.45788656  3.70964434  1.27969457 -7.21538999  7.84821922
 -6.92608478 -2.35392424  8.32831454  3.61464536 -0.28674184 -2.87386245
  4.00068395  3.69653861  9.68378531 -7.278385   -7.90873047 -2.41469174
 -3.51370882  6.49444659  7.17689955 -4.77777697 -2.97134134 -0.6338146
  0.75355793 -9.65106708 -8.89989426  6.99530143  1.12347522  0.50798599
 -9.02110719  9.84011966 -6.16446615 -2.1919846   6.83936162  7.01254857
 -5.30766773  4.18034124 -9.68850209 -5.30243924  8.6308934   6.94886296
 -3.74310371 -7.17015727  0.50585453  6.5284797  -4.27368615  0.68288962
  7.23974203  8.16376838 -9.54528187 -8.04775584 -0.09982886 -6.10924399
 -9.93495247  2.91907289 -2.923125   -3.53381594  1.95797436 -7.16019353
  2.04554788 -1.98013462  7.88027302  6.2077106   2.85884957 -1.95785965
 -7.63680301  8.4529672  -7.24006824  2.41572699  8.19679144 -8.88831243
 -8.57612565  7.081674   -5.74549524 -1.44629504]

Funkcja straty

Do określenia jak mocno przewidziane wartości y są różne zostanie użyta kwadratowa funkcja straty.

def loss_mse(ys, y_bar):
    return sum((ys - y_bar)*(ys - y_bar)) / len(ys)

Gradient

Funkcja przyjmuje współczynniki wielomianu, x, y oraz parametr prękości uczenia i wylicza za pomocą gradientu nowe wartości współczynników wielomianu.

def calc_gradient_poly_for_GD(coeffs, inputs_x, outputs_y, lr, k): 
    partial_derivatives = [[] for i in range(k+1)]
    new_coeffs = [None for i in range(k+1)]
        
    y_bars = eval_polynomial(coeffs, inputs_x, k)

    for x,y,y_bar in list(zip(inputs_x, outputs_y, y_bars)):
        for i in range(k+1):
            partial = x**(k-i) * (y - y_bar)
            partial_derivatives[i].append(partial)
    
    n = len(y_bars)
    
    for i in range(k+1):
        derivative = (-2 / n) * sum(partial_derivatives[i])
        new_coeffs[i] = coeffs[i] - lr * derivative
        
    
    #update with these new coeffs:
    new_y_bar = eval_polynomial(new_coeffs, inputs_x, k)
    
    updated_model_loss = loss_mse(outputs_y, new_y_bar)
    return updated_model_loss, new_coeffs, new_y_bar

Minimalizacja

Funkcja powtarza proces wyznaczenia nowych współczynników wielomianu zadaną ilość razy - epoch.

def gradient_descent(epochs, lr, k):
    losses = []
    rand_coeffs_to_test = rand_coeffs
    for i in range(epochs):
        loss = calc_gradient_poly_for_GD(rand_coeffs_to_test, hundred_xs, ys, lr, k)
        rand_coeffs_to_test = loss[1]
        losses.append(loss[0])
    #print(losses)
    return loss[0], loss[1], loss[2], losses  #(updated_model_loss, new_model_coeffs, new_y_bar, saved loss updates)

Uruchomienie

Wartości wyznaczonego wielomionu zostają wyświetlone na wykresie z porównaniem do zbioru danych.
UWAGA gdy otrzymamy błąd double scalars - przekroczenie wartości, należy ustawić mniejszy learning rate.

GD = gradient_descent(1500, 0.0001, k)

plt.figure(figsize=(20,10))
plt.plot(xs, ys, 'g+', label = 'original')
plt.plot(xs, GD[2], 'b.', label = 'final_prediction')
plt.title('Original vs Final prediction after Gradient Descent')
plt.legend(loc="lower right")
plt.show()

Wyznaczone współczynniki wielomianu

print("Początkowe współczynniki {}".format(rand_coeffs))
print("Wyznaczone współczynniki {}".format(GD[1]))
print("Oryginalne współczynniki {}".format(coeffs))

Początkowe współczynniki [-8, -1, -3]
Wyznaczone współczynniki [2.0470388756371727, -5.710817671350023, -1.943138081266058]
Oryginalne współczynniki [2, -5, 4]

Funkcja straty

Wykres przedstawia jak zmianiała się wartość funkcji straty w kolejnych krokach algorytmu.

plt.figure(figsize=(20,10))
plt.plot(GD[3], 'b-', label = 'loss')
plt.title('Loss over 1500 iterations')
plt.legend(loc="lower right")
plt.xlabel('Iterations')
plt.ylabel('MSE')
plt.show()