编程作业ex2：Logistic 回归

1、logistic回归

1.1 可视化数据

打开ex2data1.txt观察数据

  第一列和第二列为两次考试的成绩，第三列代表该生是否能被录取，1代表录取，0代表不录取

读取数据：

data = load('ex2data1.txt');
X = data(:, [1, 2]); y = data(:, 3); % 将数据集前两列赋给X，将数据集最后一列赋给y

将样本数据绘制在图上，坐标轴是两次考试的成绩，两种标记代表了他们能否被录取。

function plotData(X,y)
figure;
hold on;
% 用find函数寻找1和0的项，返回的位置分别存在pos和neg中
pos = find(y == 1); neg = find(y == 0);
% 画图
plot(X(pos, 1), X(pos, 2), 'k+','LineWidth', 2, 'MarkerSize', 7);
% X(pos,1),X(pos,2)为能够被录取样本的两次考试成绩，k+表示线型为加号型，黑色，线宽为2，标记大小为7
plot(X(neg, 1), X(neg, 2), 'ko', 'MarkerFaceColor', 'y','MarkerSize', 7);
% X(neg,1),X(neg,2)为不能被录取样本的两次考试成绩值，ko表示线型为圆圈，黑色，填充色为黄色，标记大小为7
hold off;
end

plot函数的使用方法

find函数的使用方法

绘制的样本如下图：

  1.2 实现

1.2.1 sigmoid 函数

逻辑回归的假设函数：

其中g是sigmoid函数：

   完成sigmoid.m的代码：

function g = sigmoid(z)
g = zeros(size(z)); % 初始化g
g = 1./ (1+exp(-z));
end

1.2.2 代价函数和梯度（完成costFunction.m中的代码）

复习一下logistic回归的代价函数：

求导得：

function [J, grad] = costFunction(theta, X, y)
m = length(y); % number of training examples
J = 0;
grad = zeros(size(theta)); % 初始化J和grad（梯度）
J = -1*sum(y.*log(sigmoid(X*theta))+(1-y).*log((1-sigmoid(X*theta))))/m;
grad = (X'*(sigmoid(X*theta)-y))/m;  % 求梯度
end

运行结果：

 1.2.3 使用fminunc函数寻找最佳参数theta

%  Set options for fminunc
options = optimset('GradObj', 'on', 'MaxIter', 400);

%  Run fminunc to obtain the optimal theta
%  This function will return theta and the cost 
[theta, cost] = ...
    fminunc(@(t)(costFunction(t, X, y)), initial_theta, options);

GradObj:On 表示使用自己的梯度公式，Maxlter:400表示最大迭代次数

运行结果：

得到的theta将被用于接下来决策边界的绘制。

  绘制决策边界：

决策边界的表达式，这里我们只考虑本题的输入，sigmoid函数的分界点是z=0，即

function plotDecisionBoundary(theta, X, y)

% Plot Data
plotData(X(:,2:3), y); % X的第2、3两列和y进行作图
hold on

if size(X, 2) <= 3  % size(X,2)返回X的列数，若列数小于等于3，则可以绘制直线作为决策边界
    % Only need 2 points to define a line, so choose two endpoints 两点确定直线
    plot_x = [min(X(:,2))-2,  max(X(:,2))+2]; %取第二列的最小值和最大值作为x_2

    % Calculate the decision boundary line 根据边界函数的方程计算出x_3的值
    plot_y = (-1./theta(3)).*(theta(2).*plot_x + theta(1));

    % Plot, and adjust axes for better viewing
    plot(plot_x, plot_y)
    
    % Legend, specific for the exercise
    legend('Admitted', 'Not admitted', 'Decision Boundary')
    axis([30, 100, 30, 100])
else   % 若X的列数大于三，则需要绘制等高线图来作为决策边界
    % Here is the grid range
    u = linspace(-1, 1.5, 50);
    v = linspace(-1, 1.5, 50);

    z = zeros(length(u), length(v));
    % Evaluate z = theta*x over the grid
    for i = 1:length(u)
        for j = 1:length(v)
            z(i,j) = mapFeature(u(i), v(j))*theta;
        end
    end
    z = z'; % important to transpose z before calling contour

    % Plot z = 0
    % Notice you need to specify the range [0, 0]
    contour(u, v, z, [0, 0], 'LineWidth', 2)
end
hold off

end

运行结果：

 1.2.4 评估模型精确度（完成predict.m中的代码）

输入X（预测集）和theta（参数），返回预测结果向量。再进行判断，每个值是 0 或 1 ；

function p = predict(theta, X)
m = size(X, 1);
%初始化p
p = zeros(m, 1);
%————————————————————————————————————————————————————————————————————————
k = find(sigmoid( X * theta) >= 0.5 ); % 寻找大于0.5的位置
p(k)= 1; 
d = find(sigmoid( X * theta) < 0.5 ); % 寻找小于0.5的位置
p(d)= 0;
end
%————————————————————————————————————————————————————————————————————————

　　p_Matric=sigmoid(X*theta);%计算sigmoid (m*1矩阵)

　　pos=find(p_Matric>0.5);%寻找>0.5的位置

　　p(pos,:)=1;%设置为1

end
%————————————————————————————————————————————————————————————————————————

调用：

prob = sigmoid([1 45 85] * theta);
fprintf(['For a student with scores 45 and 85, we predict an admission ' ...
         'probability of %f
'], prob);
fprintf('Expected value: 0.775 +/- 0.002

');

% Compute accuracy on our training set
p = predict(theta, X);

fprintf('Train Accuracy: %f
', mean(double(p == y)) * 100);
fprintf('Expected accuracy (approx): 89.0
');
fprintf('
');

运行结果：

 2.正则化逻辑回归

 假设你是芯片工厂的产品经理，要判断一批芯片是否通过工厂的质量检测（QA），现在你有一批芯片在两次测试中的结果，用这两次结果的数据集来构建回归模型。

2.1 可视化数据

仍然使用plotdata.m中的代码，结果如下：

 此数据集我们不能再用线性决策边界进行分割

2.2 特征映射

将特征映射到x1，x2的所有多项式，直到六次幂

 由于此映射，我们的两个特征的向量（两个QA测试的得分）已转换为28维向量。 在此较高维特征向量上训练的逻辑回归分类器将具有更复杂的决策边界，并且在我们的二维图中绘制时将显示非线性。

虽然特征映射能够帮助我们建立更具表现力的分类器，但是容易过度使用，因此我们可以使用正则化的方法来处理数据

2.3 代价函数和梯度

Logistic回归中的正则化代价函数：

梯度：

 

 注意：参数θ_0不需要正则化，而在MATLAB中索引是从1开始的，因此theta(1) （在这里代表了θ_0）不需要正则化

function [J, grad] = costFunctionReg(theta, X, y, lambda)
% Initialize some useful values
m = length(y); % number of training examples

% You need to return the following variables correctly 
J = 0;
grad = zeros(size(theta));

theta_1=[0;theta(2:end)];    % 先把theta(1)拿掉，不参与正则化
J= -1 * sum( y .* log( sigmoid(X*theta) ) + (1 - y ) .* log( (1 - sigmoid(X*theta)) ) ) / m  + lambda/(2*m) * theta_1' * theta_1; 
grad = ( X' * (sigmoid(X*theta) - y ) )/ m + lambda/m * theta_1; 

end

运行结果：

 2.4 绘制决策边界

仍然调用plotdecisionboundary.m中的函数

绘制结果：

 接下来尝试更改正则化参数lambda，可以看到决策边界的变化：

lambda=100

  lambda=0