编程作业（五）

正则化的线性回归以及偏差VS方差

正则化的线性回归

背景：数据集中包含水位变化的历史记录x和水坝的水量y。

任务一 可视化数据集

我们将数据集分为三部分：

• 训练集：X，y
• 交叉验证集：Xval，yval
• 测试集：Xtest，ytest

因此，本任务只需将训练集可视化即可。在ex5.m文件已将该任务代码准备好了，我们只需运行即可：

% Load from ex5data1: % You will have X, y, Xval, yval, Xtest, ytest in your environmentload ('ex5data1.mat');% m = Number of examplesm = size(X, 1);% Plot training dataplot(X, y, 'rx', 'MarkerSize', 10, 'LineWidth', 1.5);xlabel('Change in water level (x)');ylabel('Water flowing out of the dam (y)');

其运行结果为：

任务二 线性回归的正则化代价函数

我们先将正则化代价函数公式列出：

其向量化后公式为：
　J(θ) = ((Xθ - y)^T(Xθ - y) + λθ_t^Tθ_t) / 2m
其中，θ_t表示将θ的第一列替换为0，即θ₀不参与正则化操作。

因此，我们可根据上述公式在linearRegCostFunction.m文件中键入如下代码：

theta_1 = [0; theta(2:end)];J = ((X * theta - y)' * (X * theta - y)) / (2 * m) + lambda / (2 * m) * theta_1' * theta_1;

任务三 线性回归的正则化下降梯度

同样的，我们先将公式列出：

其向量化后的公式为：
　grad = (X^T(Xθ - y) + λθ_t) / m
其中θ_t同上。

因此，在linearRegCostFunction.m文件中继续键入如下代码：

grad = (X' * (X * theta - y) + lambda * theta_1) / m;

任务四 拟合线性回归

一旦你的代价函数和下降梯度运行正常，下一步就是在ex5.m文件调用并运行trainLinearReg.m文件中的代码，通过使用fmincg函数计算出使得代价函数最小化的θ。

在该任务中参数θ为2维向量，因此我们将正则化参数λ的值设为0。为何将正则化参数λ的值设为0？这是因为正则化对于低维度的θ没有太大的帮助。

trainLinearReg.m文件中的代码如下：

function [theta] = trainLinearReg(X, y, lambda)%TRAINLINEARREG Trains linear regression given a dataset (X, y) and a%regularization parameter lambda%   [theta] = TRAINLINEARREG (X, y, lambda) trains linear regression using%   the dataset (X, y) and regularization parameter lambda. Returns the%   trained parameters theta.%% Initialize Thetainitial_theta = zeros(size(X, 2), 1); % Create "short hand" for the cost function to be minimizedcostFunction = @(t) linearRegCostFunction(X, y, t, lambda);% Now, costFunction is a function that takes in only one argumentoptions = optimset('MaxIter', 200, 'GradObj', 'on');% Minimize using fmincgtheta = fmincg(costFunction, initial_theta, options);end

ex5.m文件中该部分代码如下：

%  Train linear regression with lambda = 0lambda = 0;[theta] = trainLinearReg([ones(m, 1) X], y, lambda);%  Plot fit over the dataplot(X, y, 'rx', 'MarkerSize', 10, 'LineWidth', 1.5);xlabel('Change in water level (x)');ylabel('Water flowing out of the dam (y)');hold on;plot(X, [ones(m, 1) X]*theta, '--', 'LineWidth', 2)hold off;fprintf('Program paused. Press enter to continue.\n');pause;

运行结果为：

偏差与方差

高偏差的模型通常对训练集的拟合不太好，即欠拟合问题；高方差的模型通常对训练集的拟合非常完美，但对于交叉验证集或测试集的拟合不太好，即过拟合问题。

因此，该小节将练习绘制学习曲线来诊断偏差与方差的问题。

任务一 学习曲线

为了绘制学习曲线，我们需要计算出J_train(θ)和J_CV(θ)。

其中J_train(θ)的计算公式为：

J_CV(θ)的计算公式为：

因此，我们先需要利用trainLinearReg函数计算出使得代价函数最下化的θ的值；然后在使用linearRegCostFunction函数分别计算J_train(θ)和J_CV(θ)。

注：在使用linearRegCostFunction函数时，要注意将正则化参数λ = 0。

learningCurve.m文件中的具体代码如下：

for i = 1 : m    theta = trainLinearReg(X(1:i, :), y(1:i), lambda);    error_train(i) = linearRegCostFunction(X(1:i, :), y(1:i), theta, 0);    error_val(i) = linearRegCostFunction(Xval, yval, theta, 0);end

该部分的运行结果为：

多项式回归

在之前的部分，我们的线性模型对数据的拟合不太好，即出现欠拟合问题。在本小节，我们通过增加特征变量来解决欠拟合问题。

对于多项式回归，我们的假设函数h_θ(x)为：
h_θ(x) = θ₀ + θ₁ * (waterLevel) + θ₂ * (waterLevel)² + ... + θ_p * (waterLevel)^p

现在，我们需要在数据集中增添高阶幂的特征变量。因此，我们需要在polyFeatures.m文件中键入相关代码，使得数据集X变为一个m*p的矩阵。

polyFeatures.m文件的相关代码如下：

for i = 1 : p    X_poly(:, i) = X .^ i;end

任务一 学习多项式回归

对于该部分的练习，我们使用8次幂的多项式回归模型。由于特征变量在多项式回归模型中，其取值范围各不相同。因此，我们需要对特征变量归一化。

featureNormalize.m文件中特征变量归一化代码如下：

function [X_norm, mu, sigma] = featureNormalize(X)%FEATURENORMALIZE Normalizes the features in X %   FEATURENORMALIZE(X) returns a normalized version of X where%   the mean value of each feature is 0 and the standard deviation%   is 1. This is often a good preprocessing step to do when%   working with learning algorithms.mu = mean(X);X_norm = bsxfun(@minus, X, mu);sigma = std(X_norm);X_norm = bsxfun(@rdivide, X_norm, sigma);% ============================================================end

然后，我们将正则化参数λ = 0，利用trainLinearReg函数计算出使得代价函数最下化的θ的值。

最后，我们利用linearRegCostFunction函数分别计算Jtrain(θ)和JCV(θ)，绘制出学习曲线。

该部分代码如下：

%% =========== Part 6: Feature Mapping for Polynomial Regression =============%  One solution to this is to use polynomial regression. You should now%  complete polyFeatures to map each example into its powers%p = 8;% Map X onto Polynomial Features and NormalizeX_poly = polyFeatures(X, p);[X_poly, mu, sigma] = featureNormalize(X_poly);  % NormalizeX_poly = [ones(m, 1), X_poly];                   % Add Ones% Map X_poly_test and normalize (using mu and sigma)X_poly_test = polyFeatures(Xtest, p);X_poly_test = bsxfun(@minus, X_poly_test, mu);X_poly_test = bsxfun(@rdivide, X_poly_test, sigma);X_poly_test = [ones(size(X_poly_test, 1), 1), X_poly_test];         % Add Ones% Map X_poly_val and normalize (using mu and sigma)X_poly_val = polyFeatures(Xval, p);X_poly_val = bsxfun(@minus, X_poly_val, mu);X_poly_val = bsxfun(@rdivide, X_poly_val, sigma);X_poly_val = [ones(size(X_poly_val, 1), 1), X_poly_val];           % Add Onesfprintf('Normalized Training Example 1:\n');fprintf('  %f  \n', X_poly(1, :));fprintf('\nProgram paused. Press enter to continue.\n');pause;%% =========== Part 7: Learning Curve for Polynomial Regression =============%  Now, you will get to experiment with polynomial regression with multiple%  values of lambda. The code below runs polynomial regression with %  lambda = 0. You should try running the code with different values of%  lambda to see how the fit and learning curve change.%lambda = 0;[theta] = trainLinearReg(X_poly, y, lambda);% Plot training data and fitfigure(1);plot(X, y, 'rx', 'MarkerSize', 10, 'LineWidth', 1.5);plotFit(min(X), max(X), mu, sigma, theta, p);xlabel('Change in water level (x)');ylabel('Water flowing out of the dam (y)');title (sprintf('Polynomial Regression Fit (lambda = %f)', lambda));figure(2);[error_train, error_val] = ...    learningCurve(X_poly, y, X_poly_val, yval, lambda);plot(1:m, error_train, 1:m, error_val);title(sprintf('Polynomial Regression Learning Curve (lambda = %f)', lambda));xlabel('Number of training examples')ylabel('Error')axis([0 13 0 100])legend('Train', 'Cross Validation')fprintf('Polynomial Regression (lambda = %f)\n\n', lambda);fprintf('# Training Examples\tTrain Error\tCross Validation Error\n');for i = 1:m    fprintf('  \t%d\t\t%f\t%f\n', i, error_train(i), error_val(i));endfprintf('Program paused. Press enter to continue.\n');pause;

运行结果为：

学习曲线

任务二 调整正则化参数（选做）

当正则化参数λ = 1时，其运行结果为：

学习曲线

当正则化参数λ = 100时，其运行结果为：

学习曲线

任务三 通过交叉验证集选择正则化参数

当正则化参数λ∈{0, 0.001, 0.003, 0.03, 0.1, 0.3, 1, 3, 10}时，分别计算出J_train(θ)和J_CV(θ)。

在validationCurve.m文件中键入如下代码：

for i = 1 : length(lambda_vec)    lambda = lambda_vec(i);    theta = trainLinearReg(X, y, lambda);    error_train(i) = linearRegCostFunction(X, y, theta, 0);    error_val(i) = linearRegCostFunction(Xval, yval, theta, 0);end

然后，ex5.m文件中的相关代码通过J_train(θ)和J_CV(θ)的值，绘制出相关函数图。

[lambda_vec, error_train, error_val] = ...    validationCurve(X_poly, y, X_poly_val, yval);close all;plot(lambda_vec, error_train, lambda_vec, error_val);legend('Train', 'Cross Validation');xlabel('lambda');ylabel('Error');fprintf('lambda\t\tTrain Error\tValidation Error\n');for i = 1:length(lambda_vec)    fprintf(' %f\t%f\t%f\n', ...            lambda_vec(i), error_train(i), error_val(i));endfprintf('Program paused. Press enter to continue.\n');pause;

其运行结果为：

任务四 计算测试集误差（选做）

在实际开发中，除了计算J_train(θ)和J_CV(θ)，我们还需计算J_test(θ)。

参考代码：

for i = 1 : m    theta = trainLinearReg(X(1:i, :), y(1:i), lambda);    error_test(i) = linearRegCostFunction(Xtest(1:i, :), ytest(1:i), theta, 0);end

任务五 绘制随机数据集的学习曲线（选做）

从数据集中随机抽取60%的数据作为训练集，20%的数据作为交叉验证集和20%的数据作为测试集。由于此部分答案不唯一就不过多叙述。