Uczenie maszynowe i sztuczne sieci neuronowe/Ćwiczenia 1

Z Brain-wiki

Powtórka podstawowych rachunków wektorowych i macierzowych w ptyhonie

Wektor wierszowy i kolumnowy

Natywnym typem zmiennych w numpy są tablice czyli array. Aczkolwiek są one wielowymiarowe i przy pomocy indeksowania i pobierania wycinków można się nimi sprawnie posługiwać to nie są one domyślnie macierzami w sensie matematycznym. Aby uprawiać przy ich pomocy algebrę musimy nadać im kształt :-). Służy do tego metoda reshape. Aby zbadać własności tych obiektów proszę wykonać następujący kod:

import numpy as np
x = np.array([1,2,3,4])
print len(x), x.shape
print 'x: ',x
print 'transponowany x: ',x.T

Czy tablica x i transponowana tablica x się różnią?

A teraz proszę wykonać następujący kod:

import numpy as np
x = np.array([1,2,3,4]).reshape(4,1)
print len(x), x.shape
print 'x: ',x
print 'transponowany x: ',x.T

Proszę sprawdzić kształt i transpozycję macierzy 2x2, np.

[math]A = \left[ \begin{array}{cc} 1 & 2 \\ 3 & 4 \end{array} \right] [/math]
A = np.array([[1,2],[3,4]])

Operator * służy do mnożenia macierzy element po elemencie, albo do mnożenia macierzy przez skalar:

print A*A
print 2*A
v = np.array([3])
print v*A

Do wykonywania mnożenia w sensie macierzowym służy funkcja np.dot:

print np.dot(A,A)

Co robią następujące polecenia:

import numpy as np
x = np.array([1,2,3,4]).reshape(4,1)
print np.dot(x.T,x)
print np.dot(x,x.T)

Zapoznanie się z regresją liniową

W poniższych ćwiczeniach korzystać będziemy z następującej funkcji generującej dane. Proszę ją wykonać i obejrzeć dane. (X,Y) to ciąg uczący.

# -*- coding: utf-8 -*-
import scipy.stats as st
import pylab as py
import numpy as np


def symuluj_dane(theta, ile=30, odch_std = 3):
	# symuowana zalezność ma nastepującą postać y = theta[0] + theta[1]*x 
	# wartości parametrów  
	X = np.ones((ile,2)).reshape(ile,2) # X ma być macierzą, w której wierszach zapisane są kolejne wejścia ciągu
	X[:,1] = np.linspace(0, 10,ile)
	#błąd 
	epsilon = st.norm.rvs(size = ile, loc=0, scale = odch_std).reshape(ile,1)
	#epsilon = (st.uniform.rvs(size = ile, loc=-0.5, scale = odch_std)**3).reshape(ile,1)
		
	Y =  np.dot(X,theta) + epsilon 
	return (X,Y)

theta = np.array([1.,3.]).reshape(2,1) # theta ma być wektorem kolumnowym
(X,Y) = symuluj_dane(theta,  ile = 7, odch_std = 3) #symulujemy dane
py.plot(X[:,1], Y,'bo') # przedstawiamy nasze dane na wykresie
py.show()

Algorytm równań normalnych

Proszę napisać funkcję, która:

  • na wejściu przyjmuje ciąg uczący, implementuje wzór na parametry optymalne na podstawie równań normalnych.
  • Funkcja powinna zwracać estymowane parametry theta.
  • Proszę dorysować prostą reprezentującą hipotezę do wykresu punktów ciągu uczącego.
  • dla przypomnienia: odwrotność macierzy można obliczyć w numpy funkcją: numpy.linalg.inv

Algortym gradientowy zbiorczy

Proszę napisać funkcję, która znajduje optymalne parametry theta wg. algorytmy gradientowego zbiorczego. Funkcja jako argumenty przyjmuje ciąg uczący, wartości początkowe theta i parametr szybkości zbiegania alpha. Na wyjściu funkcja powinna zwracać wyestymowane wartości parametrów.

W ramach ilustracji po każdej iteracji proszę dorysować prostą parametryzowaną przez aktualne wartości parametrów. Prostą animację można wykonać analogicznie jak tu.

Algorytm gradientowy stochastyczny

Proszę napisać funkcję, która znajduje optymalne parametry theta wg. algorytmu gradientowego stochastycznego. Funkcja jako argumenty przyjmuje ciąg uczący, wartości początkowe theta i parametr szybkości zbiegania alpha. Na wyjściu funkcja powinna zwracać wyestymowane wartości parametrów.

W ramach ilustracji po każdej iteracji proszę dorysować prostą parametryzowaną przez aktualne wartości parametrów. Prostą animację można wykonać analogicznie jak tu.

Porównanie algorytmów

  • Proszę sprawdzić zbieżność algorytmów w zależności od parametrów szybkości zbiegania.
  • Proszę sprawdzić czy algorytmy optymalizacyjne działają poprawnie dla danych gdzie błąd podlega innym rozkładom prawdopodobieństwa niż normalny. np. rozkład jednorodny lub t o 3 st. swobody.
*


Przykładowe kody

# -*- coding: utf-8 -*-
import matplotlib
matplotlib.use('TkAgg')
import scipy.stats as st
import pylab as py
import numpy as np


def symuluj_dane(theta, ile=30, odch_std = 3):
	# symuowana zalezność ma nastepującą postać y = theta[0] + theta[1]*x 
	# wartości parametrów  
	X = np.ones((ile,2)).reshape(ile,2) # X ma być macierzą, w której wierszach zapisane są kolejne wejścia ciągu
	X[:,1] = np.linspace(0, 10,ile)
	#błąd 
	epsilon = st.norm.rvs(size = ile, loc=0, scale = odch_std).reshape(ile,1)
	#epsilon = (st.uniform.rvs(size = ile, loc=-0.5, scale = odch_std)**3).reshape(ile,1)
		
	Y =  np.dot(X,theta) + epsilon 
	return (X,Y)

def licz_rownania_normalne(X,Y):
	theta = np.dot( np.dot(np.linalg.inv( np.dot(X.T,X)  ),X.T) ,Y)
	return theta

def koszt(X,Y,theta):
	return 0.5* np.dot( ( np.dot(X,theta) - Y).T, ( np.dot(X,theta) -Y) )

def licz_iteracyjnie_stoch(X,Y,theta0 = np.array([0,0]).reshape(2,1), alpha = 0.01):
	
	py.ion()	
	py.plot(X[:,1], Y,'bo')
	y_reg = np.dot(X,theta0)
	linia, = py.plot(X[:,1],y_reg,'r')
	for i in range(200):
		ind = np.random.randint(X.shape[0])		
		# mały x to wejście pojedynczego przykładu o indeksie in
		x = X[ind,:]		
		x = x.reshape(1,len(x))
		theta0 = theta0 - alpha * ( np.dot(x,theta0)- Y[ind])*x.T
		#print 'J: ', koszt(X,Y,theta)
		linia.set_ydata(np.dot(X,theta0))
		py.draw()
	return theta0

def licz_iteracyjnie_batch(X,Y,theta0 = np.array([0,0]).reshape(2,1), alpha = 0.005):
	py.ion()	
	py.plot(X[:,1], Y,'bo')
	y_reg = np.dot(X,theta0)
	linia, = py.plot(X[:,1],y_reg,'m')
	delta = np.ones(len(theta0)).reshape(len(theta0),1)
	while np.max(np.abs(delta)) > 0.001:
		delta = np.zeros(len(theta0)).reshape(len(theta0),1)
		for ind, x in enumerate(X):		
			# mały x to wejście pojedynczego przykładu o indeksie in
			x = x.reshape(1,len(x))
			delta += ( np.dot(x,theta0)- Y[ind])*x.T
		theta0 = theta0 - alpha * delta 
		#print 'J: ', koszt(X,Y,theta)
		linia.set_ydata(np.dot(X,theta0))
		py.draw()
	return theta0
	
theta = np.array([1.,3.]).reshape(2,1) # theta ma być wektorem kolumnowym
(X,Y) = symuluj_dane(theta,  ile = 7, odch_std = 3)

theta_est = licz_rownania_normalne(X,Y)
theta_est_iter_stoch = licz_iteracyjnie_stoch(X,Y)
theta_est_iter_batch = licz_iteracyjnie_batch(X,Y)

print theta
print theta_est
print theta_est_iter_stoch
print theta_est_iter_batch
py.plot(X[:,1], Y,'bo')
y_reg = np.dot(X,theta_est) 
py.plot(X[:,1],y_reg,'g')
py.show()