随机森林法

Halcom

随机森林法： 区别：本视频是算法底层代码原理实现，本帖内容为工具箱下的随机森林法实现。

百度网盘视频链接： 视频链接： 视频1： http://pan.baidu.com/s/1kVBPUsB 视频2：http://pan.baidu.com/s/1hszKNjA 录制的视频是算法底层原理讲解，底层代码实现，方便大家真正掌握算法实质，开发出更加出色的算法。录制视频的初衷是：避免读者朋友利用大把时间学习已有常见算法，本系列视频旨在让读者朋友快速深入了解这些常见算法原理，有更多的时间去研究更加高大上算法（价值）。 具体链接在halcom.cn论坛，联系人QQ：3283892722 该论坛是一个学习交流平台，我会逐一的和大家分享学习。 欢迎大家录制视频，并提交给我，我来设置视频，你可在论坛进行打赏分享。 视频专用播放器：http://halcom.cn/forum.php?mod=viewthread&tid=258&extra=page%3D1

使用环境：Win7-32bit-matlab2014a
主程序（MATLAB工具箱使用）：

1. load fisheriris
   
2. X = meas;
   
3. Y = species;
   
4. BaggedEnsemble = generic_random_forests(X,Y,60,'classification')
   
5. predict(BaggedEnsemble,[5 3 5 1.8])

复制代码

子函数程序如下：

1. BaggedEnsemble = TreeBagger(iNumBags,X,Y,'OOBPred','On','Method',str_method)
   
2. 
   
3. % plot out of bag prediction error
   
4. oobErrorBaggedEnsemble = oobError(BaggedEnsemble);
   
5. figID = figure;
   
6. plot(oobErrorBaggedEnsemble)
   
7. xlabel 'Number of grown trees';
   
8. ylabel 'Out-of-bag classification error';
   
9. print(figID, '-dpdf', sprintf('randomforest_errorplot_%s.pdf', date));
   
10. 
    
11. oobPredict(BaggedEnsemble)
    
12. 
    
13. % view trees
    
14. view(BaggedEnsemble.Trees{1}) % text description
    
15. view(BaggedEnsemble.Trees{1},'mode','graph') % graphic description
    

复制代码

参考：https://github.com/kjw0612/awesome-random-forest
http://database.51cto.com/art/201407/444788.htm

Halcom

1. % random forest随机森林 RF
   
2. clc,clear,close all
   
3. warning off
   
4. feature jit off
   
5. data = xlsread('数据3.xlsx',3,'B3:P137'); % 加载数据
   
6. rng(1945,'twister') % For reproducibility
   
7. % OOB(Out-of-Bag)样本
   
8. x = data(:,2:end);    % 各个指标
   
9. y = data(:,1);        % 输出，按照年份聚类
   
10. ntree = 300;          % CART决策树的数量
    
11. BaggedEnsemble = TreeBagger(ntree,x,y,'Method','C','oobvarimp','on','OOBPred','On','NVarToSample',2,'MinLeaf',5 );
    
12. %% 原始样本分布情况
    
13. % RF来计算样本的相似度矩阵(Proximity Matrix)
    
14. Prox = fillProximities(BaggedEnsemble);
    
15. PM = Prox.Proximity;
    
16. figure(1),
    
17. [mdsPM,eigvalue] = mdsProx(Prox,'Colors','rgbcmshpkx*.o><');
    
18. xlabel('Dim1');ylabel('Dim2');title('相似度矩阵投影图')
    
19. %% Outliers
    
20. % 异常样本也称作离群数据(Outliers)，通常定义为远离数据主体的样本。
    
21. OM = Prox.OutlierMeasure;
    
22. figure(2),bar(OM);
    
23. xlabel('样本数');ylabel('Outliers');
    
24. % 剔除异常样本
    
25. abnormal = find(OM>3);
    
26. data(abnormal,:) = [];  % 删除异常样本
    
27. x = data(:,2:end);      % 更新各个指标
    
28. y = data(:,1);          % 更新输出，按照年份聚类
    
29. ntree = 300;          % CART决策树的数量
    
30. BaggedEnsemble = TreeBagger(ntree,x,y,'Method','C','oobvarimp','on','OOBPred','On','NVarToSample',2,'MinLeaf',5 );
    
31. Prox=[];Prox = fillProximities(BaggedEnsemble);
    
32. PM=[];PM = Prox.Proximity;
    
33. figure(111),
    
34. [mdsPM,eigvalue] = mdsProx(Prox,'Colors','rgbcmshpkx*.o><');
    
35. xlabel('Dim1');ylabel('Dim2');title('去除异常样本的相似度矩阵投影图')
    
36. xx = x;
    
37. yy = y;
    
38. %% OOB估计误差
    
39. ntree = 300;          % CART决策树的数量
    
40. BaggedEnsemble=[];
    
41. BaggedEnsemble = TreeBagger(ntree,x,y,'Method','C','oobvarimp','on','OOBPred','On','NVarToSample',1,'MinLeaf',5 );
    
42. oobErrorBaggedEnsemble=[];oobErrorBaggedEnsemble = oobError(BaggedEnsemble);
    
43. figure(3),plot(oobErrorBaggedEnsemble)
    
44. xlabel(['决策树数量,ntree=',num2str(ntree)]);
    
45. ylabel('OOB分类误差');
    
46. 
    
47. ntree = 500;          % CART决策树的数量
    
48. BaggedEnsemble=[];
    
49. BaggedEnsemble = TreeBagger(ntree,x,y,'Method','C','oobvarimp','on','OOBPred','On','NVarToSample',1,'MinLeaf',5 );
    
50. oobErrorBaggedEnsemble=[];oobErrorBaggedEnsemble = oobError(BaggedEnsemble);
    
51. figure(4),plot(oobErrorBaggedEnsemble)
    
52. xlabel(['决策树数量,ntree=',num2str(ntree)]);
    
53. ylabel('OOB分类误差');
    
54. 
    
55. ntree = 1000;          % CART决策树的数量
    
56. BaggedEnsemble=[];
    
57. BaggedEnsemble = TreeBagger(ntree,x,y,'Method','C','oobvarimp','on','OOBPred','On','NVarToSample',1,'MinLeaf',5 );
    
58. oobErrorBaggedEnsemble=[];oobErrorBaggedEnsemble = oobError(BaggedEnsemble);
    
59. figure(5),plot(oobErrorBaggedEnsemble)
    
60. xlabel(['决策树数量,ntree=',num2str(ntree)]);
    
61. ylabel('OOB分类误差');
    
62. %% 确定mtry
    
63. for mtry=1:15
    
64.     ntree = 300;          % CART决策树的数量
    
65.     BaggedEnsemble=[];
    
66.     BaggedEnsemble = TreeBagger(ntree,x,y,'Method','C','OOBPred','On','NVarToSample',mtry,'MinLeaf',5 );
    
67.     oobErrorBaggedEnsemble = oobError(BaggedEnsemble);
    
68.     ROC(mtry) = mean(oobErrorBaggedEnsemble);
    
69. end
    
70. figure(6),plot(ROC,'o-');
    
71. title(['best mtry = ',num2str( find(ROC==max(ROC)) )])
    
72. xlabel('mtry'); ylabel('ROC(Repeated Cross-Validation)');
    
73. 
    
74. %% 变量的重要性求解
    
75. ntree = 300;          % CART决策树的数量
    
76. BaggedEnsemble=[];
    
77. BaggedEnsemble = TreeBagger(ntree,x,y,'Method','C','oobvarimp','on','OOBPred','On','NVarToSample',1,'MinLeaf',5 );
    
78. Prox=[];Prox = fillProximities(BaggedEnsemble);  % RF来计算样本的相似度矩阵(Proximity Matrix)
    
79. MDA = Prox.OOBPermutedVarDeltaError;        
    
80. MDG = Prox.OOBPermutedVarDeltaMeanMargin;  
    
81. % MDA、MDG画图
    
82. [wu,sortMDA] = sort(MDA,'ascend');  % 升序排列
    
83. [wu,sortMDG] = sort(MDG,'ascend');  % 升序排列
    
84. index = {'应急响应能力X21','耕地灌溉率X2','饮水保障X3','抗旱经费支出X10',...
    
85.     '水库蓄水能力X20','作物种植比例X4','农业总产值比例X18','人均粮食产量X16',...
    
86.     '旱地比重X5','人均耕地面积X1','农用设施配备比例X19','外出务工比例X12',...
    
87.     '土壤墒情X8','年平均降水量X7'}';
    
88. figure(7),
    
89. axes1 = axes('Position',[0.3 0.15 0.6 0.7]);
    
90. axis(axes1);
    
91. for i=1:length(sortMDA)
    
92.     plot(i,MDA( sortMDA(i) ),'ro');
    
93.     hold on
    
94.     text(-6.2,min(MDA)+(i-1)*(max(MDA)-min(MDA))/length(sortMDA),index( sortMDA(i) ));
    
95. end
    
96. set(gca,'ytick',[])        %  同时去掉 x 轴和 y 轴的刻度
    
97. set(gca, 'XGrid', 'on');   %  网格
    
98. set(gca, 'YGrid', 'on');   %  网格
    
99. xlabel('MDA')
    
100. hold off
     
101. 
     
102. figure(8),
     
103. axes1 = axes('Position',[0.3 0.15 0.6 0.7]);
     
104. axis(axes1);
     
105. for i=1:length(sortMDG)
     
106.     plot(i,MDG( sortMDG(i) ),'ro');
     
107.     hold on
     
108.     text(-6.2,min(MDG)+(i-1)*(max(MDG)-min(MDG))/length(sortMDG),index( sortMDG(i) ));
     
109. end
     
110. set(gca,'ytick',[])        %  同时去掉 x 轴和 y 轴的刻度
     
111. set(gca, 'XGrid', 'on');   %  网格
     
112. set(gca, 'YGrid', 'on');   %  网格
     
113. xlabel('MDG')
     
114. hold off
     
115. %% 按照MDG进行删除指标计算ROC
     
116. ntree = 300;          % CART决策树的数量
     
117. mtry = 1;             % CART决策树生长过程中分枝时随机选择的变量个数
     
118. indexP = ones(length(sortMDG),1);
     
119. for i=1:13             % 删除前8个不重要指标
     
120.     x=xx;
     
121.     x(:, sortMDG(1:i) ) = [];
     
122.     BaggedEnsemble=[];
     
123.     BaggedEnsemble = TreeBagger(ntree,x,y,'Method','C','OOBPred','On','NVarToSample',mtry,'MinLeaf',5 );
     
124.     oobErrorBaggedEnsemble=[];
     
125.     oobErrorBaggedEnsemble = oobError(BaggedEnsemble);
     
126.     AUC(i) = mean(oobErrorBaggedEnsemble);
     
127.     indexP(i) = i/13;
     
128. end
     
129. figure(9),plot(AUC,'r-.')
     
130. xlabel('删除不重要指标的个数');ylabel('OOB-AUC');
     
131. 
     
132. %% 重复交叉验证过程指标选择概率排序
     
133. figure(10),
     
134. axes2 = axes('Position',[0.3 0.15 0.6 0.7]);
     
135. axis(axes2);
     
136. for i=1:length(sortMDG)
     
137.     plot(indexP(i),i,'ro');
     
138.     hold on
     
139.     text(-0.43,i,index( sortMDG(i) ));
     
140. end
     
141. set(gca,'ytick',[])        %  同时去掉 x 轴和 y 轴的刻度
     
142. set(gca, 'XGrid', 'on');   %  网格
     
143. set(gca, 'YGrid', 'on');   %  网格
     
144. xlabel('指标被选中的概率')
     
145. hold off
     

复制代码