機器學習 - 編程練習（一）：線性回歸

編程練習（一）：線性回歸

---

文件清單

• ex1.m
• ex1_multi.m
• ex1data1.txt - ex1.m 用到的數據組
• ex1data2.txt - ex1_multi.m 用到的數據組
• submit.m - 提交代碼
• [*] warmUpExercise.m
• [*] plotData.m
• [*] computeCost.m
• [*] gradientDescent.m
• [+] computeCostMulti.m
• [+] gradientDescentMulti.m
• [+] featureNormalize.m
• [+] normalEqn.m

* 為必須要完成的
+ 為可選

1 簡單的Octave/MATLAB函數

完成 warmUpExercise.m ，使得 A 為 5x5 單位矩陣，送分的。

A = eye(5); 

輸出結果：

ans =
Diagonal Matrix
	1 	0 	0 	0 	0
	0 	1 	0 	0 	0
	0 	0 	1 	0 	0
	0 	0 	0 	1 	0
	0 	0 	0 	0 	1

1.1 提交代碼

使用 submit

解決不能上傳的問題：

Windows

進入練習的 /lib 文件夾中，編輯 submitWithConfiguration.m 修改第66行，然后重啟octave

responseBody = urlread(submissionUrl, 'post', params);
% 改為
[code, responseBody] = system(sprintf('echo jsonBody=%s | curl -k -X POST -d @- %s', body, submissionUrl));

OSX Linux

進入練習的 /lib 文件夾中，編輯 submitWithConfiguration.m 修改第66行，然后重啟octave

responseBody = urlread(submissionUrl, 'post', params);
% 改為
[code, responseBody] = system(sprintf('echo ''jsonBody=%s'' | curl -k -X POST -d @- %s', body, submissionUrl));

2 單變量線性回歸

背景：假設我們現在是個連鎖餐廳的老板，已經在很多城市開了連鎖店（提供訓練組），現在想再開新店，需要通過以前的數據預測新開的店的收益。

ex1data1.txt 提供所需要的訓練組，第一列是城市人口，第二列是對應的收益。負值代表着虧損。

ex1.m 已經加載了數據。

2.1 可視化數據

ex1.m 中加載了變量 X 和 y：

data = load('ex1data1.txt'); % read comma separated data
X = data(:, 1); y = data(:, 2);
m = length(y); % number of training examples

接着 plotData 函數將其繪制出來，我們需要完成 plotData.m 進行繪圖：

plot(x, y, 'rx', 'MarkerSize', 10); % Plot the data
ylabel('Profit in $10,000s'); % Set the y-axis label
xlabel('Population of City in 10,000s'); % Set the x-axis label

現在運行 ex1.m 就能看到類似圖1的圖像了：

* 通過看手冊學習函數是個很好的方式 help plot 。

2.2 梯度下降

在這一結中我們要通過梯度下降計算出合適的線性回歸參數theta。

2.2.1 更新等式

線性回歸的目標是將成本函數最小化：

\[J(\theta) = \frac{1}{2m}\sum_{i=1}^m (h_{\theta}(x_{i})-y_{i})^{2} \]

其中線性模型的假設 \(h_\theta(x)\) 為：

\[h(x) = \theta_{0} + \theta_{1}x \]

記住我們的變量是 \(\theta\) ，通過調整 \(\theta\) 的值來最小化成本函數。一種方式就是梯度下降，每次迭代都會更新 \(\theta\) ：

\[\theta_{j} := \theta_{j} - \alpha\frac{1}{m}\sum_{i=1}^m ((h_{\theta}(x^{(i)})-y^{(i)})x_j^{(i)}) \]

* 與 = 不同，:= 代表着 同時更新(simultaneously update)，簡單來說就是先將計算的 \(\theta\) 存入一個臨時變量，最后所有 \(\theta\) 都計算完了一起再賦值回去。

例如：

temp1 = theta1 - （theta1 - 10 * theta2 * x1） ;
temp2 = theta2 - （theta1 - 10 * theta2 * x2） ;
theta1 = temp1 ;
theta2 = temp2 ;

注意，我們使用梯度下降時在矩陣 X 第一列增加了一列 1

X = [ones(m, 1), data(:,1)]; % Add a column of ones to x

原因是成本函數的公式中 \(\theta_0\) 的系數為1.

2.2.2 工具

X = [ones(m, 1), data(:,1)]; % Add a column of ones to x
theta = zeros(2, 1); % initialize fitting parameters
iterations = 1500;
alpha = 0.01;

初始化 theta ，迭代次數 iteration ，Learning rate alpha 。

2.2.3 計算成本函數

在這一節我們需要根據成本函數完成 computeCost.m，注意 X，y是訓練組中的數據矩陣而不是非矢量數(scalar values)。

當正確完成成本函數后，ex1.m 下一步將 theta 初始化為0，運行成功后會輸出結果 32.07 。

* 完成后使用 submit 提交代碼。

參考代碼：

% ====================== YOUR CODE HERE ===================================

J = sum((X*theta-y).^2)/(2*length(y));

% =========================================================================

2.2.4 梯度下降

接下來我們需要完成 gradientDescent.m 。其中循環結構已經寫好了，僅需要寫每次更新 theta 的函數即可。

編程時要明確的知道，成本函數J是改變 theta 來最小化的，而不是 X 或者 y。

怎么樣驗證自己的梯度下降工作正常？

可以通過在循環中打印 computeCost 的結果，如果梯度下降工作正常，則打印的結果應該慢慢減少然后無變化。

完成函數后，ex1.m將會使用最后的theta繪制出一條匹配的直線。如圖：

最后得到的 \(\theta\) 我們用來預測人口 35,000 和 70,000 地區的收益。

* 完成后使用 submit 提交代碼。

*注意，Octave/MATLAB中 A*B 默認為矩陣乘法，使用 A.*B 可以用逐元素乘法。

參考代碼：

% ====================== YOUR CODE HERE ======================

theta = theta - X'*(X*theta-y)/m*alpha;

% ============================================================

2.4 可視化成本函數

為了更好的理解成本函數，在這一節我們將構建一個三維的圖形。

ex1.m 中已經寫好了這一部分的代碼，但是我們需要理解為什么代碼這么寫。

% initialize J vals to a matrix of 0's
J vals = zeros(length(theta0 vals), length(theta1 vals));
% Fill out J vals
for i = 1:length(theta0 vals)
	for j = 1:length(theta1 vals)
		t = [theta0 vals(i); theta1 vals(j)];
		J_vals(i,j) = computeCost(x, y, t);
	end
end

通過 computeCost 函數，我們計算出對應的值並儲存到一個以 theta0 和 theta1 為坐標的矩陣 J_vals 中

通過使用 surf 和 contour 命令，我們可以構建出圖形：

這些圖像的目的是為了展示 \(J(\theta)\) 隨着 \(\theta_0\) 和 \(\theta_1\) 的改變的值。（在2D輪廓圖中比3D的更直觀）。最小點是 \(\theta_0\) 和 \(\theta_1\) 的最適點， 每一步梯度下降都會更靠近這個點。

我們可以通過旋轉看到為什么叫“輪廓圖”：

附加練習

---

3 多元線性回歸

在這一節中，我們將會利用多元線性回歸預測房子的價格。假設我們要出售我們的房子，現在想知道在房市中能賣一個什么價格。一種方式是收集最近出售的房價數據，構建一個房子價格的模型。

ex1data2.txt中包含了房子價格的訓練組。第一列是房子的尺寸（平方英尺），第二列是卧室的數量，第三列是房子的價格。

我們通過 ex1_multi.m 完成本次任務。

3.1 特性一般化(Feature Normalization)

首先我們需要用特性放縮讓數據的范圍縮小，使得梯度下降計算的更快：

我們的任務是完成 featureNormalize.m：

• 計算每個特性的平均值(mean)
• 計算標准差(standard deviations)
• 特性放縮(feature scaling)

* 我們這里利用的是標准差(standard deviation)，也可以使用差值(max - min)。

* 完成后使用 submit 提交代碼。

參考代碼：

% ====================== YOUR CODE HERE ======================   

mu = mean(X);
sigma = std(X);
X_norm = (X - mu) ./ sigma;

% ============================================================

3.2 梯度下降

我們之前用梯度下降解決單變量線性回歸問題，和多元線性回歸問題的唯一區別就是矩陣 X 的特性不止一個，假設函數的公式沒有變化。

我們需要完成 computeCostMulti.m 和 gradientDescentMulti.m 完成多元線性回歸的成本函數和梯度下降。如果之前的梯度下降適用於多變量，這里也可以使用。

* 確保代碼支持多個特性，而且滿足向量運算條件。可以使用 size(X, 2) 找出訓練組中有多少特性。

* 完成后使用 submit 提交代碼。

* 多變量中， 成本函數也可以寫成向量計算形式：

\(J(\theta) = \frac{1}{2m} (X\theta - \vec{y})^{T} (X\theta - \vec{y})\）

參考代碼：

% ====================== YOUR CODE HERE ======================

theta = theta - X'*(X*theta-y)/m*alpha;

% ============================================================

3.2.1 附加練習：學習速率的選擇（Learning rate）

對於不同的訓練組，選擇合適的學習速率可以使得結果計算的更快更准確。我們可以在 ex1_multi.m 中修改 alpha 值。

實際經驗我們可以一步一步的改變值，例如 0.3 -> 0.1 -> 0.03 -> 0.01 等等。

不合適學習速率

過大：導致找不到最適點，Octave/MATLAB會返回一個 NaNs 代表 not a number ，也就是正無窮或者負無窮。

成本函數和迭代次數的圖形如下：

（！圖片）

過小：導致梯度下降計算耗時嚴重。

合適的學習速率

接下來，我們需要用訓練好的 theta 值來預測 1650 平方英尺 3 個卧室的房子的價格。

注意不要忘記放縮！

參考代碼：

% ====================== YOUR CODE HERE ======================
P = [1 (([1650 3] - mu) ./ sigama)];
price = P*theta;
% ============================================================

3.3 一般等式

這一結我們需要完成 normalEqn.m 利用公式：

\(\theta = (X^T X)^{-1}X^T y\)

完成之后我們可以預測 1650 平方英尺 3 個卧室的房子的價格。

參考代碼：

% ====================== YOUR CODE HERE ======================
price = 1650; % You should change this
price = [1,1650,3]*theta;
% ============================================================