经常被np 的axis 搞的乱七八糟。。

在numpy中数组都有着[]标记，则axis=0对应着最外层的[]，axis=1对应第二外层的[]，以此类推，axis=n对应第n外层的[]。计算完后，移除最外层[]。脑补链接

前向传播

import matplotlib.pyplot as plt
import numpy as np
import scipy.io as sio
import matplotlib
import scipy.optimize as opt
from sklearn.metrics import classification_report#这个包是评价报告def load_data(path, transpose=True):data = sio.loadmat(path)y = data.get('y')  # (5000,1)y = y.reshape(y.shape[0])  X = data.get('X')  # (5000,400)if transpose:X = np.array([im.reshape((20, 20)).T for im in X])X = np.array([im.reshape(400) for im in X])  #  (5000,400) 应该跟上面那个X一样的。。return X, yX_raw, y_raw = load_data('ex4data1.mat', transpose=False)
X = np.insert(X_raw, 0, np.ones(X_raw.shape[0]), axis=1) #增加全部为1的一列
print(X.shape,y_raw.shape)
print(np.unique(y_raw))
>>>
(5000, 401) (5000,)
[ 1  2  3  4  5  6  7  8  9 10]     # 十个标签# from sklearn.preprocessing import OneHotEncoder
# encoder = OneHotEncoder(sparse=False)
# y_onehot = encoder.fit_transform(y)
# y_onehot.shapedef expand_y(y):
#     where y=10 -> [0, 0, 0, 0, 0, 0, 0, 0, 0, 1]: ndarray,就这个效果res = []for i in y:y_array = np.zeros(10)y_array[i - 1] = 1res.append(y_array)return np.array(res)
# 跟前面的 OneHotEncoder 一个效果。。
y = expand_y(y_raw)
y.shape >>>(5000, 10)def load_weight(path):           # 读取权重data = sio.loadmat(path)return data['Theta1'], data['Theta2']
t1, t2 = load_weight('ex4weights.mat')
t1.shape, t2.shape   >>> ((25, 401), (10, 26))def serialize(a, b):  # 序列化，，，，return np.concatenate((np.ravel(a), np.ravel(b)))   # 在这个nn架构中，我们有theta1（25,401），theta2（10,26），它们的梯度是delta1，delta2  
theta = serialize(t1, t2)  # 扁平化参数，25*401+10*26=10285，所有的theta。
theta.shape  >>>(10285,)def sigmoid(z):return 1 / (1 + np.exp(-z))def deserialize(seq):   # 这个也要通过一维的 theta 在把两个权重矩阵搞回来。。return seq[:25 * 401].reshape(25, 401), seq[25 * 401:].reshape(10, 26)def feed_forward(theta, X):          # 前向传播t1, t2 = deserialize(theta)  # t1: (25,401) t2: (10,26)m = X.shape[0]a1 = X          # 5000 * 401z2 = a1 @ t1.T  # 5000 * 25a2 = np.insert(sigmoid(z2), 0, np.ones(m), axis=1)  # 5000*26z3 = a2 @ t2.T  # 5000 * 10h = sigmoid(z3)  # 5000*10, this is h_theta(X)return a1, z2, a2, z3, h    # 就用这些来反向
_, _, _, _, h = feed_forward(theta, X)   
np.unique(h.argmax(axis=1)) # 5000*10
>>> array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9], dtype=int64)  # 从最大值的索引来看。。def cost(theta, X, y):      # 代价函数m = X.shape[0]  _, _, _, _, h = feed_forward(theta, X) # h 5000,10  X (5000, 401) y (5000, 10)pair_computation = -np.multiply(y, np.log(h)) - np.multiply((1 - y), np.log(1 - h))# pari  5000,10 return pair_computation.sum() / m  # 进行 m 5000 and K 10 的求和cost(theta, X, y) >>>0.28762916516131892

def regularized_cost(theta, X, y, l=1):           # 正则化代价函数
# t1: (25,401) t2: (10,26)  t1,t2 的第一列是截距项，正则化是忽略t1, t2 = deserialize(theta)  # t1: (25,401) t2: (10,26)m = X.shape[0]reg_t1 = (l / (2 * m)) * np.power(t1[:, 1:], 2).sum()  # this is how you ignore first colreg_t2 = (l / (2 * m)) * np.power(t2[:, 1:], 2).sum()return cost(theta, X, y) + reg_t1 + reg_t2   # 就是加了两项。。regularized_cost(theta, X, y)   >>>0.38376985909092365   # 代价确实上升了。。

反向传播

X.shape,y.shape，t1.shape, t2.shape， theta.shape
>>> ((5000, 401), (5000, 10), (25, 401), (10, 26), (10285,))def sigmoid_gradient(z):return np.multiply(sigmoid(z), 1 - sigmoid(z))def gradient(theta, X, y):    # 注意各种维度的转换。。。吐了吐了t1, t2 = deserialize(theta)  # t1: (25,401) t2: (10,26)m = X.shape[0]delta1 = np.zeros(t1.shape)  # (25, 401)  # 后一项都是加了偏置单元delta2 = np.zeros(t2.shape)  # (10, 26)a1, z2, a2, z3, h = feed_forward(theta, X)for i in range(m):a1i = a1[i, :]  # (1, 401)z2i = z2[i, :]  # (1, 25)a2i = a2[i, :]  # (1, 26)hi = h[i, :]    # (1, 10)yi = y[i, :]    # (1, 10)d3i = hi - yi   # (1, 10)z2i = np.insert(z2i, 0, np.ones(1))  # make it (1, 26) to compute d2id2i = np.multiply(t2.T @ d3i, sigmoid_gradient(z2i))  # (1, 26) 一个@ 一个点乘# a2i 的导# careful with np vector transposedelta2 += np.matrix(d3i).T @ np.matrix(a2i)  # (1, 10).T @ (1, 26) -> (10, 26)delta1 += np.matrix(d2i[1:]).T @ np.matrix(a1i)  # (1, 25).T @ (1, 401) -> (25, 401)# 这个维度确实乱，，不过确实补充了没有讲到的东西。。。累加偏导，就是算得了所有theta的偏导。维度都是对应的嘛。。delta1 = delta1 / mdelta2 = delta2 / mreturn serialize(delta1, delta2)    d1, d2 = deserialize(gradient(theta, X, y))
d1.shape, d2.shape
>>> ((25, 401), (10, 26))def regularized_gradient(theta, X, y, l=1):        # 加上正则项"""don't regularize theta of bias terms"""m = X.shape[0]delta1, delta2 = deserialize(gradient(theta, X, y))t1, t2 = deserialize(theta)   # ((25, 401), (10, 26))t1[:, 0] = 0reg_term_d1 = (l / m) * t1      # t1, t2 对应两层，delta1 = delta1 + reg_term_d1t2[:, 0] = 0reg_term_d2 = (l / m) * t2delta2 = delta2 + reg_term_d2return serialize(delta1, delta2)

训练我们的模型

def random_init(size):     # 随机初始化，破坏对称性return np.random.uniform(-0.12, 0.12, size)  # 从一个均匀分布中采样def nn_training(X, y):init_theta = random_init(10285)        # 25*401 + 10*26res = opt.minimize(fun=regularized_cost,x0=init_theta,args=(X, y, 1),method='TNC',jac=regularized_gradient,options={'maxiter': 400})return resres = nn_training(X, y)   # 慢，十分钟，，终于有机器学习的感觉了，，pytorch 要赶快安排上了。。
res
>>>fun: 0.312513531441238jac: array([ 3.78273626e-04, -1.92256503e-07,  2.48737126e-07, ...,-1.32757264e-04, -1.18317799e-04, -1.63580336e-04])message: 'Max. number of function evaluations reached'nfev: 402nit: 25status: 3success: Falsex: array([-2.76786178e+00, -9.61282515e-04,  1.24368563e-03, ...,   # 这个就是我们的theta9.54831564e-01,  9.78685642e-01,  8.95277691e-01])_, y_answer = load_data('ex4data1.mat')
final_theta = res.x
def show_accuracy(theta, X, y):_, _, _, _, h = feed_forward(theta, X)y_pred = np.argmax(h, axis=1) + 1print(y_pred)print(classification_report(y, y_pred))show_accuracy(final_theta,X,y_answer)              # 正确率好高啊。。。precision    recall  f1-score   support1       0.99      1.00      0.99       5002       0.99      1.00      0.99       5003       0.99      0.99      0.99       5004       1.00      1.00      1.00       5005       1.00      1.00      1.00       5006       1.00      1.00      1.00       5007       0.99      1.00      0.99       5008       1.00      1.00      1.00       5009       1.00      0.99      0.99       50010       1.00      1.00      1.00       500accuracy                           1.00      5000macro avg       1.00      1.00      1.00      5000
weighted avg       1.00      1.00      1.00      5000

不自量力的用一下。。。。

def sigmoid(inx):a=inx.shapeinx = np.ravel(inx)for i in range(inx.shape[0]):if inx[i]>=0:     inx[i] = 1.0/(1+np.exp(-inx[i]))else:inx[i] = np.exp(inx[i])/(1+np.exp(inx[i]))return inx.reshape((a))import cv2
p1 = cv2.imread(r'D:\jupyter-notebook\opencv\num\p_0.jpg',0)
p2 = cv2.imread(r'D:\jupyter-notebook\opencv\num\p_1.bmp',0)
p3 = cv2.imread(r'D:\jupyter-notebook\opencv\num\p_2.bmp',0)
p_0 = np.insert(p1.reshape(1,400),0,np.ones(1)).reshape(1,401)
p_1 = np.insert(p2.reshape(1,400),0,np.ones(1)).reshape(1,401)
p_2 = np.insert(p3.reshape(1,400),0,np.ones(1)).reshape(1,401)
print(p_0.shape) 
>>>(1, 401)p = np.sum((p_0,p_1,p_2),axis=1)
p.shape   # 3, 401)_, _, _, _, h = feed_forward(res.x, p)
np.argmax(h,axis=1) 
>>> array([2, 2, 2], dtype=int64)  # 试了好多个全都是2 .。。。

汗颜啊，，到底怎么使啊。。有没有大佬能告诉我问题出在哪呢