線性回歸與梯度下降法[二]——優化與比較

接着上文——機器學習基礎——梯度下降法（Gradient Descent）往下講。這次我們主要用matlab來實現更一般化的梯度下降法。由上文中的幾個變量到多個變量。改變算法的思路，使用矩陣來進行計算。同時對算法的優化和調參進行總結。即特征縮放（feature scaling）問題和學習速率\(\alpha\)的取值問題。還有在擬合線性模型時，如何選擇正確的算法，梯度下降 or 最小二乘法？

matlab的基本用法，已經總結在matlab基礎教程——根據Andrew Ng的machine learning整理上。相對於python來說，的確matlab的語法更簡單，上手更快，只要簡單學習下就能上手編寫算法，或者用於計算。

matlab實現一般梯度下降法

針對於Andrew Ng的machine learning課當中的數據，我們這里給出了matlab的一般梯度下降法的代碼。
在了解了算法原理后，實現起來就非常簡單了，主要分為以下幾步:

1. 數據加載和預處理
   從文件中讀入數據，並分塊放到相應的矩陣中去。這里我們選取的模型是\(y = \theta_0+\theta_1x_1\)。顯然\(x_0\)取1，所以需要在原有的X矩陣中添加一列全1的列。

data = load ('data1.txt');
X = data(:, 1);
y = data(:, 2);
X = [ones(length(data),1),data(:,1)];

2. 設置算法的學習參數

X = [ones(length(data),1),data(:,1)];
theta = zeros(2,1);%theta初始值為0
alpha = 0.01;
max_iter = 5000;
m = length(y);%數據組數
J_history = zeros(max_iter,1);%初始化迭代誤差變量

3. 執行算法。

for iter = 1:max_iter
	theta = theta - alpha / m * X' * (X * theta - y);	%梯度下降
    J_history(iter) = sum((X * theta-y) .^ 2) / (2*m);	%記錄每次 迭代后的全局誤差
    fprintf('iter:%d ------ Error:%f\n',iter,J_history(iter));
end

下面是完整的代碼：https://github.com/maoqyhz/machine_learning_practice

%加載數據和數據預處理
data = load ('data1.txt');
fprintf('Running Gradient Descent ...\n');

X = data(:, 1); 	%取X集合
y = data(:, 2);		%取y集合

%畫出數據的散點圖
figure;
plot(X, y, 'rx', 'MarkerSize', 10);
ylabel('Profit in $10,000s'); 
xlabel('Population of City in 10,000s'); 
hold on; 

%設置學習參數
X = [ones(length(data),1),data(:,1)]; 	% y = theta0*x0 + theta1*x1 默認x0為1
theta = zeros(2,1);	%theta初始值為0
alpha = 0.01;
max_iter = 5000;

m = length(y);	%數據組數
J_history = zeros(max_iter,1);	%初始化迭代誤差變量
iter = 1;
for iter = 1:max_iter
	%每迭代100次畫一條曲線
    if  mod(iter,100) == 0
        plot(X(:,2), X*theta, 'g')
    end

	theta = theta - alpha / m * X' * (X * theta - y);	%梯度下降
    J_history(iter) = sum((X * theta-y) .^ 2) / (2*m);	%記錄每次迭代后的全局誤差
    fprintf('iter:%d ------ Error:%f\n',iter,J_history(iter));
end

%輸出最終的theta值和最終的擬合曲線
disp('Theta found by gradient descent:');
disp(theta);
plot(X(:,2), X*theta, 'k')
legend('Training data', 'Linear regression')
hold off 

運行結果如下：

• 圖1為得到的結果，可以看到在現有的學習速率下，基本上4600次左右就能收斂了。得出了最終的\(\theta\)值。
  

• 圖2為擬合的過程，每迭代100次會畫出一條綠色曲線，不斷迭代最終收斂的曲線是黑色的那條，可以清楚看見擬合的過程。
  

特征縮放和數據的歸一化

下圖為上例中的數據部分，可以看到第一列和第二列的數據差距並不大，所以不需要很多的迭代次數就能收斂。但是如果數據之間差距很大，在計算過程中會導致數據超出數據類型的最大存儲值，或者是需要額外的迭代次數來到達收斂的情況，所以此時我們對數據進行處理，也就是所說的特征縮放。

數據歸一化

數據歸一化是非常有用的特征縮放方法，可以把數據縮放到\([-1,1]\)中。常見的歸一化方法有兩種。

min-max標准化（Min-Max Normalization）

可以使數據進行線性變換，使結果值映射到\([0,1]\)之間。設原區間為\([a,b]\)，轉換函數如下：

\[x^* = \frac{x-a}{b-a} \]

Z-score標准化方法（Z-socre Normalization）

這種方法給予原始數據的均值（mean）和標准差（standard deviation）進行數據的標准化。經過處理的數據符合標准正態分布，即均值為0，標准差為1，轉化函數為：

\[x^* = \frac{x-\mu}{\sigma} \]

其中\(\mu\)為數據的均值，\(\sigma\)為數據的標准差。

matlab來處理數據歸一化問題

新的數據集如下，顯然，不同的特征之前的差距比較大，此時就需要進行歸一化處理。

1. 數據歸一化

%Z-score Normalization
function [ X_norm,avg,sigma ] = normalize(X)
	avg = mean(X,1); %均值
	sigma = std(X,1);%標准差
	%repmat函數用來數據填充
	X_norm  = (X - repmat(avg,size(X,1),1)) ./  repmat(sigma,size(X,1),1);
end

2. 測試樣本需要進行歸一化。

% 使用1 1650 3對數據進行估計
x = [1650 3];
price = [1 (([1650 3] - avg) ./ sigma)] * theta ;

學習速率的取值

可以看到學習速率太小，需要進行多次迭代才能收斂；學習速率太大，可能會miss最小值兒無法收斂。因此我們需要選取合適的學習速率Andrew Ng已經告訴我們方法，如圖所示，選擇學習速率的區間，不斷地去調整以滿足自己的數據和算法。說到底機器學習就是個調參的過程，啊哈哈~~

梯度下降法VS最小二乘法

對於線性回歸問題，我們既可以用梯度下降法也可以用最小二乘法來解決。那么這兩種算法有和不同？簡而言之，前者是需要多次迭代來收斂到全局最小值，后者則是提供一種解析解法，直接一次性求得\(\theta\)的最優值。

對於最小二乘法，我們這里只給出具體的結論。證明牽涉到矩陣的求導，可以直接看西瓜書。

那么，當\(X^TX\)為滿秩矩陣或正定矩陣時，我們可以一步得到\(\theta\)的最優值。

\[\theta = (X^TX)^{-1}X^Ty \]

而在matlab中，最小二乘法的求解過程也的確只需要一行代碼就可以完成。把下面的代碼替換掉梯度下降的代碼即可。

theta = pinv( X' * X ) * X' * y;

How to choose？

梯度下降法

• 需要設置學習速率。
• 需要進行多次迭代。
• 數據比較大時（百萬？），可以很好的工作。

最小二乘法

• 不需要設置學習速率。
• 不需要迭代。
• 需要矩陣計算，數據大時，計算速度慢。

參考文獻

• Coursera機器學習課
• 數據歸一化和兩種常用的歸一化方法
• 在博客中使用LaTeX插入數學公式
• 機器學習基礎——梯度下降法（Gradient Descent）
• matlab基礎教程——根據Andrew Ng的machine learning整理