一、Adaboost算法
1.1 Adaboost分类算法
adaboost_discrete_c.py
import numpy as np
import copy
from ch4.decision_tree_C import DecisionTreeClassifierclass AdaBoostClassifier:"""adaboost分类算法:既可以做二分类、也可以做多分类,取决于基本分类器1. 同质学习器:非列表形式,按同种类型的基学习器构造2. 异质学习器:列表传递[logisticsregression, svm, cart, ...]"""def __init__(self, base_estimator=None, n_estimators=10, learning_rate=1.0):""":param base_estimator: 基学习器:param n_estimators: 基学习器的个数T:param learning_rate: 学习率,降低后续训练的基分类器的权重,避免过拟合"""self.base_estimator = base_estimatorself.n_estimators = n_estimatorsself.learning_rate = learning_rate# 如果不提供基学习器,则默认按照深度为2的决策树作为基分类器if self.base_estimator is None:self.base_estimator = DecisionTreeClassifier(max_depth=2)if type(base_estimator) != list:# 同质(同种类型)的分类器,深拷贝self.base_estimator = [copy.deepcopy(self.base_estimator)for _ in range(self.n_estimators)]else:# 异质(不同种类型)的分类器self.n_estimators = len(self.base_estimator)self.estimator_weights = [] # 每个基学习器的权重系数def fit(self, x_train, y_train):"""训练AdaBoost每个基学习器,计算权重分布,每个基学习器的误差率和权重系数α,:param x_train: 训练集,二维数组:m * k:param y_train: 目标集:return:"""x_train, y_train = np.asarray(x_train), np.asarray(y_train)n_samples, n_class = x_train.shape[0], len(set(y_train)) # 样本量,类别数sample_weight = np.ones(n_samples) # 为适应自写的基学习器,设置样本均匀权重为1.0# 针对每一个基学习器,根据带有权重分布的训练集训练基学习器,计算相关参数for idx in range(self.n_estimators):# 1. 使用具有权重分布的Dm的训练数据集学习,并预测self.base_estimator[idx].fit(x_train, y_train, sample_weight=sample_weight)# 只关心分类错误的,如果分类错误,则为0,正确则为1y_hat_0 = (self.base_estimator[idx].predict(x_train) == y_train).astype(int)# 2. 计算分类误差率error_rate = sample_weight.dot(1.0 - y_hat_0) / n_samplesif error_rate > 0.5:self.estimator_weights.append(0) # 当前基分类器不起作用continue# 3. 计算基分类器的权重系数,考虑溢出alpha_rate = 0.5 * np.log((1 - error_rate) / error_rate + 1e-8) + np.log(n_class - 1)alpha_rate = min(10.0, alpha_rate) # 避免权重系数过大self.estimator_weights.append(alpha_rate)# 4. 更新样本权重,为了适应多分类,yi*Gm(xi)计算np.power(-1.0, 1 - y_hat_0)sample_weight *= np.exp(-1.0 * alpha_rate * np.power(-1.0, 1 - y_hat_0))sample_weight = sample_weight / np.sum(sample_weight) * n_samples# 5. 更新estimator的权重系数,按照学习率for i in range(self.n_estimators):self.estimator_weights[i] *= np.power(self.learning_rate, i)def predict_proba(self, x_test):"""预测测试样本所属类别概率,软投票:param x_test: 测试样本集:return:"""x_test = np.asarray(x_test)# 按照加法模型,线性组合基学习器# 每个测试样本,每个基学习器预测概率(10,[(0.68, 0.32),(0.55, 0.45)]...)y_hat_prob = np.sum([self.base_estimator[i].predict_proba(x_test) *self.estimator_weights[i] for i in range(self.n_estimators)], axis=0)return y_hat_prob / y_hat_prob.sum(axis=1, keepdims=True)def predict(self, x_test):"""预测测试样本所属类别:param x_test: 测试样本集:return:"""return np.argmax(self.predict_proba(x_test), axis=1)
1.2 Adaboost回归算法
adaboost_regressor.py
import numpy as np
import copy
from ch4.decision_tree_R import DecisionTreeRegression # CARTclass AdaBoostRegressior:"""adaboost回归算法:结合(集成)策略:加权中位数、预测值的平均加权1. 同质学习器,异质学习器2. 回归误差率依赖于相对误差:平方误差、线性误差、指数误差"""def __init__(self, base_estimator=None, n_estimators=10, learning_rate=1.0,loss="square", comb_strategy="weight_median"):""":param base_estimator: 基学习器:param n_estimators: 基学习器的个数T:param learning_rate: 学习率,降低后续训练的基分类器的权重,避免过拟合:param loss: 损失函数:linear、square、exp:param comb_strategy: weight_median、weight_mean"""self.base_estimator = base_estimatorself.n_estimators = n_estimatorsself.learning_rate = learning_rateself.loss = loss # 相对误差的损失函数self.comb_strategy = comb_strategy # 结合策略# 如果不提供基学习器,则默认按照深度为2的决策树作为基分类器if self.base_estimator is None:self.base_estimator = DecisionTreeRegression(max_depth=2)if type(base_estimator) != list:# 同质(同种类型)的分类器,深拷贝self.base_estimator = [copy.deepcopy(self.base_estimator)for _ in range(self.n_estimators)]else:# 异质(不同种类型)的分类器self.n_estimators = len(self.base_estimator)self.estimator_weights = [] # 每个基学习器的权重系数def _cal_loss(self, y_true, y_hat):"""根据损失函数计算相对误差:param y_true: 真值:param y_hat: 预测值:return:"""errors = np.abs(y_true - y_hat) # 绝对值误差if self.loss.lower() == "linear": # 线性return errors / np.max(errors)elif self.loss.lower() == "square": # 平方errors_s = (y_true - y_hat) ** 2return errors_s / np.max(errors) ** 2elif self.loss.lower() == "exp": # 指数return 1 - np.exp(-errors / np.max(errors))else:raise ValueError("仅支持linear、square和exp...")def fit(self, x_train, y_train):"""Adaboost回归算法,T个基学习器的训练:1. 基学习器基于权重分布Dt的训练集训练2. 计算最大绝对误差、相对误差、回归误差率3. 计算当前ht的置信度4. 更新下一轮的权重分布:param x_train::param y_train::return:"""x_train, y_train = np.asarray(x_train), np.asarray(y_train)n_samples, n_class = x_train.shape[0], len(set(y_train)) # 样本量,类别数sample_weight = np.ones(n_samples) # 为适应自写的基学习器,设置样本均匀权重为1.0for idx in range(self.n_estimators):# 1. 基学习器基于权重分布Dt的训练集训练以及预测self.base_estimator[idx].fit(x_train, y_train, sample_weight=sample_weight)y_hat = self.base_estimator[idx].predict(x_train) # 当前训练集的预测值# 2. 计算最大绝对误差、相对误差、回归误差率errors = self._cal_loss(y_train, y_hat) # 相对误差error_rate = np.dot(errors, sample_weight / n_samples) # 回归误差率# 3. 计算当前ht的置信度,基学习器的权重参数alpha_rate = error_rate / (1 - error_rate)self.estimator_weights.append(alpha_rate)# 4. 更新下一轮的权重分布sample_weight *= np.power(alpha_rate, 1 - errors)sample_weight = sample_weight / np.sum(sample_weight) * n_samples# 5. 计算基学习器的权重系数以及考虑学习率self.estimator_weights = np.log(1 / np.asarray(self.estimator_weights))for i in range(self.n_estimators):self.estimator_weights[i] *= np.power(self.learning_rate, i)def predict(self, x_test):"""Adaboost回归算法预测,按照加权中位数以及加权平均两种结合策略:param x_test: 测试样本集:return:"""x_test = np.asarray(x_test)if self.comb_strategy == "weight_mean": # 加权平均self.estimator_weights /= np.sum(self.estimator_weights)# n * Ty_hat_mat = np.array([self.estimator_weights[i] *self.base_estimator[i].predict(x_test)for i in range(self.n_estimators)])# print(y_hat_mat.shape) (10, 5160)return np.sum(y_hat_mat, axis=0)elif self.comb_strategy == "weight_median": # 加权中位数# T个基学习器的预测结果构成一个二维数组(10, 5160)y_hat_mat = np.array([self.estimator_weights[i] *self.base_estimator[i].predict(x_test)for i in range(self.n_estimators)]).Tsorted_idx = np.argsort(y_hat_mat, axis=1) # 二维数组# 按照每个样本预测值的升序排列序号,排序权重系数,然后累加计算weight_cdf = np.cumsum(self.estimator_weights[sorted_idx], axis=1)# 选择最小的t,如下代码产生二维bool数组median_or_above = weight_cdf >= 0.5 * weight_cdf[:, -1][:, np.newaxis]# print(median_idx)median_idx = np.argmax(median_or_above, axis=1) # 返回每个样本的t索引值median_estimators = sorted_idx[np.arange(x_test.shape[0]), median_idx]return y_hat_mat[np.arange(x_test.shape[0]), median_estimators]
1.3 SAMME算法
samme_r_muti_classifier.py
import numpy as np
import copy
from ch4.decision_tree_C import DecisionTreeClassifierclass SAMMERClassifier:"""SAMME.R算法是将SAMME拓展到连续数值型的范畴。基学习器的输出为连续型,一般为类别概率的预测值。"""def __init__(self, base_estimator=None, n_estimators=10):""":param base_estimator: 基学习器:param n_estimators: 基学习器的个数T"""self.base_estimator = base_estimatorself.n_estimators = n_estimators# 如果不提供基学习器,则默认按照深度为2的决策树作为基分类器if self.base_estimator is None:self.base_estimator = DecisionTreeClassifier(max_depth=2)if type(base_estimator) != list:# 同质(同种类型)的分类器,深拷贝self.base_estimator = [copy.deepcopy(self.base_estimator)for _ in range(self.n_estimators)]else:# 异质(不同种类型)的分类器self.n_estimators = len(self.base_estimator)self.estimator_weights = [] # 每个基学习器的权重系数self.n_samples, self.n_class = None, None # 样本量和类别数def _target_encoding(self, y_train):"""对目标值进行编码:param y_train: 训练目标集:return:"""self.n_samples, self.n_class = len(y_train), len(set(y_train))target = -1 / (self.n_class - 1) * np.ones((self.n_samples, self.n_class))for i in range(self.n_samples):target[i, y_train[i]] = 1 # 对应该样本的类别所在编码中的列改为1return targetdef fit(self, x_train, y_train):"""训练SAMME.Rt每个基学习器,根据预测类别概率计算权重分布:param x_train: 训练集,二维数组:m * k:param y_train: 目标集:return:"""x_train, y_train = np.asarray(x_train), np.asarray(y_train)target = self._target_encoding(y_train) # 编码sample_weight = np.ones(self.n_samples) # 为适应自写的基学习器,设置样本均匀权重为1.0# 针对每一个基学习器,根据带有权重分布的训练集训练基学习器,计算相关参数c = (self.n_class - 1) / self.n_classfor idx in range(self.n_estimators):# 1. 使用具有权重分布的Dm的训练数据集学习,并预测self.base_estimator[idx].fit(x_train, y_train, sample_weight=sample_weight)# 根据训练的基学习器,获得其样本的预测类别概率pred_p = self.base_estimator[idx].predict_proba(x_train)# 针对预测概率,小于eps的值替换为eps,避免log函数溢出np.clip(pred_p, np.finfo(pred_p.dtype).eps, None, out=pred_p)# 2. 更新样本权重sample_weight *= np.exp(-c * (target * np.log(pred_p)).sum(axis=1))sample_weight = sample_weight / np.sum(sample_weight) * self.n_samples@staticmethoddef softmax_func(x):"""softmax函数,为避免上溢或下溢,对参数x做限制:param x: 数组: batch_size * n_classes:return: 1 * n_classes"""exps = np.exp(x - np.max(x)) # 避免溢出,每个数减去其最大值exp_sum = np.sum(exps, axis=1, keepdims=True)return exps / exp_sumdef predict_proba(self, x_test):"""预测测试样本所属类别概率,软投票:param x_test: 测试样本集:return:"""x_test = np.asarray(x_test)C_x = np.zeros((x_test.shape[0], self.n_class))for i in range(self.n_estimators):y_prob = self.base_estimator[i].predict_proba(x_test)np.clip(y_prob, np.finfo(y_prob.dtype).eps, None, out=y_prob)y_ln = np.log(y_prob)C_x += (self.n_class - 1) * (y_ln - np.sum(y_ln, axis=1, keepdims=True) / self.n_class)return C_xdef predict(self, x_test):"""预测测试样本所属类别:param x_test: 测试样本集:return:"""return np.argmax(self.predict_proba(x_test), axis=1)
1.4 Adaboost分类算法测试
test_adaboost_c.py
from sklearn.datasets import make_classification
from sklearn.metrics import classification_reportfrom ch4.decision_tree_C import DecisionTreeClassifier # 基学习器,决策树
from ch3.logistic_regression_2class import LogisticRegression # 逻辑回归
from ch6.svm_smo_classifier import SVMClassifier # 支持向量机
from adaboost_discrete_c import AdaBoostClassifier
from ch8.plt_decision_function import plot_decision_functionX, y = make_classification(n_samples=300, n_features=2, n_informative=1, n_redundant=0, n_repeated=0, n_classes=2,n_clusters_per_class=1, class_sep=1, random_state=42)
# 同质:同种类型的基学习器
base_tree = DecisionTreeClassifier(max_depth=3, is_feature_all_R=True, max_bins=20)
ada_bc = AdaBoostClassifier(base_estimator=base_tree, n_estimators=10, learning_rate=1.0)
ada_bc.fit(X, y) # adaboost训练
print("基学习器的权重系数:\n", ada_bc.estimator_weights)
y_pred = ada_bc.predict(X) # 预测类别
print(classification_report(y, y_pred))
plot_decision_function(X, y, ada_bc)# 异质:不同类型的基学习器
log_reg = LogisticRegression(batch_size=20, max_epochs=5)
cart = DecisionTreeClassifier(max_depth=4, is_feature_all_R=True)
svm = SVMClassifier(C=5.0, max_epochs=20)
ada_bc2 = AdaBoostClassifier(base_estimator=[log_reg, cart, svm], learning_rate=1.0)
ada_bc2.fit(X, y) # adaboost训练
print("异质基学习器的权重系数:", ada_bc2.estimator_weights)
y_pred = ada_bc2.predict(X) # 预测类别
print(classification_report(y, y_pred))
plot_decision_function(X, y, ada_bc2)
test_adaboost_c2.py
import numpy as np
import matplotlib.pyplot as plt
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import KFold
from sklearn.metrics import accuracy_score
from sklearn.datasets import make_blobs
from ch4.decision_tree_C import DecisionTreeClassifier
from ch8.adaboost_discrete_c import AdaBoostClassifierX, y = make_blobs(n_samples=1000, n_features=10, centers=5, cluster_std=[1.5, 2, 0.9, 3, 2.8], random_state=0)
X = StandardScaler().fit_transform(X)base_em = DecisionTreeClassifier(max_depth=4, is_feature_all_R=True, max_bins=10)
acc_scores = [] # 存储每次交叉验证的均分
# 用10折交叉验证评估不同基学习器个数T下的分类正确率
for n in range(1, 21):scores = [] # 一次交叉验证的acc均值k_fold = KFold(n_splits=10)for idx_train, idx_test in k_fold.split(X, y):classifier = AdaBoostClassifier(base_estimator=base_em, n_estimators=n, learning_rate=1)classifier.fit(X[idx_train, :], y[idx_train])y_test_pred = classifier.predict(X[idx_test, :])scores.append(accuracy_score(y[idx_test], y_test_pred))acc_scores.append(np.mean(scores))print(n, ":", acc_scores[-1])plt.figure(figsize=(7, 5))
plt.plot(range(1, 21), acc_scores, "ko-", lw=1)
plt.xlabel("Number of Estimations", fontdict={"fontsize": 12})
plt.ylabel("Accuracy Score", fontdict={"fontsize": 12})
plt.title("Cross Validation Scores of Different Number of Base Learners", fontdict={"fontsize": 14})
plt.grid(ls=":")
plt.show()
1.5 Adaboost回归算法测试
test_adaboost_regressor.py
import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import fetch_california_housing
from sklearn.metrics import r2_score
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import train_test_split
from ch4.decision_tree_R import DecisionTreeRegression
from ch8.adaboost_regressor import AdaBoostRegressiorhousing = fetch_california_housing()
X, y = housing.data, housing.target
# print(X.shape)
# print(y.shape)
X = StandardScaler().fit_transform(X)
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.25, random_state=42)base_ht = DecisionTreeRegression(max_bins=50, max_depth=5)
plt.figure(figsize=(14, 15))def train_plot(cs, loss, i):abr = AdaBoostRegressior(base_estimator=base_ht, n_estimators=30,comb_strategy=cs, loss=loss)abr.fit(X_train, y_train)y_hat = abr.predict(X_test)# print(r2_score(y_test, y_hat))plt.subplots(231 + i)idx = np.argsort(y_test) # 对真值排序plt.plot(y_test[idx], "k-", lw=1.5, label="Test True")plt.plot(y_hat[idx], "r-", lw=1, label="Predict")plt.legend(frameon=False)plt.title("%s, %s, R2 = %.5f, MSE = %.5f" %(cs, loss, r2_score(y_test, y_hat), ((y_test - y_hat) ** 2).mean()))plt.xlabel("Test Samples Serial Number", fontdict={"fontsize": 12})plt.ylabel("True VS Predict", fontdict={"fontsize": 12})plt.grid(ls=":")print(cs, loss)loss_func = ["linear", "square", "exp"]
comb_strategy = ["weight_mean", "weight_median"]
i = 0
for loss in loss_func:for cs in comb_strategy:train_plot(cs, loss, i)i += 1
plt.show()
1.6 SAMME算法测试
test_samme_r_c.py
import numpy as np
import matplotlib.pyplot as plt
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import KFold
from sklearn.metrics import accuracy_score
from sklearn.datasets import make_blobs
# from ch4.decision_tree_C import DecisionTreeClassifierR
from sklearn.tree import DecisionTreeClassifier
from ch8.samme_r_muti_classifier import SAMMERClassifierX, y = make_blobs(n_samples=1000, n_features=10, centers=5, cluster_std=[1.5, 2, 0.9, 3, 2.8], random_state=0)
X = StandardScaler().fit_transform(X)base_em = DecisionTreeClassifier(max_depth=4)
acc_scores = [] # 存储每次交叉验证的均分
# 用10折交叉验证评估不同基学习器个数T下的分类正确率
for n in range(1, 21):scores = [] # 一次交叉验证的acc均值k_fold = KFold(n_splits=10)for idx_train, idx_test in k_fold.split(X, y):classifier = SAMMERClassifier(base_estimator=base_em, n_estimators=n)classifier.fit(X[idx_train, :], y[idx_train])y_test_pred = classifier.predict(X[idx_test, :])scores.append(accuracy_score(y[idx_test], y_test_pred))acc_scores.append(np.mean(scores))print(n, ":", acc_scores[-1])plt.figure(figsize=(7, 5))
plt.plot(range(1, 21), acc_scores, "ko-", lw=1)
plt.xlabel("Number of Estimations", fontdict={"fontsize": 12})
plt.ylabel("Accuracy Score", fontdict={"fontsize": 12})
plt.title("Cross Validation Scores of Different Number of Base Learners", fontdict={"fontsize": 14})
plt.grid(ls=":")
plt.show()
1.7 可视化分类边界函数
plt_decision_function.py
import matplotlib.pyplot as plt
import numpy as npdef plot_decision_function(X, y, clf, is_show=True):"""可视化分类边界函数:param X: 测试样本:param y: 测试样本的类别:param clf: 分类模型:param is_show: 是否在当前显示图像,用于父函数绘制子图:return:"""if is_show:plt.figure(figsize=(7, 5))x_min, x_max = X[:, 0].min() - 1, X[:, 0].max() + 1y_min, y_max = X[:, 1].min() - 1, X[:, 1].max() + 1xi, yi = np.meshgrid(np.linspace(x_min, x_max, 100),np.linspace(y_min, y_max, 100))y_pred = clf.predict(np.c_[xi.ravel(), yi.ravel()]) # 模型预测值y_pred = y_pred.reshape(xi.shape)plt.contourf(xi, yi, y_pred, cmap="winter", alpha=0.4)plt.scatter(X[:, 0], X[:, 1], c=y, edgecolors="k")plt.xlabel("Feature 1", fontdict={"fontsize": 12})plt.ylabel("Feature 2", fontdict={"fontsize": 12})plt.title("Model Classification Boundary", fontdict={"fontsize": 14})if is_show:plt.show()
二、提升树算法boosting tree
2.1 提升树回归算法
boostingtree_r.py
import numpy as np
import copy
from ch4.decision_tree_R import DecisionTreeRegression # CARTclass BoostTreeRegressor:"""提升树回归算法,采用平方误差损失"""def __init__(self, base_estimator=None, n_estimators=10, learning_rate=1.0):""":param base_estimator: 基学习器:param n_estimators: 基学习器的个数T:param learning_rate: 学习率,降低后续训练的基分类器的权重,避免过拟合"""self.base_estimator = base_estimatorself.n_estimators = n_estimatorsself.learning_rate = learning_rate# 如果不提供基学习器,则默认按照深度为2的决策树作为基分类器if self.base_estimator is None:self.base_estimator = DecisionTreeRegression(max_depth=2)if type(base_estimator) != list:# 同质(同种类型)的分类器,深拷贝self.base_estimator = [copy.deepcopy(self.base_estimator)for _ in range(self.n_estimators)]else:# 异质(不同种类型)的分类器self.n_estimators = len(self.base_estimator)def fit(self, x_train, y_train):"""提升树的训练,针对每个基决策树算法,拟合上一轮的残差:param x_train: 训练集:param y_train: 目标集:return:"""x_train, y_train = np.asarray(x_train), np.asarray(y_train)# 1. 训练第一棵回归决策树,并预测self.base_estimator[0].fit(x_train, y_train)y_hat = self.base_estimator[0].predict(x_train)y_residual = y_train - y_hat # 残差,MSE的负梯度# 2. 从第二棵树开始,每一轮拟合上一轮的残差for idx in range(1, self.n_estimators):self.base_estimator[idx].fit(x_train, y_residual) # 拟合残差# 累加第m-1棵树的预测值,当前模型是f_(m-1)y_hat += self.base_estimator[idx].predict(x_train) * self.learning_ratey_residual = y_train - y_hat # 当前模型的残差def predict(self, x_test):"""回归提升树的预测:param x_test: 测试样本集:return:"""x_test = np.asarray(x_test)y_hat_mat = np.sum([self.base_estimator[0].predict(x_test)] +[np.power(self.learning_rate, i) * self.base_estimator[i].predict(x_test)for i in range(1, self.n_estimators - 1)] +[self.base_estimator[-1].predict(x_test)], axis=0)return y_hat_mat
2.2 梯度提升树分类算法
gradientboosting_c
import numpy as np
import copy
from ch4.decision_tree_R import DecisionTreeRegression # CARTclass GradientBoostClassifier:"""梯度提升树多分类算法:多分类也可用回归树来做,即训练与类别数相同的几组回归树,每一组代表一个类别,然后对所有组的输出进行softmax操作将其转换为概率分布,再通过交叉熵或者KL一类的损失函数求每棵树相应的负梯度,指导下一轮的训练。"""def __init__(self, base_estimator=None, n_estimators=10, learning_rate=1.0):""":param base_estimator: 基学习器:param n_estimators: 基学习器的个数T:param learning_rate: 学习率,降低后续训练的基分类器的权重,避免过拟合"""self.base_estimator = base_estimatorself.n_estimators = n_estimatorsself.learning_rate = learning_rate# 如果不提供基学习器,则默认按照深度为2的决策树作为基分类器if self.base_estimator is None:self.base_estimator = DecisionTreeRegression(max_depth=2)if type(base_estimator) != list:# 同质(同种类型)的分类器,深拷贝self.base_estimator = [copy.deepcopy(self.base_estimator)for _ in range(self.n_estimators)]else:# 异质(不同种类型)的分类器self.n_estimators = len(self.base_estimator)self.base_estimators = [] # 扩展到class_num组分类器@staticmethoddef one_hot_encoding(target):class_labels = np.unique(target)target_y = np.zeros((len(target), len(class_labels)), dtype=np.int32)for i, label in enumerate(target):target_y[i, label] = 1 # 对应类别所在的列为1return target_y@staticmethoddef softmax_func(x):exps = np.exp(x - np.max(x))return exps / np.sum(exps, axis=1, keepdims=True)def fit(self, x_train, y_train):"""梯度提升分类算法的训练,共训练M * K个基学习器:param x_train: 训练集:param y_train: 目标集:return:"""x_train, y_train = np.asarray(x_train), np.asarray(y_train)class_num = len(np.unique(y_train)) # 类别数y_encoded = self.one_hot_encoding(y_train)# 深拷贝class_num组分类器,每组(每个类别)n_estimators个基学习器# 假设是三分类:[[0, 1, 2, ..., 9], [10], [10]]self.base_estimators = [copy.deepcopy(self.base_estimator) for _ in range(class_num)]# 初始化第一轮基学习器,针对每个类别,分别训练一个基学习器y_hat_scores = [] # 用于存储每个类别的预测值for c_idx in range(class_num):self.base_estimators[c_idx][0].fit(x_train, y_encoded[:, c_idx])y_hat_scores.append(self.base_estimators[c_idx][0].predict(x_train))y_hat_scores = np.c_[y_hat_scores].T # 把每个类别的预测值构成一列,(120, 3) (n_samples, class_num)# print(np.asarray(y_hat_vals).shape)grad_y = y_encoded - self.softmax_func(y_hat_scores) # 按类别计算样本的负梯度值# 训练后续基学习器,共M - 1轮,每轮针对每个类别,分别训练一个基学习器for idx in range(1, self.n_estimators):y_hat_values = [] # 用于存储每个类别的预测值for c_idx in range(class_num):self.base_estimators[c_idx][idx].fit(x_train, grad_y[:, c_idx])y_hat_values.append(self.base_estimators[c_idx][idx].predict(x_train))y_hat_scores += np.c_[y_hat_values].T * self.learning_rate# print(np.asarray(y_hat_vals).shape)grad_y = y_encoded - self.softmax_func(y_hat_scores) # 按类别计算样本的负梯度值def predict_proba(self, x_test):"""预测测试样本所属类别的概率:param x_test: 测试样本集:return:"""x_test = np.asarray(x_test)y_hat_scores = []for c_idx in range(len(self.base_estimators)):# 取当前类别的M个基学习器estimator = self.base_estimators[c_idx]y_hat_scores.append(np.sum([estimator[0].predict(x_test)] +[self.learning_rate * estimator[i].predict(x_test)for i in range(1, self.n_estimators - 1)] +[estimator[-1].predict(x_test)], axis=0))# y_hat_scores的维度(3 * 30)return self.softmax_func(np.c_[y_hat_scores].T)def predict(self, x_test):"""预测测试样本所属类别,概率大的idx标记为类别:param x_test: 测试样本集:return:"""print(self.predict_proba(x_test))return np.argmax(self.predict_proba(x_test), axis=1)
2.3 梯度提升树回归算法
gradientboosting_r
import numpy as np
import copy
from ch4.decision_tree_R import DecisionTreeRegression # CARTclass GradientBoostRegressor:"""梯度提升树回归算法,损失函数:五个,以损失函数在当前模型的负梯度近似为残差"""def __init__(self, base_estimator=None, n_estimators=10, learning_rate=1.0,loss="ls", huber_threshold=0.1, quantile_threshold=0.5):""":param base_estimator: 基学习器:param n_estimators: 基学习器的个数T:param learning_rate: 学习率,降低后续训练的基分类器的权重,避免过拟合"""self.base_estimator = base_estimatorself.n_estimators = n_estimatorsself.learning_rate = learning_rate# 如果不提供基学习器,则默认按照深度为2的决策树作为基分类器if self.base_estimator is None:self.base_estimator = DecisionTreeRegression(max_depth=2)if type(base_estimator) != list:# 同质(同种类型)的分类器,深拷贝self.base_estimator = [copy.deepcopy(self.base_estimator)for _ in range(self.n_estimators)]else:# 异质(不同种类型)的分类器self.n_estimators = len(self.base_estimator)self.loss = loss # 损失函数的类型self.huber_threshold = huber_threshold # 仅对Huber损失有效self.quantile_threshold = quantile_threshold # 仅对分位数损失函数有效def _cal_negative_gradient(self, y_true, y_pred):"""计算负梯度值:param y_true: 真值:param y_pred: 预测值:return:"""if self.loss.lower() == "ls": # MSEreturn y_true - y_predelif self.loss.lower() == "lae": # MAEreturn np.sign(y_true - y_pred)elif self.loss.lower() == "huber": # 平滑平均绝对损失return np.where(np.abs(y_true - y_pred) > self.huber_threshold,self.huber_threshold * np.sign(y_true - y_pred),y_true - y_pred)elif self.loss.lower() == "quantile": # 分位数损失return np.where(y_true > y_pred, self.quantile_threshold,self.quantile_threshold - 1)elif self.loss.lower() == "logcosh": # 双曲余弦的对数的负梯度return -np.tanh(y_pred - y_true)else:raise ValueError("仅限于ls、lae、huber、quantile和logcosh,选择有误...")def fit(self, x_train, y_train):"""提升树的训练,针对每个基决策树算法,拟合上一轮的残差1. 假设回归决策树以mse构建的,针对不同的损失函数,计算不同的基尼指数划分标准2. 预测,集成,也根据不同的损失函数,预测叶子结点的输出...:param x_train: 训练集:param y_train: 目标集:return:"""x_train, y_train = np.asarray(x_train), np.asarray(y_train)# 1. 训练第一棵回归决策树,并预测self.base_estimator[0].fit(x_train, y_train)y_hat = self.base_estimator[0].predict(x_train)y_residual = self._cal_negative_gradient(y_train, y_hat) # 负梯度# 2. 从第二棵树开始,每一轮拟合上一轮的残差for idx in range(1, self.n_estimators):self.base_estimator[idx].fit(x_train, y_residual) # 拟合残差# 累加第m-1棵树的预测值,当前模型是f_(m-1)y_hat += self.base_estimator[idx].predict(x_train) * self.learning_ratey_residual = self._cal_negative_gradient(y_train, y_hat) # 负梯度def predict(self, x_test):"""回归提升树的预测:param x_test: 测试样本集:return:"""x_test = np.asarray(x_test)y_hat_mat = np.sum([self.base_estimator[0].predict(x_test)] +[np.power(self.learning_rate, i) * self.base_estimator[i].predict(x_test)for i in range(1, self.n_estimators - 1)] +[self.base_estimator[-1].predict(x_test)], axis=0)return y_hat_mat
2.4 提升树算法测试
test_boosting_tree_r.py
import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import fetch_california_housing
from sklearn.metrics import r2_score
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import train_test_split
from ch4.decision_tree_R import DecisionTreeRegression
from ch8.boostingtree_r import BoostTreeRegressor
from sklearn.tree import DecisionTreeRegressor# housing = fetch_california_housing()
# X, y = housing.data[0:20000:100, :], housing.target[0:20000:100]
# print(X.shape)
# print(y.shape)
# X = StandardScaler().fit_transform(X)
# X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.25, random_state=42)X = np.linspace(1, 10, 10).reshape(-1, 1)
y = np.array([5.56, 5.70, 5.91, 6.40, 6.80, 7.05, 8.90, 8.70, 9.00, 9.05])# base_ht = DecisionTreeRegression(max_bins=10, max_depth=1)
base_ht = DecisionTreeRegressor(max_depth=1)
# n_estimators = np.linspace(2, 31, 29, dtype=np.int32)
# r2_scores = []
# for n in n_estimators:
# btr = BoostTreeRegressior(base_estimator=base_ht, n_estimators=n)
# btr.fit(X_train, y_train)
# y_hat = btr.predict(X_test)
# # print(r2_score(y_test, y_hat))
# r2_scores.append(r2_score(y_test, y_hat))
# print(n, ":", r2_scores[-1])r2_scores = []
for n in range(1, 7):btr = BoostTreeRegressor(base_estimator=base_ht, n_estimators=n)btr.fit(X, y)y_hat = btr.predict(X)# print(r2_score(y_test, y_hat))r2_scores.append(r2_score(y, y_hat))print(n, ":", r2_scores[-1], np.sum((y - y_hat) ** 2))# plt.figure(figsize=(7, 5))
# plt.plot(n_estimators, r2_scores, "ko-", lw=1)
# plt.show()# idx = np.argsort(y_test) # 对真值排序
#
# plt.figure(figsize=(7, 5))
# plt.plot(y_test[idx], "k-", lw=1.5, label="Test True")
# plt.plot(y_hat[idx], "r-", lw=1, label="Predict")
# plt.legend(frameon=False)
# plt.title("Regression Boosting Tree, R2 = %.5f, MSE = %.5f" %
# (r2_score(y_test, y_hat), ((y_test - y_hat) ** 2).mean()))
# plt.xlabel("Test Samples Serial Number", fontdict={"fontsize": 12})
# plt.ylabel("True VS Predict", fontdict={"fontsize": 12})
# plt.grid(ls=":")
#
# plt.show()
2.5 梯度提升树算法测试
test_gradboost_c1.py
from ch8.gradientboosting_c import GradientBoostClassifier
from sklearn.datasets import load_iris, load_digits, load_breast_cancer
from sklearn.decomposition import PCA
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import train_test_split
from sklearn.metrics import classification_report
from ch4.decision_tree_R import DecisionTreeRegression
from sklearn.tree import DecisionTreeRegressor# iris = load_iris()
# X, y = iris.data, iris.targetdigits = load_digits()
X, y = digits.data, digits.target# bc_data = load_breast_cancer()
# X, y = bc_data.data, bc_data.targetX = PCA(n_components=10).fit_transform(X)
X = StandardScaler().fit_transform(X)X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, shuffle=True, random_state=42)# base_es = DecisionTreeRegression(max_bins=50, max_depth=3)
base_es = DecisionTreeRegressor(max_depth=3)gbc = GradientBoostClassifier(base_estimator=base_es, n_estimators=50)
gbc.fit(X_train, y_train)
y_hat = gbc.predict(X_test)
print(classification_report(y_test, y_hat))
三、Bagging算法
3.1 Bagging算法
bagging_c_r.py
import numpy as np
import copy
from ch4.decision_tree_R import DecisionTreeRegression # CART
from ch4.decision_tree_C import DecisionTreeClassifier
from sklearn.metrics import accuracy_score, r2_scoreclass BaggingClassifierRegressor:"""1. Bagging的基本流程:采样出T个含m个训练样本的采样集,然后基于每个采样集训练出一个基学习器,再集成。2. 预测输出进行结合:Bagging通常对分类任务采用简单投票法,对回归任务使用简单平均法。3. 把回归任务与分类任务集成到一个算法中,右参数task来控制,包外估计OOB控制"""def __init__(self, base_estimator=None, n_estimators=10, task="C", OOB=False):""":param base_estimator: 基学习器:param n_estimators: 基学习器的个数T:param task: 任务:C代表分类任务,R代表回归任务:param OOB: 布尔变量,True表示进行包外估计"""self.base_estimator = base_estimatorself.n_estimators = n_estimatorsself.task = taskif task.lower() not in ["c", "r"]:raise ValueError("Bagging任务仅限分类(C/c)、回归(R/r)")# 如果不提供基学习器,则默认按照深度为2的决策树作为基分类器if self.base_estimator is None:if self.task.lower() == "c":self.base_estimator = DecisionTreeClassifier()elif self.task.lower() == "r":self.base_estimator = DecisionTreeRegression()if type(base_estimator) != list:# 同质(同种类型)的分类器,深拷贝self.base_estimator = [copy.deepcopy(self.base_estimator)for _ in range(self.n_estimators)]else:# 异质(不同种类型)的分类器self.n_estimators = len(self.base_estimator)self.OOB = OOB # 是否进行包外估计self.oob_indices = [] # 保存每次有放回采样未被使用的样本索引self.y_oob_hat = None # 包括估计样本预测值(回归)或预测类别概率(分类)self.oob_score = None # 包外估计的评分,分类和回归def fit(self, x_train, y_train):"""Bagging算法(包含分类和回归)的训练:param x_train: 训练集:param y_train: 目标集:return:"""x_train, y_train = np.asarray(x_train), np.asarray(y_train)n_samples = x_train.shape[0]for estimator in self.base_estimator:# 1. 有放回的随机重采样训练集indices = np.random.choice(n_samples, n_samples, replace=True) # 采样样本索引indices = np.unique(indices)x_bootstrap, y_bootstrap = x_train[indices, :], y_train[indices]# 2. 基于采样数据,训练基学习器estimator.fit(x_bootstrap, y_bootstrap)# 存储每个基学习器未使用的样本索引n_indices = set(np.arange(n_samples)).difference(set(indices))self.oob_indices.append(list(n_indices)) # 每个基学习器未参与训练的样本索引# 3. 包外估计if self.OOB:if self.task.lower() == "c":self._oob_score_classifier(x_train, y_train)else:self._oob_score_regressor(x_train, y_train)def _oob_score_classifier(self, x_train, y_train):"""分类任务的包外估计:param x_train::param y_train::return:"""self.y_oob_hat, y_true = [], []for i in range(x_train.shape[0]): # 针对每个训练样本y_hat_i = [] # 当前样本在每个基学习器下的预测概率,个数未必等于self.n_estimatorsfor idx in range(self.n_estimators): # 针对每个基学习器if i in self.oob_indices[idx]: # 如果该样本属于包外估计y_hat = self.base_estimator[idx].predict_proba(x_train[i, np.newaxis])y_hat_i.append(y_hat[0])# print(y_hat_i)if y_hat_i: # 非空,计算各基学习器预测类别概率的均值self.y_oob_hat.append(np.mean(np.c_[y_hat_i], axis=0))y_true.append(y_train[i]) # 存储对应的真值self.y_oob_hat = np.asarray(self.y_oob_hat)self.oob_score = accuracy_score(y_true, np.argmax(self.y_oob_hat, axis=1))def _oob_score_regressor(self, x_train, y_train):"""回归任务的包外估计:param x_train::param y_train::return:"""self.y_oob_hat, y_true = [], []for i in range(x_train.shape[0]): # 针对每个训练样本y_hat_i = [] # 当前样本在每个基学习器下的预测概率,个数未必等于self.n_estimatorsfor idx in range(self.n_estimators): # 针对每个基学习器if i in self.oob_indices[idx]: # 如果该样本属于包外估计y_hat = self.base_estimator[idx].predict(x_train[i, np.newaxis])y_hat_i.append(y_hat[0])# print(y_hat_i)if y_hat_i: # 非空,计算各基学习器预测类别概率的均值self.y_oob_hat.append(np.mean(y_hat_i))y_true.append(y_train[i]) # 存储对应的真值self.y_oob_hat = np.asarray(self.y_oob_hat)self.oob_score = r2_score(y_true, self.y_oob_hat)def predict_proba(self, x_test):"""分类任务中测试样本所属类别的概率预测:param x_test::return:"""if self.task.lower() != "c":raise ValueError("predict_proba()仅适用于分类任务。")x_test = np.asarray(x_test)y_test_hat = [] # 用于存储测试样本所属类别概率for estimator in self.base_estimator:y_test_hat.append(estimator.predict_proba(x_test))# print(y_test_hat)return np.mean(y_test_hat, axis=0)def predict(self, x_test):"""分类任务:预测测试样本所属类别,类别概率大者索引为所属类别回归任务:预测测试样本,对每个基学习器预测值简单平均:param x_test::return:"""if self.task.lower() == "c":return np.argmax(self.predict_proba(x_test), axis=1)elif self.task.lower() == "r":y_hat = [] # 预测值for estimator in self.base_estimator:y_hat.append(estimator.predict(x_test))return np.mean(y_hat, axis=0)
3.2 Bagging算法测试
test_bagging_c1.py
from sklearn.datasets import load_iris
from ch8.bagging_c_r import BaggingClassifierRegressor
from ch4.decision_tree_C import DecisionTreeClassifier
from sklearn.metrics import classification_report
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaleriris = load_iris()
X, y = iris.data, iris.target
X = StandardScaler().fit_transform(X)
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.25, shuffle=True, random_state=42)base_es = DecisionTreeClassifier(max_depth=10, max_bins=50, is_feature_all_R=True)
bagcr = BaggingClassifierRegressor(base_estimator=base_es, n_estimators=20, task="c", OOB=True)
bagcr.fit(X_train, y_train)
y_hat = bagcr.predict(X_test)
print(classification_report(y_test, y_hat))
print("包外估计的精度:", bagcr.oob_score)
test_bagging_c2.py
from sklearn.datasets import load_iris
from ch8.bagging_c_r import BaggingClassifierRegressor
from ch4.decision_tree_C import DecisionTreeClassifier
from sklearn.metrics import classification_report, accuracy_score
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler, LabelEncoder
import pandas as pd
import numpy as np
import matplotlib.pyplot as pltnursery = pd.read_csv("../ch4/data/nursery.csv").dropna()
X, y = np.asarray(nursery.iloc[:, :-1]), np.asarray(nursery.iloc[:, -1])
y = LabelEncoder().fit_transform(y)X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.30, shuffle=True, random_state=42)base_es = DecisionTreeClassifier(max_depth=10)
bagcr = BaggingClassifierRegressor(base_estimator=base_es, n_estimators=30, task="c")
bagcr.fit(X_train, y_train)
y_hat = bagcr.predict(X_test)
print(classification_report(y_test, y_hat))
# print("包外估计的精度:", bagcr.oob_score)y_test_scores = []
for i in range(30):bagcr = BaggingClassifierRegressor(base_estimator=base_es, n_estimators=1, task="c")bagcr.fit(X_train, y_train)y_hat = bagcr.predict(X_test)y_test_scores.append(accuracy_score(y_test, y_hat))plt.figure(figsize=(7, 5))
plt.plot(range(1, 31), y_test_scores, "ko-", lw=1.5)
plt.xlabel("Training Times", fontsize=12)
plt.ylabel("Test Accuracy", fontsize=12)
plt.grid(ls=":")
plt.show()
test_bagging_r.py
import numpy as np
import matplotlib.pyplot as plt
from ch4.decision_tree_R import DecisionTreeRegression
from ch8.bagging_c_r import BaggingClassifierRegressor
from sklearn.metrics import r2_scoref = lambda x: 0.5 * np.exp(-(x + 3) ** 2) + np.exp(-x ** 2) + 1.5 * np.exp(-(x - 3) ** 2)np.random.seed(0)
N = 200
X = np.random.rand(N) * 10 - 5
X = np.sort(X)
y = f(X) + 0.05 * np.random.randn(N)
X = X.reshape(-1, 1)
# print(X)base_estimator = DecisionTreeRegression(max_bins=30, max_depth=8)
model = BaggingClassifierRegressor(base_estimator=base_estimator, n_estimators=100, task="r")
model.fit(X, y)X_test = np.linspace(1.1 * X.min(axis=0), 1.1 * X.max(axis=0), 1000).reshape(-1, 1)y_bagging_hat = model.predict(X_test)base_estimator.fit(X, y)
y_cart_hat = base_estimator.predict(X_test)plt.figure(figsize=(7, 5))
plt.scatter(X, y, s=10, c="k", label="Raw Data")
plt.plot(X_test, f(X_test), "k-", lw=1.5, label="True F(x)")
plt.plot(X_test, y_bagging_hat, "r-", label="Bagging(R2 = %.5f)" % r2_score(f(X_test), y_bagging_hat))
plt.plot(X_test, y_cart_hat, "b-", label="CART(R2 = %.5f)" % r2_score(f(X_test), y_cart_hat))
plt.legend(frameon=False)
plt.xlabel("X", fontsize=12)
plt.ylabel("Y", fontsize=12)
plt.grid(ls=":")
plt.title("Bagging(100 estimators) VS CART Regression", fontsize=14)
plt.show()
四、随机森林算法
4.1 随机森林算法
rf_classifier_regressor.py
import numpy as np
import copy
from ch4.decision_tree_R import DecisionTreeRegression # CART
from ch4.decision_tree_C import DecisionTreeClassifier
from sklearn.metrics import accuracy_score, r2_scoreclass RandomForestClassifierRegressor:"""随机森林RF是Bagging的一个扩展变体。 RF在以决策树为基学习器构建Bagging集成的基础上,进一步在决策树的训练过程中引入了随机属性选择, 即对训练样本和输入变量增加随机扰动。"""def __init__(self, base_estimator=None, n_estimators=10, feature_sampling_rate=0.5,task="C", OOB=False, feature_importance=False):""":param base_estimator: 基学习器:param n_estimators: 基学习器的个数T:param task: 任务:C代表分类任务,R代表回归任务:param OOB: 布尔变量,True表示进行包外估计:param feature_sampling_rate: 特征变量的抽样率:param feature_importance: 布尔变量,表示是否进行特征重要性的评估"""self.base_estimator = base_estimatorself.n_estimators = n_estimatorsself.feature_sampling_rate = feature_sampling_rateif task.lower() not in ["c", "r"]:raise ValueError("Bagging任务仅限分类(C/c)、回归(R/r)")self.task = task# 如果不提供基学习器,则默认决策树作为基分类器if self.base_estimator is None:if self.task.lower() == "c":base_estimator = DecisionTreeClassifier()elif self.task.lower() == "r":base_estimator = DecisionTreeRegression()self.base_estimator = [copy.deepcopy(base_estimator)for _ in range(self.n_estimators)]self.OOB = OOB # 是否进行包外估计self.oob_indices = [] # 保存每次有放回采样未被使用的样本索引self.y_oob_hat = None # 包括估计样本预测值(回归)或预测类别概率(分类)self.oob_score = None # 包外估计的评分,分类和回归self.feature_importance = feature_importanceself.feature_importance_scores = None # 特征变量的重要性评分self.feature_importance_indices = [] # 针对每个基学习器,存储特征变量的抽样索引def fit(self, x_train, y_train):"""随机森林算法(包含分类和回归)的训练:param x_train: 训练集:param y_train: 目标集:return:"""x_train, y_train = np.asarray(x_train), np.asarray(y_train)n_samples, n_features = x_train.shapefor estimator in self.base_estimator:# 1. 有放回的随机重采样训练集indices = np.random.choice(n_samples, n_samples, replace=True) # 采样样本索引indices = np.unique(indices)x_bootstrap, y_bootstrap = x_train[indices, :], y_train[indices]# 2. 对特征属性变量进行抽样fb_num = int(self.feature_sampling_rate * n_features) # 抽样特征数feature_idx = np.random.choice(n_features, fb_num, replace=False) # 不放回self.feature_importance_indices.append(feature_idx)x_bootstrap = x_bootstrap[:, feature_idx] # 获取特征变量抽样后的训练样本# 3. 基于采样数据,训练基学习器estimator.fit(x_bootstrap, y_bootstrap)# 存储每个基学习器未使用的样本索引n_indices = set(np.arange(n_samples)).difference(set(indices))self.oob_indices.append(list(n_indices)) # 每个基学习器未参与训练的样本索引# 4. 包外估计if self.OOB:if self.task.lower() == "c":self._oob_score_classifier(x_train, y_train)else:self._oob_score_regressor(x_train, y_train)# 5. 特征重要性估计if self.feature_importance:if self.task.lower() == "c":self._feature_importance_score_classifier(x_train, y_train)else:self._feature_importance_score_regressor(x_train, y_train)def _oob_score_classifier(self, x_train, y_train):"""分类任务的包外估计:param x_train::param y_train::return:"""self.y_oob_hat, y_true = [], []for i in range(x_train.shape[0]): # 针对每个训练样本y_hat_i = [] # 当前样本在每个基学习器下的预测概率,个数未必等于self.n_estimatorsfor idx in range(self.n_estimators): # 针对每个基学习器if i in self.oob_indices[idx]: # 如果该样本属于包外估计x_sample = x_train[i, self.feature_importance_indices[idx]]y_hat = self.base_estimator[idx].predict_proba(x_sample.reshape(1, -1))y_hat_i.append(y_hat[0])# print(y_hat_i)if y_hat_i: # 非空,计算各基学习器预测类别概率的均值self.y_oob_hat.append(np.mean(np.c_[y_hat_i], axis=0))y_true.append(y_train[i]) # 存储对应的真值self.y_oob_hat = np.asarray(self.y_oob_hat)self.oob_score = accuracy_score(y_true, np.argmax(self.y_oob_hat, axis=1))def _oob_score_regressor(self, x_train, y_train):"""回归任务的包外估计:param x_train::param y_train::return:"""self.y_oob_hat, y_true = [], []for i in range(x_train.shape[0]): # 针对每个训练样本y_hat_i = [] # 当前样本在每个基学习器下的预测概率,个数未必等于self.n_estimatorsfor idx in range(self.n_estimators): # 针对每个基学习器if i in self.oob_indices[idx]: # 如果该样本属于包外估计x_sample = x_train[i, self.feature_importance_indices[idx]]y_hat = self.base_estimator[idx].predict(x_sample.reshape(1, -1))y_hat_i.append(y_hat[0])# print(y_hat_i)if y_hat_i: # 非空,计算各基学习器预测类别概率的均值self.y_oob_hat.append(np.mean(y_hat_i))y_true.append(y_train[i]) # 存储对应的真值self.y_oob_hat = np.asarray(self.y_oob_hat)self.oob_score = r2_score(y_true, self.y_oob_hat)def _feature_importance_score_classifier(self, x_train, y_train):"""分类问题的特征变量重要性评估计算:param x_train::param y_train::return:"""n_feature = x_train.shape[1]self.feature_importance_scores = np.zeros(n_feature) # 特征变量重要性评分for f_j in range(n_feature): # 针对每个特征变量f_j_scores = [] # 当前第j个特征变量在所有基学习器预测的OOB误差变化for idx, estimator in enumerate(self.base_estimator):f_s_indices = list(self.feature_importance_indices[idx]) # 获取当前基学习器的特征变量索引if f_j in f_s_indices: # 表示当前基学习器中存在第j个特征变量# 1. 计算基于OOB的测试误差errorx_samples = x_train[self.oob_indices[idx], :][:, f_s_indices] # OOB样本以及特征抽样y_hat = estimator.predict(x_samples)error = 1 - accuracy_score(y_train[self.oob_indices[idx]], y_hat)# 2. 计算第j个特征随机打乱顺序后的测试误差np.random.shuffle(x_samples[:, f_s_indices.index(f_j)]) # 原地打乱第j个特征变量取值,其他特征取值不变y_hat_j = estimator.predict(x_samples)error_j = 1 - accuracy_score(y_train[self.oob_indices[idx]], y_hat_j)f_j_scores.append(error_j - error)# 3. 计算所有基学习器对当前第j个特征评分的均值self.feature_importance_scores[f_j] = np.mean(f_j_scores)return self.feature_importance_scoresdef _feature_importance_score_regressor(self, x_train, y_train):"""回归任务的特征变量重要性评估计算:param x_train::param y_train::return:"""n_feature = x_train.shape[1]self.feature_importance_scores = np.zeros(n_feature) # 特征变量重要性评分for f_j in range(n_feature): # 针对每个特征变量f_j_scores = [] # 当前第j个特征变量在所有基学习器预测的OOB误差变化for idx, estimator in enumerate(self.base_estimator):f_s_indices = list(self.feature_importance_indices[idx]) # 获取当前基学习器的特征变量索引if f_j in f_s_indices: # 表示当前基学习器中存在第j个特征变量# 1. 计算基于OOB的测试误差errorx_samples = x_train[self.oob_indices[idx], :][:, f_s_indices] # OOB样本以及特征抽样y_hat = estimator.predict[x_samples]error = 1 - r2_score(y_train[self.oob_indices[idx]], y_hat)# 2. 计算第j个特征随机打乱顺序后的测试误差np.random.shuffle(x_samples[:, f_s_indices.index(f_j)]) # 原地打乱第j个特征变量取值,其他特征取值不变y_hat_j = estimator.predict[x_samples]error_j = 1 - r2_score(y_train[self.oob_indices[idx]], y_hat_j)f_j_scores.append(error_j - error)# 3. 计算所有基学习器对当前第j个特征评分的均值self.feature_importance_scores[f_j] = np.mean(f_j_scores)return self.feature_importance_scoresdef predict_proba(self, x_test):"""分类任务中测试样本所属类别的概率预测:param x_test::return:"""if self.task.lower() != "c":raise ValueError("predict_proba()仅适用于分类任务。")x_test = np.asarray(x_test)y_test_hat = [] # 用于存储测试样本所属类别概率for idx, estimator in enumerate(self.base_estimator):x_test_bootstrap = x_test[:, self.feature_importance_indices[idx]]y_test_hat.append(estimator.predict_proba(x_test_bootstrap))# print(y_test_hat)return np.mean(y_test_hat, axis=0)def predict(self, x_test):"""分类任务:预测测试样本所属类别,类别概率大者索引为所属类别回归任务:预测测试样本,对每个基学习器预测值简单平均:param x_test::return:"""if self.task.lower() == "c":return np.argmax(self.predict_proba(x_test), axis=1)elif self.task.lower() == "r":y_hat = [] # 预测值for idx, estimator in enumerate(self.base_estimator):x_test_bootstrap = x_test[:, self.feature_importance_indices[idx]]y_hat.append(estimator.predict(x_test_bootstrap))return np.mean(y_hat, axis=0)
4.2 随机森林算法测试
test_rf_c1.py
from sklearn.datasets import load_iris, load_wine, load_digits
from ch8.randomforest.rf_classifier_regressor import RandomForestClassifierRegressor
# from ch4.decision_tree_C import DecisionTreeClassifier
from sklearn.tree import DecisionTreeClassifier
from sklearn.metrics import classification_report
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler
import matplotlib.pyplot as plt
import seaborn as sns
import pandas as pdiris = load_iris()
X, y = iris.data, iris.target
X = StandardScaler().fit_transform(X)
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.25, shuffle=True, random_state=42)# base_es = DecisionTreeClassifier(max_depth=10, max_bins=50, is_feature_all_R=True)
base_es = DecisionTreeClassifier(max_depth=10)
rf_model = RandomForestClassifierRegressor(base_estimator=base_es, n_estimators=30,task="c", OOB=True, feature_importance=True)
rf_model.fit(X_train, y_train)
y_hat = rf_model.predict(X_test)
print(classification_report(y_test, y_hat))
print("包外估计的精度:", rf_model.oob_score)
print("特征重要性评分:", rf_model.feature_importance_scores)plt.figure(figsize=(9, 5))
data_pd = pd.DataFrame([iris.feature_names, rf_model.feature_importance_scores]).T
data_pd.columns = ["Feature Names", "Importance"]
sns.barplot(x="Importance", y="Feature Names", data=data_pd)
plt.title("Iris DataSet Feature Importance Scores", fontdict={"fontsize": 14})
plt.grid(ls=":")
print(data_pd)
plt.show()