[机器学习实验3]Logistic Regression and Newton Method

第三篇实验记录，logic回归的，一个二分类问题，输入变量的维度比较低（二维），所以可以在图像画出边界线来直观感受的。

题目链接：
http://openclassroom.stanford.edu/MainFolder/DocumentPage.php?course=MachineLearning&doc=exercises/ex4/ex4.html

这个是测试的数据集直接的显示，数据可以在这里下载
http://openclassroom.stanford.edu/MainFolder/courses/MachineLearning/exercises/ex4materials/ex4Data.zip
我们的假设模型的函数如下：

我们要进行求最优解的代价函数如下：

这个公式是由假设这些训练集的每个数据是独立的情况下推导出来的


我们要求的就是J(θ)的最大值（极大似然估计），我们可以选用之前实验使用的梯度下降法，但是该方法的迭代次数较多，所以本次实验中使用的是牛顿迭代法，这个方法网上有很多其他的资料，这里就不细讲了，推导也不难。

我们求出代价函数一阶导的零点即可（极大极小都这么求）

用Hessian矩阵表示就是:

在贴代码前，先分析下各个计算公式中变量的维度（矩阵行列数），这样方便写代码，θ的维度是3维的，所以H矩阵的肯定要是3X3的，其次，对于θ 的迭代公式中的各项要注意计算H阵时是x(i)各项与h(x(i))各项相乘得到一个新的矩阵后再乘x矩阵，这个计算H矩阵的地方要转换成对角阵来求才能达到公式中的效果，在代码里面会有说明

下面贴上代码

function Logistic_Regression_and_Newton_Method()
x = load('ex4x.dat');
y = load('ex4y.dat');
[m, n] = size(x);
x = [ones(m, 1), x];

% find returns the indices of the
% rows meeting the specified condition
pos = find(y == 1); neg = find(y == 0);

% Assume the features are in the 2nd and 3rd
% columns of x
plot(x(pos, 2), x(pos,3), '+'); hold on
plot(x(neg, 2), x(neg, 3), 'o')

g = inline('1.0 ./ (1.0 + exp(-z))'); 
% Usage: To find the value of the sigmoid 
% evaluated at 2, call g(2)
% Initialize fitting parameters
theta = zeros(n+1, 1);
iter = 10;%迭代次数
J = zeros(iter, 1);
    for i =1:iter
        z = x*theta;%列向量 80x1
        %h=g(z)
%            temp1 = diag(g(z));
%            temp2 = diag(1-g(z));
        p=(1/m).*x' * (g(z)-y);
        %这里对公式里面的h阵转化为对角阵（除主对角线外都是0)才能进行计算（要满足矩阵乘法的规则还要满足公式中每项的对应的元素相乘，而不是乘积和）
        H = (1/m).*x'*(diag(g(z))*diag(1-g(z))*x);%3x3
        theta = theta-inv(H)*p;
        temp = y.*log(g(z));
        J(i) =(1/m)*sum(-y.*log(g(z)) - (1-y).*log(1-g(z)));
    end

theta

% Plot Newton's method result
% Only need 2 points to define a line, so choose two endpoints
plot_x = [min(x(:,2))-2,  max(x(:,2))+2];
% Calculate the decision boundary line
plot_y = (-1./theta(3)).*(theta(2).*plot_x +theta(1));
plot(plot_x, plot_y)
legend('Admitted', 'Not admitted', 'Decision Boundary')
hold off

% Plot J
figure
plot(0:iter-1, J, 'o--', 'MarkerFaceColor', 'r', 'MarkerSize', 8)
xlabel('Iteration'); ylabel('J')
% Display J
J

prob = 1 - g([1, 20, 80]*theta)

end

得到的分界线如下：

分界线的公式:θ0+θ1*x+θ2*y=0;

最后的运行结果回答几个问题：

1、
2、
下面是代价函数的收敛过程的示意，可以看出，在迭代到第四次的时候基本就已经收敛了，所以是远远比梯度下降的方法是要快的。

但是牛顿迭代法在H矩阵为非正定阵时是无法收敛的，需要使用拟牛顿法，这方面还没有研究到，后续遇到再了解吧。