Matlab实现线性回归和逻辑回归: Linear Regression & Logistic Regression

超、凢脫俗 2022-08-05 12:28 444阅读 0赞

本文为Maching Learning 栏目补充内容，为上几章中所提到[单参数线性回归][Link 1]、[多参数线性回归][Link 2]和 [逻辑回归][Link 3]的总结版。旨在帮助大家更好地理解回归，所以我在Matlab中分别对他们予以实现，在本文中由易到难地逐个介绍。

**本讲内容：**

**Matlab 实现各种回归函数**

**=========================**

**基本模型**

**Y=θ0\+θ1X1型---线性回归（直线拟合）**

**解决过拟合问题---Regularization**

**Y=1/(1+e^X)型---逻辑回归（sigmod 函数拟合）**

**=========================**

**第一部分：基本模型**

在解决拟合问题的解决之前，我们首先回忆一下线性回归和逻辑回归的基本模型。

设待拟合参数 θn\*1 和输入参数\[ xm\*n, ym\*1 \] 。

对于各类拟合我们都要根据梯度下降的算法，给出两部分：

① cost function（指出真实值y与拟合值h<hypothesis>之间的距离）：给出cost function 的表达式，每次迭代保证cost function的量减小；给出梯度gradient，即cost function对每一个参数θ的求导结果。

function \[ jVal,gradient \] = costFunction ( theta )

② Gradient\_descent（主函数）：用来运行梯度下降算法，调用上面的cost function进行不断迭代，直到最大迭代次数达到给定标准或者cost function返回值不再减小。

function \[optTheta,functionVal,exitFlag\]=Gradient\_descent( )

线性回归：拟合方程为hθ(x)=θ0x0\+θ1x1\+…+θnxn，当然也可以有xn的幂次方作为线性回归项（如![1341472718_8627.jpg][]），这与普通意义上的线性不同，而是类似多项式的概念。

其cost function 为：![1341901237_1654.jpg][]

逻辑回归：拟合方程为hθ(x)=1/(1+e^(θTx))，其cost function 为：![1341658176_1292.jpg][]

cost function对各θj的求导请自行求取，看[第三章][Link 3]最后一图，或者参见后文代码。

后面，我们分别对几个模型方程进行拟合，给出代码，并用matlab中的fit函数进行验证。

**第二部分：Y=θ0\+θ1X1型---线性回归（直线拟合）**

在[Matlab 线性拟合 & 非线性拟合][Matlab _ _]中我们已经讲过如何用matlab自带函数fit进行直线和曲线的拟合，非常实用。而这里我们是进行ML课程的学习，因此研究如何利用前面讲到的梯度下降法（gradient descent）进行拟合。

cost function：

**\[cpp\]** [view plain][] [copy][view plain]

1.  function \[ jVal,gradient \] = costFunction2( theta ) 
2.  %COSTFUNCTION2 Summary of **this** function goes here 
3.  % linear regression -> y=theta0 + theta1\*x 
4.  % parameter: x:m\*n theta:n\*1 y:m\*1 (m=4,n=1) 
5.  % 
6.  
7.  %Data 
8.  x=\[1;2;3;4\]; 
9.  y=\[1.1;2.2;2.7;3.8\]; 
10. m=size(x,1); 
11. 
12. hypothesis = h\_func(x,theta); 
13. delta = hypothesis - y; 
14. jVal=sum(delta.^2); 
15. 
16. gradient(1)=sum(delta)/m; 
17. gradient(2)=sum(delta.\*x)/m; 
18. 
19. end

其中，h\_func是hypothesis的结果：

**\[cpp\]** [view plain][] [copy][view plain]

1.  function \[res\] = h\_func(inputx,theta) 
2.  %H\_FUNC Summary of **this** function goes here 
3.  % Detailed explanation goes here 
4.  
5.  
6.  %cost function 2 
7.  res= theta(1)+theta(2)\*inputx;function \[res\] = h\_func(inputx,theta) 
8.  end

Gradient\_descent：

**\[cpp\]** [view plain][] [copy][view plain]

1.  function \[optTheta,functionVal,exitFlag\]=Gradient\_descent( ) 
2.  %GRADIENT\_DESCENT Summary of **this** function goes here 
3.  % Detailed explanation goes here 
4.  
5.   options = optimset('GradObj','on','MaxIter',100); 
6.   initialTheta = zeros(2,1); 
7.   \[optTheta,functionVal,exitFlag\] = fminunc(@costFunction2,initialTheta,options); 
8.  
9.  end

result：

**\[cpp\]** [view plain][] [copy][view plain]

1.  >> \[optTheta,functionVal,exitFlag\] = Gradient\_descent() 
2.  
3.  Local minimum found. 
4.  
5.  Optimization completed because the size of the gradient is less than 
6.  the **default** value of the function tolerance. 
7.  
8.  <stopping criteria details> 
9.  
10. 
11. optTheta = 
12. 
13.  0.3000 
14.  0.8600 
15. 
16. 
17. functionVal = 
18. 
19.  0.0720 
20. 
21. 
22. exitFlag = 
23. 
24.  1

即得y=0.3+0.86x;

验证：

**\[cpp\]** [view plain][] [copy][view plain]

1.  function \[ parameter \] = checkcostfunc( ) 
2.  %CHECKC2 Summary of **this** function goes here 
3.  % check **if** the cost function works well 
4.  % check with the matlab fit function as standard 
5.  
6.  %check cost function 2 
7.  x=\[1;2;3;4\]; 
8.  y=\[1.1;2.2;2.7;3.8\]; 
9.  
10. EXPR= \{ 'x','1'\}; 
11. p=fittype(EXPR); 
12. parameter=fit(x,y,p); 
13. 
14. end

运行结果：

**\[cpp\]** [view plain][] [copy][view plain]

1.  >> checkcostfunc() 
2.  
3.  ans = 
4.  
5.   Linear model: 
6.   ans(x) = a\*x + b 
7.   Coefficients (with 95% confidence bounds): 
8.   a = 0.86 (0.4949, 1.225) 
9.   b = 0.3 (-0.6998, 1.3)

和我们的结果一样。下面画图：

**\[cpp\]** [view plain][] [copy][view plain]

1.  function PlotFunc( xstart,xend ) 
2.  %PLOTFUNC Summary of **this** function goes here 
3.  % draw original data and the fitted 
4.  
5.  
6.  
7.  %===================cost function 2====linear regression 
8.  %original data 
9.  x1=\[1;2;3;4\]; 
10. y1=\[1.1;2.2;2.7;3.8\]; 
11. %plot(x1,y1,'ro-','MarkerSize',10); 
12. plot(x1,y1,'rx','MarkerSize',10); 
13. hold on; 
14. 
15. %fitted line - 拟合曲线 
16. x\_co=xstart:0.1:xend; 
17. y\_co=0.3+0.86\*x\_co; 
18. %plot(x\_co,y\_co,'g'); 
19. plot(x\_co,y\_co); 
20. 
21. hold off; 
22. end

![1341902243_7379.jpg][]

**第三部分：** **解决过拟合问题---Regularization**

过拟合问题解决方法我们已在第三章中讲过，利用Regularization的方法就是在cost function中加入关于 θ的项，使得部分θ的值偏小，从而达到fit效果。

例如定义 costfunction J(θ)： jVal=(theta(1)-5)^2+(theta(2)-5)^2;

在每次迭代中，按照gradient descent的方法更新参数θ：θ(i)-=gradient(i),其中gradient(i)是J(θ)对θi求导的函数式，在此例中就有gradient(1)=2\*(theta(1)-5), gradient(2)=2\*(theta(2)-5)。

函数costFunction, 定义jVal=J(θ)和对两个θ的gradient：

**\[cpp\]** [view plain][] [copy][view plain]

1.  function \[ jVal,gradient \] = costFunction( theta ) 
2.  %COSTFUNCTION Summary of **this** function goes here 
3.  % Detailed explanation goes here 
4.  
5.  jVal= (theta(1)-5)^2+(theta(2)-5)^2; 
6.  
7.  gradient = zeros(2,1); 
8.  %code to compute derivative to theta 
9.  gradient(1) = 2 \* (theta(1)-5); 
10. gradient(2) = 2 \* (theta(2)-5); 
11. 
12. end

Gradient\_descent，进行参数优化

**\[cpp\]** [view plain][] [copy][view plain]

1.  function \[optTheta,functionVal,exitFlag\]=Gradient\_descent( ) 
2.  %GRADIENT\_DESCENT Summary of **this** function goes here 
3.  % Detailed explanation goes here 
4.  
5.   options = optimset('GradObj','on','MaxIter',100); 
6.   initialTheta = zeros(2,1) 
7.   \[optTheta,functionVal,exitFlag\] = fminunc(@costFunction,initialTheta,options); 
8.  
9.  end

matlab主窗口中调用，得到优化厚的参数(θ1,θ2)=(5,5)

**\[cpp\]** [view plain][] [copy][view plain]

1.   \[optTheta,functionVal,exitFlag\] = Gradient\_descent() 
2.  
3.  initialTheta = 
4.  
5.   0 
6.   0 
7.  
8.  
9.  Local minimum found. 
10. 
11. Optimization completed because the size of the gradient is less than 
12. the **default** value of the function tolerance. 
13. 
14. <stopping criteria details> 
15. 
16. 
17. optTheta = 
18. 
19.  5 
20.  5 
21. 
22. 
23. functionVal = 
24. 
25.  0 
26. 
27. 
28. exitFlag = 
29. 
30.  1

**第四部分：Y=1/(1+e^X)型---逻辑回归（sigmod 函数拟合）**

hypothesis function:

**\[cpp\]** [view plain][] [copy][view plain]

1.  function \[res\] = h\_func(inputx,theta) 
2.  
3.  %cost function 3 
4.  tmp=theta(1)+theta(2)\*inputx;%m\*1 
5.  res=1./(1+exp(-tmp));%m\*1 
6.  
7.  end

cost function:

**\[cpp\]** [view plain][] [copy][view plain]

1.  function \[ jVal,gradient \] = costFunction3( theta ) 
2.  %COSTFUNCTION3 Summary of **this** function goes here 
3.  % Logistic Regression 
4.  
5.  x=\[-3; -2; -1; 0; 1; 2; 3\]; 
6.  y=\[0.01; 0.05; 0.3; 0.45; 0.8; 1.1; 0.99\]; 
7.  m=size(x,1); 
8.  
9.  %hypothesis data 
10. hypothesis = h\_func(x,theta); 
11. 
12. %jVal-cost function & gradient updating 
13. jVal=-sum(log(hypothesis+0.01).\*y + (1-y).\*log(1-hypothesis+0.01))/m; 
14. gradient(1)=sum(hypothesis-y)/m; %reflect to theta1 
15. gradient(2)=sum((hypothesis-y).\*x)/m; %reflect to theta 2 
16. 
17. end

Gradient\_descent:

**\[cpp\]** [view plain][] [copy][view plain]

1.  function \[optTheta,functionVal,exitFlag\]=Gradient\_descent( ) 
2.  
3.   options = optimset('GradObj','on','MaxIter',100); 
4.   initialTheta = \[0;0\]; 
5.   \[optTheta,functionVal,exitFlag\] = fminunc(@costFunction3,initialTheta,options); 
6.  
7.  end

运行结果：

**\[cpp\]** [view plain][] [copy][view plain]

1.   \[optTheta,functionVal,exitFlag\] = Gradient\_descent() 
2.  
3.  Local minimum found. 
4.  
5.  Optimization completed because the size of the gradient is less than 
6.  the **default** value of the function tolerance. 
7.  
8.  <stopping criteria details> 
9.  
10. 
11. optTheta = 
12. 
13.  0.3526 
14.  1.7573 
15. 
16. 
17. functionVal = 
18. 
19.  0.2498 
20. 
21. 
22. exitFlag = 
23. 
24.  1

画图验证：

**\[cpp\]** [view plain][] [copy][view plain]

1.  function PlotFunc( xstart,xend ) 
2.  %PLOTFUNC Summary of **this** function goes here 
3.  % draw original data and the fitted 
4.  
5.  %===================cost function 3=====logistic regression 
6.  
7.  %original data 
8.  x=\[-3; -2; -1; 0; 1; 2; 3\]; 
9.  y=\[0.01; 0.05; 0.3; 0.45; 0.8; 1.1; 0.99\]; 
10. plot(x,y,'rx','MarkerSize',10); 
11. hold on 
12. 
13. %fitted line 
14. x\_co=xstart:0.1:xend; 
15. theta = \[0.3526,1.7573\]; 
16. y\_co=h\_func(x\_co,theta); 
17. plot(x\_co,y\_co); 
18. hold off 
19. 
20. end

![1341902876_4751.jpg][]

有朋友问，这里就补充一下logistic regression中gradient的推导：

令

[![z = \\frac\{1\}\{1+e^\{-\\theta x\}\}][z _ _frac_1_1_e_-_theta x]][z _ _frac_1_1_e_-_theta x_z _ _frac_1_1_e_-_theta x]

则有

[![z'\_\{\\theta\}=\\frac\{e^\{-\\theta x\}\}\{(1+e^\{-\\theta x\})^2\} \\cdot (-x) = z(z-1)(-x)\\\\][z_theta_frac_e_-_theta x_1_e_-_theta x_2_ _cdot _-x_ _ z_z-1_-x]][z_theta_frac_e_-_theta x_1_e_-_theta x_2_ _cdot _-x_ _ z_z-1_-x_z_theta_frac_e_-_theta x_1_e_-_theta x_2_ _cdot _-x_ _ z_z-1_-x]

由于cost function

[![J = ylnz+(1-y)ln(1-z)][J _ ylnz_1-y_ln_1-z]][J _ ylnz_1-y_ln_1-z_J _ ylnz_1-y_ln_1-z]

可得

[![J'\_\{\\theta\} = y\\frac\{1\}\{z\}z'\_\{\\theta\}+(1-y)\\frac\{-z'\_\\theta\}\{1-z\}\\\\ J'\_\{\\theta\}=z'\_\\theta(\\frac\{y\}\{z\}-\\frac\{1-y\}\{1-z\}) = z(z-1)(-x)\\frac\{y-yz-z+yz\}\{z(1-z)\} = (y-z)x][J_theta_ _ y_frac_1_z_z_theta_1-y_frac_-z_theta_1-z_ J_theta_z_theta_frac_y_z_-_frac_1-y_1-z_ _ z_z-1_-x_frac_y-yz-z_yz_z_1-z_ _ _y-z_x]][J_theta_ _ y_frac_1_z_z_theta_1-y_frac_-z_theta_1-z_ J_theta_z_theta_frac_y_z_-_frac_1-y_1-z_ _ z_z-1_-x_frac_y-yz-z_yz_z_1-z_ _ _y-z_x_J_theta_ _ y_frac_1_z_z_theta_1-y_frac_-z_theta_1-z_ J_theta_z_theta_frac_y_z_-_frac_1-y_1-z_ _ z_z-1_-x_frac_y-yz-z_yz_z_1-z_ _ _y-z_x]

所以gradient = -J'(theta) = (z-y)x

关于Machine Learning更多的学习资料将继续更新，敬请关注本博客和新浪微博[Sophia\_qing][Sophia_qing]。

**==============================**

Reference:

1. [http://blog.csdn.net/abcjennifer/article/details/7691571][Link 1]

2.[http://blog.csdn.net/abcjennifer/article/details/7700772][Link 2]

3.[http://blog.csdn.net/abcjennifer/article/details/7716281][Link 3]

4.[http://blog.csdn.net/abcjennifer/article/details/7684836][Matlab _ _]

[Link 1]: http://blog.csdn.net/abcjennifer/article/details/7691571
[Link 2]: http://blog.csdn.net/abcjennifer/article/details/7700772
[Link 3]: http://blog.csdn.net/abcjennifer/article/details/7716281
[1341472718_8627.jpg]: https://img-my.csdn.net/uploads/201207/05/1341472718_8627.jpg
[1341901237_1654.jpg]: https://img-my.csdn.net/uploads/201207/10/1341901237_1654.jpg
[1341658176_1292.jpg]: https://img-my.csdn.net/uploads/201207/07/1341658176_1292.jpg
[Matlab _ _]: http://blog.csdn.net/abcjennifer/article/details/7684836
[view plain]: http://blog.csdn.net/abcjennifer/article/details/7732417#
[1341902243_7379.jpg]: https://img-my.csdn.net/uploads/201207/10/1341902243_7379.jpg
[1341902876_4751.jpg]: https://img-my.csdn.net/uploads/201207/10/1341902876_4751.jpg
[z _ _frac_1_1_e_-_theta x]: http://latex.codecogs.com/gif.latex?z%20=%20\frac%7b1%7d%7b1+e^%7b-\theta%20x%7d%7d
[z _ _frac_1_1_e_-_theta x_z _ _frac_1_1_e_-_theta x]: http://private.codecogs.com/eqnedit.php?latex=z%20=%20\frac%7b1%7d%7b1@plus;e^%7b-\theta%20x%7d%7d
[z_theta_frac_e_-_theta x_1_e_-_theta x_2_ _cdot _-x_ _ z_z-1_-x]: http://latex.codecogs.com/gif.latex?z%27_%7b\theta%7d=\frac%7be^%7b-\theta%20x%7d%7d%7b%281+e^%7b-\theta%20x%7d%29^2%7d%20\cdot%20%28-x%29%20=%20z%28z-1%29%28-x%29\\
[z_theta_frac_e_-_theta x_1_e_-_theta x_2_ _cdot _-x_ _ z_z-1_-x_z_theta_frac_e_-_theta x_1_e_-_theta x_2_ _cdot _-x_ _ z_z-1_-x]: http://private.codecogs.com/eqnedit.php?latex=z%27_%7b\theta%7d=\frac%7be^%7b-\theta%20x%7d%7d%7b%281@plus;e^%7b-\theta%20x%7d%29^2%7d%20\cdot%20%28-x%29%20=%20z%28z-1%29%28-x%29\\
[J _ ylnz_1-y_ln_1-z]: http://latex.codecogs.com/gif.latex?J%20=%20ylnz+%281-y%29ln%281-z%29\\
[J _ ylnz_1-y_ln_1-z_J _ ylnz_1-y_ln_1-z]: /images/20220805/4421b81ede6140eb9ee9cafb7dd6d401.png\
[J_theta_ _ y_frac_1_z_z_theta_1-y_frac_-z_theta_1-z_ J_theta_z_theta_frac_y_z_-_frac_1-y_1-z_ _ z_z-1_-x_frac_y-yz-z_yz_z_1-z_ _ _y-z_x]: http://latex.codecogs.com/gif.latex?J%27_%7b\theta%7d%20=%20y\frac%7b1%7d%7bz%7dz%27_%7b\theta%7d+%281-y%29\frac%7b-z%27_\theta%7d%7b1-z%7d\\%20J%27_%7b\theta%7d=z%27_\theta%28\frac%7by%7d%7bz%7d-\frac%7b1-y%7d%7b1-z%7d%29%20=%20z%28z-1%29%28-x%29\frac%7by-yz-z+yz%7d%7bz%281-z%29%7d%20=%20%28y-z%29x
[J_theta_ _ y_frac_1_z_z_theta_1-y_frac_-z_theta_1-z_ J_theta_z_theta_frac_y_z_-_frac_1-y_1-z_ _ z_z-1_-x_frac_y-yz-z_yz_z_1-z_ _ _y-z_x_J_theta_ _ y_frac_1_z_z_theta_1-y_frac_-z_theta_1-z_ J_theta_z_theta_frac_y_z_-_frac_1-y_1-z_ _ z_z-1_-x_frac_y-yz-z_yz_z_1-z_ _ _y-z_x]: http://private.codecogs.com/eqnedit.php?latex=J%27_%7b\theta%7d%20=%20y\frac%7b1%7d%7bz%7dz%27_%7b\theta%7d@plus;%281-y%29\frac%7b-z%27_\theta%7d%7b1-z%7d\\%20J%27_%7b\theta%7d=z%27_\theta%28\frac%7by%7d%7bz%7d-\frac%7b1-y%7d%7b1-z%7d%29%20=%20z%28z-1%29%28-x%29\frac%7by-yz-z@plus;yz%7d%7bz%281-z%29%7d%20=%20%28y-z%29x
[Sophia_qing]: http://weibo.com/u/2607574543