吴恩达老师机器学习-ex4

梯度检测没有实现。有借鉴网上的部分

导入相关库，读取数据

因为这次的数据是mat文件，需要使用scipy库中的loadmat进行读取数据。

通过对数据类型的分析，发现是字典类型，查看该字典的键，可以发现又X，y等关键字。

import numpy as np
import scipy.io as sio
from scipy.optimize import minimize
from sklearn.preprocessing import OneHotEncoder

#读取数据
path = "./ex4data1.mat"
data = sio.loadmat(path)
# print(type(data))
# print(data.keys())
X = data.get("X")
Y = data.get("y").flatten()
# X = np.insert(X,0,values=1,axis=1)
# print(Y.shape)
# print(Y)

one-hot编码

在之前没有涉及神经网络的方向传播等时，可以直接将Y降维。

而在神经网络中需要对Y进行编码，使每一个y值都是一个10维的向量

#one-hot编码
encoder = OneHotEncoder(sparse=False)
Y_onehot = encoder.fit_transform(Y.reshape(-1,1))

随机初始化

当使用高级优化算法或者梯度下降算法时，需要对向量进行赋值。

第一种想法是全部设为0，在逻辑回归中，是完全被允许的，但在训练神经网络时，这种初始值起不到任何作用。因为这会导致第二层激活项的值完全相同（因为前一层的权重相同），这个问题称为对称权重问题。

为了解决这个问题，我们采用随机初始化的方法，也就是将每一个都初始化为一个范围在中的一个随机数。

#初始化参数
input_size = 400
hidden_size = 25
num_labels = 10
lamda = 1

size = (input_size+1)*hidden_size+(hidden_size+1)*num_labels
params = np.random.uniform(-1.2,1.2,size)

序列化

因为到后面优化的时候fun和jac参数要求代价函数和梯度的第一个参数必须为一维向量，所以这里需要序列化。

#序列化
def serialize(theta1,theta2):
    return np.r_[theta1.flatten(),theta2.flatten()]
def deserialize(params):
    return params[:(input_size+1)*hidden_size].reshape(hidden_size,input_size+1),params[(input_size+1)*hidden_size:].reshape(num_labels,hidden_size+1)

theta1,theta2 = deserialize(params)

前向传播

同上一题一致的前向传播

#前向传播
def feed_forward(params,X):
    theta1, theta2 = deserialize(params)
    a1 = np.insert(X, 0, values=np.ones(X.shape[0]), axis=1)
    z2 = a1 @ theta1.T
    a2 = 1/(1+np.exp(-z2))
    a2 = np.insert(a2,0,values=1,axis=1)
    z3 = a2 @ theta2.T
    h = 1/(1+np.exp(-z3))
    return a1,z2,a2,z3,h

代价函数

不需要正则化

公式如下：

#代价函数
def cost_func(params,X,Y,lamda):
    theta1, theta2 = deserialize(params)
    a1,z2, a2, z3, h = feed_forward(params,X)
    m = len(X)
    cost = -np.sum(Y * np.log(h) + (1-Y) * np.log(1-h))/m
    reg = (np.sum(np.power(theta1[:,1:],2)) + np.sum(np.power(theta2[:,1:],2)))* lamda/(2*m)
    return reg+cost

反向传播

相关公式如下：

#反向传播
def sigmoid_gradient(z):
    return 1/(1+np.exp(-z))*(1-1/(1+np.exp(-z)))

def gradientReg(params,X,Y,lamda=1):
    theta1, theta2 = deserialize(params)
    m = len(X)
    a1, z2, a2, z3, h = feed_forward(params,X)
    delta3 = h - Y
    delta2 = delta3 @ theta2[:,1:] * sigmoid_gradient(z2)
    Delta2 = delta3.T @ a2 / m
    Delta1 = delta2.T @ a1 / m
    theta1[:, 0] = 0
    theta2[:, 0] = 0
    regDelta1 = Delta1 + (lamda / m) * theta1
    regDelta2 = Delta2 + (lamda / m) * theta2
    return serialize(regDelta1,regDelta2)

优化

fmin = minimize(fun=cost_func,
                x0=params,
                args=(X,Y_onehot,lamda),
                method="TNC",
                jac=gradientReg,
                options={'maxiter': 400})

预测

a1,z2,a2,z3,h = feed_forward(fmin.x,X)
y_pred = np.argmax(h,axis=1)
y_pred = y_pred + 1
acc = np.mean(y_pred==Y)
print(acc)

准确率达到了98.64%

---------------------------------------------------------------------------------------------------------------------------------

可以先将数据可视化，展示手写数字。

def one_image(X):
    pick_one = np.random.randint(5000)
    image = X[pick_one,:]
    fig,ax = plt.subplots()
    ax.imshow(image.reshape(20,20).T,cmap="gray_r")
    plt.xticks([])
    plt.yticks([])
    plt.show()

one_image(X)
def more_image(X):
    pick_more = np.random.choice(len(X),100)
    images = X[pick_more,:]
    fig,ax = plt.subplots(nrows=10,ncols=10, figsize=(8, 8), sharex=True, sharey=True)
    for i in range(10):
        for j in range(10):
            ax[i,j].imshow(images[10*i+j].reshape(20,20).T,cmap="gray_r")
    plt.xticks([])
    plt.yticks([])
    plt.show()
more_image(X)

最后可以使用上面写出的函数将隐藏层可视化。

thetafinal1, thetafinal2 = deserialize(fmin.x)
hidden_layer = thetafinal1[:, 1:] 
more_image(hidden_layer)