本文已参加「新人创造礼」活动，一起敞开创造之路。

  本文介绍根据Python的随机森林（Random Forest，RF）回归代码，以及模型超参数（包含决策树个数与最大深度、最小别离样本数、最小叶子节点样本数、最大别离特征数等）主动优化的代码。

  本文是在上一篇Python完结随机森林RF回归、自变量重要性排序的基础上完结的，因而本次仅对随机森林模型超参数主动择优部分的代码加以详细解说；而数据预备、模型树立、精度鉴定等其他部分的代码详细解说，咱们直接点击上述文章Python完结随机森林RF回归、自变量重要性排序检查即可。

  其中，关于根据MATLAB完结相同进程的代码与实战，咱们能够点击检查MATLAB完结随机森林回归预测并对变量的影响程度进行排序。

  本文分为两部分，第一部分为代码的分段讲解，第二部分为完好代码。

1 代码分段讲解

1.1 数据与模型预备

  本部分是对随机森林算法的数据与模型预备，因为在之前的博客中现已详细介绍过了，本文就不再赘述~咱们直接检查Python完结随机森林RF回归、自变量重要性排序即可。

import pydot
import numpy as np
import pandas as pd
import scipy.stats as stats
import matplotlib.pyplot as plt
from pprint import pprint
from sklearn import metrics
from openpyxl import load_workbook
from sklearn.tree import export_graphviz
from sklearn.ensemble import RandomForestRegressor
from sklearn.model_selection import GridSearchCV
from sklearn.model_selection import RandomizedSearchCV
# Attention! Data Partition
# Attention! One-Hot Encoding
train_data_path='G:/CropYield/03_DL/00_Data/AllDataAll_Train.csv'
test_data_path='G:/CropYield/03_DL/00_Data/AllDataAll_Test.csv'
write_excel_path='G:/CropYield/03_DL/05_NewML/ParameterResult_ML.xlsx'
tree_graph_dot_path='G:/CropYield/03_DL/05_NewML/tree.dot'
tree_graph_png_path='G:/CropYield/03_DL/05_NewML/tree.png'
random_seed=44
random_forest_seed=np.random.randint(low=1,high=230)
# Data import
train_data=pd.read_csv(train_data_path,header=0)
test_data=pd.read_csv(test_data_path,header=0)
# Separate independent and dependent variables
train_Y=np.array(train_data['Yield'])
train_X=train_data.drop(['ID','Yield'],axis=1)
train_X_column_name=list(train_X.columns)
train_X=np.array(train_X)
test_Y=np.array(test_data['Yield'])
test_X=test_data.drop(['ID','Yield'],axis=1)
test_X=np.array(test_X)

1.2 超参数规模给定

  首要，咱们需求对随机森林模型超参数各自的规模加以确认，之后咱们将在这些规模内确认各个超参数的终究最优取值。换句话说，咱们现在先给每一个需求择优的超参数划定一个很大很大的规模（例如对于“决策树个数”这个超参数，咱们能够将其规模划定在10到5000这样一个很大的规模），然后后期将用择优算法在每一个超参数的这个规模内进行查找。

  在此，咱们先要确认对哪些超参数进行择优。本文挑选在随机森林算法中比较重要的几个超参数进行调优，别离是：决策树个数n_estimators，决策树最大深度max_depth，最小别离样本数（即拆分决策树节点所需的最小样本数）min_samples_split，最小叶子节点样本数（即一个叶节点所需包含的最小样本数）min_samples_leaf，最大别离特征数（即寻找最佳节点分割时要考虑的特征变量数量）max_features，以及是否进行随机抽样bootstrap等六种。关于上述超参数假如咱们不是太了解详细的意义，能够检查Python完结随机森林RF回归、自变量重要性排序的1.5部分，或许就会比较好理解了（不过其实不理解也不影响接下来的操作）。

  这儿提一句，其实随机森林的超参数并不止上述这些，我这儿也是结合数据状况与终究的精度需求，挑选了相比照较常用的几个超参数；咱们根据各自实际需求，挑选需求调整的超参数，并用相同的代码思路履行即可。

# Search optimal hyperparameter
n_estimators_range=[int(x) for x in np.linspace(start=50,stop=3000,num=60)]
max_features_range=['auto','sqrt']
max_depth_range=[int(x) for x in np.linspace(10,500,num=50)]
max_depth_range.append(None)
min_samples_split_range=[2,5,10]
min_samples_leaf_range=[1,2,4,8]
bootstrap_range=[True,False]
random_forest_hp_range={'n_estimators':n_estimators_range,
                        'max_features':max_features_range,
                        'max_depth':max_depth_range,
                        'min_samples_split':min_samples_split_range,
                        'min_samples_leaf':min_samples_leaf_range
                        # 'bootstrap':bootstrap_range
                        }
pprint(random_forest_hp_range)

  能够看到，上述代码首要是对六种超参数划定了一个规模，然后将其别离存入了一个超参数规模字典random_forest_hp_range。在这儿咱们能够看到，我在存入字典时，将bootstrap的规模这一句注释掉了，这是因为其时运转后我发现bootstrap还是处于True这个状态比较好（也便是得到的成果精度比较高），因而就取消了这一超参数的择优；咱们根据个人数据与模型的实际状况来即可~

  咱们能够看一下random_forest_hp_range变量的取值状况：

  没错，它是一个字典，键便是超参数的称号，值便是超参数的规模。因为我将bootstrap注释掉了，因而这个字典里就没有bootstrap这一项了~

1.3 超参数随机匹配择优

  上面咱们确认了每一种超参数各自的规模，那么接下来咱们就将他们别离组合，比照每一个超参数取值组合所得到的模型成果，然后确认最优超参数组合。

  其实咱们会发现，咱们上面划定六种超参数（除掉我后来删去的bootstrap的话是五种），假如依照排列组合来核算的话，会有很多很多种组合方法，假如要一一测验未免也太麻烦了。因而，咱们用到RandomizedSearchCV这一功能——其将随机匹配每一种超参数组合，并输出最优的组合。换句话说，咱们用RandomizedSearchCV来进行随机的排列，而不是对所有的超参数排列组合方法进行遍历。这姿态的确能够节约很多时刻。

random_forest_model_test_base=RandomForestRegressor()
random_forest_model_test_random=RandomizedSearchCV(estimator=random_forest_model_test_base,
                                                   param_distributions=random_forest_hp_range,
                                                   n_iter=200,
                                                   n_jobs=-1,
                                                   cv=3,
                                                   verbose=1,
                                                   random_state=random_forest_seed
                                                   )
random_forest_model_test_random.fit(train_X,train_Y)
best_hp_now=random_forest_model_test_random.best_params_
pprint(best_hp_now)

  由代码能够看到，咱们首要树立一个随机森林模型random_forest_model_test_base，并将其带入到RandomizedSearchCV中；其中，RandomizedSearchCV的参数组合便是刚刚咱们看的random_forest_hp_range，n_iter便是详细随机调配超参数组合的次数（这个次数因而肯定是越大包括的组合数越多，作用越好，可是也越费时刻），cv是穿插验证的折数（RandomizedSearchCV衡量每一种组合方法的作用便是用穿插验证来进行的），n_jobs与verbose是关于模型线程、日志相关的信息，咱们不用太介意，random_state是随机森林中随机抽样的随机数种子。

  之后，咱们对random_forest_model_test_random加以训练，并获取其所得到的最优超参数匹配组合best_hp_now。在这儿，模型的训练次数便是n_iter与cv的乘积（因为穿插验证有几折，那么就需求运转几次；而一共有n_iter个参数匹配组合，因而总次数便是二者相乘）。例如，用上述代码那么一共就需求运转600次。运转进程在程序中将主动显现，如下图。

  能够看到，一共有600次fit，我这儿共花了11.7min完结。详细速度和电脑配置、自变量与因变量数据量大小，以及电脑此时内存等等都有关。

  运转结束，咱们来看看找到的最有超参数组合best_hp_now。

  能够看到，通过200种组合匹配方法的核算，现在五种超参数最优的组合调配方法现已得到了。其实每一次得到的超参数最优组合成果距离也是蛮大的——例如同一批数据，有的时分我得到的n_estimators最优值是如图所示的100，有的时分也会是2350；所以咱们根据实际状况来判别即可~

  那么接下来，咱们就持续对这一best_hp_now所示的成果进行更进一步的择优。

1.4 超参数遍历匹配择优

  刚刚咱们根据RandomizedSearchCV，完结了200次的超参数随机匹配与择优；可是此时的成果是一个随机不完全遍历后所得的成果，因而其最优组合或许并不是全局最优的，而仅仅一个大概的最优规模。因而接下来，咱们需求根据上述所得到的随机最优匹配成果，进行遍历全部组合的匹配择优。

  遍历匹配即在随机匹配最优成果的基础上，在其接近规模内选取几个数值，并通过GridSearchCV对每一种匹配都遍历，然后选出比较好的超参数终究取值成果。

# Grid Search
random_forest_hp_range_2={'n_estimators':[60,100,200],
                          'max_features':[12,13],
                          'max_depth':[350,400,450],
                          'min_samples_split':[2,3] # Greater than 1
                          # 'min_samples_leaf':[1,2]
                          # 'bootstrap':bootstrap_range
                          }
random_forest_model_test_2_base=RandomForestRegressor()
random_forest_model_test_2_random=GridSearchCV(estimator=random_forest_model_test_2_base,
                                               param_grid=random_forest_hp_range_2,
                                               cv=3,
                                               verbose=1,
                                               n_jobs=-1)
random_forest_model_test_2_random.fit(train_X,train_Y)
best_hp_now_2=random_forest_model_test_2_random.best_params_
pprint(best_hp_now_2)

  咱们能够看到，本部分代码其实和1.3部分比较相似。咱们侧重讲解random_forest_hp_range_2。其中，n_estimators设定为了[60,100,200]，这是因为咱们刚刚得到的best_hp_now中n_estimators为100，那么咱们就在100附近选取几个值，作为新的n_estimators规模；max_features也是相同的，因为best_hp_now中max_features为'sqrt'，也便是输入数据特征（自变量）的个数的平方根，而我这儿自变量个数大概是150多个，因而其开平方之后便是12.24左右，那么就挑选其附近的两个数（需求为整数），因而就挑选了[12,13]。其他的超参数取值也是相似的。这儿我将'min_samples_leaf'也给注释掉了是因为我跑了很多次发现，'min_samples_leaf'还是取1最好，那么就直接挑选为默许1（'min_samples_leaf'在不指定的状况下默许为1）即可，因为超参数规模越小，程序跑的就越快。

  这儿程序运转的次数便是每一种超参数取值个数的排列组合次数乘以穿插验证的折数，也便是(2*3*2*3)*3=108次，咱们来看看是不是108次：

  很明显，没有问题，便是108个fit。和前面的600次fit比起来，这样就快很多了（这也是为什么我直接将'min_samples_leaf'与'bootstrap'注释掉的原因；要是这两个超参数也参加的话，那么假设他们两个各有2个取值的话，总时刻至少就要翻2*2=4倍）。

  再来看看通过遍历择优后的最优超参数匹配取值best_hp_now_2。

  以上便是咱们通过一次随机择优、一次遍历择优之后的超参数成果（不要忘记了'min_samples_leaf'与'bootstrap'还要别离取1和True，也便是默许值）。假如咱们感觉这个组合调配还不是很好，那么能够持续履行本文“1.4 超参数遍历匹配择优”部分1到2次，精度或许会有更进一步的提高。

1.5 模型运转与精度鉴定

  结束了上述超参数择优进程，咱们就能够进行模型运转、精度鉴定与成果输出等操作。本部分内容除了第一句代码（将最优超参数组合分配给模型）之外，其余部分因为在之前的博客中现已详细介绍过了，本文就不再赘述~咱们直接检查Python完结随机森林RF回归、自变量重要性排序即可。

# Build RF regression model with optimal hyperparameters
random_forest_model_final=random_forest_model_test_2_random.best_estimator_
# Predict test set data
random_forest_predict=random_forest_model_test_2_random.predict(test_X)
random_forest_error=random_forest_predict-test_Y
# Draw test plot
plt.figure(1)
plt.clf()
ax=plt.axes(aspect='equal')
plt.scatter(test_Y,random_forest_predict)
plt.xlabel('True Values')
plt.ylabel('Predictions')
Lims=[0,10000]
plt.xlim(Lims)
plt.ylim(Lims)
plt.plot(Lims,Lims)
plt.grid(False)
plt.figure(2)
plt.clf()
plt.hist(random_forest_error,bins=30)
plt.xlabel('Prediction Error')
plt.ylabel('Count')
plt.grid(False)
# Verify the accuracy
random_forest_pearson_r=stats.pearsonr(test_Y,random_forest_predict)
random_forest_R2=metrics.r2_score(test_Y,random_forest_predict)
random_forest_RMSE=metrics.mean_squared_error(test_Y,random_forest_predict)**0.5
print('Pearson correlation coefficient is {0}, and RMSE is {1}.'.format(random_forest_pearson_r[0],
                                                                        random_forest_RMSE))
# Save key parameters
excel_file=load_workbook(write_excel_path)
excel_all_sheet=excel_file.sheetnames
excel_write_sheet=excel_file[excel_all_sheet[0]]
excel_write_sheet=excel_file.active
max_row=excel_write_sheet.max_row
excel_write_content=[random_forest_pearson_r[0],random_forest_R2,random_forest_RMSE,
                     random_seed,random_forest_seed]
for i in range(len(excel_write_content)):
        exec("excel_write_sheet.cell(max_row+1,i+1).value=excel_write_content[i]")
excel_file.save(write_excel_path)
# Draw decision tree visualizing plot
random_forest_tree=random_forest_model_final.estimators_[5]
export_graphviz(random_forest_tree,out_file=tree_graph_dot_path,
                feature_names=train_X_column_name,rounded=True,precision=1)
(random_forest_graph,)=pydot.graph_from_dot_file(tree_graph_dot_path)
random_forest_graph.write_png(tree_graph_png_path)
# Calculate the importance of variables
random_forest_importance=list(random_forest_model_final.feature_importances_)
random_forest_feature_importance=[(feature,round(importance,8)) 
                                  for feature, importance in zip(train_X_column_name,
                                                                 random_forest_importance)]
random_forest_feature_importance=sorted(random_forest_feature_importance,key=lambda x:x[1],reverse=True)
plt.figure(3)
plt.clf()
importance_plot_x_values=list(range(len(random_forest_importance)))
plt.bar(importance_plot_x_values,random_forest_importance,orientation='vertical')
plt.xticks(importance_plot_x_values,train_X_column_name,rotation='vertical')
plt.xlabel('Variable')
plt.ylabel('Importance')
plt.title('Variable Importances')

2 完好代码

  本文所用完好代码如下。

# -*- coding: utf-8 -*-
"""
Created on Sun Mar 21 22:05:37 2021
@author: fkxxgis
"""
import pydot
import numpy as np
import pandas as pd
import scipy.stats as stats
import matplotlib.pyplot as plt
from pprint import pprint
from sklearn import metrics
from openpyxl import load_workbook
from sklearn.tree import export_graphviz
from sklearn.ensemble import RandomForestRegressor
from sklearn.model_selection import GridSearchCV
from sklearn.model_selection import RandomizedSearchCV
# Attention! Data Partition
# Attention! One-Hot Encoding
train_data_path='G:/CropYield/03_DL/00_Data/AllDataAll_Train.csv'
test_data_path='G:/CropYield/03_DL/00_Data/AllDataAll_Test.csv'
write_excel_path='G:/CropYield/03_DL/05_NewML/ParameterResult_ML.xlsx'
tree_graph_dot_path='G:/CropYield/03_DL/05_NewML/tree.dot'
tree_graph_png_path='G:/CropYield/03_DL/05_NewML/tree.png'
random_seed=44
random_forest_seed=np.random.randint(low=1,high=230)
# Data import
train_data=pd.read_csv(train_data_path,header=0)
test_data=pd.read_csv(test_data_path,header=0)
# Separate independent and dependent variables
train_Y=np.array(train_data['Yield'])
train_X=train_data.drop(['ID','Yield'],axis=1)
train_X_column_name=list(train_X.columns)
train_X=np.array(train_X)
test_Y=np.array(test_data['Yield'])
test_X=test_data.drop(['ID','Yield'],axis=1)
test_X=np.array(test_X)
# Search optimal hyperparameter
n_estimators_range=[int(x) for x in np.linspace(start=50,stop=3000,num=60)]
max_features_range=['auto','sqrt']
max_depth_range=[int(x) for x in np.linspace(10,500,num=50)]
max_depth_range.append(None)
min_samples_split_range=[2,5,10]
min_samples_leaf_range=[1,2,4,8]
bootstrap_range=[True,False]
random_forest_hp_range={'n_estimators':n_estimators_range,
                        'max_features':max_features_range,
                        'max_depth':max_depth_range,
                        'min_samples_split':min_samples_split_range,
                        'min_samples_leaf':min_samples_leaf_range
                        # 'bootstrap':bootstrap_range
                        }
pprint(random_forest_hp_range)
random_forest_model_test_base=RandomForestRegressor()
random_forest_model_test_random=RandomizedSearchCV(estimator=random_forest_model_test_base,
                                                   param_distributions=random_forest_hp_range,
                                                   n_iter=200,
                                                   n_jobs=-1,
                                                   cv=3,
                                                   verbose=1,
                                                   random_state=random_forest_seed
                                                   )
random_forest_model_test_random.fit(train_X,train_Y)
best_hp_now=random_forest_model_test_random.best_params_
pprint(best_hp_now)
# Grid Search
random_forest_hp_range_2={'n_estimators':[60,100,200],
                          'max_features':[12,13],
                          'max_depth':[350,400,450],
                          'min_samples_split':[2,3] # Greater than 1
                          # 'min_samples_leaf':[1,2]
                          # 'bootstrap':bootstrap_range
                          }
random_forest_model_test_2_base=RandomForestRegressor()
random_forest_model_test_2_random=GridSearchCV(estimator=random_forest_model_test_2_base,
                                               param_grid=random_forest_hp_range_2,
                                               cv=3,
                                               verbose=1,
                                               n_jobs=-1)
random_forest_model_test_2_random.fit(train_X,train_Y)
best_hp_now_2=random_forest_model_test_2_random.best_params_
pprint(best_hp_now_2)
# Build RF regression model with optimal hyperparameters
random_forest_model_final=random_forest_model_test_2_random.best_estimator_
# Predict test set data
random_forest_predict=random_forest_model_test_2_random.predict(test_X)
random_forest_error=random_forest_predict-test_Y
# Draw test plot
plt.figure(1)
plt.clf()
ax=plt.axes(aspect='equal')
plt.scatter(test_Y,random_forest_predict)
plt.xlabel('True Values')
plt.ylabel('Predictions')
Lims=[0,10000]
plt.xlim(Lims)
plt.ylim(Lims)
plt.plot(Lims,Lims)
plt.grid(False)
plt.figure(2)
plt.clf()
plt.hist(random_forest_error,bins=30)
plt.xlabel('Prediction Error')
plt.ylabel('Count')
plt.grid(False)
# Verify the accuracy
random_forest_pearson_r=stats.pearsonr(test_Y,random_forest_predict)
random_forest_R2=metrics.r2_score(test_Y,random_forest_predict)
random_forest_RMSE=metrics.mean_squared_error(test_Y,random_forest_predict)**0.5
print('Pearson correlation coefficient is {0}, and RMSE is {1}.'.format(random_forest_pearson_r[0],
                                                                        random_forest_RMSE))
# Save key parameters
excel_file=load_workbook(write_excel_path)
excel_all_sheet=excel_file.sheetnames
excel_write_sheet=excel_file[excel_all_sheet[0]]
excel_write_sheet=excel_file.active
max_row=excel_write_sheet.max_row
excel_write_content=[random_forest_pearson_r[0],random_forest_R2,random_forest_RMSE,
                     random_seed,random_forest_seed]
for i in range(len(excel_write_content)):
        exec("excel_write_sheet.cell(max_row+1,i+1).value=excel_write_content[i]")
excel_file.save(write_excel_path)
# Draw decision tree visualizing plot
random_forest_tree=random_forest_model_final.estimators_[5]
export_graphviz(random_forest_tree,out_file=tree_graph_dot_path,
                feature_names=train_X_column_name,rounded=True,precision=1)
(random_forest_graph,)=pydot.graph_from_dot_file(tree_graph_dot_path)
random_forest_graph.write_png(tree_graph_png_path)
# Calculate the importance of variables
random_forest_importance=list(random_forest_model_final.feature_importances_)
random_forest_feature_importance=[(feature,round(importance,8)) 
                                  for feature, importance in zip(train_X_column_name,
                                                                 random_forest_importance)]
random_forest_feature_importance=sorted(random_forest_feature_importance,key=lambda x:x[1],reverse=True)
plt.figure(3)
plt.clf()
importance_plot_x_values=list(range(len(random_forest_importance)))
plt.bar(importance_plot_x_values,random_forest_importance,orientation='vertical')
plt.xticks(importance_plot_x_values,train_X_column_name,rotation='vertical')
plt.xlabel('Variable')
plt.ylabel('Importance')
plt.title('Variable Importances')