梯度下降找方程参数?

机器学习里边,梯度下降用来寻找模型参数.其实它就是在寻找方程的系数,什么方程?能量函数.今天使用形式简单的函数来验证一下,是不是所有形式的方程都能顺利地用梯度下降求解.

线性函数

这个其实都不用验证,这就是线性回归嘛.不过还是按照从简单到复杂一步步来吧.
有这么一个线性函数y=Ax+B,参数A和B就用theta表示了.要用梯度下降拟合参数需要知道两个东西,能量函数和能量函数对于theta的偏导,下面两个函数就做的这个事.

1
2
3
4
5
6
7
% linearFunc.m
function y = linearFunc(x, theta)
% y = Ax+B
%%
    y = theta(1)*x + theta(2);
    
end

1
2
3
4
5
6
7
8
9
% linearFuncCost.m 返回能量以及梯度
function [J, grad] = linearFuncCost(theta, x, y)
    J = sum( (theta(1)*x + theta(2) - y).^2); %能量函数,这里为了简单问题没有对m个样本进行归一化,
					%因为我们只是寻找方程参数,并不是线性回归要考虑到样本个数.
    grad = zeros(size(theta));
    grad(1) = sum(2*(theta(1)*x + theta(2) - y) .* x); %两个偏导
    grad(2) = sum(2*(theta(1)*x + theta(2) - y));

end

然后呢?我们给定一组x,再用公式来产生正确的y,试试能否通过梯度下降来找到正确的参数theta.

1
2
3
4
5
6
7
8
9
10
% LinearTest.m
x = [-10:0.1:10];
realTheta = [5, 2];
y = linearFunc(x, realTheta);

options = optimset('GradObj', 'on', 'MaxIter', 50);
theta = [1;1];
theta = fmincg (@(t)(linearFuncCost(t, x, y)), ...
                [1;1], options);
fprintf(theta)

这里用到一个优化函数fmincg,其实比简单的梯度下降稍微高端一点,用了共厄梯度,是Ng公开课里推荐用的.用fmincg找参数,需要把linearFuncCost当作参数传进去,还有一组初始theta值,随便给个[1;1]就行,还有一个options,比如这里的options规定最多进行50次迭代.当然还有输入x和y,可以看成是样本.结果当然是能找到了,只经过13次迭代.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
Iteration     1 | Cost: 3.776981e+03
Iteration     2 | Cost: 9.778718e+02
Iteration     3 | Cost: 1.037721e+02
Iteration     4 | Cost: 5.217549e+00
Iteration     5 | Cost: 1.193253e-02
Iteration     6 | Cost: 1.155018e-02
Iteration     7 | Cost: 4.344808e-03
Iteration     8 | Cost: 2.320043e-03
Iteration     9 | Cost: 3.970677e-04
Iteration    10 | Cost: 7.156246e-06
Iteration    11 | Cost: 5.619422e-08
Iteration    12 | Cost: 3.901228e-09
Iteration    13 | Cost: 1.045008e-09
Iteration    14 | Cost: 5.914386e-27
Iteration    15 | Cost: 3.352659e-30

    5.0000
    2.0000

给y加点噪声呢?

1
y = linearFunc(x, realTheta) + rand(1, length(x))*0.1;

再迭代,效果还是很好.

1
2
4.9993
2.0481

复杂一点的函数

稍微复杂一点的,什么交叉熵之类的就不验证了,它们是真能找到参数.下面是一个我这几天实际碰到的问题中的一个函数.具体问题不谈,问题的求解转化成了用给定数据拟合这样一个形式的函数的问题.
y=theta(3)√theta(2)+x−theta(4)√theta(1)+x形式略复杂,不过能量函数和偏导还是能写出来的.这个函数由于有4个参数,所以有4个偏导.

1
2
3
4
% f_x.m
function y = f_x(x, theta)
    y = theta(3)./sqrt(theta(2)+x) - theta(4)./sqrt(theta(1)+x);
end

1
2
3
4
5
6
7
8
9
10
11
function [J, grad] = f_xcost(theta, x, y)
    J = sum((y - theta(3)./sqrt(theta(2)+x) + theta(4)./sqrt(theta(1)+x)).^2);   
    J_a_grad = sum(    2*( y - theta(3)./sqrt(theta(2) + x) + theta(4)./sqrt(theta(1)+x) ).* (-0.5*theta(4)*(theta(1)+x).^(-1.5))    );
    J_b_grad = sum(    2*( y - theta(3)./sqrt(theta(2) + x) + theta(4)./sqrt(theta(1)+x) ).* (0.5*theta(3)*(theta(2)+x).^(-1.5))    );    

    J_c_grad = sum(    2*( y - theta(3)./sqrt(theta(2) + x) + theta(4)./sqrt(theta(1)+x) ).* ( -1./sqrt(theta(2)+x) )    );
    J_d_grad = sum(    2*( y - theta(3)./sqrt(theta(2) + x) + theta(4)./sqrt(theta(1)+x) ).* ( 1./sqrt(theta(1)+x) )    );

    grad = [J_a_grad; J_b_grad;J_c_grad; J_d_grad];

end

再试试.

1
2
3
4
5
6
7
8
% f_xTest.m
x = (5:0.1:9);
realTheta = [2.4;2.9;0.8;0.2];
y = f_x(x, realTheta);
initTheta = [1;1;    1;1];
options = optimset('GradObj', 'on', 'MaxIter', 4000);
estimateTheta = fmincg (@(t)(f_xcost(t, x, y)), ...
                initTheta, options)

1
2
3
4
5
6
7
Iteration  3699 | Cost: 2.437550e-16
estimateTheta =

    2.5463
    2.8467
    1.0387
    0.4387

经过了3699次迭代收敛,但是估计出来的参数却和实际参数偏差很大.只能解释为遇到了局部最小值,掉坑里了.
改进的想法是,一次用梯度下降更新4个参数可能有点过快了,尝试一次更新3个参数,留一个不更新.这样每次估计出3个新参和一个不变的参数,用这4个来进行第二次更新,依然更新其中3个,不过不是第一次那3个了.这样经过4次就能更新全了.

结果

使用部分更新的办法有时能正确的找到参数,有时也不行,不知道为什么,希望了解数值优化的同学帮帮忙啊.