Machine-Learning-Mastery-机器学习算法教程-三-

Machine Learning Mastery 机器学习算法教程（三）

原文：Machine Learning Mastery

协议：CC BY-NC-SA 4.0

如何用 Python 从零开始实现装袋

原文： machinelearningmastery.com/implement-bagging-scratch-python/

决策树是一种简单而强大的预测性建模技术，但它们存在高度差异。

这意味着在给定不同的训练数据的情况下，树可以得到非常不同

使决策树更加健壮并实现更好表现的技术称为引导程序聚合或简称包装。

在本教程中，您将了解如何使用 Python 从零开始使用决策树实现装袋过程。

完成本教程后，您将了解：

• 如何创建数据集的引导样本。
• 如何使用自举模型做出预测。
• 如何将装袋应用于您自己的预测性建模问题。

让我们开始吧。

• 2017 年 1 月更新：将 cross_validation_split（）中的 fold_size 计算更改为始终为整数。修复了 Python 3 的问题。
• 2017 年 2 月更新：修复了 build_tree 中的错误。
• 2017 年 8 月更新：修正了基尼计算中的一个错误，根据组大小添加了组基尼评分缺失的权重（感谢迈克尔！）。
• 更新 Aug / 2018 ：经过测试和更新，可与 Python 3.6 配合使用。

如何用 Python 实现套装 照片由 Michael Cory 拍摄，保留一些权利。

说明

本节简要介绍 Bootstrap Aggregation 和将在本教程中使用的 Sonar 数据集。

Bootstrap 聚合算法

引导程序是具有替换的数据集的样本。

这意味着从现有数据集的随机样本创建新数据集，其中可以选择给定行并将其多次添加到样本中。

当您只有可用的有限数据集时，在估算诸如更广泛数据集的均值等值时使用它是一种有用的方法。通过创建数据集的样本并估算这些样本的均值，您可以获取这些估计的平均值，并更好地了解潜在问题的真实均值。

这种方法可以与具有高方差的机器学习算法一起使用，例如决策树。针对每个数据引导样本以及用于做出预测的那些模型的平均输出，训练单独的模型。这种技术简称为 bootstrap 聚合或装袋。

方差意味着算法的表现对训练数据敏感，高方差表明训练数据的变化越多，算法的表现就越差。

通过训练许多树并取其预测的平均值，可以改善诸如未修剪的决策树之类的高方差机器学习算法的表现。结果通常优于单个决策树。

除了提高表现之外，装袋的另一个好处是袋装决策树不能过度配合问题。可以继续添加树木，直到达到最大表现。

声纳数据集

我们将在本教程中使用的数据集是 Sonar 数据集。

这是一个描述声纳啁啾返回从不同表面反弹的数据集。 60 个输入变量是不同角度的回报强度。这是一个二分类问题，需要一个模型来区分岩石和金属圆柱。共有 208 个观测结果。

这是一个众所周知的数据集。所有变量都是连续的，通常在 0 到 1 的范围内。输出变量是我的字符串“M”和摇滚的“R”，需要将其转换为整数 1 和 0。

通过预测数据集（M 或矿）中具有最多观测值的类，零规则算法可以实现 53％的准确度。

您可以在 UCI 机器学习库中了解有关此数据集的更多信息。

免费下载数据集并将其放在工作目录中，文件名为 sonar.all-data.csv 。

教程

本教程分为两部分：

1. Bootstrap Resample。
2. 声纳数据集案例研究。

这些步骤提供了实现和将决策树的引导聚合应用于您自己的预测性建模问题所需的基础。

1. Bootstrap Resample

让我们首先深入了解自举法的工作原理。

我们可以通过从数据集中随机选择行并将它们添加到新列表来创建数据集的新样本。我们可以针对固定数量的行重复此操作，或者直到新数据集的大小与原始数据集的大小的比率匹配为止。

我们可以通过不删除已选择的行来允许替换采样，以便将来可以选择。

下面是一个名为 **subsample（）**的函数，它实现了这个过程。来自随机模块的 **randrange（）**函数用于选择随机行索引以在循环的每次迭代中添加到样本。样本的默认大小是原始数据集的大小。

# Create a random subsample from the dataset with replacement
def subsample(dataset, ratio=1.0):
	sample = list()
	n_sample = round(len(dataset) * ratio)
	while len(sample) < n_sample:
		index = randrange(len(dataset))
		sample.append(dataset[index])
	return sample

我们可以使用此函数来估计人为数据集的平均值。

首先，我们可以创建一个包含 20 行和 0 到 9 之间的单列随机数的数据集，并计算平均值。

然后，我们可以制作原始数据集的引导样本，计算平均值，并重复此过程，直到我们有一个均值列表。取这些样本均值的平均值可以给出我们对整个数据集平均值的可靠估计。

下面列出了完整的示例。

每个 bootstrap 样本创建为原始 20 个观察数据集的 10％样本（或 2 个观察值）。然后，我们通过创建原始数据集的 1,10,100 个引导样本，计算它们的平均值，然后平均所有这些估计的平均值来进行实验。

from random import seed
from random import random
from random import randrange

# Create a random subsample from the dataset with replacement
def subsample(dataset, ratio=1.0):
	sample = list()
	n_sample = round(len(dataset) * ratio)
	while len(sample) < n_sample:
		index = randrange(len(dataset))
		sample.append(dataset[index])
	return sample

# Calculate the mean of a list of numbers
def mean(numbers):
	return sum(numbers) / float(len(numbers))

seed(1)
# True mean
dataset = [[randrange(10)] for i in range(20)]
print('True Mean: %.3f' % mean([row[0] for row in dataset]))
# Estimated means
ratio = 0.10
for size in [1, 10, 100]:
	sample_means = list()
	for i in range(size):
		sample = subsample(dataset, ratio)
		sample_mean = mean([row[0] for row in sample])
		sample_means.append(sample_mean)
	print('Samples=%d, Estimated Mean: %.3f' % (size, mean(sample_means)))

运行该示例将打印我们要估计的原始平均值。

然后我们可以从各种不同数量的自举样本中看到估计的平均值。我们可以看到，通过 100 个样本，我们可以很好地估计平均值。

True Mean: 4.450
Samples=1, Estimated Mean: 4.500
Samples=10, Estimated Mean: 3.300
Samples=100, Estimated Mean: 4.480

我们可以从每个子样本创建一个模型，而不是计算平均值。

接下来，让我们看看如何组合多个 bootstrap 模型的预测。

2.声纳数据集案例研究

在本节中，我们将随机森林算法应用于 Sonar 数据集。

该示例假定数据集的 CSV 副本位于当前工作目录中，文件名为 sonar.all-data.csv 。

首先加载数据集，将字符串值转换为数字，并将输出列从字符串转换为 0 到 1 的整数值。这可以通过辅助函数 load_csv（）， **str_column_to_float（ ）**和 **str_column_to_int（）**加载和准备数据集。

我们将使用 k-fold 交叉验证来估计学习模型在看不见的数据上的表现。这意味着我们将构建和评估 k 模型并将表现估计为平均模型误差。分类精度将用于评估每个模型。这些行为在 cross_validation_split（）， **accuracy_metric（）**和 **evaluate_algorithm（）**辅助函数中提供。

我们还将使用适用于装袋的分类和回归树（CART）算法的实现，包括辅助函数 **test_split（）**将数据集分成组， **gini_index（）**来评估分裂点， **get_split（）**找到最佳分裂点， to_terminal（）， **split（）**和 **build_tree（）**使用创建单个决策树，**预测（）用决策树和上一步骤中描述的子样本（）**函数做出预测，以制作训练数据集的子样本

开发了一个名为 **bagging_predict（）**的新函数，负责使用每个决策树做出预测并将预测组合成单个返回值。这是通过从袋装树所做的预测列表中选择最常见的预测来实现的。

最后，开发了一个名为 **bagging（）**的新功能，负责创建训练数据集的样本，在每个样本上训练决策树，然后使用袋装树列表对测试数据集做出预测。

The complete example is listed below.

# Bagging Algorithm on the Sonar dataset
from random import seed
from random import randrange
from csv import reader

# Load a CSV file
def load_csv(filename):
	dataset = list()
	with open(filename, 'r') as file:
		csv_reader = reader(file)
		for row in csv_reader:
			if not row:
				continue
			dataset.append(row)
	return dataset

# Convert string column to float
def str_column_to_float(dataset, column):
	for row in dataset:
		row[column] = float(row[column].strip())

# Convert string column to integer
def str_column_to_int(dataset, column):
	class_values = [row[column] for row in dataset]
	unique = set(class_values)
	lookup = dict()
	for i, value in enumerate(unique):
		lookup[value] = i
	for row in dataset:
		row[column] = lookup[row[column]]
	return lookup

# Split a dataset into k folds
def cross_validation_split(dataset, n_folds):
	dataset_split = list()
	dataset_copy = list(dataset)
	fold_size = int(len(dataset) / n_folds)
	for i in range(n_folds):
		fold = list()
		while len(fold) < fold_size:
			index = randrange(len(dataset_copy))
			fold.append(dataset_copy.pop(index))
		dataset_split.append(fold)
	return dataset_split

# Calculate accuracy percentage
def accuracy_metric(actual, predicted):
	correct = 0
	for i in range(len(actual)):
		if actual[i] == predicted[i]:
			correct += 1
	return correct / float(len(actual)) * 100.0

# Evaluate an algorithm using a cross validation split
def evaluate_algorithm(dataset, algorithm, n_folds, *args):
	folds = cross_validation_split(dataset, n_folds)
	scores = list()
	for fold in folds:
		train_set = list(folds)
		train_set.remove(fold)
		train_set = sum(train_set, [])
		test_set = list()
		for row in fold:
			row_copy = list(row)
			test_set.append(row_copy)
			row_copy[-1] = None
		predicted = algorithm(train_set, test_set, *args)
		actual = [row[-1] for row in fold]
		accuracy = accuracy_metric(actual, predicted)
		scores.append(accuracy)
	return scores

# Split a dataset based on an attribute and an attribute value
def test_split(index, value, dataset):
	left, right = list(), list()
	for row in dataset:
		if row[index] < value:
			left.append(row)
		else:
			right.append(row)
	return left, right

# Calculate the Gini index for a split dataset
def gini_index(groups, classes):
	# count all samples at split point
	n_instances = float(sum([len(group) for group in groups]))
	# sum weighted Gini index for each group
	gini = 0.0
	for group in groups:
		size = float(len(group))
		# avoid divide by zero
		if size == 0:
			continue
		score = 0.0
		# score the group based on the score for each class
		for class_val in classes:
			p = [row[-1] for row in group].count(class_val) / size
			score += p * p
		# weight the group score by its relative size
		gini += (1.0 - score) * (size / n_instances)
	return gini

# Select the best split point for a dataset
def get_split(dataset):
	class_values = list(set(row[-1] for row in dataset))
	b_index, b_value, b_score, b_groups = 999, 999, 999, None
	for index in range(len(dataset[0])-1):
		for row in dataset:
		# for i in range(len(dataset)):
		# 	row = dataset[randrange(len(dataset))]
			groups = test_split(index, row[index], dataset)
			gini = gini_index(groups, class_values)
			if gini < b_score:
				b_index, b_value, b_score, b_groups = index, row[index], gini, groups
	return {'index':b_index, 'value':b_value, 'groups':b_groups}

# Create a terminal node value
def to_terminal(group):
	outcomes = [row[-1] for row in group]
	return max(set(outcomes), key=outcomes.count)

# Create child splits for a node or make terminal
def split(node, max_depth, min_size, depth):
	left, right = node['groups']
	del(node['groups'])
	# check for a no split
	if not left or not right:
		node['left'] = node['right'] = to_terminal(left + right)
		return
	# check for max depth
	if depth >= max_depth:
		node['left'], node['right'] = to_terminal(left), to_terminal(right)
		return
	# process left child
	if len(left) <= min_size:
		node['left'] = to_terminal(left)
	else:
		node['left'] = get_split(left)
		split(node['left'], max_depth, min_size, depth+1)
	# process right child
	if len(right) <= min_size:
		node['right'] = to_terminal(right)
	else:
		node['right'] = get_split(right)
		split(node['right'], max_depth, min_size, depth+1)

# Build a decision tree
def build_tree(train, max_depth, min_size):
	root = get_split(train)
	split(root, max_depth, min_size, 1)
	return root

# Make a prediction with a decision tree
def predict(node, row):
	if row[node['index']] < node['value']:
		if isinstance(node['left'], dict):
			return predict(node['left'], row)
		else:
			return node['left']
	else:
		if isinstance(node['right'], dict):
			return predict(node['right'], row)
		else:
			return node['right']

# Create a random subsample from the dataset with replacement
def subsample(dataset, ratio):
	sample = list()
	n_sample = round(len(dataset) * ratio)
	while len(sample) < n_sample:
		index = randrange(len(dataset))
		sample.append(dataset[index])
	return sample

# Make a prediction with a list of bagged trees
def bagging_predict(trees, row):
	predictions = [predict(tree, row) for tree in trees]
	return max(set(predictions), key=predictions.count)

# Bootstrap Aggregation Algorithm
def bagging(train, test, max_depth, min_size, sample_size, n_trees):
	trees = list()
	for i in range(n_trees):
		sample = subsample(train, sample_size)
		tree = build_tree(sample, max_depth, min_size)
		trees.append(tree)
	predictions = [bagging_predict(trees, row) for row in test]
	return(predictions)

# Test bagging on the sonar dataset
seed(1)
# load and prepare data
filename = 'sonar.all-data.csv'
dataset = load_csv(filename)
# convert string attributes to integers
for i in range(len(dataset[0])-1):
	str_column_to_float(dataset, i)
# convert class column to integers
str_column_to_int(dataset, len(dataset[0])-1)
# evaluate algorithm
n_folds = 5
max_depth = 6
min_size = 2
sample_size = 0.50
for n_trees in [1, 5, 10, 50]:
	scores = evaluate_algorithm(dataset, bagging, n_folds, max_depth, min_size, sample_size, n_trees)
	print('Trees: %d' % n_trees)
	print('Scores: %s' % scores)
	print('Mean Accuracy: %.3f%%' % (sum(scores)/float(len(scores))))

k 值为 5 用于交叉验证，每次迭代时评估每个折叠 208/5 = 41.6 或仅超过 40 个记录。

构建深度树，最大深度为 6，每个节点为 2 的最小训练行数。训练数据集的样本创建为原始数据集大小的 50％。这是为了强制用于训练每棵树的数据集子样本中的某些变体。装袋的默认设置是使样本数据集的大小与原始训练数据集的大小相匹配。

评估了一系列 4 种不同数量的树以显示算法的行为。

打印每个折叠的精度和每个配置的平均精度。随着树木数量的增加，我们可以看到表现略有提升的趋势。

Trees: 1
Scores: [87.8048780487805, 65.85365853658537, 65.85365853658537, 65.85365853658537, 73.17073170731707]
Mean Accuracy: 71.707%

Trees: 5
Scores: [60.97560975609756, 80.48780487804879, 78.04878048780488, 82.92682926829268, 63.41463414634146]
Mean Accuracy: 73.171%

Trees: 10
Scores: [60.97560975609756, 73.17073170731707, 82.92682926829268, 80.48780487804879, 68.29268292682927]
Mean Accuracy: 73.171%

Trees: 50
Scores: [63.41463414634146, 75.60975609756098, 80.48780487804879, 75.60975609756098, 85.36585365853658]
Mean Accuracy: 76.098%

这种方法的一个难点是，即使构建了深树，创建的袋装树也非常相似。反过来，这些树的预测也是相似的，并且我们希望在训练数据集的不同样本上训练的树之间的高方差减小。

这是因为在构造选择相同或相似分裂点的树时使用的贪婪算法。

本教程试图通过约束用于训练每棵树的样本大小来重新注入此方差。更强大的技术是约束在创建每个分割点时可以评估的特征。这是随机森林算法中使用的方法。

扩展

• 调整示例。探索树木数量甚至单个树配置的不同配置，以了解您是否可以进一步改善结果。
• Bag 另一种算法。其他算法可与套袋一起使用。例如，具有低 k 值的 k-最近邻算法将具有高方差并且是用于装袋的良好候选者。
• 回归问题。套袋可以与回归树一起使用。您可以从袋装树中返回预测的平均值，而不是从预测集中预测最常见的类值。回归问题的实验。

你有没有试过这些扩展？ 在下面的评论中分享您的经验。

评论

在本教程中，您了解了如何使用 Python 从零开始实现引导程序聚合。

具体来说，你学到了：

• 如何创建子样本并估计引导数量。
• 如何创建决策树集合并使用它们做出预测。
• 如何将装袋应用于现实世界的预测性建模问题。

你有什么问题吗？ 在下面的评论中提出您的问题，我会尽力回答。

如何用 Python 从零开始实现基线机器学习算法

原文： machinelearningmastery.com/implement-baseline-machine-learning-algorithms-scratch-python/

在预测性建模问题上建立基线表现非常重要。

基线为您稍后评估的更高级方法提供了比较点。

在本教程中，您将了解如何在 Python 中从零开始实现基线机器学习算法。

完成本教程后，您将了解：

• 如何实现随机预测算法。
• 如何实现零规则预测算法。

让我们开始吧。

• 更新 Aug / 2018 ：经过测试和更新，可与 Python 3.6 配合使用。

如何使用 Python 从零开始实现基线机器学习算法 照片由 Vanesser III ，保留一些权利。

描述

有许多机器学习算法可供选择。事实上数以百计。

您必须知道给定算法的预测是否良好。但你怎么知道的？

答案是使用基线预测算法。基线预测算法提供了一组预测，您可以像对问题的任何预测一样进行评估，例如分类准确度或 RMSE。

在评估问题的所有其他机器学习算法时，这些算法的分数提供了所需的比较点。

一旦建立，您可以评论给定算法与朴素基线算法相比有多好，提供给定方法实际有多好的背景。

两种最常用的基线算法是：

• 随机预测算法。
• 零规则算法。

当开始一个比传统分类或回归问题更具粘性的新问题时，首先设计一个特定于您的预测问题的随机预测算法是个好主意。稍后您可以对此进行改进并设计零规则算法。

让我们实现这些算法，看看它们是如何工作的。

教程

本教程分为两部分：

1. 随机预测算法。
2. 零规则算法。

这些步骤将为您的机器学习算法实现和计算基线表现提供所需的基础。

1.随机预测算法

随机预测算法预测在训练数据中观察到的随机结果。

它可能是最简单的算法。

它要求您将所有不同的结果值存储在训练数据中，这可能对具有许多不同值的回归问题很大。

因为随机数用于做出决策，所以在使用算法之前修复随机数种子是个好主意。这是为了确保我们获得相同的随机数集，并且每次运行算法时都会得到相同的决策。

下面是名为 **random_algorithm（）**的函数中的随机预测算法的实现。

该函数既包含包含输出值的训练数据集，也包含必须预测输出值的测试数据集。

该功能适用​​于分类和回归问题。它假定训练数据中的输出值是每行的最后一列。

首先，从训练数据中收集该组唯一输出值。然后，为测试集中的每一行选择随机选择的输出值。

# Generate random predictions
def random_algorithm(train, test):
	output_values = [row[-1] for row in train]
	unique = list(set(output_values))
	predicted = list()
	for row in test:
		index = randrange(len(unique))
		predicted.append(unique[index])
	return predicted

为简单起见，我们可以使用仅包含输出列的小数据集来测试此函数。

训练数据集中的输出值为“0”或“1”，表示算法将从中选择的预测集是{0,1}。测试集还包含单个列，没有数据，因为预测未知。

from random import seed
from random import randrange

# Generate random predictions
def random_algorithm(train, test):
	output_values = [row[-1] for row in train]
	unique = list(set(output_values))
	predicted = list()
	for row in test:
		index = randrange(len(unique))
		predicted.append(unique[index])
	return predicted

seed(1)
train = [[0], [1], [0], [1], [0], [1]]
test = [[None], [None], [None], [None]]
predictions = random_algorithm(train, test)
print(predictions)

运行该示例计算测试数据集的随机预测并打印这些预测。

[0, 0, 1, 0]

随机预测算法易于实现且运行速度快，但我们可以做得更好作为基线。

2.零规则算法

零规则算法是比随机算法更好的基线。

它使用有关给定问题的更多信息来创建一个规则以做出预测。此规则因问题类型而异。

让我们从分类问题开始，预测一个类标签。

分类

对于分类问题，一条规则是预测训练数据集中最常见的类值。这意味着如果训练数据集具有 90 个类“0”的实例和 10 个类“1”的实例，则它将预测“0”并且实现 90/100 或 90％的基线准确度。

这比随机预测算法要好得多，后者平均只能达到 82％的准确率。有关如何计算随机搜索估计值的详细信息，请参阅以下内容：

= ((0.9 * 0.9) + (0.1 * 0.1)) * 100
= 82%

下面是一个名为 **zero_rule_algorithm_classification（）**的函数，它为分类情况实现了这个功能。

# zero rule algorithm for classification
def zero_rule_algorithm_classification(train, test):
	output_values = [row[-1] for row in train]
	prediction = max(set(output_values), key=output_values.count)
	predicted = [prediction for i in range(len(test))]
	return predicted

该函数利用 **max（）**函数和 key 属性，这有点聪明。

给定训练数据中观察到的类值列表， **max（）**函数采用一组唯一的类值，并调用集合中每个类值的类值列表上的计数。

结果是它返回在训练数据集中观察到的类值列表中具有最高观察值计数的类值。

如果所有类值具有相同的计数，那么我们将选择在数据集中观察到的第一个类值。

一旦我们选择了一个类值，它就会用于对测试数据集中的每一行做出预测。

下面是一个包含设计数据集的工作示例，其中包含 4 个类“0”的示例和 2 个类“1”的示例。我们希望算法选择类值“0”作为测试数据集中每行的预测。

from random import seed
from random import randrange

# zero rule algorithm for classification
def zero_rule_algorithm_classification(train, test):
	output_values = [row[-1] for row in train]
	prediction = max(set(output_values), key=output_values.count)
	predicted = [prediction for i in range(len(train))]
	return predicted

seed(1)
train = [['0'], ['0'], ['0'], ['0'], ['1'], ['1']]
test = [[None], [None], [None], [None]]
predictions = zero_rule_algorithm_classification(train, test)
print(predictions)

运行此示例做出预测并将其打印到屏幕。正如所料，选择并预测了类值“0”。

['0', '0', '0', '0', '0', '0']

现在，让我们看一下回归问题的零规则算法。

回归

回归问题需要预测实际价值。

对实际值的良好默认预测是预测集中趋势。这可以是均值或中位数。

一个好的默认值是使用训练数据中观察到的输出值的平均值（也称为平均值）。

这可能比随机预测具有更低的误差，随机预测将返回任何观察到的输出值。

下面是一个名为 **zero_rule_algorithm_regression（）**的函数。它的工作原理是计算观察到的输出值的平均值。

mean = sum(value) / total values

一旦计算出，则对训练数据中的每一行预测平均值。

from random import randrange

# zero rule algorithm for regression
def zero_rule_algorithm_regression(train, test):
	output_values = [row[-1] for row in train]
	prediction = sum(output_values) / float(len(output_values))
	predicted = [prediction for i in range(len(test))]
	return predicted

可以使用一个简单的示例测试此功能。

我们可以设计一个小数据集，其中已知平均值为 15。

10
15
12
15
18
20

mean = (10 + 15 + 12 + 15 + 18 + 20) / 6
mean = 90 / 6
mean = 15

以下是完整的示例。我们期望测试数据集中的 4 行中的每一行都预测平均值 15。

from random import seed
from random import randrange

# zero rule algorithm for regression
def zero_rule_algorithm_regression(train, test):
	output_values = [row[-1] for row in train]
	prediction = sum(output_values) / float(len(output_values))
	predicted = [prediction for i in range(len(test))]
	return predicted

seed(1)
train = [[10], [15], [12], [15], [18], [20]]
test = [[None], [None], [None], [None]]
predictions = zero_rule_algorithm_regression(train, test)
print(predictions)

运行该示例计算打印的预测输出值。正如所料，测试数据集中每行的平均值为 15。

[15.0, 15.0, 15.0, 15.0, 15.0, 15.0]

扩展

以下是基线算法的一些扩展，您可能希望将工具作为本教程的扩展进行研究。

• 替代中心趋势，其中预测中位数，模式或其他集中趋势计算而不是平均值。
• 移动平均值用于预测最后 n 个记录的平均值的时间序列问题。

评论

在本教程中，您发现了计算机器学习问题的表现基线的重要性。

你现在知道了：

• 如何实现分类和回归问题的随机预测算法。
• 如何实现分类和回归问题的零规则算法。

你有什么问题吗？ 在评论中提出您的问题，我会尽力回答。

如何在 Python 中从零开始实现决策树算法

原文： machinelearningmastery.com/implement-decision-tree-algorithm-scratch-python/

决策树是一种强大的预测方法，非常受欢迎。

它们很受欢迎，因为最终模型很容易被从业者和领域专家所理解。最终决策树可以准确解释为什么进行特定预测，使其对操作使用非常有吸引力。

决策树还为更先进的集合方法提供了基础，例如装袋，随机森林和梯度增强。

在本教程中，您将了解如何使用 Python 从零开始实现分类和回归树算法。

完成本教程后，您将了解：

• 如何计算和评估数据中的候选分裂点。
• 如何安排拆分为决策树结构。
• 如何将分类和回归树算法应用于实际问题。

让我们开始吧。

• 2017 年 1 月更新：将 cross_validation_split（）中的 fold_size 计算更改为始终为整数。修复了 Python 3 的问题。
• 2017 年 2 月更新：修复了 build_tree 中的错误。
• 2017 年 8 月更新：修正了基尼计算中的一个错误，根据组大小添加了组基尼评分缺失的权重（感谢迈克尔！）。
• 更新 Aug / 2018 ：经过测试和更新，可与 Python 3.6 配合使用。

如何在 Python 中从零开始实现决策树算法 Martin Cathrae 的照片，保留一些权利。

说明

本节简要介绍了本教程中使用的分类和回归树算法以及 Banknote 数据集。

分类和回归树

分类和回归树或简称 CART 是 Leo Breiman 引用的首字母缩略词，用于指代可用于分类或回归预测性建模问题的决策树算法。

我们将在本教程中专注于使用 CART 进行分类。

CART 模型的表示是二叉树。这是来自算法和数据结构的相同二叉树，没什么太花哨的（每个节点可以有零个，一个或两个子节点）。

假设变量是数字，节点表示单个输入变量（X）和该变量上的分割点。树的叶节点（也称为终端节点）包含用于做出预测的输出变量（y）。

一旦创建，就可以使用拆分在每个分支之后使用新的数据行来导航树，直到进行最终预测。

创建二元决策树实际上是划分输入空间的过程。贪婪的方法用于划分称为递归二进制分裂的空间。这是一个数值程序，其中所有值都排成一行，并使用成本函数尝试和测试不同的分裂点。

选择具有最佳成本（最低成本，因为我们最小化成本）的分割。基于成本函数，以贪婪的方式评估和选择所有输入变量和所有可能的分裂点。

• 回归：为选择分割点而最小化的成本函数是落在矩形内的所有训练样本的总和平方误差。
• 分类：使用基尼成本函数，其表示节点的纯度，其中节点纯度是指分配给每个节点的训练数据的混合程度。

拆分继续，直到节点包含最少数量的训练示例或达到最大树深度。

钞票数据集

钞票数据集涉及根据从照片中采取的若干措施来预测给定钞票是否是真实的。

数据集包含 1,372 行，包含 5 个数字变量。这是两个类的分类问题（二分类）。

下面提供了数据集中五个变量的列表。

1. 小波变换图像的方差（连续）。
2. 小波变换图像的偏度（连续）。
3. 小波峰度变换图像（连续）。
4. 图像熵（连续）。
5. class（整数）。

下面是数据集的前 5 行的示例

3.6216,8.6661,-2.8073,-0.44699,0
4.5459,8.1674,-2.4586,-1.4621,0
3.866,-2.6383,1.9242,0.10645,0
3.4566,9.5228,-4.0112,-3.5944,0
0.32924,-4.4552,4.5718,-0.9888,0
4.3684,9.6718,-3.9606,-3.1625,0

使用零规则算法预测最常见的类值，问题的基线准确率约为 50％。

您可以从 UCI 机器学习库了解更多信息并下载数据集。

下载数据集并将其放在当前工作目录中，文件名为 data_banknote_authentication.csv 。

教程

本教程分为 5 个部分：

1. 基尼指数。
2. 创建拆分。
3. 建树。
4. 做一个预测。
5. 钞票案例研究。

这些步骤将为您提供从零开始实现 CART 算法所需的基础，并将其应用于您自己的预测性建模问题。

基尼系数

Gini 索引是用于评估数据集中拆分的成本函数的名称。

数据集中的拆分涉及一个输入属性和该属性的一个值。它可用于将训练模式划分为两组行。

基尼分数通过分割创建的两个组中的类的混合程度，可以了解分割的好坏程度。完美分离导致基尼评分为 0，而最差情况分裂导致每组 50/50 分类导致基尼评分为 0.5（对于 2 类问题）。

通过示例可以最好地演示计算基尼系数。

我们有两组数据，每组有 2 行。第一组中的行都属于类 0，第二组中的行属于类 1，因此它是完美的分割。

我们首先需要计算每组中班级的比例。

proportion = count(class_value) / count(rows)

这个例子的比例是：

group_1_class_0 = 2 / 2 = 1
group_1_class_1 = 0 / 2 = 0
group_2_class_0 = 0 / 2 = 0
group_2_class_1 = 2 / 2 = 1

然后为每个子节点计算 Gini，如下所示：

gini_index = sum(proportion * (1.0 - proportion))
gini_index = 1.0 - sum(proportion * proportion)

然后，必须根据组的大小，相对于父组中的所有样本，对每组的基尼系数进行加权。当前正在分组的所有样本。我们可以将此权重添加到组的 Gini 计算中，如下所示：

gini_index = (1.0 - sum(proportion * proportion)) * (group_size/total_samples)

在此示例中，每组的基尼评分计算如下：

Gini(group_1) = (1 - (1*1 + 0*0)) * 2/4
Gini(group_1) = 0.0 * 0.5 
Gini(group_1) = 0.0 
Gini(group_2) = (1 - (0*0 + 1*1)) * 2/4
Gini(group_2) = 0.0 * 0.5 
Gini(group_2) = 0.0

然后在分割点处的每个子节点上添加分数，以给出可以与其他候选分割点进行比较的分割点的最终基尼分数。

然后，此分裂点的基尼计算为 0.0 + 0.0 或完美基尼得分为 0.0。

下面是一个名为 **gini_index（）**的函数，它计算组列表的 Gini 索引和已知类值的列表。

你可以看到那里有一些安全检查，以避免空组的除以零。

# Calculate the Gini index for a split dataset
def gini_index(groups, classes):
	# count all samples at split point
	n_instances = float(sum([len(group) for group in groups]))
	# sum weighted Gini index for each group
	gini = 0.0
	for group in groups:
		size = float(len(group))
		# avoid divide by zero
		if size == 0:
			continue
		score = 0.0
		# score the group based on the score for each class
		for class_val in classes:
			p = [row[-1] for row in group].count(class_val) / size
			score += p * p
		# weight the group score by its relative size
		gini += (1.0 - score) * (size / n_instances)
	return gini

我们可以使用上面的工作示例测试此函数。我们还可以测试每组中 50/50 分裂的最坏情况。下面列出了完整的示例。

# Calculate the Gini index for a split dataset
def gini_index(groups, classes):
	# count all samples at split point
	n_instances = float(sum([len(group) for group in groups]))
	# sum weighted Gini index for each group
	gini = 0.0
	for group in groups:
		size = float(len(group))
		# avoid divide by zero
		if size == 0:
			continue
		score = 0.0
		# score the group based on the score for each class
		for class_val in classes:
			p = [row[-1] for row in group].count(class_val) / size
			score += p * p
		# weight the group score by its relative size
		gini += (1.0 - score) * (size / n_instances)
	return gini

# test Gini values
print(gini_index([[[1, 1], [1, 0]], [[1, 1], [1, 0]]], [0, 1]))
print(gini_index([[[1, 0], [1, 0]], [[1, 1], [1, 1]]], [0, 1]))

运行该示例打印两个 Gini 分数，首先是最差情况的分数为 0.5，然后是最佳情况的分数为 0.0。

0.5
0.0

现在我们知道如何评估拆分的结果，让我们看一下创建拆分。

2.创建拆分

拆分由数据集中的属性和值组成。

我们可以将其概括为要拆分的属性的索引以及在该属性上拆分行的值。这只是索引数据行的有用简写。

创建拆分涉及三个部分，第一部分我们已经看过计算基尼评分。剩下的两部分是：

1. 拆分数据集。
2. 评估所有拆分。

我们来看看每一个。

2.1。拆分数据集

拆分数据集意味着在给定属性索引和该属性的拆分值的情况下将数据集分成两个行列表。

一旦我们拥有这两个组，我们就可以使用上面的 Gini 评分来评估拆分的成本。

拆分数据集涉及迭代每一行，检查属性值是低于还是高于拆分值，并分别将其分配给左侧或右侧组。

下面是一个名为 **test_split（）**的函数，它实现了这个过程。

# Split a dataset based on an attribute and an attribute value
def test_split(index, value, dataset):
	left, right = list(), list()
	for row in dataset:
		if row[index] < value:
			left.append(row)
		else:
			right.append(row)
	return left, right

不是很多。

请注意，右侧组包含索引值大于或等于拆分值的所有行。

2.2。评估所有拆分

通过上面的 Gini 函数和测试分割函数，我们现在拥有评估分割所需的一切。

给定一个数据集，我们必须检查每个属性的每个值作为候选分割，评估分割的成本并找到我们可以做出的最佳分割。

找到最佳拆分后，我们可以将其用作决策树中的节点。

这是一个详尽而贪婪的算法。

我们将使用字典来表示决策树中的节点，因为我们可以按名称存储数据。当选择最佳分割并将其用作树的新节点时，我们将存储所选属性的索引，要分割的属性的值以及由所选分割点分割的两组数据。

每组数据都是其自己的小数据集，只有那些通过拆分过程分配给左或右组的行。您可以想象我们如何在构建决策树时递归地再次拆分每个组。

下面是一个名为 **get_split（）**的函数，它实现了这个过程。您可以看到它迭代每个属性（类值除外），然后遍历该属性的每个值，分割和评估拆分。

记录最佳分割，然后在所有检查完成后返回。

# Select the best split point for a dataset
def get_split(dataset):
	class_values = list(set(row[-1] for row in dataset))
	b_index, b_value, b_score, b_groups = 999, 999, 999, None
	for index in range(len(dataset[0])-1):
		for row in dataset:
			groups = test_split(index, row[index], dataset)
			gini = gini_index(groups, class_values)
			if gini < b_score:
				b_index, b_value, b_score, b_groups = index, row[index], gini, groups
	return {'index':b_index, 'value':b_value, 'groups':b_groups}

我们可以设计一个小数据集来测试这个函数和我们的整个数据集拆分过程。

X1			X2			Y
2.771244718		1.784783929		0
1.728571309		1.169761413		0
3.678319846		2.81281357		0
3.961043357		2.61995032		0
2.999208922		2.209014212		0
7.497545867		3.162953546		1
9.00220326		3.339047188		1
7.444542326		0.476683375		1
10.12493903		3.234550982		1
6.642287351		3.319983761		1

我们可以为每个类使用单独的颜色绘制此数据集。您可以看到，手动选择 X1 的值（图中的 x 轴）来分割此数据集并不困难。

CART Contrived Dataset

下面的例子将所有这些放在一起。

# Split a dataset based on an attribute and an attribute value
def test_split(index, value, dataset):
	left, right = list(), list()
	for row in dataset:
		if row[index] < value:
			left.append(row)
		else:
			right.append(row)
	return left, right

# Calculate the Gini index for a split dataset
def gini_index(groups, classes):
	# count all samples at split point
	n_instances = float(sum([len(group) for group in groups]))
	# sum weighted Gini index for each group
	gini = 0.0
	for group in groups:
		size = float(len(group))
		# avoid divide by zero
		if size == 0:
			continue
		score = 0.0
		# score the group based on the score for each class
		for class_val in classes:
			p = [row[-1] for row in group].count(class_val) / size
			score += p * p
		# weight the group score by its relative size
		gini += (1.0 - score) * (size / n_instances)
	return gini

# Select the best split point for a dataset
def get_split(dataset):
	class_values = list(set(row[-1] for row in dataset))
	b_index, b_value, b_score, b_groups = 999, 999, 999, None
	for index in range(len(dataset[0])-1):
		for row in dataset:
			groups = test_split(index, row[index], dataset)
			gini = gini_index(groups, class_values)
			print('X%d < %.3f Gini=%.3f' % ((index+1), row[index], gini))
			if gini < b_score:
				b_index, b_value, b_score, b_groups = index, row[index], gini, groups
	return {'index':b_index, 'value':b_value, 'groups':b_groups}

dataset = [[2.771244718,1.784783929,0],
	[1.728571309,1.169761413,0],
	[3.678319846,2.81281357,0],
	[3.961043357,2.61995032,0],
	[2.999208922,2.209014212,0],
	[7.497545867,3.162953546,1],
	[9.00220326,3.339047188,1],
	[7.444542326,0.476683375,1],
	[10.12493903,3.234550982,1],
	[6.642287351,3.319983761,1]]
split = get_split(dataset)
print('Split: [X%d < %.3f]' % ((split['index']+1), split['value']))

**get_split（）**功能被修改为打印出每个分割点，并且在评估时它是基尼指数。

运行该示例打印所有 Gini 分数，然后在 X1 的数据集中打印最佳分割的分数＆lt; 6.642，基尼指数为 0.0 或完美分裂。

X1 < 2.771 Gini=0.444
X1 < 1.729 Gini=0.500
X1 < 3.678 Gini=0.286
X1 < 3.961 Gini=0.167
X1 < 2.999 Gini=0.375
X1 < 7.498 Gini=0.286
X1 < 9.002 Gini=0.375
X1 < 7.445 Gini=0.167
X1 < 10.125 Gini=0.444
X1 < 6.642 Gini=0.000
X2 < 1.785 Gini=0.500
X2 < 1.170 Gini=0.444
X2 < 2.813 Gini=0.320
X2 < 2.620 Gini=0.417
X2 < 2.209 Gini=0.476
X2 < 3.163 Gini=0.167
X2 < 3.339 Gini=0.444
X2 < 0.477 Gini=0.500
X2 < 3.235 Gini=0.286
X2 < 3.320 Gini=0.375
Split: [X1 < 6.642]

现在我们知道如何在数据集或行列表中找到最佳分割点，让我们看看如何使用它来构建决策树。

3.建造一棵树

创建树的根节点很简单。

我们使用整个数据集调用上面的 **get_split（）**函数。

向树中添加更多节点更有趣。

构建树可以分为 3 个主要部分：

1. 终端节点。
2. 递归拆分。
3. 建造一棵树。

3.1。终端节点

我们需要决定何时停止种树。

我们可以使用节点在训练数据集中负责的行数和行数来实现。

• 最大树深。这是树的根节点的最大节点数。一旦满足树的最大深度，我们必须停止拆分添加新节点。更深的树木更复杂，更有可能过拟合训练数据。
• 最小节点记录。这是给定节点负责的最小训练模式数。一旦达到或低于此最小值，我们必须停止拆分和添加新节点。预计训练模式太少的节点过于具体，可能会过度训练训练数据。

这两种方法将是用户指定的树构建过程参数。

还有一个条件。可以选择所有行属于一个组的拆分。在这种情况下，我们将无法继续拆分和添加子节点，因为我们将无法在一侧或另一侧拆分记录。

现在我们有一些关于何时停止种植树木的想法。当我们在给定点停止增长时，该节点被称为终端节点并用于进行最终预测。

这是通过获取分配给该节点的行组并选择组中最常见的类值来完成的。这将用于做出预测。

下面是一个名为 **to_terminal（）**的函数，它将为一组行选择一个类值。它返回行列表中最常见的输出值。

# Create a terminal node value
def to_terminal(group):
	outcomes = [row[-1] for row in group]
	return max(set(outcomes), key=outcomes.count)

3.2。递归拆分

我们知道如何以及何时创建终端节点，现在我们可以构建我们的树。

构建决策树涉及在为每个节点创建的组上反复调用上面开发的 **get_split（）**函数。

添加到现有节点的新节点称为子节点。节点可以具有零个子节点（终端节点），一个子节点（一侧直接做出预测）或两个子节点。我们将在给定节点的字典表示中将子节点称为左和右。

创建节点后，我们可以通过再次调用相同的函数，对拆分中的每组数据递归创建子节点。

下面是一个实现此递归过程的函数。它将节点作为参数以及节点中的最大深度，最小模式数和节点的当前深度。

您可以想象这可能首先如何在根节点中传递，并且深度为 1.此函数最好用以下步骤解释：

1. 首先，提取节点分割的两组数据以供使用并从节点中删除。当我们处理这些组时，节点不再需要访问这些数据。
2. 接下来，我们检查左侧或右侧行组是否为空，如果是，我们使用我们拥有的记录创建终端节点。
3. 然后我们检查是否已达到最大深度，如果是，我们创建一个终端节点。
4. 然后我们处理左子节点，如果行组太小则创建终端节点，否则以深度优先方式创建和添加左节点，直到在该分支上到达树的底部。
5. 然后以相同的方式处理右侧，因为我们将构造的树恢复到根。

# Create child splits for a node or make terminal
def split(node, max_depth, min_size, depth):
	left, right = node['groups']
	del(node['groups'])
	# check for a no split
	if not left or not right:
		node['left'] = node['right'] = to_terminal(left + right)
		return
	# check for max depth
	if depth >= max_depth:
		node['left'], node['right'] = to_terminal(left), to_terminal(right)
		return
	# process left child
	if len(left) <= min_size:
		node['left'] = to_terminal(left)
	else:
		node['left'] = get_split(left)
		split(node['left'], max_depth, min_size, depth+1)
	# process right child
	if len(right) <= min_size:
		node['right'] = to_terminal(right)
	else:
		node['right'] = get_split(right)
		split(node['right'], max_depth, min_size, depth+1)

3.3。建造一棵树

我们现在可以将所有部分组合在一起。

构建树包括创建根节点并调用 **split（）**函数，然后递归调用自身以构建整个树。

下面是实现此过程的小 **build_tree（）**函数。

# Build a decision tree
def build_tree(train, max_depth, min_size):
	root = get_split(train)
	split(root, max_depth, min_size, 1)
	return root

我们可以使用上面设计的小数据集测试整个过程。

以下是完整的示例。

还包括一个小的 **print_tree（）**函数，它递归地打印出决策树的节点，每个节点一行。虽然没有真正的决策树图那么引人注目，但它给出了树结构和决策的概念。

# Split a dataset based on an attribute and an attribute value
def test_split(index, value, dataset):
	left, right = list(), list()
	for row in dataset:
		if row[index] < value:
			left.append(row)
		else:
			right.append(row)
	return left, right

# Calculate the Gini index for a split dataset
def gini_index(groups, classes):
	# count all samples at split point
	n_instances = float(sum([len(group) for group in groups]))
	# sum weighted Gini index for each group
	gini = 0.0
	for group in groups:
		size = float(len(group))
		# avoid divide by zero
		if size == 0:
			continue
		score = 0.0
		# score the group based on the score for each class
		for class_val in classes:
			p = [row[-1] for row in group].count(class_val) / size
			score += p * p
		# weight the group score by its relative size
		gini += (1.0 - score) * (size / n_instances)
	return gini

# Select the best split point for a dataset
def get_split(dataset):
	class_values = list(set(row[-1] for row in dataset))
	b_index, b_value, b_score, b_groups = 999, 999, 999, None
	for index in range(len(dataset[0])-1):
		for row in dataset:
			groups = test_split(index, row[index], dataset)
			gini = gini_index(groups, class_values)
			if gini < b_score:
				b_index, b_value, b_score, b_groups = index, row[index], gini, groups
	return {'index':b_index, 'value':b_value, 'groups':b_groups}

# Create a terminal node value
def to_terminal(group):
	outcomes = [row[-1] for row in group]
	return max(set(outcomes), key=outcomes.count)

# Create child splits for a node or make terminal
def split(node, max_depth, min_size, depth):
	left, right = node['groups']
	del(node['groups'])
	# check for a no split
	if not left or not right:
		node['left'] = node['right'] = to_terminal(left + right)
		return
	# check for max depth
	if depth >= max_depth:
		node['left'], node['right'] = to_terminal(left), to_terminal(right)
		return
	# process left child
	if len(left) <= min_size:
		node['left'] = to_terminal(left)
	else:
		node['left'] = get_split(left)
		split(node['left'], max_depth, min_size, depth+1)
	# process right child
	if len(right) <= min_size:
		node['right'] = to_terminal(right)
	else:
		node['right'] = get_split(right)
		split(node['right'], max_depth, min_size, depth+1)

# Build a decision tree
def build_tree(train, max_depth, min_size):
	root = get_split(train)
	split(root, max_depth, min_size, 1)
	return root

# Print a decision tree
def print_tree(node, depth=0):
	if isinstance(node, dict):
		print('%s[X%d < %.3f]' % ((depth*' ', (node['index']+1), node['value'])))
		print_tree(node['left'], depth+1)
		print_tree(node['right'], depth+1)
	else:
		print('%s[%s]' % ((depth*' ', node)))

dataset = [[2.771244718,1.784783929,0],
	[1.728571309,1.169761413,0],
	[3.678319846,2.81281357,0],
	[3.961043357,2.61995032,0],
	[2.999208922,2.209014212,0],
	[7.497545867,3.162953546,1],
	[9.00220326,3.339047188,1],
	[7.444542326,0.476683375,1],
	[10.12493903,3.234550982,1],
	[6.642287351,3.319983761,1]]
tree = build_tree(dataset, 1, 1)
print_tree(tree)

我们可以在运行此示例时更改最大深度参数，并查看对打印树的影响。

最大深度为 1（调用 **build_tree（）**函数时的第二个参数），我们可以看到树使用了我们在上一节中发现的完美分割。这是一个具有一个节点的树，也称为决策树桩。

[X1 < 6.642]
 [0]
 [1]

将最大深度增加到 2，即使不需要，我们也会强制树进行拆分。然后，根节点的左右子节点再次使用 X1 属性来拆分已经完美的类混合。

[X1 < 6.642]
 [X1 < 2.771]
  [0]
  [0]
 [X1 < 7.498]
  [1]
  [1]

最后，反过来说，我们可以强制一个更高级别的分裂，最大深度为 3。

[X1 < 6.642]
 [X1 < 2.771]
  [0]
  [X1 < 2.771]
   [0]
   [0]
 [X1 < 7.498]
  [X1 < 7.445]
   [1]
   [1]
  [X1 < 7.498]
   [1]
   [1]

这些测试表明，很有可能优化实现以避免不必要的拆分。这是一个扩展。

现在我们可以创建一个决策树，让我们看看如何使用它来对新数据做出预测。

4.做出预测

使用决策树做出预测涉及使用专门提供的数据行导航树。

同样，我们可以使用递归函数实现此功能，其中使用左子节点或右子节点再次调用相同的预测例程，具体取决于拆分如何影响提供的数据。

我们必须检查子节点是否是要作为预测返回的终端值，或者它是否是包含要考虑的另一级树的字典节点。

下面是实现此过程的 **predict（）**函数。您可以看到给定节点中的索引和值

您可以看到给定节点中的索引和值如何用于评估提供的数据行是否位于拆分的左侧或右侧。

# Make a prediction with a decision tree
def predict(node, row):
	if row[node['index']] < node['value']:
		if isinstance(node['left'], dict):
			return predict(node['left'], row)
		else:
			return node['left']
	else:
		if isinstance(node['right'], dict):
			return predict(node['right'], row)
		else:
			return node['right']

我们可以使用我们设计的数据集来测试这个功能。下面是一个使用硬编码决策树的示例，该决策树具有最佳分割数据的单个节点（决策树桩）。

该示例对数据集中的每一行做出预测。

# Make a prediction with a decision tree
def predict(node, row):
	if row[node['index']] < node['value']:
		if isinstance(node['left'], dict):
			return predict(node['left'], row)
		else:
			return node['left']
	else:
		if isinstance(node['right'], dict):
			return predict(node['right'], row)
		else:
			return node['right']

dataset = [[2.771244718,1.784783929,0],
	[1.728571309,1.169761413,0],
	[3.678319846,2.81281357,0],
	[3.961043357,2.61995032,0],
	[2.999208922,2.209014212,0],
	[7.497545867,3.162953546,1],
	[9.00220326,3.339047188,1],
	[7.444542326,0.476683375,1],
	[10.12493903,3.234550982,1],
	[6.642287351,3.319983761,1]]

#  predict with a stump
stump = {'index': 0, 'right': 1, 'value': 6.642287351, 'left': 0}
for row in dataset:
	prediction = predict(stump, row)
	print('Expected=%d, Got=%d' % (row[-1], prediction))

运行该示例将按预期为每行打印正确的预测。

Expected=0, Got=0
Expected=0, Got=0
Expected=0, Got=0
Expected=0, Got=0
Expected=0, Got=0
Expected=1, Got=1
Expected=1, Got=1
Expected=1, Got=1
Expected=1, Got=1
Expected=1, Got=1

我们现在知道如何创建决策树并使用它来做出预测。现在，让我们将它应用于真实的数据集。

5.钞票案例研究

本节将 CART 算法应用于 Bank Note 数据集。

第一步是加载数据集并将加载的数据转换为可用于计算分割点的数字。为此，我们将使用辅助函数 **load_csv（）**来加载文件，使用 **str_column_to_float（）**将字符串数转换为浮点数。

我们将使用 5 倍折叠交叉验证来评估算法。这意味着每个折叠中将使用 1372/5 = 274.4 或仅超过 270 个记录。我们将使用辅助函数 **evaluate_algorithm（）**来评估具有交叉验证的算法和 **accuracy_metric（）**来计算预测的准确率。

开发了一个名为 **decision_tree（）**的新函数来管理 CART 算法的应用，首先从训练数据集创建树，然后使用树对测试数据集做出预测。

下面列出了完整的示例。

# CART on the Bank Note dataset
from random import seed
from random import randrange
from csv import reader

# Load a CSV file
def load_csv(filename):
	file = open(filename, "rb")
	lines = reader(file)
	dataset = list(lines)
	return dataset

# Convert string column to float
def str_column_to_float(dataset, column):
	for row in dataset:
		row[column] = float(row[column].strip())

# Split a dataset into k folds
def cross_validation_split(dataset, n_folds):
	dataset_split = list()
	dataset_copy = list(dataset)
	fold_size = int(len(dataset) / n_folds)
	for i in range(n_folds):
		fold = list()
		while len(fold) < fold_size:
			index = randrange(len(dataset_copy))
			fold.append(dataset_copy.pop(index))
		dataset_split.append(fold)
	return dataset_split

# Calculate accuracy percentage
def accuracy_metric(actual, predicted):
	correct = 0
	for i in range(len(actual)):
		if actual[i] == predicted[i]:
			correct += 1
	return correct / float(len(actual)) * 100.0

# Evaluate an algorithm using a cross validation split
def evaluate_algorithm(dataset, algorithm, n_folds, *args):
	folds = cross_validation_split(dataset, n_folds)
	scores = list()
	for fold in folds:
		train_set = list(folds)
		train_set.remove(fold)
		train_set = sum(train_set, [])
		test_set = list()
		for row in fold:
			row_copy = list(row)
			test_set.append(row_copy)
			row_copy[-1] = None
		predicted = algorithm(train_set, test_set, *args)
		actual = [row[-1] for row in fold]
		accuracy = accuracy_metric(actual, predicted)
		scores.append(accuracy)
	return scores

# Split a dataset based on an attribute and an attribute value
def test_split(index, value, dataset):
	left, right = list(), list()
	for row in dataset:
		if row[index] < value:
			left.append(row)
		else:
			right.append(row)
	return left, right

# Calculate the Gini index for a split dataset
def gini_index(groups, classes):
	# count all samples at split point
	n_instances = float(sum([len(group) for group in groups]))
	# sum weighted Gini index for each group
	gini = 0.0
	for group in groups:
		size = float(len(group))
		# avoid divide by zero
		if size == 0:
			continue
		score = 0.0
		# score the group based on the score for each class
		for class_val in classes:
			p = [row[-1] for row in group].count(class_val) / size
			score += p * p
		# weight the group score by its relative size
		gini += (1.0 - score) * (size / n_instances)
	return gini

# Select the best split point for a dataset
def get_split(dataset):
	class_values = list(set(row[-1] for row in dataset))
	b_index, b_value, b_score, b_groups = 999, 999, 999, None
	for index in range(len(dataset[0])-1):
		for row in dataset:
			groups = test_split(index, row[index], dataset)
			gini = gini_index(groups, class_values)
			if gini < b_score:
				b_index, b_value, b_score, b_groups = index, row[index], gini, groups
	return {'index':b_index, 'value':b_value, 'groups':b_groups}

# Create a terminal node value
def to_terminal(group):
	outcomes = [row[-1] for row in group]
	return max(set(outcomes), key=outcomes.count)

# Create child splits for a node or make terminal
def split(node, max_depth, min_size, depth):
	left, right = node['groups']
	del(node['groups'])
	# check for a no split
	if not left or not right:
		node['left'] = node['right'] = to_terminal(left + right)
		return
	# check for max depth
	if depth >= max_depth:
		node['left'], node['right'] = to_terminal(left), to_terminal(right)
		return
	# process left child
	if len(left) <= min_size:
		node['left'] = to_terminal(left)
	else:
		node['left'] = get_split(left)
		split(node['left'], max_depth, min_size, depth+1)
	# process right child
	if len(right) <= min_size:
		node['right'] = to_terminal(right)
	else:
		node['right'] = get_split(right)
		split(node['right'], max_depth, min_size, depth+1)

# Build a decision tree
def build_tree(train, max_depth, min_size):
	root = get_split(train)
	split(root, max_depth, min_size, 1)
	return root

# Make a prediction with a decision tree
def predict(node, row):
	if row[node['index']] < node['value']:
		if isinstance(node['left'], dict):
			return predict(node['left'], row)
		else:
			return node['left']
	else:
		if isinstance(node['right'], dict):
			return predict(node['right'], row)
		else:
			return node['right']

# Classification and Regression Tree Algorithm
def decision_tree(train, test, max_depth, min_size):
	tree = build_tree(train, max_depth, min_size)
	predictions = list()
	for row in test:
		prediction = predict(tree, row)
		predictions.append(prediction)
	return(predictions)

# Test CART on Bank Note dataset
seed(1)
# load and prepare data
filename = 'data_banknote_authentication.csv'
dataset = load_csv(filename)
# convert string attributes to integers
for i in range(len(dataset[0])):
	str_column_to_float(dataset, i)
# evaluate algorithm
n_folds = 5
max_depth = 5
min_size = 10
scores = evaluate_algorithm(dataset, decision_tree, n_folds, max_depth, min_size)
print('Scores: %s' % scores)
print('Mean Accuracy: %.3f%%' % (sum(scores)/float(len(scores))))

该示例使用 5 层的最大树深度和每个节点的最小行数为 10.这些 CART 参数通过一些实验选择，但绝不是最佳的。

运行该示例打印每个折叠的平均分类准确度以及所有折叠的平均表现。

您可以看到 CART 和所选配置的平均分类精度达到了约 97％，这明显优于达到 50％精度的零规则算法。

Scores: [96.35036496350365, 97.08029197080292, 97.44525547445255, 98.17518248175182, 97.44525547445255]
Mean Accuracy: 97.299%

扩展

本节列出了您可能希望探索的本教程的扩展。

• 算法调整。未调整 CART 在 Bank Note 数据集中的应用。尝试使用不同的参数值，看看是否可以获得更好的表现。
• 交叉熵。用于评估分裂的另一个成本函数是交叉熵（logloss）。您可以实现和试验此替代成本函数。
• 树修剪。减少训练数据集过拟合的一项重要技术是修剪树木。调查并实现树修剪方法。
• 分类数据集。该示例设计用于具有数字或序数输入属性的输入数据，尝试分类输入数据和可能使用相等而不是排名的拆分。
• 回归。使用不同的成本函数和方法调整树以进行回归以创建终端节点。
• 更多数据集。将算法应用于 UCI 机器学习库中的更多数据集。

你有没有探索过这些扩展？ 在下面的评论中分享您的经验。

评论

在本教程中，您了解了如何使用 Python 从零开始实现决策树算法。

具体来说，你学到了：

• 如何选择和评估训练数据集中的分割点。
• 如何从多个拆分中递归构建决策树。
• 如何将 CART 算法应用于现实世界的分类预测性建模问题。

你有什么问题吗？ 在下面的评论中提出您的问题，我会尽力回答。

如何用 Python 从零开始实现学习向量量化

原文： machinelearningmastery.com/implement-learning-vector-quantization-scratch-python/

k-最近邻居的限制是您必须保留一个大型训练样例数据库才能做出预测。

学习向量量化算法通过学习最能代表训练数据的更小的模式子集来解决这个问题。

在本教程中，您将了解如何使用 Python 从零开始实现学习向量量化算法。

完成本教程后，您将了解：

• 如何从训练数据集中学习一组码本向量。
• 如何使用学习的码本向量做出预测。
• 如何将学习向量量化应用于实际预测性建模问题。

让我们开始吧。

• 2017 年 1 月更新：将 cross_validation_split（）中的 fold_size 计算更改为始终为整数。修复了 Python 3 的问题。
• 更新 Aug / 2018 ：经过测试和更新，可与 Python 3.6 配合使用。

如何用 Python 从零开始实现学习向量量化 照片由 Tony Faiola ，保留一些权利。

描述

本节简要介绍了学习向量量化算法和我们将在本教程中使用的电离层分类问题

学习向量量化

学习向量量化（LVQ）算法很像 K 最近邻。

通过在模式库中找到最佳匹配来做出预测。不同之处在于，模式库是从训练数据中学习的，而不是使用训练模式本身。

模式库称为码本向量，每个模式称为码本。将码本向量初始化为来自训练数据集的随机选择的值。然后，在许多时期，它们适于使用学习算法最佳地总结训练数据。

学习算法一次显示一个训练记录，在码本向量中找到最佳匹配单元，如果它们具有相同的类，则将其移动到更接近训练记录，或者如果它们具有不同的类，则更远离训练记录。

一旦准备好，码本向量用于使用 K 最近邻算法做出预测，其中 k = 1。

该算法是为分类预测性建模问题而开发的，但可以适用于回归问题。

电离层数据集

电离层数据集根据雷达返回数据预测电离层的结构。

每个实例都描述了大气层雷达回波的特性，任务是预测电离层中是否存在结构。

共有 351 个实例和 34 个数字输入变量，每对雷达脉冲有 17 对 2，通常具有 0-1 的相同比例。类值是一个字符串，其值为“g”表示良好返回，“b”表示不良返回。

使用零规则算法预测具有最多观测值的类，可以实现 64.286％的基线准确度。

您可以从 UCI 机器学习库了解更多信息并下载数据集。

下载数据集并将其放在当前工作目录中，名称为 ionosphere.csv 。

教程

本教程分为 4 个部分：

1. 欧几里德距离。
2. 最佳匹配单位。
3. 训练码本向量。
4. 电离层案例研究。

这些步骤将为实现 LVQ 算法并将其应用于您自己的预测性建模问题奠定基础。

欧几里德距离

需要的第一步是计算数据集中两行之间的距离。

数据行主要由数字组成，计算两行或数字向量之间的距离的简单方法是绘制一条直线。这在 2D 或 3D 中是有意义的，并且可以很好地扩展到更高的尺寸。

我们可以使用欧几里德距离测量来计算两个向量之间的直线距离。它被计算为两个向量之间的平方差之和的平方根。

distance = sqrt( sum( (x1_i - x2_i)² )

其中 x1 是第一行数据， x2 是第二行数据， i 是特定列的索引，因为我们对所有列求和。

对于欧几里德距离，值越小，两个记录就越相似。值为 0 表示两个记录之间没有差异。

下面是一个名为 **euclidean_distance（）**的函数，它在 Python 中实现了这一功能。

# calculate the Euclidean distance between two vectors
def euclidean_distance(row1, row2):
	distance = 0.0
	for i in range(len(row1)-1):
		distance += (row1[i] - row2[i])**2
	return sqrt(distance)

您可以看到该函数假定每行中的最后一列是从距离计算中忽略的输出值。

我们可以用一个小的人为分类数据集测试这个距离函数。当我们构造 LVQ 算法所需的元素时，我们将使用该数据集几次。

X1			X2			Y
2.7810836		2.550537003		0
1.465489372		2.362125076		0
3.396561688		4.400293529		0
1.38807019		1.850220317		0
3.06407232		3.005305973		0
7.627531214		2.759262235		1
5.332441248		2.088626775		1
6.922596716		1.77106367		1
8.675418651		-0.242068655		1
7.673756466		3.508563011		1

综上所述，我们可以编写一个小例子，通过打印第一行和所有其他行之间的距离来测试我们的距离函数。我们希望第一行和它自己之间的距离为 0，这是一个值得注意的好事。

下面列出了完整的示例。

from math import sqrt

# calculate the Euclidean distance between two vectors
def euclidean_distance(row1, row2):
	distance = 0.0
	for i in range(len(row1)-1):
		distance += (row1[i] - row2[i])**2
	return sqrt(distance)

# Test distance function
dataset = [[2.7810836,2.550537003,0],
	[1.465489372,2.362125076,0],
	[3.396561688,4.400293529,0],
	[1.38807019,1.850220317,0],
	[3.06407232,3.005305973,0],
	[7.627531214,2.759262235,1],
	[5.332441248,2.088626775,1],
	[6.922596716,1.77106367,1],
	[8.675418651,-0.242068655,1],
	[7.673756466,3.508563011,1]]

row0 = dataset[0]
for row in dataset:
	distance = euclidean_distance(row0, row)
	print(distance)

运行此示例将打印数据集中第一行和每一行之间的距离，包括其自身。

0.0
1.32901739153
1.94946466557
1.55914393855
0.535628072194
4.85094018699
2.59283375995
4.21422704263
6.52240998823
4.98558538245

现在是时候使用距离计算来定位数据集中的最佳匹配单位。

2.最佳匹配单位

最佳匹配单元或 BMU 是与新数据最相似的码本向量。

要在数据集中找到 BMU 以获取新的数据，我们必须首先计算每个码本与新数据之间的距离。我们可以使用上面的距离函数来做到这一点。

计算距离后，我们必须按照与新数据的距离对所有码本进行排序。然后我们可以返回第一个或最相似的码本向量。

我们可以通过跟踪数据集中每个记录的距离作为元组来进行此操作，按距离（按降序排序）对元组列表进行排序，然后检索 BMU。

下面是一个名为 **get_best_matching_unit（）**的函数，它实现了这个功能。

# Locate the best matching unit
def get_best_matching_unit(codebooks, test_row):
	distances = list()
	for codebook in codebooks:
		dist = euclidean_distance(codebook, test_row)
		distances.append((codebook, dist))
	distances.sort(key=lambda tup: tup[1])
	return distances[0][0]

您可以看到上一步中开发的 **euclidean_distance（）**函数用于计算每个码本与新 test_row 之间的距离。

在使用自定义键的情况下对码本和距离元组的列表进行排序，以确保在排序操作中使用元组中的第二项（ tup [1] ）。

最后，返回顶部或最相似的码本向量作为 BMU。

我们可以使用上一节中准备的小型人为数据集来测试此功能。

下面列出了完整的示例。

from math import sqrt

# calculate the Euclidean distance between two vectors
def euclidean_distance(row1, row2):
	distance = 0.0
	for i in range(len(row1)-1):
		distance += (row1[i] - row2[i])**2
	return sqrt(distance)

# Locate the best matching unit
def get_best_matching_unit(codebooks, test_row):
	distances = list()
	for codebook in codebooks:
		dist = euclidean_distance(codebook, test_row)
		distances.append((codebook, dist))
	distances.sort(key=lambda tup: tup[1])
	return distances[0][0]

# Test best matching unit function
dataset = [[2.7810836,2.550537003,0],
	[1.465489372,2.362125076,0],
	[3.396561688,4.400293529,0],
	[1.38807019,1.850220317,0],
	[3.06407232,3.005305973,0],
	[7.627531214,2.759262235,1],
	[5.332441248,2.088626775,1],
	[6.922596716,1.77106367,1],
	[8.675418651,-0.242068655,1],
	[7.673756466,3.508563011,1]]
test_row = dataset[0]
bmu = get_best_matching_unit(dataset, test_row)
print(bmu)

运行此示例将数据集中的 BMU 打印到第一个记录。正如预期的那样，第一条记录与自身最相似，位于列表的顶部。

[2.7810836, 2.550537003, 0]

使用一组码本向量做出预测是一回事。

我们使用 1 最近邻居算法。也就是说，对于我们希望做出预测的每个新模式，我们在集合中找到最相似的码本向量并返回其关联的类值。

现在我们知道如何从一组码本向量中获得最佳匹配单元，我们需要学习如何训练它们。

3.训练码本向量

训练一组码本向量的第一步是初始化该集合。

我们可以使用训练数据集中随机特征构建的模式对其进行初始化。

下面是一个名为 **random_codebook（）**的函数，它实现了这个功能。从训练数据中选择随机输入和输出特征。

# Create a random codebook vector
def random_codebook(train):
	n_records = len(train)
	n_features = len(train[0])
	codebook = [train[randrange(n_records)][i] for i in range(n_features)]
	return codebook

在将码本向量初始化为随机集之后，必须调整它们以最好地总结训练数据。

这是迭代完成的。

1. 时期：在顶层，对于固定数量的时期或训练数据的曝光重复该过程。
2. 训练数据集：在一个时期内，每次使用一个训练模式来更新该码本向量集。
3. 模式特征：对于给定的训练模式，更新最佳匹配码本向量的每个特征以使其移近或远离。

为每个训练模式找到最佳匹配单元，并且仅更新该最佳匹配单元。训练模式和 BMU 之间的差异被计算为误差。比较类值（假定为列表中的最后一个值）。如果它们匹配，则将错误添加到 BMU 以使其更接近训练模式，否则，将其减去以将其推得更远。

调整 BMU 的量由学习率控制。这是对所有 BMU 所做更改量的加权。例如，学习率为 0.3 意味着 BMU 仅移动了训练模式和 BMU 之间的误差或差异的 30％。

此外，调整学习率以使其在第一时期具有最大效果并且随着训练继续进行直到其在最后时期中具有最小效果的效果较小。这称为线性衰减学习率计划，也可用于人工神经网络。

我们可以按时期总结学习率的衰减如下：

rate = learning_rate * (1.0 - (epoch/total_epochs))

我们可以通过假设学习率为 0.3 和 10 个时期来测试这个等式。每个时期的学习率如下：

Epoch		Effective Learning Rate
0		0.3
1		0.27
2		0.24
3		0.21
4		0.18
5		0.15
6		0.12
7		0.09
8		0.06
9		0.03

我们可以把所有这些放在一起。下面是一个名为 **train_codebooks（）**的函数，它实现了在给定训练数据集的情况下训练一组码本向量的过程。

该函数对训练数据集，创建和训练的码本向量的数量，初始学习率和训练码本向量的时期数量采用 3 个附加参数。

您还可以看到该函数记录每个时期的总和平方误差，并打印一条消息，显示时期编号，有效学习率和总和平方误差分数。在调试训练函数或给定预测问题的特定配置时，这很有用。

您可以看到使用 **random_codebook（）**初始化码本向量和 **get_best_matching_unit（）**函数来查找一个迭代内每个训练模式的 BMU。

# Train a set of codebook vectors
def train_codebooks(train, n_codebooks, lrate, epochs):
	codebooks = [random_codebook(train) for i in range(n_codebooks)]
	for epoch in range(epochs):
		rate = lrate * (1.0-(epoch/float(epochs)))
		sum_error = 0.0
		for row in train:
			bmu = get_best_matching_unit(codebooks, row)
			for i in range(len(row)-1):
				error = row[i] - bmu[i]
				sum_error += error**2
				if bmu[-1] == row[-1]:
					bmu[i] += rate * error
				else:
					bmu[i] -= rate * error
		print('>epoch=%d, lrate=%.3f, error=%.3f' % (epoch, rate, sum_error))
	return codebooks

我们可以将它与上面的示例结合起来，为我们设计的数据集学习一组代码簿向量。

以下是完整的示例。

from math import sqrt
from random import randrange
from random import seed

# calculate the Euclidean distance between two vectors
def euclidean_distance(row1, row2):
	distance = 0.0
	for i in range(len(row1)-1):
		distance += (row1[i] - row2[i])**2
	return sqrt(distance)

# Locate the best matching unit
def get_best_matching_unit(codebooks, test_row):
	distances = list()
	for codebook in codebooks:
		dist = euclidean_distance(codebook, test_row)
		distances.append((codebook, dist))
	distances.sort(key=lambda tup: tup[1])
	return distances[0][0]

# Create a random codebook vector
def random_codebook(train):
	n_records = len(train)
	n_features = len(train[0])
	codebook = [train[randrange(n_records)][i] for i in range(n_features)]
	return codebook

# Train a set of codebook vectors
def train_codebooks(train, n_codebooks, lrate, epochs):
	codebooks = [random_codebook(train) for i in range(n_codebooks)]
	for epoch in range(epochs):
		rate = lrate * (1.0-(epoch/float(epochs)))
		sum_error = 0.0
		for row in train:
			bmu = get_best_matching_unit(codebooks, row)
			for i in range(len(row)-1):
				error = row[i] - bmu[i]
				sum_error += error**2
				if bmu[-1] == row[-1]:
					bmu[i] += rate * error
				else:
					bmu[i] -= rate * error
		print('>epoch=%d, lrate=%.3f, error=%.3f' % (epoch, rate, sum_error))
	return codebooks

# Test the training function
seed(1)
dataset = [[2.7810836,2.550537003,0],
	[1.465489372,2.362125076,0],
	[3.396561688,4.400293529,0],
	[1.38807019,1.850220317,0],
	[3.06407232,3.005305973,0],
	[7.627531214,2.759262235,1],
	[5.332441248,2.088626775,1],
	[6.922596716,1.77106367,1],
	[8.675418651,-0.242068655,1],
	[7.673756466,3.508563011,1]]
learn_rate = 0.3
n_epochs = 10
n_codebooks = 2
codebooks = train_codebooks(dataset, n_codebooks, learn_rate, n_epochs)
print('Codebooks: %s' % codebooks)

运行该示例训练一组 2 个码本向量用于 10 个时期，初始学习率为 0.3。每个时期打印细节，并显示从训练数据中学习的一组 2​​个码本向量。

我们可以看到学习率的变化符合我们上面针对每个时期探讨的期望。我们还可以看到，每个时期的总和平方误差在训练结束时继续下降，并且可能有机会进一步调整示例以实现更少的错误。

>epoch=0, lrate=0.300, error=43.270
>epoch=1, lrate=0.270, error=30.403
>epoch=2, lrate=0.240, error=27.146
>epoch=3, lrate=0.210, error=26.301
>epoch=4, lrate=0.180, error=25.537
>epoch=5, lrate=0.150, error=24.789
>epoch=6, lrate=0.120, error=24.058
>epoch=7, lrate=0.090, error=23.346
>epoch=8, lrate=0.060, error=22.654
>epoch=9, lrate=0.030, error=21.982
Codebooks: [[2.432316086217663, 2.839821664184211, 0], [7.319592257892681, 1.97013382654341, 1]]

现在我们知道如何训练一组码本向量，让我们看看如何在真实数据集上使用这个算法。

4.电离层案例研究

在本节中，我们将学习向量量化算法应用于电离层数据集。

第一步是加载数据集并将加载的数据转换为我们可以与欧氏距离计算一起使用的数字。为此我们将使用辅助函数 **load_csv（）**来加载文件， **str_column_to_float（）**将字符串数转换为浮点数， **str_column_to_int（）**转换​​为 class 列到整数值。

我们将使用 5 倍折叠交叉验证来评估算法。这意味着每个折叠中将有 351/5 = 70.2 或仅超过 70 个记录。我们将使用辅助函数 **evaluate_algorithm（）**来评估具有交叉验证的算法和 **accuracy_metric（）**来计算预测的准确率。

The complete example is listed below.

# LVQ for the Ionosphere Dataset
from random import seed
from random import randrange
from csv import reader
from math import sqrt

# Load a CSV file
def load_csv(filename):
	dataset = list()
	with open(filename, 'r') as file:
		csv_reader = reader(file)
		for row in csv_reader:
			if not row:
				continue
			dataset.append(row)
	return dataset

# Convert string column to float
def str_column_to_float(dataset, column):
	for row in dataset:
		row[column] = float(row[column].strip())

# Convert string column to integer
def str_column_to_int(dataset, column):
	class_values = [row[column] for row in dataset]
	unique = set(class_values)
	lookup = dict()
	for i, value in enumerate(unique):
		lookup[value] = i
	for row in dataset:
		row[column] = lookup[row[column]]
	return lookup

# Split a dataset into k folds
def cross_validation_split(dataset, n_folds):
	dataset_split = list()
	dataset_copy = list(dataset)
	fold_size = int(len(dataset) / n_folds)
	for i in range(n_folds):
		fold = list()
		while len(fold) < fold_size:
			index = randrange(len(dataset_copy))
			fold.append(dataset_copy.pop(index))
		dataset_split.append(fold)
	return dataset_split

# Calculate accuracy percentage
def accuracy_metric(actual, predicted):
	correct = 0
	for i in range(len(actual)):
		if actual[i] == predicted[i]:
			correct += 1
	return correct / float(len(actual)) * 100.0

# Evaluate an algorithm using a cross validation split
def evaluate_algorithm(dataset, algorithm, n_folds, *args):
	folds = cross_validation_split(dataset, n_folds)
	scores = list()
	for fold in folds:
		train_set = list(folds)
		train_set.remove(fold)
		train_set = sum(train_set, [])
		test_set = list()
		for row in fold:
			row_copy = list(row)
			test_set.append(row_copy)
			row_copy[-1] = None
		predicted = algorithm(train_set, test_set, *args)
		actual = [row[-1] for row in fold]
		accuracy = accuracy_metric(actual, predicted)
		scores.append(accuracy)
	return scores

# calculate the Euclidean distance between two vectors
def euclidean_distance(row1, row2):
	distance = 0.0
	for i in range(len(row1)-1):
		distance += (row1[i] - row2[i])**2
	return sqrt(distance)

# Locate the best matching unit
def get_best_matching_unit(codebooks, test_row):
	distances = list()
	for codebook in codebooks:
		dist = euclidean_distance(codebook, test_row)
		distances.append((codebook, dist))
	distances.sort(key=lambda tup: tup[1])
	return distances[0][0]

# Make a prediction with codebook vectors
def predict(codebooks, test_row):
	bmu = get_best_matching_unit(codebooks, test_row)
	return bmu[-1]

# Create a random codebook vector
def random_codebook(train):
	n_records = len(train)
	n_features = len(train[0])
	codebook = [train[randrange(n_records)][i] for i in range(n_features)]
	return codebook

# Train a set of codebook vectors
def train_codebooks(train, n_codebooks, lrate, epochs):
	codebooks = [random_codebook(train) for i in range(n_codebooks)]
	for epoch in range(epochs):
		rate = lrate * (1.0-(epoch/float(epochs)))
		for row in train:
			bmu = get_best_matching_unit(codebooks, row)
			for i in range(len(row)-1):
				error = row[i] - bmu[i]
				if bmu[-1] == row[-1]:
					bmu[i] += rate * error
				else:
					bmu[i] -= rate * error
	return codebooks

# LVQ Algorithm
def learning_vector_quantization(train, test, n_codebooks, lrate, epochs):
	codebooks = train_codebooks(train, n_codebooks, lrate, epochs)
	predictions = list()
	for row in test:
		output = predict(codebooks, row)
		predictions.append(output)
	return(predictions)

# Test LVQ on Ionosphere dataset
seed(1)
# load and prepare data
filename = 'ionosphere.csv'
dataset = load_csv(filename)
for i in range(len(dataset[0])-1):
	str_column_to_float(dataset, i)
# convert class column to integers
str_column_to_int(dataset, len(dataset[0])-1)
# evaluate algorithm
n_folds = 5
learn_rate = 0.3
n_epochs = 50
n_codebooks = 20
scores = evaluate_algorithm(dataset, learning_vector_quantization, n_folds, n_codebooks, learn_rate, n_epochs)
print('Scores: %s' % scores)
print('Mean Accuracy: %.3f%%' % (sum(scores)/float(len(scores))))

运行此示例将打印每个折叠的分类准确度以及所有折叠的平均分类精度。

我们可以看出，87.143％的准确率优于 64.286％的基线。我们还可以看到，我们的 20 个码本向量库远远少于保存整个训练数据集。

Scores: [90.0, 88.57142857142857, 84.28571428571429, 87.14285714285714, 85.71428571428571]
Mean Accuracy: 87.143%

扩展

本节列出了您可能希望探索的教程的扩展。

• 调谐参数。上述示例中的参数未进行调整，请尝试使用不同的值来提高分类准确度。
• 不同的距离测量。尝试不同的距离测量，如曼哈顿距离和闵可夫斯基距离。
• 多次通过 LVQ 。可以通过多次训练运行来更新码本向量。通过大学习率的训练进行实验，接着是大量具有较小学习率的时期来微调码本。
• 更新更多 BMU 。尝试在训练时选择多个 BMU，并将其从训练数据中拉出。
• 更多问题。将 LVQ 应用于 UCI 机器学习存储库中的更多分类问题。

你有没有探索过这些扩展？ 在下面的评论中分享您的经验。

评论

在本教程中，您了解了如何在 Python 中从零开始实现学习向量量化算法。

具体来说，你学到了：

• 如何计算模式之间的距离并找到最佳匹配单元。
• 如何训练一组码本向量以最好地总结训练数据集。
• 如何将学习向量量化算法应用于实际预测性建模问题。

你有什么问题吗？ 在下面的评论中提出您的问题，我会尽力回答。

如何利用 Python 从零开始实现线性回归

原文： machinelearningmastery.com/implement-linear-regression-stochastic-gradient-descent-scratch-python/

许多机器学习算法的核心是优化。

机器学习算法使用优化算法在给定训练数据集的情况下找到一组好的模型参数。

机器学习中最常用的优化算法是随机梯度下降。

在本教程中，您将了解如何实现随机梯度下降，以便从零开始使用 Python 优化线性回归算法。

完成本教程后，您将了解：

• 如何使用随机梯度下降估计线性回归系数。
• 如何预测多元线性回归。
• 如何利用随机梯度下降实现线性回归来预测新数据。

让我们开始吧。

• 2017 年 1 月更新：将 cross_validation_split（）中的 fold_size 计算更改为始终为整数。修复了 Python 3 的问题。
• 更新 Aug / 2018 ：经过测试和更新，可与 Python 3.6 配合使用。

如何使用 Python 随机梯度下降进行线性回归 照片 star5112 ，保留一些权利。

描述

在本节中，我们将描述线性回归，随机梯度下降技术和本教程中使用的葡萄酒质量数据集。

多元线性回归

线性回归是一种预测实际价值的技术。

令人困惑的是，这些预测真实价值的问题被称为回归问题。

线性回归是一种使用直线来模拟输入和输出值之间关系的技术。在两个以上的维度中，该直线可以被认为是平面或超平面。

预测作为输入值的组合来预测输出值。

使用系数（b）对每个输入属性（x）进行加权，并且学习算法的目标是发现导致良好预测（y）的一组系数。

y = b0 + b1 * x1 + b2 * x2 + ...

使用随机梯度下降可以找到系数。

随机梯度下降

梯度下降是通过遵循成本函数的梯度来最小化函数的过程。

这包括了解成本的形式以及衍生物，以便从给定的点知道梯度并且可以在该方向上移动，例如，向下走向最小值。

在机器学习中，我们可以使用一种技术来评估和更新称为随机梯度下降的每次迭代的系数，以最小化模型对我们的训练数据的误差。

此优化算法的工作方式是每个训练实例一次显示给模型一个。该模型对训练实例做出预测，计算误差并更新模型以减少下一次预测的误差。该过程重复固定次数的迭代。

该过程可用于在模型中找到导致训练数据上模型的最小误差的系数集。每次迭代，机器学习语言中的系数（b）使用以下等式更新：

b = b - learning_rate * error * x

其中 b 是被优化的系数或权重， learning_rate 是您必须配置的学习率（例如 0.01），错误是模型的预测误差关于归因于重量的训练数据， x 是输入值。

葡萄酒质量数据集

在我们开发具有随机梯度下降的线性回归算法之后，我们将使用它来模拟葡萄酒质量数据集。

该数据集由 4,898 种白葡萄酒的细节组成，包括酸度和 pH 值等测量值。目标是使用这些客观测量来预测 0 到 10 之间的葡萄酒质量。

以下是此数据集中前 5 个记录的示例。

7,0.27,0.36,20.7,0.045,45,170,1.001,3,0.45,8.8,6
6.3,0.3,0.34,1.6,0.049,14,132,0.994,3.3,0.49,9.5,6
8.1,0.28,0.4,6.9,0.05,30,97,0.9951,3.26,0.44,10.1,6
7.2,0.23,0.32,8.5,0.058,47,186,0.9956,3.19,0.4,9.9,6
7.2,0.23,0.32,8.5,0.058,47,186,0.9956,3.19,0.4,9.9,6

必须将数据集标准化为 0 到 1 之间的值，因为每个属性具有不同的单位，并且反过来具有不同的比例。

通过预测归一化数据集上的平均值（零规则算法），可以实现 0.148 的基线均方根误差（RMSE）。

您可以在 UCI 机器学习库上了解有关数据集的更多信息。

您可以下载数据集并将其保存在当前工作目录中，名称为 winequality-white.csv 。您必须从文件的开头删除标头信息，并将“;”值分隔符转换为“，”以符合 CSV 格式。

教程

本教程分为 3 个部分：

1. 做出预测。
2. 估计系数。
3. 葡萄酒质量预测。

这将为您自己的预测性建模问题提供实现和应用具有随机梯度下降的线性回归所需的基础。

1.做出预测

第一步是开发一个可以做出预测的功能。

在随机梯度下降中的候选系数值的评估中以及在模型完成之后，我们希望开始对测试数据或新数据做出预测。

下面是一个名为 **predict（）**的函数，它预测给定一组系数的行的输出值。

第一个系数 in 总是截距，也称为偏差或 b0，因为它是独立的，不负责特定的输入值。

# Make a prediction with coefficients
def predict(row, coefficients):
	yhat = coefficients[0]
	for i in range(len(row)-1):
		yhat += coefficients[i + 1] * row[i]
	return yhat

我们可以设计一个小数据集来测试我们的预测函数。

x, y
1, 1
2, 3
4, 3
3, 2
5, 5

下面是该数据集的图表。

线性回归的小受控数据集

我们还可以使用先前准备的系数来对该数据集做出预测。

综合这些，我们可以测试下面的 **predict（）**函数。

# Make a prediction with coefficients
def predict(row, coefficients):
	yhat = coefficients[0]
	for i in range(len(row)-1):
		yhat += coefficients[i + 1] * row[i]
	return yhat

dataset = [[1, 1], [2, 3], [4, 3], [3, 2], [5, 5]]
coef = [0.4, 0.8]
for row in dataset:
	yhat = predict(row, coef)
	print("Expected=%.3f, Predicted=%.3f" % (row[-1], yhat))

有一个输入值（x）和两个系数值（b0 和 b1）。我们为这个问题建模的预测方程是：

y = b0 + b1 * x

或者，我们手动选择的具体系数值为：

y = 0.4 + 0.8 * x

运行此函数，我们得到的结果与预期的输出（y）值相当接近。

Expected=1.000, Predicted=1.200
Expected=3.000, Predicted=2.000
Expected=3.000, Predicted=3.600
Expected=2.000, Predicted=2.800
Expected=5.000, Predicted=4.400

现在我们准备实现随机梯度下降来优化我们的系数值。

2.估计系数

我们可以使用随机梯度下降来估计训练数据的系数值。

随机梯度下降需要两个参数：

• 学习率：用于限制每次更新时每个系数的校正量。
• 时期：更新系数时运行训练数据的次数。

这些以及训练数据将是该函数的参数。

我们需要在函数中执行 3 个循环：

1. 循环每个时代。
2. 循环遍历训练数据中的每一行以获得一个迭代。
3. 循环遍历每个系数并将其更新为一个迭代中的一行。

如您所见，我们更新训练数据中每一行的每个系数，每个时期。

系数根据模型产生的误差进行更新。该误差被计算为用候选系数进行的预测与预期输出值之间的差异。

error = prediction - expected

有一个系数可以对每个输入属性进行加权，并且这些系数以一致的方式更新，例如：

b1(t+1) = b1(t) - learning_rate * error(t) * x1(t)

列表开头的特殊系数（也称为截距或偏差）以类似的方式更新，除非没有输入，因为它与特定输入值无关：

b0(t+1) = b0(t) - learning_rate * error(t)

现在我们可以将所有这些放在一起。下面是一个名为 **coefficients_sgd（）**的函数，它使用随机梯度下降计算训练数据集的系数值。

# Estimate linear regression coefficients using stochastic gradient descent
def coefficients_sgd(train, l_rate, n_epoch):
	coef = [0.0 for i in range(len(train[0]))]
	for epoch in range(n_epoch):
		sum_error = 0
		for row in train:
			yhat = predict(row, coef)
			error = yhat - row[-1]
			sum_error += error**2
			coef[0] = coef[0] - l_rate * error
			for i in range(len(row)-1):
				coef[i + 1] = coef[i + 1] - l_rate * error * row[i]
		print('>epoch=%d, lrate=%.3f, error=%.3f' % (epoch, l_rate, sum_error))
	return coef

您可以看到，此外，我们会跟踪每个时期的平方误差（正值）的总和，以便我们可以在外部循环中打印出一条好消息。

我们可以在上面的同样小的人为数据集上测试这个函数。

# Make a prediction with coefficients
def predict(row, coefficients):
	yhat = coefficients[0]
	for i in range(len(row)-1):
		yhat += coefficients[i + 1] * row[i]
	return yhat

# Estimate linear regression coefficients using stochastic gradient descent
def coefficients_sgd(train, l_rate, n_epoch):
	coef = [0.0 for i in range(len(train[0]))]
	for epoch in range(n_epoch):
		sum_error = 0
		for row in train:
			yhat = predict(row, coef)
			error = yhat - row[-1]
			sum_error += error**2
			coef[0] = coef[0] - l_rate * error
			for i in range(len(row)-1):
				coef[i + 1] = coef[i + 1] - l_rate * error * row[i]
		print('>epoch=%d, lrate=%.3f, error=%.3f' % (epoch, l_rate, sum_error))
	return coef

# Calculate coefficients
dataset = [[1, 1], [2, 3], [4, 3], [3, 2], [5, 5]]
l_rate = 0.001
n_epoch = 50
coef = coefficients_sgd(dataset, l_rate, n_epoch)
print(coef)

我们使用 0.001 的小学习率并且将模型训练 50 个时期，或者将系数的 50 次曝光训练到整个训练数据集。

运行该示例在每个时期打印一条消息，该消息包含该时期的总和平方误差和最后一组系数。

>epoch=45, lrate=0.001, error=2.650
>epoch=46, lrate=0.001, error=2.627
>epoch=47, lrate=0.001, error=2.607
>epoch=48, lrate=0.001, error=2.589
>epoch=49, lrate=0.001, error=2.573
[0.22998234937311363, 0.8017220304137576]

你可以看到即使在最后一个时代，错误仍会继续下降。我们可以训练更长时间（更多迭代）或增加每个时期更新系数的量（更高的学习率）。

试验并看看你提出了什么。

现在，让我们将这个算法应用于真实数据集。

3.葡萄酒质量预测

在本节中，我们将使用随机梯度下降在葡萄酒质量数据集上训练线性回归模型。

该示例假定数据集的 CSV 副本位于当前工作目录中，文件名为 winequality-white.csv 。

首先加载数据集，将字符串值转换为数字，并将每列标准化为 0 到 1 范围内的值。这可以通过辅助函数 **load_csv（）**和 **str_column_to_float（）[来实现]。 HTG3]加载并准备数据集和 **dataset_minmax（）**和 **normalize_dataset（）来规范化它。

我们将使用 k-fold 交叉验证来估计学习模型在看不见的数据上的表现。这意味着我们将构建和评估 k 模型并将表现估计为平均模型误差。均方根误差将用于评估每个模型。这些行为在 cross_validation_split（）， **rmse_metric（）**和 **evaluate_algorithm（）**辅助函数中提供。

我们将使用上面创建的 predict（）， **coefficient_sgd（）**和 **linear_regression_sgd（）**函数来训练模型。

以下是完整的示例。

# Linear Regression With Stochastic Gradient Descent for Wine Quality
from random import seed
from random import randrange
from csv import reader
from math import sqrt

# Load a CSV file
def load_csv(filename):
	dataset = list()
	with open(filename, 'r') as file:
		csv_reader = reader(file)
		for row in csv_reader:
			if not row:
				continue
			dataset.append(row)
	return dataset

# Convert string column to float
def str_column_to_float(dataset, column):
	for row in dataset:
		row[column] = float(row[column].strip())

# Find the min and max values for each column
def dataset_minmax(dataset):
	minmax = list()
	for i in range(len(dataset[0])):
		col_values = [row[i] for row in dataset]
		value_min = min(col_values)
		value_max = max(col_values)
		minmax.append([value_min, value_max])
	return minmax

# Rescale dataset columns to the range 0-1
def normalize_dataset(dataset, minmax):
	for row in dataset:
		for i in range(len(row)):
			row[i] = (row[i] - minmax[i][0]) / (minmax[i][1] - minmax[i][0])

# Split a dataset into k folds
def cross_validation_split(dataset, n_folds):
	dataset_split = list()
	dataset_copy = list(dataset)
	fold_size = int(len(dataset) / n_folds)
	for i in range(n_folds):
		fold = list()
		while len(fold) < fold_size:
			index = randrange(len(dataset_copy))
			fold.append(dataset_copy.pop(index))
		dataset_split.append(fold)
	return dataset_split

# Calculate root mean squared error
def rmse_metric(actual, predicted):
	sum_error = 0.0
	for i in range(len(actual)):
		prediction_error = predicted[i] - actual[i]
		sum_error += (prediction_error ** 2)
	mean_error = sum_error / float(len(actual))
	return sqrt(mean_error)

# Evaluate an algorithm using a cross validation split
def evaluate_algorithm(dataset, algorithm, n_folds, *args):
	folds = cross_validation_split(dataset, n_folds)
	scores = list()
	for fold in folds:
		train_set = list(folds)
		train_set.remove(fold)
		train_set = sum(train_set, [])
		test_set = list()
		for row in fold:
			row_copy = list(row)
			test_set.append(row_copy)
			row_copy[-1] = None
		predicted = algorithm(train_set, test_set, *args)
		actual = [row[-1] for row in fold]
		rmse = rmse_metric(actual, predicted)
		scores.append(rmse)
	return scores

# Make a prediction with coefficients
def predict(row, coefficients):
	yhat = coefficients[0]
	for i in range(len(row)-1):
		yhat += coefficients[i + 1] * row[i]
	return yhat

# Estimate linear regression coefficients using stochastic gradient descent
def coefficients_sgd(train, l_rate, n_epoch):
	coef = [0.0 for i in range(len(train[0]))]
	for epoch in range(n_epoch):
		for row in train:
			yhat = predict(row, coef)
			error = yhat - row[-1]
			coef[0] = coef[0] - l_rate * error
			for i in range(len(row)-1):
				coef[i + 1] = coef[i + 1] - l_rate * error * row[i]
			# print(l_rate, n_epoch, error)
	return coef

# Linear Regression Algorithm With Stochastic Gradient Descent
def linear_regression_sgd(train, test, l_rate, n_epoch):
	predictions = list()
	coef = coefficients_sgd(train, l_rate, n_epoch)
	for row in test:
		yhat = predict(row, coef)
		predictions.append(yhat)
	return(predictions)

# Linear Regression on wine quality dataset
seed(1)
# load and prepare data
filename = 'winequality-white.csv'
dataset = load_csv(filename)
for i in range(len(dataset[0])):
	str_column_to_float(dataset, i)
# normalize
minmax = dataset_minmax(dataset)
normalize_dataset(dataset, minmax)
# evaluate algorithm
n_folds = 5
l_rate = 0.01
n_epoch = 50
scores = evaluate_algorithm(dataset, linear_regression_sgd, n_folds, l_rate, n_epoch)
print('Scores: %s' % scores)
print('Mean RMSE: %.3f' % (sum(scores)/float(len(scores))))

k 值为 5 用于交叉验证，每次迭代时评估每个折叠 4,898 / 5 = 979.6 或略低于 1000 个记录。通过一些实验选择了 0.01 和 50 个训练时期的学习率。

您可以尝试自己的配置，看看是否可以打败我的分数。

运行此示例将打印 5 个交叉验证折叠中每个折叠的分数，然后打印平均 RMSE。

我们可以看到 RMSE（在标准化数据集上）是 0.126，如果我们只是预测平均值（使用零规则算法），则低于 0.148 的基线值。

Scores: [0.12248058224159092, 0.13034017509167112, 0.12620370547483578, 0.12897687952843237, 0.12446990678682233]
Mean RMSE: 0.126

扩展

本节列出了本教程的一些扩展，您可能希望考虑这些扩展。

• 调整示例。调整学习率，时期数，甚至数据准备方法，以获得葡萄酒质量数据集的改进分数。
• 批随机梯度下降。改变随机梯度下降算法以在每个时期累积更新，并且仅在时期结束时批量更新系数。
• 其他回归问题。将该技术应用于 UCI 机器学习库中的其他回归问题。

你有没有探索过这些扩展？ 请在下面的评论中告诉我。

评论

在本教程中，您了解了如何使用 Python 从零开始使用随机梯度下降来实现线性回归。

你学到了

• 如何预测多元线性回归问题。
• 如何使用随机梯度下降来优化一组系数。
• 如何将该技术应用于真实的回归预测性建模问题。

你有什么问题吗？ 在下面的评论中提出您的问题，我会尽力回答。