梯度检测没有实现。有借鉴网上的部分
导入相关库,读取数据
因为这次的数据是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)