読者です 読者をやめる 読者になる 読者になる

けつあご日記

こんちは.菊地です.なんか浅ーいことをいろいろ書きます

Coursera ML ex1をpythonで(バカ正直に)やってみた【なれない日記20160728】

なれない日記 機械学習 Python

Coursera Machine Learningの課題ex1を,特にライブラリなどを使わずにバカ正直にやってみた版.バカ正直にやり過ぎてほぼOctaveのときと変わらない.Coursera Honor Codeに触れそう...大丈夫かな...

ライブラリを駆使したりするのは今度やる.その場合ココらへんが参考になりそうだ.
Coursera Machine Learningの課題をPythonで: ex1(線形回帰) - Qiita
Coursera / Machine Learningの教材を2度楽しむ - Qiita

とりあえずこんな感じ.

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

data = np.genfromtxt("ex1data1.txt", delimiter=',')
x = data[:, 0]
y = data[:, 1]

X = np.transpose(np.array([np.ones_like(x), x]))
theta = np.zeros(2)
m = len(y)

iterations = 1500
alpha = 0.01 # learning rate

def compute_cost(X, y, theta):
    J = sum((np.dot(X, theta) - y) ** 2) /(2 * m)
    return J

# Update theta and record history of cost function
def gradient_descent(X, y, theta, alpha, num_iters):
    J_history = np.zeros(num_iters)

    for i in range(num_iters):
        theta = theta - alpha / m * np.dot(np.transpose(X), (np.dot(X, theta) - y))
        J_history[i] = compute_cost(X, y, theta)
    
    return(theta, J_history)

theta, J_history = gradient_descent(X, y, theta, alpha, iterations)


# Prediction examples
predict1 = np.dot([1, 3.5], theta)
predict2 = np.dot([1, 7], theta)
print('For population = 35,000, we predict a profit of\n' + '...' + str(predict1*10000))
print('For population = 70,000, we predict a profit of\n' + '...' + str(predict2*10000))


# Calculate cost function about ALL theta
theta0_vals = np.arange(-10, 10, 0.2)
theta1_vals = np.arange(-1, 4, 0.05)
J_vals = np.zeros((len(theta0_vals), len(theta1_vals)))

for i in range(len(theta0_vals)):
  for j in range(len(theta1_vals)):
    t = np.array([theta0_vals[i], theta1_vals[j]])
    J_vals[i, j] = compute_cost(X, y, t)


fig = plt.figure(figsize=(10,8))

# Data plot and Regression line
ax = fig.add_subplot(2, 2, 1)
ax.plot(x, y, 'x', label='Traning data')
ax.plot(X[:, 1], np.dot(X, theta), label='Linear regression')
ax.legend(loc='lower right', shadow=True, fontsize='small')
ax.set_title('Regression with Gradient Descent')
ax.set_xlim(4, 25)
ax.set_ylim(-5, 30)
ax.set_xlabel('Population of City in 10,000s')
ax.set_ylabel('Profit in $10,000s')

# Cost function v.s. Iterations
ax = fig.add_subplot(2, 2, 2)
ax.plot(np.arange(iterations), J_history, color='green')
ax.set_title('Convergence of gradient descent')
ax.set_xlabel('No. of iterations')
ax.set_ylabel('Cost function')

# Cost function v.s. theta (3D)
J_vals = np.transpose(J_vals) # NEED TRANSPOSE!!!
ax = fig.add_subplot(2, 2, 3, projection='3d')
theta0_grid, theta1_grid = np.meshgrid(theta0_vals, theta1_vals)
ax.plot_wireframe(theta0_grid, theta1_grid, J_vals, rstride=5, cstride=5)
ax.set_title('3D Cost function')
ax.set_xlabel(r'$\theta_0$')
ax.set_ylabel(r'$\theta_1$')
ax.set_zlabel('Cost function')

# Cost function v.s. theta (contour)
ax = fig.add_subplot(2, 2, 4)
CS = plt.contour(theta0_grid, theta1_grid, J_vals, 100)
ax.set_title('Contour of cost function')
ax.set_xlabel(r'$\theta_0$')
ax.set_ylabel(r'$\theta_1$')

fig.subplots_adjust(hspace=0.3, wspace=0.3)

plt.show()

図はこんな感じ.

f:id:kichiku_kikuchi:20160728005745p:plain

うん.すべてがそのままだ.まぁライブラリを使ったときの恩恵が強く感じられるだろうから良いことにしよう.



↓アフィカスリンク

ITエンジニアのための機械学習理論入門

ITエンジニアのための機械学習理論入門

Python機械学習プログラミング 達人データサイエンティストによる理論と実践 (impress top gear)

Python機械学習プログラミング 達人データサイエンティストによる理論と実践 (impress top gear)