吴恩达机器学习系列内容的学习目录 → \rightarrow →吴恩达机器学习系列内容汇总。
本次练习对应的基础知识总结 → \rightarrow → 神经网络:Representation。
本次练习对应的文档说明和提供的MATLAB代码 → \rightarrow → 提取码:mwww 。
本次练习对应的完整代码实现(MATLAB + Python版本) → \rightarrow →Github链接。
一、多类别分类
在本练习中,我们将使用Logistic回归和神经网络识别手写数字(从0到9)。如今自动手写数字识别已广泛使用,从识别邮件信封上的邮政编码到识别银行支票上的金额。本练习将展示如何将我们之前学到的方法用于此分类任务。
在练习的第一部分,我们将扩展先前的Logistic回归的实现,并将其应用于一对多分类。
1.1 数据集
在ex3data1.mat中为我们提供了一个数据集,其中包含5000个手写数字的训练样本。.mat格式表示该数据已保存为MATLAB矩阵格式,而不是文本(ASCII)格式。使用load命令可以将这些矩阵直接读取到程序中,加载后我们程序的内存中将会出现有着正确尺寸和值的矩阵。该矩阵已经被命名,因此我们无需为其分配名称。
% Load saved matrices from file
load('ex3data1.mat');
% The matrices X and y will now be in your Octave environment
ex3data1.mat中有5000个训练样本,其中每个训练样本都是像素 20 × 20 20\times20 20×20的数字灰度图像。每个像素由一个浮点数表示,该数字指示该位置的灰度强度。将 20 × 20 20\times20 20×20像素网格“展开”为400维向量。这些训练样本中的每一个都成为数据矩阵 X X X中的一行。这为我们提供了一个 5000 × 400 5000\times400 5000×400的矩阵 X X X,其中每一行都是一个手写数字图像的训练样本。
X = [ − ( x ( 1 ) ) T − − ( x ( 2 ) ) T − ⋮ − ( x ( m ) ) T − ] X=\begin{bmatrix} -\left ( x^{\left ( 1 \right )} \right )^{T}-\\ -\left ( x^{\left ( 2 \right )} \right )^{T}-\\ \vdots \\ -\left ( x^{\left ( m \right )} \right )^{T}- \end{bmatrix} X=⎣⎢⎢⎢⎢⎡−(x(1))T−−(x(2))T−⋮−(x(m))T−⎦⎥⎥⎥⎥⎤
训练集的第二部分是一个5000维向量
y
y
y,其中包含训练集的标签。因为 MATLAB中没有零索引,所以我们将数字零映射到值10。因此,将“ 0”数字标记为“ 10”,而将数字“ 1”至“ 9”按照其自然顺序标记为“ 1”至“ 9”。
1.2 可视化数据
我们将从可视化训练集的子集开始。在ex3.m的第1部分中,代码从X中随机选择100行,并将这些行传递给displayData函数。displayData函数将每一行映射到像素
20
×
20
20\times20
20×20 的灰度图像,并将这些图像一起显示。代码中已经提供了displayData函数,执行完此步骤后,我们应该看到一个如图1所示的图像。
1.3 向量化Logistic回归
我们将使用多个一对一的Logistic回归模型来构建多类别分类器。由于有10个类别,我们将需要训练10个单独的Logistic回归分类器。为了提高训练效率,确保代码被很好地向量化是很重要的。在本节中,我们将实现不使用任何for循环的Logistic回归向量化版本。我们可以使用ex2练习中的代码作为此练习的起点。
1.3.1 向量化代价函数
我们将从编写代价函数的向量化版本开始。回想一下,在(非正则化的)Logistic回归中,代价函数为
J
(
θ
)
=
1
m
∑
i
=
1
m
[
−
y
(
i
)
×
l
o
g
(
h
θ
(
x
(
i
)
)
)
−
(
1
−
y
(
i
)
)
×
l
o
g
(
1
−
h
θ
(
x
(
i
)
)
)
]
J(\theta )=\frac{1}{m}\sum_{i=1}^{m}[-y^{(i)}\times log(h_{\theta}(x^{(i)}))-(1-y^{(i)})\times log(1-h_{\theta}(x^{(i)}))]
J(θ)=m1i=1∑m[−y(i)×log(hθ(x(i)))−(1−y(i))×log(1−hθ(x(i)))]
要计算求和中的每个元素,我们必须为每个样本
i
i
i计算
h
θ
(
x
(
i
)
)
h_{θ}(x^{(i)})
hθ(x(i)),其中
h
θ
(
x
(
i
)
)
=
g
(
θ
T
x
(
i
)
)
h_{θ}(x^{(i)})=g(θ^{T}x^{(i)})
hθ(x(i))=g(θTx(i))和
g
(
z
)
=
1
1
+
e
−
z
g(z)= \frac{1}{1+e^{-z}}
g(z)=1+e−z1。事实证明,我们可以使用矩阵乘法快速计算出所有样本。
X
X
X和
θ
θ
θ定义为
X
=
[
−
(
x
(
1
)
)
T
−
−
(
x
(
2
)
)
T
−
⋮
−
(
x
(
m
)
)
T
−
]
a
n
d
θ
=
[
θ
0
θ
1
⋮
θ
n
]
\begin{matrix} X=\begin{bmatrix} -\left ( x^{\left ( 1 \right )} \right )^{T}-\\ -\left ( x^{\left ( 2 \right )} \right )^{T}-\\ \vdots \\ -\left ( x^{\left ( m \right )} \right )^{T}- \end{bmatrix} & and & \theta =\begin{bmatrix} \theta _{0}\\ \theta _{1}\\ \vdots \\ \theta _{n} \end{bmatrix} \end{matrix}
X=⎣⎢⎢⎢⎢⎡−(x(1))T−−(x(2))T−⋮−(x(m))T−⎦⎥⎥⎥⎥⎤andθ=⎣⎢⎢⎢⎡θ0θ1⋮θn⎦⎥⎥⎥⎤
然后,通过计算矩阵乘积
X
θ
X θ
Xθ,我们可以得到
X
θ
=
[
−
(
x
(
1
)
)
T
θ
−
−
(
x
(
2
)
)
T
θ
−
⋮
−
(
x
(
m
)
)
T
θ
−
]
=
[
−
θ
T
(
x
(
1
)
)
−
−
θ
T
(
x
(
2
)
)
−
⋮
−
θ
T
(
x
(
m
)
)
−
]
X\theta =\begin{bmatrix} -\left ( x^{\left ( 1 \right )} \right )^{T}\theta -\\ -\left ( x^{\left ( 2 \right )} \right )^{T}\theta -\\ \vdots \\ -\left ( x^{\left ( m \right )} \right )^{T}\theta - \end{bmatrix}=\begin{bmatrix} -\theta^{T}\left ( x^{\left ( 1 \right )} \right ) -\\ -\theta^{T}\left ( x^{\left ( 2 \right )} \right ) -\\ \vdots\\ -\theta^{T}\left ( x^{\left ( m \right )} \right ) - \end{bmatrix}
Xθ=⎣⎢⎢⎢⎢⎡−(x(1))Tθ−−(x(2))Tθ−⋮−(x(m))Tθ−⎦⎥⎥⎥⎥⎤=⎣⎢⎢⎢⎡−θT(x(1))−−θT(x(2))−⋮−θT(x(m))−⎦⎥⎥⎥⎤
在最后一个等式中,如果a和b是向量,则有
a
T
b
=
b
T
a
a^{T}b=b^{T}a
aTb=bTa,这使得我们能够在一行代码中计算出所有样本
i
i
i的乘积
θ
T
x
(
i
)
θ^{T}x^{(i)}
θTx(i)。
我们要完成的是在文件lrCostFunction.m中编写非正则化的代价函数,实现时应使用我们上面介绍的策略来计算
θ
T
x
(
i
)
θ^{T}x^{(i)}
θTx(i)。除此之外,我们还应该对代价函数使用向量化方法,完全向量化的lrCostFunction.m版本不应包含任何循环。
1.3.2 向量化梯度
回想一下,(非正则化的)Logistic回归代价函数的梯度是一个向量,其中第
j
j
j个元素的定义为
∂
J
(
θ
)
∂
θ
j
=
1
m
∑
i
=
1
m
(
(
h
θ
(
x
(
i
)
)
−
y
(
i
)
)
x
j
(
i
)
)
\frac{\partial J(\theta )}{\partial \theta _{j}}=\frac{1}{m}\sum_{i=1}^{m} \left ( (h _{\theta}(x^{(i)})-y^{(i)})x_{j}^{(i)}\right )
∂θj∂J(θ)=m1i=1∑m((hθ(x(i))−y(i))xj(i))
为了对数据集进行向量化处理,我们首先针对所有
θ
j
θ_{j}
θj明确写出其偏导数
[
∂
J
∂
θ
0
∂
J
∂
θ
1
∂
J
∂
θ
2
⋮
∂
J
∂
θ
n
]
=
1
m
[
∑
i
=
1
m
(
(
h
θ
(
x
(
i
)
)
−
y
(
i
)
)
x
0
(
i
)
)
∑
i
=
1
m
(
(
h
θ
(
x
(
i
)
)
−
y
(
i
)
)
x
1
(
i
)
)
∑
i
=
1
m
(
(
h
θ
(
x
(
i
)
)
−
y
(
i
)
)
x
2
(
i
)
)
⋮
∑
i
=
1
m
(
(
h
θ
(
x
(
i
)
)
−
y
(
i
)
)
x
n
(
i
)
)
]
=
1
m
∑
i
=
1
m
(
(
h
θ
(
x
(
i
)
)
−
y
(
i
)
)
x
(
i
)
)
=
1
m
X
T
(
(
h
θ
(
x
)
−
y
)
\begin{matrix} \begin{bmatrix} \frac{\partial J}{\partial \theta _{0}}\\ \frac{\partial J}{\partial \theta _{1}}\\ \frac{\partial J}{\partial \theta _{2}}\\ \vdots \\ \frac{\partial J}{\partial \theta _{n}} \end{bmatrix} & =\frac{1}{m}\begin{bmatrix} \sum_{i=1}^{m}\left ( (h _{\theta}(x^{(i)})-y^{(i)})x_{0}^{(i)} \right )\\ \sum_{i=1}^{m}\left ( (h _{\theta}(x^{(i)})-y^{(i)})x_{1}^{(i)} \right )\\ \sum_{i=1}^{m}\left ( (h _{\theta}(x^{(i)})-y^{(i)})x_{2}^{(i)} \right )\\ \vdots \\ \sum_{i=1}^{m}\left ( (h _{\theta}(x^{(i)})-y^{(i)})x_{n}^{(i)} \right ) \end{bmatrix}\\ &\\ & =\frac{1}{m}\sum_{i=1}^{m}\left ( (h _{\theta}(x^{(i)})-y^{(i)})x^{(i)} \right )\ _{}\ _{}\ _{}\ _{}\\ &\\ & =\frac{1}{m}X^{T}\left ((h _{\theta}(x) -y \right )\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{} \end{matrix}
⎣⎢⎢⎢⎢⎢⎢⎡∂θ0∂J∂θ1∂J∂θ2∂J⋮∂θn∂J⎦⎥⎥⎥⎥⎥⎥⎤=m1⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎡∑i=1m((hθ(x(i))−y(i))x0(i))∑i=1m((hθ(x(i))−y(i))x1(i))∑i=1m((hθ(x(i))−y(i))x2(i))⋮∑i=1m((hθ(x(i))−y(i))xn(i))⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎤=m1∑i=1m((hθ(x(i))−y(i))x(i)) =m1XT((hθ(x)−y)
其中, h θ ( x ) − y = [ h θ ( x ( 1 ) ) − y ( 1 ) h θ ( x ( 2 ) ) − y ( 2 ) ⋮ h θ ( x ( m ) ) − y ( m ) ] h _{\theta}(x) -y=\begin{bmatrix} h _{\theta}(x^{(1)}) -y^{(1)}\\ h _{\theta}(x^{(2)}) -y^{(2)}\\ \vdots \\ h _{\theta}(x^{(m)}) -y^{(m)} \end{bmatrix} hθ(x)−y=⎣⎢⎢⎢⎡hθ(x(1))−y(1)hθ(x(2))−y(2)⋮hθ(x(m))−y(m)⎦⎥⎥⎥⎤。
需要注意的是,
x
(
i
)
x^{(i)}
x(i)是向量,而
(
h
θ
(
x
(
i
)
)
−
y
(
i
)
)
\left (h _{\theta}(x^{(i)}) -y^{(i)}\right )
(hθ(x(i))−y(i))是标量。为了理解推导的最后一步,令
β
i
=
(
h
θ
(
x
(
i
)
)
−
y
(
i
)
)
β_{i}=\left (h _{\theta}(x^{(i)}) -y^{(i)}\right )
βi=(hθ(x(i))−y(i)),我们可以得到
∑
i
β
i
x
(
i
)
=
[
∣
x
(
1
)
∣
∣
x
(
2
)
∣
⋯
∣
x
(
m
)
∣
]
[
β
1
β
2
⋮
β
m
]
=
X
T
β
\sum_{i}^{ } \beta _{i}x^{(i)}=\begin{bmatrix} \begin{matrix} |\\ x^{(1)}\\ | \end{matrix} &\begin{matrix} |\\ x^{(2)}\\ | \end{matrix} & \cdots &\begin{matrix} |\\ x^{(m)}\\ | \end{matrix} \end{bmatrix}\begin{bmatrix} \beta _{1}\\ \beta _{2}\\ \vdots \\ \beta _{m} \end{bmatrix}=X^{T}\beta
i∑βix(i)=⎣⎡∣x(1)∣∣x(2)∣⋯∣x(m)∣⎦⎤⎣⎢⎢⎢⎡β1β2⋮βm⎦⎥⎥⎥⎤=XTβ
上面的表达式使我们无需循环即可计算所有偏导数。现在,我们应该通过编写代码实现公式来计算正确的向量化梯度。
完成lrCostFunction.m时需要填写以下代码:
h = sigmoid(X * theta);
J = 1/m .* ( -y' * log(h)- (1-y)' * log(1-h)) ;
grad = 1/m .* X' * (h - y);
1.3.3 向量化Logistic回归的正则化
在Logistic回归实现向量化之后,现在将向代价函数添加正则化项。回想一下,对于正则化Logistic回归,代价函数定义为
J
(
θ
)
=
1
m
∑
i
=
1
m
[
−
y
(
i
)
×
l
o
g
(
h
θ
(
x
(
i
)
)
)
−
(
1
−
y
(
i
)
)
×
l
o
g
(
1
−
h
θ
(
x
(
i
)
)
)
]
+
λ
2
m
∑
j
=
1
n
θ
j
2
J(\theta)=\frac{1}{m}\sum_{i=1}^{m} [-y^{(i)}\times log(h_{\theta}(x^{(i)}))-(1-y^{(i)})\times log(1-h_{\theta}(x^{(i)}))]+\frac{\lambda }{2m} \sum_{j=1}^{n} \theta _{j}^{2}
J(θ)=m1i=1∑m[−y(i)×log(hθ(x(i)))−(1−y(i))×log(1−hθ(x(i)))]+2mλj=1∑nθj2
要注意的是,我们不应对用于偏差项的
θ
0
θ_{0}
θ0进行正则化。
相应地,正则化Logistic回归代价函数对
θ
j
θ_{j}
θj的偏导数定义为
∂ J ( θ ) ∂ θ 0 = 1 m ∑ i = 1 m ( h θ ( x ( i ) ) − y ( i ) ) x j ( i ) f o r j = 0 ∂ J ( θ ) ∂ θ j = ( 1 m ∑ i = 1 m ( h θ ( x ( i ) ) − y ( i ) ) x j ( i ) ) + λ m θ j f o r j ⩾ 1 \begin{matrix} \frac{\partial J(\theta )}{\partial \theta _{0}}=\frac{1}{m}\sum_{i=1}^{m} (h _{\theta}(x^{(i)})-y^{(i)})x_{j}^{(i)}\ _{}\ _{}\ _{}\ _{}for\ _{}\ _{}j=0\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\ _{}\\ \\ \frac{\partial J(\theta )}{\partial \theta _{j}}=\left (\frac{1}{m}\sum_{i=1}^{m} (h _{\theta}(x^{(i)})-y^{(i)})x_{j}^{(i)} \right )+\frac{\lambda }{m}\theta _{j}\ _{}\ _{}\ _{}\ _{}for\ _{}\ _{}j\geqslant 1 \end{matrix} ∂θ0∂J(θ)=m1∑i=1m(hθ(x(i))−y(i))xj(i) for j=0 ∂θj∂J(θ)=(m1∑i=1m(hθ(x(i))−y(i))xj(i))+mλθj for j⩾1
现在,我们需要在lrCostFunction中对之前编写的代码进行修改,在代码中加入正则化部分。
修改之前编写的lrCostFunction.m为如下代码:
temp=[0;theta(2:end)]; % 先把theta(1)拿掉,不参与正则化
h = sigmoid(X * theta);
J = 1/m .* ( -y' * log(h)- (1-y)' * log(1-h)) + lambda /(2*m) * temp' * temp ;
grad = 1/m .* X' * (h - y) + lambda/m * temp;
完成lrCostFunction.m中的代码后,运行ex3.m脚本将可以得到如下结果:
Testing lrCostFunction() with regularization
Cost: 2.534819
Expected cost: 2.534819
Gradients:
0.146561
-0.548558
0.724722
1.398003
Expected gradients:
0.146561
-0.548558
0.724722
1.398003
1.4 一对多分类
在本部分练习中,对于我们数据集中的每个K类(如图1所示),我们将通过训练多个正则化Logistic回归分类器来实现一对多分类。在手写数字数据集中,K=10,但是我们编写的代码应适用于任何K值。
现在,我们应该在oneVsAll.m中完成代码来为每个类别训练一个分类器。我们的代码应以矩阵
Θ
∈
R
K
×
(
N
+
1
)
\Theta \in R^{K\times (N+1) }
Θ∈RK×(N+1)的形式返回所有分类器参数,其中
Θ
\Theta
Θ每一行对应于一类的Logistic回归学习参数。我们可以使用从1到K的“ for”循环进行此操作,独立训练每个分类器。
要注意的是,这个函数的y参数是一个从1到10的标签向量,我们将数字“0”映射到标签10(以避免与索引混淆)。
在训练类别
k
∈
{
1
,
.
.
.
,
K
}
k\in \left \{ 1,...,K \right \}
k∈{1,...,K}的分类器时,我们将需要一个标注y的多维向量,其中
y
j
∈
0
,
1
y_{j}\in0,1
yj∈0,1,表示第
j
个
j个
j个训练样本是属于K类(
y
j
=
1
y_{j}= 1
yj=1),还是属于另一个类(
y
j
=
0
y_{j}= 0
yj=0)。
此外,我们将在本练习中使用fmincg(而不是fminunc)。 fmincg的工作方式与fminunc相似,但处理大量参数的效率更高。正确完成oneVsAll.m的代码后,脚本ex3.m将继续使用oneVsAll函数来训练多类分类器。
完成oneVsAll.m时需要填写以下代码:
initial_theta = zeros(n+1 , 1);
options = optimset('Gradobj', 'on', 'MaxIter', 50);
for c = 1: num_labels
[all_theta(c,:)] = fmincg(@(t)lrCostFunction(t, X, (y==c), lambda) ,initial_theta, options);%使选项中指定的优化选项最小化。
end
%当我们调用这个fminunc函数时,它会自动的从众多高级优化算法中挑选一个来使用(你也可以把它当做一个可以自动选择合适的学习速率的梯度下降算法)。
对上述代码有如下解释:
- Logistic回归分类的目标是区分两种类型。写好代价函数之后,调用一次fminunc函数得到theta,就可以画出boundary。但是多分类中需要训练K个分类器,所以多了一个oneVsAll.m。其实就是循环K次,得到一个theta矩阵(row为 K,column为特征值个数)。
- 多分类使用的Octave函数是fmincg不是fminunc,fmincg更适合参数多的情况。
- 这里用到了y == c,这个操作将y中每个数据和c进行比较。得到的矩阵中(这里是向量)相同的为1,不同的为0。表示第j个训练实例是属于K类(yj = 1),还是属于另一个类(yj = 0)。
完成oneVsAll.m中的代码后,运行ex3.m脚本将可以得到如下结果:
1.4.1 一对多预测
训练完所有函数后,我们现在可以使用它来预测给定图像中包含的数字。对于每个输入,我们应该使用训练有素的Logistic回归分类器计算它属于每个类别的“概率”。一对多预测函数将会选择相应的Logistic回归分类器所输出概率最高的类别,并返回类标签(1,2,…,K)作为输入样本的预测。
现在,我们应该在predictOneVsAll.m中完成代码以使用“一对多”分类器进行预测。
完成predictOneVsAll.m时需要填写以下代码:
temp = X * all_theta'; %5000X10
[x, p] = max(temp,[],2); % 返回每行最大值的索引位置,也就是预测的数字
对上述代码有如下解释:
- 将X与theta相乘得到预测结果,得到一个 5000 × 10 5000\times10 5000×10的矩阵,每行对应一个分类结果(只有一个1,其余为0)。题目要求返回的矩阵p是一个5000*1的矩阵,每行是1到K中的数字。
- 使用max函数,得到两个返回值。第一个x是一个全1的向量,没用。第二个p是1所在的index,也就是预测的类别。
完成后,ex3.m将使用学习到的
Θ
\Theta
Θ值调用prepareOneVsAll函数。我们应该看到训练集的准确性约为94.9%,即在训练集中将94.9%的样本进行了正确地分类。
我的运行结果如下:
Training Set Accuracy: 95.020000
二、神经网络
在本练习的上一部分中,我们实现了多类别Logistic回归来识别手写数字。但是,Logistic回归不能构成更加复杂的假设,因为它只是线性分类器。
在本部分的练习中,我们将实现一个神经网络,使用与以前相同的训练集来识别手写数字,神经网络将能够代表形成非线性假设的复杂模型。接下来我们将使用吴恩达老师已经训练过的神经网络参数,我们的目标是实现前馈传播算法,并使用我们的权重进行预测。在下次练习中,我们将编写用于学习神经网络参数的反向传播算法。
提供的脚本ex3_nn.m将帮助我们逐步完成此练习。
2.1 模型表示
我们的神经网络如图2所示。它有3层——输入层,隐藏层和输出层。回想一下,我们的输入是数字图像的像素值。由于图像的尺寸为
20
×
20
20×20
20×20,因此这给了我们400个输入单元(不包括输出始终加一的额外偏置单元)。像以前一样,训练数据将被加载到变量
X
X
X和
y
y
y中。
本练习的资料中已经为我们提供了一组网络参数 ( Θ ( 1 ) , Θ ( 2 ) ) (\Theta^{(1)},\Theta^{(2)}) (Θ(1),Θ(2))。这些参数存储在ex3weights.mat中,并将由ex3_nn.m加载到Theta1和Theta2中,实现代码如下所示。参数的大小适用于神经网络,中间层有25个单元,输出层有10个单元(对应于10位数字类别) 。
% Load saved matrices from file
load('ex3weights.mat');
% The matrices Theta1 and Theta2 will now be in your Octave
% environment
% Theta1 has size 25 x 401
% Theta2 has size 10 x 26
2.2 前馈传播和预测
现在,我们将为神经网络实现前馈传播,需要完成predict.m中的代码以返回神经网络的预测。
我们应该实现前馈计算,该计算为每个样本
i
i
i计算
h
θ
(
x
(
i
)
)
h_{\theta}(x^{(i)})
hθ(x(i))并返回相关预测。与一对多分类策略相似,神经网络的预测结果将是输出
(
h
θ
(
x
)
)
k
\left ( h_{\theta}(x)\right )_{k}
(hθ(x))k最大时的标签。
完成predict.m时需要填写以下代码:
a1 = [ones(m, 1) X]; %5000x401
a2 = sigmoid(a1 * Theta1'); %5000x401乘以401x25得到5000x25。即把401个feature映射到25
a2 = [ones(m, 1) a2]; %5000x26
a3 = sigmoid(a2 * Theta2'); %5000x26乘以26x10得到5000x10。即把26个feature映射到10
[x,p] = max(a3,[],2);%和上面逻辑回归多分类一样,最后使用max函数获得分类结果
%[C, I] = max(a, [], dim)表示a是个二维矩阵,dim = 2表示比较的是行,返回size(a, 1) 行,每行元素是a该行最大的元素;dim = 1 表示比较的是列。
完成后,ex3_nn.m将使用已加载的Theta1和Theta2参数集调用预测函数,我们应该得到准确度约为97.5%。之后,交互序列将一次一个地从训练集启动显示图像,而控制台将打印出所显示图像的预测标签。要停止图像序列,请按Ctrl+C。
准确度结果如下:
Training Set Accuracy: 97.520000
显示的样本图像(部分)和预测结果如下:
- Neural Network Prediction: 8 (digit 8)
- Neural Network Prediction: 2 (digit 2)
- Neural Network Prediction: 10 (digit 0)
- Neural Network Prediction: 6 (digit 6)
三、MATLAB实现
3.1 ex3.m
%% Machine Learning Online Class - Exercise 3 | Part 1: One-vs-all
% Instructions
% ------------
%
% This file contains code that helps you get started on the
% linear exercise. You will need to complete the following functions
% in this exericse:
%
% lrCostFunction.m (logistic regression cost function)
% oneVsAll.m
% predictOneVsAll.m
% predict.m
%
% For this exercise, you will not need to change any code in this file,
% or any other files other than those mentioned above.
%
%% Initialization
clear ; close all; clc
%% Setup the parameters you will use for this part of the exercise
input_layer_size = 400; % 20x20 Input Images of Digits
num_labels = 10; % 10 labels, from 1 to 10
% (note that we have mapped "0" to label 10)
%% =========== Part 1: Loading and Visualizing Data =============
% We start the exercise by first loading and visualizing the dataset.
% You will be working with a dataset that contains handwritten digits.
%
% Load Training Data
fprintf('Loading and Visualizing Data ...\n')
load('ex3data1.mat'); % training data stored in arrays X, y
m = size(X, 1); %返回矩阵X的行数。(若size(X,2)返回X的列数)
% Randomly select 100 data points to display
rand_indices = randperm(m); %将1到5000整数随机打乱
sel = X(rand_indices(1:100), :); %先打乱再选取100行向量,即100个数字。display显示10*10数字矩阵
displayData(sel);
fprintf('Program paused. Press enter to continue.\n');
pause;
%% ============ Part 2a: Vectorize Logistic Regression ============
% In this part of the exercise, you will reuse your logistic regression
% code from the last exercise. You task here is to make sure that your
% regularized logistic regression implementation is vectorized. After
% that, you will implement one-vs-all classification for the handwritten
% digit dataset.
%
% Test case for lrCostFunction
fprintf('\nTesting lrCostFunction() with regularization');
theta_t = [-2; -1; 1; 2]; %theta_t:4x1
X_t = [ones(5,1) reshape(1:15,5,3)/10]; %X_t:5x4
y_t = ([1;0;1;0;1] >= 0.5); %y_t:5x1
lambda_t = 3;
[J grad] = lrCostFunction(theta_t, X_t, y_t, lambda_t);
fprintf('\nCost: %f\n', J);
fprintf('Expected cost: 2.534819\n');
fprintf('Gradients:\n');
fprintf(' %f \n', grad);
fprintf('Expected gradients:\n');
fprintf(' 0.146561\n -0.548558\n 0.724722\n 1.398003\n');
fprintf('Program paused. Press enter to continue.\n');
pause;
%% ============ Part 2b: One-vs-All Training ============
fprintf('\nTraining One-vs-All Logistic Regression...\n')
lambda = 0.1;
[all_theta] = oneVsAll(X, y, num_labels, lambda);
fprintf('Program paused. Press enter to continue.\n');
pause;
%% ================ Part 3: Predict for One-Vs-All ================
pred = predictOneVsAll(all_theta, X);
fprintf('\nTraining Set Accuracy: %f\n', mean(double(pred == y)) * 100);
3.2 ex3_nn.m
%% Machine Learning Online Class - Exercise 3 | Part 2: Neural Networks
% Instructions
% ------------
%
% This file contains code that helps you get started on the
% linear exercise. You will need to complete the following functions
% in this exericse:
%
% lrCostFunction.m (logistic regression cost function)
% oneVsAll.m
% predictOneVsAll.m
% predict.m
%
% For this exercise, you will not need to change any code in this file,
% or any other files other than those mentioned above.
%
%% Initialization
clear ; close all; clc
%% Setup the parameters you will use for this exercise
input_layer_size = 400; % 20x20 Input Images of Digits
hidden_layer_size = 25; % 25 hidden units
num_labels = 10; % 10 labels, from 1 to 10
% (note that we have mapped "0" to label 10)
%% =========== Part 1: Loading and Visualizing Data =============
% We start the exercise by first loading and visualizing the dataset.
% You will be working with a dataset that contains handwritten digits.
%
% Load Training Data
fprintf('Loading and Visualizing Data ...\n')
load('ex3data1.mat');
m = size(X, 1);
% Randomly select 100 data points to display
sel = randperm(size(X, 1));
sel = sel(1:100);
displayData(X(sel, :));
fprintf('Program paused. Press enter to continue.\n');
pause;
%% ================ Part 2: Loading Pameters ================
% In this part of the exercise, we load some pre-initialized
% neural network parameters.
fprintf('\nLoading Saved Neural Network Parameters ...\n')
% Load the weights into variables Theta1 and Theta2
load('ex3weights.mat');
%% ================= Part 3: Implement Predict =================
% After training the neural network, we would like to use it to predict
% the labels. You will now implement the "predict" function to use the
% neural network to predict the labels of the training set. This lets
% you compute the training set accuracy.
pred = predict(Theta1, Theta2, X);
fprintf('\nTraining Set Accuracy: %f\n', mean(double(pred == y)) * 100);
fprintf('Program paused. Press enter to continue.\n');
pause;
% To give you an idea of the network's output, you can also run
% through the examples one at the a time to see what it is predicting.
% Randomly permute examples
rp = randperm(m);
for i = 1:m
% Display
fprintf('\nDisplaying Example Image\n');
displayData(X(rp(i), :));
pred = predict(Theta1, Theta2, X(rp(i),:));
fprintf('\nNeural Network Prediction: %d (digit %d)\n', pred, mod(pred, 10));
% Pause with quit option
s = input('Paused - press enter to continue, q to exit:','s');
if s == 'q'
break
end
end
四、Python实现
4.1 ex3.py
import numpy as np
import matplotlib.pylab as plt
import scipy.io as sio
import math
import scipy.optimize as op
input_layer_size = 400
num_labels = 10
# =========== Part 1: Loading and Visualizing Data =============
print('Loading and Visualizing Data ...')
matinfo = sio.loadmat('ex3data1.mat')
X = matinfo['X']
Y = matinfo['y'][:, 0]
m = np.size(X, 0)
rand_indices = np.random.permutation(m) #将1到5000整数随机打乱
sel = X[rand_indices[0:100], :]#%先打乱再选取100行向量,即100个数字。display显示10*10数字矩阵
# 可视化数据集
# 显示随机100个图像, 疑问:最后的数组需要转置才会显示正的图像
def displayData(x):
width = round(math.sqrt(np.size(x, 1)))#单张图片宽度
m, n = np.shape(x)
height = int(n/width)#单张图片高度
# 显示图像的数量
drows = math.floor(math.sqrt(m))#显示图中,一行多少张图
dcols = math.ceil(m/drows)#显示图中,一列多少张图片
pad = 1#图片间的间隔
# 建立一个空白“背景布”
darray = -1*np.ones((pad+drows*(height+pad), pad+dcols*(width+pad)))#初始化图片矩阵
curr_ex = 0#当前的图片计数
#将每张小图插入图片数组中
for j in range(drows):
for i in range(dcols):
if curr_ex >= m:
break
max_val = np.max(np.abs(X[curr_ex, :]))
darray[pad+j*(height+pad):pad+j*(height+pad)+height, pad+i*(width+pad):pad+i*(width+pad)+width]\
= x[curr_ex, :].reshape((height, width))/max_val
curr_ex += 1
if curr_ex >= m:
break
plt.imshow(darray.T, cmap='gray')
plt.show()
displayData(sel)
_ = input('Press [Enter] to continue.')
# ============ Part 2: Vectorize Logistic Regression ============
print('Training One-vs-All Logistic Regression...')
# sigmoid函数
def sigmoid(z):
g = 1/(1+np.exp(-1*z))
return g
# 损失函数
def lrCostFunc(theta, x, y, lam):
m = np.size(y, 0)
h = sigmoid(x.dot(theta))
j = -1/m*(y.dot(np.log(h))+(1-y).dot(np.log(1-h)))+lamb*(theta[1:].dot(theta[1:]))/(2*m)
return j
# 梯度函数
def lrGradFunc(theta, x, y, lam):
m = np.size(y, 0)
h = sigmoid(x.dot(theta))
grad = np.zeros(np.size(theta))
grad[0] = 1 / m * (x[:, 0].dot(h-y))
grad[1:] = 1/m*(x[:, 1:].T.dot(h-y))+lam/m*theta[1:]
return grad
# 获取多个分类器的theta值
def oneVsAll(x, y, num_labels, lam):
m, n = np.shape(x)
all_theta = np.zeros((num_labels, n+1))
x = np.concatenate((np.ones((m, 1)), x), axis=1)
for i in range(num_labels):
num = 10 if i == 0 else i
init_theta = np.zeros((n+1,))
result = op.minimize(lrCostFunc, init_theta, method='BFGS', jac=lrGradFunc, args=(x, 1*(y == num), lam), options={'maxiter': 50})
all_theta[i, :] = result.x
return all_theta
lamb = 0.1
all_theta = oneVsAll(X, Y, num_labels, lamb)
_ = input('Press [Enter] to continue.')
# ================ Part 3: Predict for One-Vs-All ================
# 预测值函数
def predictOneVsAll(all_theta, x):
m = np.size(x, 0)
x = np.concatenate((np.ones((m, 1)), x), axis=1)# 这里的axis=1表示按照列进行合并
p = np.argmax(x.dot(all_theta.T), axis=1)#np.argmax(a)取出a中元素最大值所对应的索引(索引值默认从0开始),axis=1即按照行方向搜索最大值
return p
pred = predictOneVsAll(all_theta, X)
print('Training Set Accuracy: ', np.sum(pred == (Y % 10))/np.size(Y, 0))
4.2 ex3_nn.py
import numpy as np
import matplotlib.pylab as plt
import scipy.io as sio
import math
# 显示随机100个图像, 疑问:最后的数组需要转置才会显示正的图像
def displayData(x):
width = round(math.sqrt(np.size(x, 1)))
m, n = np.shape(x)
height = int(n/width)
# 显示图像的数量
drows = math.floor(math.sqrt(m))
dcols = math.ceil(m/drows)
pad = 1
# 建立一个空白“背景布”
darray = -1*np.ones((pad+drows*(height+pad), pad+dcols*(width+pad)))
curr_ex = 0
for j in range(drows):
for i in range(dcols):
if curr_ex >= m:
break
max_val = np.max(np.abs(X[curr_ex, :]))
darray[pad+j*(height+pad):pad+j*(height+pad)+height, pad+i*(width+pad):pad+i*(width+pad)+width]\
= x[curr_ex, :].reshape((height, width))/max_val
curr_ex += 1
if curr_ex >= m:
break
plt.imshow(darray.T, cmap='gray')
plt.show()
input_layer_size = 400
hidden_layer_size = 25
num_labels = 10
# =========== Part 1: Loading and Visualizing Data =============
# LoadData
print('Loading and Visualizing Data ...')
datainfo = sio.loadmat('ex3data1.mat')
X = datainfo['X']
Y = datainfo['y'][:, 0]
m = np.size(X, 0)
rand_indices = np.random.permutation(m)
sel = X[rand_indices[0:100], :]
displayData(sel)
_ = input('Press [Enter] to continue.')
# ================ Part 2: Loading Pameters ================
print('Loading Saved Neural Network Parameters ...')
weightinfo = sio.loadmat('ex3weights.mat')
theta1 = weightinfo['Theta1']
theta2 = weightinfo['Theta2']
# ================= Part 3: Implement Predict =================
# sigmoid函数
def sigmoid(z):
g = 1/(1+np.exp(-1*z))
return g
# 预测函数
def predict(t1, t2, x):
m = np.size(x, 0)
x = np.concatenate((np.ones((m, 1)), x), axis=1)
temp1 = sigmoid(x.dot(theta1.T))
temp = np.concatenate((np.ones((m, 1)), temp1), axis=1)
temp2 = sigmoid(temp.dot(theta2.T))
p = np.argmax(temp2, axis=1)+1
return p
pred = predict(theta1, theta2, X)
print('Training Set Accuracy: ', np.sum(pred == Y)/np.size(Y, 0))
_ = input('Press [Enter] to continue.')
# 随机展示图像
num = 10
rindex = np.random.permutation(m) #numpy.random.permutation(x)随机排列一个序列或者数组。
for i in range(num):
print('Displaying Example Image')
displayData(X[rindex[i]:rindex[i]+1, :])#显示一个数字
pred = predict(theta1, theta2, X[rindex[i]:rindex[i]+1, :])
print('Neural Network Prediction: %d (digit %d)' % (pred, pred % 10))