Mathematica中文论坛-非官方

标题: 最小二乘法和主成分分析的比较 matlab  儿子的papa [打印本页]

作者: meatball1982    时间: 2016-5-4 21:36
标题: 最小二乘法和主成分分析的比较 matlab  儿子的papa
比较ls和pca,实际上,当点比较少的时候,最小二乘法是个不错的选择,但要观察数据的大体方向,pca是个更好的选择。下面说明。

生成在一定范围内的矩形数据,下图中的蓝色的点(我暂时只生成5000个点。),将这些点进行旋转(这时是-30度),得到图中的红色点。分别用ls和pca分析得到红色点的主要方向。

  1. clear all
  2. clc
  3. clf
  4. %% compare lsfit and  pca
  5. % outline:
  6. % 01 gen random points in rectangle
  7. % 02 rotate and shift rectangle
  8. % 03 ls fit
  9. % 04 pca
  10. % 05 plot result

  11. %% main
  12. % parameters ---------------------------------------------
  13. rng default
  14. n=5000;          % points number
  15. the=-pi/6;       % rotate angle

  16. %% 01 gen random points in rectangle ----------------------
  17. x=rand(n,1)*4-2;
  18. y=rand(n,1)*3-1.5;

  19. %% 02 rotate and shift rectangle --------------------------
  20. x_new= cos(the)*x+sin(the)*y+3;
  21. y_new=-sin(the)*x+cos(the)*y+5;

  22. %% 03 ls fit ----------------------------------------------
  23. p=polyfit(x_new,y_new,1);
  24. y_p=p(1)*x_new+p(2);

  25. %% 04 pca -------------------------------------------------
  26. dat=[x_new';y_new']';
  27. [coeff,score,latent] = princomp(dat); % new version of matlab, use pca

  28. % rotate first PC back
  29. eye_mat=eye(2);
  30. eye_new=eye_mat*(coeff)^(-1);

  31. % gen pca line
  32. x_cen=mean(x_new);
  33. y_cen=mean(y_new);
  34. y_p_pca=eye_new(1,2)/eye_new(1,1)*(x_new-x_cen)+y_cen;

  35. % first PC angle
  36. atan(eye_new(1,2)/eye_new(1,1))*180/pi

  37. % second PC angle
  38. atan(eye_new(2,2)/eye_new(2,1))*180/pi

  39. %% 05 plot result------------------------------------------
  40. hold on
  41. plot(x,y,'b.')                          % plot ori points
  42. plot(x_new,y_new,'r.')                  % plot rotate points
  43. plot(x_new,y_p,'k-','linewidth',2)      % ls fit line
  44. plot(x_new,y_p_pca,'g-','linewidth',2)  % pca line(first PC)

  45. grid on
  46. axis equal
  47. legend('ori point','rot point','ls fit','PCA','Location','southeast')

  48. % option: output fig
  49. % h=gcf;
  50. % fig_na=['../imgs/fig_01_comp_lsfit_pca'];
  51. % fun_work_li_035_myfig_out(h,fig_na,3);

  52. %% logs
  53. % mod :
复制代码


fig_01_comp_lsfit_pca.jpg (234.86 KB, 下载次数: 930)

fig_01_comp_lsfit_pca.jpg





欢迎光临 Mathematica中文论坛-非官方 (http://ilovemathematica.com/) Powered by Discuz! X3.2