Machine-Learning-Mastery-深度学习时间序列教程-十一-

布客飞龙
56 阅读1小时+

Machine Learning Mastery 深度学习时间序列教程（十一）

原文：Machine Learning Mastery

协议：CC BY-NC-SA 4.0

如何加载，可视化和探索复杂的多变量多步时间序列预测数据集

原文： machinelearningmastery.com/how-to-load-visualize-and-explore-a-complex-multivariate-multistep-time-series-forecasting-dataset/

实时世界时间序列预测具有挑战性，其原因不仅限于问题特征，例如具有多个输入变量，需要预测多个时间步骤，以及需要对多个物理站点执行相同类型的预测。

EMC Data Science Global Hackathon 数据集或简称“_ 空气质量预测 _”数据集描述了多个站点的天气状况，需要预测随后三天的空气质量测量结果。

在本教程中，您将发现并探索空气质量预测数据集，该数据集代表具有挑战性的多变量，多站点和多步骤时间序列预测问题。

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

  • 如何加载和探索数据集的块结构。
  • 如何探索和可视化数据集的输入和目标变量。
  • 如何使用新的理解来概述一系列方法来构建问题，准备数据和建模数据集。

让我们开始吧。

  • Update Apr / 2019 ：修正了计算总缺失值的错误（感谢 zhangzhe）。

How to Load, Visualize, and Explore a Complex Multivariate Multistep Time Series Forecasting Dataset

如何加载，可视化和探索复杂的多变量多步时间序列预测数据集 照片由 H Matthew Howarth ，保留一些权利。

教程概述

本教程分为七个部分;他们是：

  1. 问题描述
  2. 加载数据集
  3. 块数据结构
  4. 输入变量
  5. 目标变量
  6. 与目标变量的皱痕
  7. 关于建模的思考

问题描述

EMC Data Science Global Hackathon 数据集或简称“_ 空气质量预测 _”数据集描述了多个站点的天气状况，需要预测随后三天的空气质量测量结果。

具体而言，对于多个站点，每小时提供 8 天的温度，压力，风速和风向等天气观测。目标是预测未来三天在多个地点的空气质量测量。预测的提前期不是连续的;相反，必须在 72 小时预测期内预测具体的提前期;他们是：

+1, +2, +3, +4, +5, +10, +17, +24, +48, +72

此外，数据集被划分为不相交但连续的数据块，其中 8 天的数据随后是需要预测的 3 天。

并非所有站点或块都可以获得所有观察结果，并且并非所有站点和块都可以使用所有输出变量。必须解决大部分缺失数据。

该数据集被用作 2012 年 Kaggle 网站上短期机器学习竞赛（或黑客马拉松）的基础。

根据从参与者中扣留的真实观察结果评估竞赛的提交，并使用平均绝对误差（MAE）进行评分。提交要求在由于缺少数据而无法预测的情况下指定-1,000,000 的值。实际上，提供了一个插入缺失值的模板，并且要求所有提交都采用（模糊的是什么）。

获胜者在滞留测试集（私人排行榜）上使用随机森林在滞后观察中获得了 0.21058 的 MAE。该帖子中提供了此解决方案的说明：

在本教程中，我们将探索此数据集，以便更好地了解预测问题的性质，并提出如何建模的方法。

加载数据集

第一步是下载数据集并将其加载到内存中。

数据集可以从 Kaggle 网站免费下载。您可能必须创建一个帐户并登录才能下载数据集。

下载整个数据集，例如“_ 将所有 _”下载到您的工作站并使用名为'AirQualityPrediction'的文件夹解压缩当前工作目录中的存档

你应该在 AirQualityPrediction / 文件夹中有五个文件;他们是：

  • SiteLocations.csv
  • SiteLocations_with_more_sites.csv
  • SubmissionZerosExceptNAs.csv
  • TrainingData.csv
  • sample_code.r

我们的重点将是包含训练数据集的'TrainingData.csv'，特别是块中的数据，其中每个块是八个连续的观察日和目标变量。

在撰写本文时，测试数据集（每个块的剩余三天）不可用于此数据集。

打开'TrainingData.csv'文件并查看内容。解压缩的数据文件相对较小（21 兆字节），很容易适应 RAM。

查看文件的内容，我们可以看到数据文件包含标题行。

我们还可以看到缺失数据标有'NA'值，Pandas 将自动转换为NumPy.NaN。

我们可以看到'_ 工作日 _'列包含作为字符串的日期，而所有其他数据都是数字。

以下是数据文件的前几行供参考。

"rowID","chunkID","position_within_chunk","month_most_common","weekday","hour","Solar.radiation_64","WindDirection..Resultant_1","WindDirection..Resultant_1018","WindSpeed..Resultant_1","WindSpeed..Resultant_1018","Ambient.Max.Temperature_14","Ambient.Max.Temperature_22","Ambient.Max.Temperature_50","Ambient.Max.Temperature_52","Ambient.Max.Temperature_57","Ambient.Max.Temperature_76","Ambient.Max.Temperature_2001","Ambient.Max.Temperature_3301","Ambient.Max.Temperature_6005","Ambient.Min.Temperature_14","Ambient.Min.Temperature_22","Ambient.Min.Temperature_50","Ambient.Min.Temperature_52","Ambient.Min.Temperature_57","Ambient.Min.Temperature_76","Ambient.Min.Temperature_2001","Ambient.Min.Temperature_3301","Ambient.Min.Temperature_6005","Sample.Baro.Pressure_14","Sample.Baro.Pressure_22","Sample.Baro.Pressure_50","Sample.Baro.Pressure_52","Sample.Baro.Pressure_57","Sample.Baro.Pressure_76","Sample.Baro.Pressure_2001","Sample.Baro.Pressure_3301","Sample.Baro.Pressure_6005","Sample.Max.Baro.Pressure_14","Sample.Max.Baro.Pressure_22","Sample.Max.Baro.Pressure_50","Sample.Max.Baro.Pressure_52","Sample.Max.Baro.Pressure_57","Sample.Max.Baro.Pressure_76","Sample.Max.Baro.Pressure_2001","Sample.Max.Baro.Pressure_3301","Sample.Max.Baro.Pressure_6005","Sample.Min.Baro.Pressure_14","Sample.Min.Baro.Pressure_22","Sample.Min.Baro.Pressure_50","Sample.Min.Baro.Pressure_52","Sample.Min.Baro.Pressure_57","Sample.Min.Baro.Pressure_76","Sample.Min.Baro.Pressure_2001","Sample.Min.Baro.Pressure_3301","Sample.Min.Baro.Pressure_6005","target_1_57","target_10_4002","target_10_8003","target_11_1","target_11_32","target_11_50","target_11_64","target_11_1003","target_11_1601","target_11_4002","target_11_8003","target_14_4002","target_14_8003","target_15_57","target_2_57","target_3_1","target_3_50","target_3_57","target_3_1601","target_3_4002","target_3_6006","target_4_1","target_4_50","target_4_57","target_4_1018","target_4_1601","target_4_2001","target_4_4002","target_4_4101","target_4_6006","target_4_8003","target_5_6006","target_7_57","target_8_57","target_8_4002","target_8_6004","target_8_8003","target_9_4002","target_9_8003"
1,1,1,10,"Saturday",21,0.01,117,187,0.3,0.3,NA,NA,NA,14.9,NA,NA,NA,NA,NA,NA,NA,NA,5.8,NA,NA,NA,NA,NA,NA,NA,NA,747,NA,NA,NA,NA,NA,NA,NA,NA,750,NA,NA,NA,NA,NA,NA,NA,NA,743,NA,NA,NA,NA,NA,2.67923294292042,6.1816228132982,NA,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,NA,2.38965627997991,NA,5.56815355612325,0.690015329704154,NA,NA,NA,NA,NA,NA,2.84349016287551,0.0920223353681394,1.69321097077376,0.368089341472558,0.184044670736279,0.368089341472558,0.276067006104418,0.892616653070952,1.74842437199465,NA,NA,5.1306307034019,1.34160578423204,2.13879182993514,3.01375212399952,NA,5.67928016629218,NA
2,1,2,10,"Saturday",22,0.01,231,202,0.5,0.6,NA,NA,NA,14.9,NA,NA,NA,NA,NA,NA,NA,NA,5.8,NA,NA,NA,NA,NA,NA,NA,NA,747,NA,NA,NA,NA,NA,NA,NA,NA,750,NA,NA,NA,NA,NA,NA,NA,NA,743,NA,NA,NA,NA,NA,2.67923294292042,8.47583334194495,NA,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,NA,1.99138023331659,NA,5.56815355612325,0.923259948195698,NA,NA,NA,NA,NA,NA,3.1011527019063,0.0920223353681394,1.94167127626774,0.368089341472558,0.184044670736279,0.368089341472558,0.368089341472558,1.73922213845783,2.14412041407765,NA,NA,5.1306307034019,1.19577906855465,2.72209869264472,3.88871241806389,NA,7.42675098668978,NA
3,1,3,10,"Saturday",23,0.01,247,227,0.5,1.5,NA,NA,NA,14.9,NA,NA,NA,NA,NA,NA,NA,NA,5.8,NA,NA,NA,NA,NA,NA,NA,NA,747,NA,NA,NA,NA,NA,NA,NA,NA,750,NA,NA,NA,NA,NA,NA,NA,NA,743,NA,NA,NA,NA,NA,2.67923294292042,8.92192983362627,NA,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,NA,1.7524146053186,NA,5.56815355612325,0.680296803933673,NA,NA,NA,NA,NA,NA,3.06434376775904,0.0920223353681394,2.52141198908702,0.460111676840697,0.184044670736279,0.368089341472558,0.368089341472558,1.7852333061419,1.93246904273093,NA,NA,5.13639545700122,1.40965825154816,3.11096993445111,3.88871241806389,NA,7.68373198968942,NA
4,1,4,10,"Sunday",0,0.01,219,218,0.2,1.2,NA,NA,NA,14,NA,NA,NA,NA,NA,NA,NA,NA,4.8,NA,NA,NA,NA,NA,NA,NA,NA,751,NA,NA,NA,NA,NA,NA,NA,NA,754,NA,NA,NA,NA,NA,NA,NA,NA,748,NA,NA,NA,NA,NA,2.67923294292042,5.09824561921501,NA,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,0.114975168664303,NA,2.38965627997991,NA,5.6776192223642,0.612267123540305,NA,NA,NA,NA,NA,NA,3.21157950434806,0.184044670736279,2.374176252498,0.460111676840697,0.184044670736279,0.368089341472558,0.276067006104418,1.86805340797323,2.08890701285676,NA,NA,5.21710200739181,1.47771071886428,2.04157401948354,3.20818774490271,NA,4.83124285639335,NA
...

我们可以使用 Pandas read_csv（）函数将数据文件加载到内存中，并在第 0 行指定标题行。

# load dataset
dataset = read_csv('AirQualityPrediction/TrainingData.csv', header=0)

我们还可以快速了解数据集中有多少缺失数据。我们可以通过首先修剪前几列来删除字符串工作日数据并将剩余列转换为浮点值来实现。

# trim and transform to floats
values = dataset.values
data = values[:, 6:].astype('float32')

然后，我们可以计算缺失观测的总数和缺失的值的百分比。

# summarize amount of missing data
total_missing = count_nonzero(isnan(data))
percent_missing = total_missing / data.size * 100
print('Total Missing: %d/%d (%.1f%%)' % (total_missing, data.size, percent_missing))

下面列出了完整的示例。

# load dataset
from numpy import isnan
from numpy import count_nonzero
from pandas import read_csv
# load dataset
dataset = read_csv('AirQualityPrediction/TrainingData.csv', header=0)
# summarize
print(dataset.shape)
# trim and transform to floats
values = dataset.values
data = values[:, 6:].astype('float32')
# summarize amount of missing data
total_missing = count_nonzero(isnan(data))
percent_missing = total_missing / data.size * 100
print('Total Missing: %d/%d (%.1f%%)' % (total_missing, data.size, percent_missing))

首先运行该示例将打印已加载数据集的形状。

我们可以看到我们有大约 37,000 行和 95 列。我们知道这些数字是误导的，因为数据实际上被分成块，并且列被分成不同站点的相同观察。

我们还可以看到有超过 40％的数据丢失。这很多。数据非常不完整，在建模问题之前我们必须很好地理解这一点。

(37821, 95)
Total Missing: 1922092/3366069 (57.1%)

块数据结构

一个很好的起点是根据块来查看数据。

块持续时间

我们可以通过'chunkID'变量（列索引 1）对数据进行分组。

如果每个块是 8 天并且观察是每小时的，那么我们期望每块（8 * 24）或 192 行数据。

如果有 37,821 行数据，则必须有大于或小于 192 小时的块，因为 37,821 / 192 大约是 196.9 块。

我们首先将数据拆分成块。我们可以先获得唯一的块标识符列表。

chunk_ids = unique(values[:, 1])

然后，我们可以收集每个块标识符的所有行，并将它们存储在字典中以便于访问。

chunks = dict()
# sort rows by chunk id
for chunk_id in chunk_ids:
	selection = values[:, chunk_ix] == chunk_id
	chunks[chunk_id] = values[selection, :]

下面定义了一个名为to_chunks()的函数，它接受加载数据的 NumPy 数组，并将chunk_id的字典返回到块的行。

# split the dataset by 'chunkID', return a dict of id to rows
def to_chunks(values, chunk_ix=1):
	chunks = dict()
	# get the unique chunk ids
	chunk_ids = unique(values[:, chunk_ix])
	# group rows by chunk id
	for chunk_id in chunk_ids:
		selection = values[:, chunk_ix] == chunk_id
		chunks[chunk_id] = values[selection, :]
	return chunks

数据文件中的'position_within_chunk'表示块内行的顺序。在这个阶段，我们假设行已经排序，不需要排序。原始数据文件的略读似乎证实了这一假设。

一旦数据被分类成块，我们就可以计算每个块中的行数并查看分布，例如使用框和胡须图。

# plot distribution of chunk durations
def plot_chunk_durations(chunks):
	# chunk durations in hours
	chunk_durations = [len(v) for k,v in chunks.items()]
	# boxplot
	pyplot.subplot(2, 1, 1)
	pyplot.boxplot(chunk_durations)
	# histogram
	pyplot.subplot(2, 1, 2)
	pyplot.hist(chunk_durations)
	# histogram
	pyplot.show()

下面列出了将所有这些联系在一起的完整示例

# split data into chunks
from numpy import unique
from pandas import read_csv
from matplotlib import pyplot

# split the dataset by 'chunkID', return a dict of id to rows
def to_chunks(values, chunk_ix=1):
	chunks = dict()
	# get the unique chunk ids
	chunk_ids = unique(values[:, chunk_ix])
	# group rows by chunk id
	for chunk_id in chunk_ids:
		selection = values[:, chunk_ix] == chunk_id
		chunks[chunk_id] = values[selection, :]
	return chunks

# plot distribution of chunk durations
def plot_chunk_durations(chunks):
	# chunk durations in hours
	chunk_durations = [len(v) for k,v in chunks.items()]
	# boxplot
	pyplot.subplot(2, 1, 1)
	pyplot.boxplot(chunk_durations)
	# histogram
	pyplot.subplot(2, 1, 2)
	pyplot.hist(chunk_durations)
	# histogram
	pyplot.show()

# load dataset
dataset = read_csv('AirQualityPrediction/TrainingData.csv', header=0)
# group data by chunks
values = dataset.values
chunks = to_chunks(values)
print('Total Chunks: %d' % len(chunks))
# plot chunk durations
plot_chunk_durations(chunks)

首先运行该示例将打印数据集中的块数。

我们可以看到有 208，这表明每小时观察的数量确实必须在各个块之间变化。

Total Chunks: 208

创建框和胡须图以及块持续时间的直方图。我们可以看到中位数确实是 192，这意味着大多数块有八天的观察或接近它。

我们还可以看到持续时间长达约 25 行的长尾。虽然这些案例并不多，但鉴于缺乏数据，我们预计这将是一个挑战。

该分布还提出了关于每个块内的观察结果可能是多么连续的问题。

Box and whisker plot and Histogram plot of chunk duration in hours

盒子和须状图以及以小时为单位的块长度的直方图

大块连续性

了解在没有完整八天数据的那些块中观察是否连续（或不连续）可能会有所帮助。

考虑这一点的一种方法是为每个不连续的块创建线图并显示观察中的间隙。

我们可以在一个地块上做到这一点。每个块都有一个唯一的标识符，从 1 到 208，我们可以使用它作为序列的值，并通过NaN值在 8 天间隔内标记缺失的观察值，这些值不会出现在图上。

反过来说，我们可以假设我们对一个块中的所有时间步都有 NaN 值，然后使用'position_within_chunk'列（索引 2）来确定具有值的时间步长并用它们标记它们。块 ID。

下面的plot_discontinuous_chunks()实现了这种行为，在同一个图上为每个缺少行的块创建一个系列或行。期望的是，突破线将帮助我们看到这些不完整的块是多么连续或不连续。

# plot chunks that do not have all data
def plot_discontiguous_chunks(chunks, row_in_chunk_ix=2):
	n_steps = 8 * 24
	for c_id,rows in chunks.items():
		# skip chunks with all data
		if rows.shape[0] == n_steps:
			continue
		# create empty series
		series = [nan for _ in range(n_steps)]
		# mark all rows with data
		for row in rows:
			# convert to zero offset
			r_id = row[row_in_chunk_ix] - 1
			# mark value
			series[r_id] = c_id
		# plot
		pyplot.plot(series)
	pyplot.show()

下面列出了完整的示例。

# plot discontiguous chunks
from numpy import nan
from numpy import unique
from pandas import read_csv
from matplotlib import pyplot

# split the dataset by 'chunkID', return a dict of id to rows
def to_chunks(values, chunk_ix=1):
	chunks = dict()
	# get the unique chunk ids
	chunk_ids = unique(values[:, chunk_ix])
	# group rows by chunk id
	for chunk_id in chunk_ids:
		selection = values[:, chunk_ix] == chunk_id
		chunks[chunk_id] = values[selection, :]
	return chunks

# plot chunks that do not have all data
def plot_discontiguous_chunks(chunks, row_in_chunk_ix=2):
	n_steps = 8 * 24
	for c_id,rows in chunks.items():
		# skip chunks with all data
		if rows.shape[0] == n_steps:
			continue
		# create empty series
		series = [nan for _ in range(n_steps)]
		# mark all rows with data
		for row in rows:
			# convert to zero offset
			r_id = row[row_in_chunk_ix] - 1
			# mark value
			series[r_id] = c_id
		# plot
		pyplot.plot(series)
	pyplot.show()

# load dataset
dataset = read_csv('AirQualityPrediction/TrainingData.csv', header=0)
# group data by chunks
values = dataset.values
chunks = to_chunks(values)
# plot discontiguous chunks
plot_discontiguous_chunks(chunks)

运行该示例会为每个缺少数据的块创建一个带有一行的图形。

每个块的行中断的数量和长度给出了每个块中的观察结果是多么不连续的概念。

许多块确实有很长的连续数据，这是建模的一个好兆头。

在某些情况下，块的观察结果非常少，而且存在的观察结果是小的连续斑块。模型可能具有挑战性。

此外，并非所有这些块都在块的末尾有观察结果：需要预测之前的时间段。对于那些寻求坚持最近观察的模型而言，这些将是一个挑战。

块内系列数据的不连续性也将使评估模型具有挑战性。例如，人们不能简单地将块数据分成两半，在前半部分进行训练，在观察结果不完整时对第二部分进行测试。至少，当考虑不完整的块数据时。

Line plots of chunks with discontinuous observations

具有不连续观察的块的线图

块内的每日覆盖率

块的不连续性也表明，查看每个块所覆盖的小时数可能很重要。

一天中的时间在环境数据中很重要，并且假设每个块包含相同的每日或每周周期的模型可能会发生绊倒，如果一天的开始和结束时间不同。

我们可以通过绘制每个块的第一个小时（每天 24 小时）的分布来快速检查这一点。

柱状图中的区间数设置为 24，因此我们可以清楚地看到 24 小时内每天每小时的分布。

此外，当收集块的第一个小时时，我们小心只从那些具有所有八天数据的块中收集它，以防丢失数据的块在块的开头没有观察，我们知道发生。

# plot distribution of chunk start hour
def plot_chunk_start_hour(chunks, hour_in_chunk_ix=5):
	# chunk start hour
	chunk_start_hours = [v[0, hour_in_chunk_ix] for k,v in chunks.items() if len(v)==192]
	# boxplot
	pyplot.subplot(2, 1, 1)
	pyplot.boxplot(chunk_start_hours)
	# histogram
	pyplot.subplot(2, 1, 2)
	pyplot.hist(chunk_start_hours, bins=24)
	# histogram
	pyplot.show()

下面列出了完整的示例。

# plot distribution of chunk start hour
from numpy import nan
from numpy import unique
from pandas import read_csv
from matplotlib import pyplot

# split the dataset by 'chunkID', return a dict of id to rows
def to_chunks(values, chunk_ix=1):
	chunks = dict()
	# get the unique chunk ids
	chunk_ids = unique(values[:, chunk_ix])
	# group rows by chunk id
	for chunk_id in chunk_ids:
		selection = values[:, chunk_ix] == chunk_id
		chunks[chunk_id] = values[selection, :]
	return chunks

# plot distribution of chunk start hour
def plot_chunk_start_hour(chunks, hour_in_chunk_ix=5):
	# chunk start hour
	chunk_start_hours = [v[0, hour_in_chunk_ix] for k,v in chunks.items() if len(v)==192]
	# boxplot
	pyplot.subplot(2, 1, 1)
	pyplot.boxplot(chunk_start_hours)
	# histogram
	pyplot.subplot(2, 1, 2)
	pyplot.hist(chunk_start_hours, bins=24)
	# histogram
	pyplot.show()

# load dataset
dataset = read_csv('AirQualityPrediction/TrainingData.csv', header=0)
# group data by chunks
values = dataset.values
chunks = to_chunks(values)
# plot distribution of chunk start hour
plot_chunk_start_hour(chunks)

运行该示例将创建一个框和胡须图以及每个块中第一个小时的直方图。

我们可以看到当天 24 小时内合理均匀的开始时间分布。

此外，这意味着每个块的预测间隔也将在 24 小时内变化。这为可能期望标准三天预测期（午夜至午夜）的模型增加了皱纹。

Distribution of first hour in observations within each chunk

每个块内观察的第一个小时的分布

现在我们已经了解了数据的块结构，让我们仔细研究描述气象观测的输入变量。

输入变量

有 56 个输入变量。

前六个（索引 0 到 5）是块的元数据信息和观察的时间。他们是：

rowID
chunkID
position_within_chunk
month_most_common
weekday
hour

其余 50 个描述了特定地点的气象信息;他们是：

Solar.radiation_64
WindDirection..Resultant_1
WindDirection..Resultant_1018
WindSpeed..Resultant_1
WindSpeed..Resultant_1018
Ambient.Max.Temperature_14
Ambient.Max.Temperature_22
Ambient.Max.Temperature_50
Ambient.Max.Temperature_52
Ambient.Max.Temperature_57
Ambient.Max.Temperature_76
Ambient.Max.Temperature_2001
Ambient.Max.Temperature_3301
Ambient.Max.Temperature_6005
Ambient.Min.Temperature_14
Ambient.Min.Temperature_22
Ambient.Min.Temperature_50
Ambient.Min.Temperature_52
Ambient.Min.Temperature_57
Ambient.Min.Temperature_76
Ambient.Min.Temperature_2001
Ambient.Min.Temperature_3301
Ambient.Min.Temperature_6005
Sample.Baro.Pressure_14
Sample.Baro.Pressure_22
Sample.Baro.Pressure_50
Sample.Baro.Pressure_52
Sample.Baro.Pressure_57
Sample.Baro.Pressure_76
Sample.Baro.Pressure_2001
Sample.Baro.Pressure_3301
Sample.Baro.Pressure_6005
Sample.Max.Baro.Pressure_14
Sample.Max.Baro.Pressure_22
Sample.Max.Baro.Pressure_50
Sample.Max.Baro.Pressure_52
Sample.Max.Baro.Pressure_57
Sample.Max.Baro.Pressure_76
Sample.Max.Baro.Pressure_2001
Sample.Max.Baro.Pressure_3301
Sample.Max.Baro.Pressure_6005
Sample.Min.Baro.Pressure_14
Sample.Min.Baro.Pressure_22
Sample.Min.Baro.Pressure_50
Sample.Min.Baro.Pressure_52
Sample.Min.Baro.Pressure_57
Sample.Min.Baro.Pressure_76
Sample.Min.Baro.Pressure_2001
Sample.Min.Baro.Pressure_3301
Sample.Min.Baro.Pressure_6005

真的，只有八个气象输入变量：

Solar.radiation
WindDirection..Resultant
WindSpeed..Resultant
Ambient.Max.Temperature
Ambient.Min.Temperature
Sample.Baro.Pressure
Sample.Max.Baro.Pressure
Sample.Min.Baro.Pressure

这些变量记录在 23 个独特的站点中;他们是：

1, 14, 22, 50, 52, 57, 64, 76, 1018, 2001, 3301, 6005

数据非常复杂。

并非所有变量都记录在所有站点上。

目标变量中使用的站点标识符有一些重叠，例如 1,50,64 等。

目标变量中使用的站点标识符未在输入变量中使用，例如 4002.还有在输入中使用的站点标识符，这些站点标识符未在目标标识符中使用，例如 15。

这表明，至少并非所有变量都记录在所有位置。这些录制站在各个站点之间是异构的。此外，对于仅收集给定类型的度量或收集所有度量的站点，可能存在一些特殊情况。

让我们仔细看看输入变量的数据。

块的输入的时间结构

我们可以从查看每个块的输入结构和分布开始。

所有八天观察的前几个块的块大小为 1,3 和 5。

我们可以枚举所有输入列并为每个输入列创建一个线图。这将为每个输入变量创建一个时间序列线图，以便大致了解每个输入变量的变化情况。

我们可以针对几个块重复这一点，以了解时间结构如何在块之间有所不同。

下面名为plot_chunk_inputs()的函数获取块格式的数据和要绘制的块 ID 列表。它将创建一个包含 50 个线图的图形，每个输入变量一个，每个图块 n 行，每个块一个。

# plot all inputs for one or more chunk ids
def plot_chunk_inputs(chunks, c_ids):
	pyplot.figure()
	inputs = range(6, 56)
	for i in range(len(inputs)):
		ax = pyplot.subplot(len(inputs), 1, i+1)
		ax.set_xticklabels([])
		ax.set_yticklabels([])
		column = inputs[i]
		for chunk_id in c_ids:
			rows = chunks[chunk_id]
			pyplot.plot(rows[:,column])
	pyplot.show()

下面列出了完整的示例。

# plot inputs for a chunk
from numpy import unique
from pandas import read_csv
from matplotlib import pyplot

# split the dataset by 'chunkID', return a dict of id to rows
def to_chunks(values, chunk_ix=1):
	chunks = dict()
	# get the unique chunk ids
	chunk_ids = unique(values[:, chunk_ix])
	# group rows by chunk id
	for chunk_id in chunk_ids:
		selection = values[:, chunk_ix] == chunk_id
		chunks[chunk_id] = values[selection, :]
	return chunks

# plot all inputs for one or more chunk ids
def plot_chunk_inputs(chunks, c_ids):
	pyplot.figure()
	inputs = range(6, 56)
	for i in range(len(inputs)):
		ax = pyplot.subplot(len(inputs), 1, i+1)
		ax.set_xticklabels([])
		ax.set_yticklabels([])
		column = inputs[i]
		for chunk_id in c_ids:
			rows = chunks[chunk_id]
			pyplot.plot(rows[:,column])
	pyplot.show()

# load data
dataset = read_csv('AirQualityPrediction/TrainingData.csv', header=0)
# group data by chunks
values = dataset.values
chunks = to_chunks(values)
# plot inputs for some chunks
plot_chunk_inputs(chunks, [1])

运行该示例将创建一个包含 50 个线图的单个图形，每个图形用于每个气象输入变量。

这些图很难看，因此您可能希望增加所创建图形的大小。

我们可以看到前五个变量的观察看起来非常完整;这些是太阳辐射，风速和风向。其余变量看起来非常不完整，至少对于这个块而言。

Parallel Time Series Line Plots For All Input Variables for 1 Chunks

1 个块的所有输入变量的并行时间序列线图

我们可以更新示例并绘制前三个块的输入变量，并进行完整的八天观察。

plot_chunk_inputs(chunks, [1, 3 ,5])

运行该示例将创建相同的 50 个线图，每个图每个图有三个系列或线，每个块一个。

同样，该图使单个图很难看到，因此您可能需要增加图的大小以更好地查看模式。

我们可以看到这三个图在每个线图中都显示出类似的结构。这是有帮助的发现，因为它表明跨多个块建模相同的变量可能是有用的。

Parallel Time Series Line Plots For All Input Variables for 3 Chunks

3 个块的所有输入变量的并行时间序列线图

它确实提出了关于变量的分布是否在不同站点之间差异很大的问题。

输入数据分布

我们可以使用 box 和 whisker 图粗略地查看输入变量的分布。

下面的plot_chunk_input_boxplots()将为每个输入要素创建一个盒子和胡须，用于一个块的数据。

# boxplot for input variables for a chuck
def plot_chunk_input_boxplots(chunks, c_id):
	rows = chunks[c_id]
	pyplot.boxplot(rows[:,6:56])
	pyplot.show()

下面列出了完整的示例。

# boxplots of inputs for a chunk
from numpy import unique
from numpy import isnan
from numpy import count_nonzero
from pandas import read_csv
from matplotlib import pyplot

# split the dataset by 'chunkID', return a dict of id to rows
def to_chunks(values, chunk_ix=1):
	chunks = dict()
	# get the unique chunk ids
	chunk_ids = unique(values[:, chunk_ix])
	# group rows by chunk id
	for chunk_id in chunk_ids:
		selection = values[:, chunk_ix] == chunk_id
		chunks[chunk_id] = values[selection, :]
	return chunks

# boxplot for input variables for a chuck
def plot_chunk_input_boxplots(chunks, c_id):
	rows = chunks[c_id]
	pyplot.boxplot(rows[:,6:56])
	pyplot.show()

# load data
dataset = read_csv('TrainingData.csv', header=0)
# group data by chunks
values = dataset.values
chunks = to_chunks(values)
# boxplot for input variables
plot_chunk_input_boxplots(chunks, 1)

运行该示例将创建 50 个箱图，每个输入变量用于训练数据集中第一个块中的观察。

我们可以看到相同类型的变量可能具有相同的观察范围，并且每组变量似乎具有不同的单位。也许是风向的度数，压力的百帕斯卡，温度的摄氏度等等。

Box and whisker plots of input variables for one chunk

一个块的输入变量的框和胡须图

对于八种变量类型中的每一种，进一步研究观察的分布和扩散可能是有趣的。这是一个进一步的练习。

我们对输入变量有一些粗略的想法，也许它们可能对预测目标变量很有用。我们无法确定。

我们现在可以将注意力转向目标变量。

目标变量

预测问题的目标是预测多个站点的多个变量，为期三天。

有 39 个时间序列变量可供预测。

从列标题中，它们是：

"target_1_57","target_10_4002","target_10_8003","target_11_1","target_11_32","target_11_50","target_11_64","target_11_1003","target_11_1601","target_11_4002","target_11_8003","target_14_4002","target_14_8003","target_15_57","target_2_57","target_3_1","target_3_50","target_3_57","target_3_1601","target_3_4002","target_3_6006","target_4_1","target_4_50","target_4_57","target_4_1018","target_4_1601","target_4_2001","target_4_4002","target_4_4101","target_4_6006","target_4_8003","target_5_6006","target_7_57","target_8_57","target_8_4002","target_8_6004","target_8_8003","target_9_4002","target_9_8003"

这些列标题的命名约定是：

target_[variable identifier]_[site identifier]]

我们可以使用一点 regex 将这些列标题转换为变量 id 和 site id 的小数据集。

结果如下：

var, site
1,57
10,4002
10,8003
11,1
11,32
11,50
11,64
11,1003
11,1601
11,4002
11,8003
14,4002
14,8003
15,57
2,57
3,1
3,50
3,57
3,1601
3,4002
3,6006
4,1
4,50
4,57
4,1018
4,1601
4,2001
4,4002
4,4101
4,6006
4,8003
5,6006
7,57
8,57
8,4002
8,6004
8,8003
9,4002
9,8003

有用的是，目标按变量 id 分组。

我们可以看到，可能需要跨多个站点预测一个变量;例如，在站点 1,32,50 等处预测的变量 11，等等：

var, site
11,1
11,32
11,50
11,64
11,1003
11,1601
11,4002
11,8003

我们可以看到，对于给定的站点，可能需要预测不同的变量。例如，站点 50 需要变量 11,3 和 4：

var, site
11,50
3,50
4,50

我们可以将目标的小数据集保存到名为“targets.txt”的文件中并加载它以进行快速分析。

# summarize targets
from numpy import unique
from pandas import read_csv
# load dataset
dataset = read_csv('targets.txt', header=0)
values = dataset.values
# summarize unique
print('Unique Variables: %d' % len(unique(values[:, 0])))
print('Unique Sites: %d' % len(unique(values[:, 1])))

运行该示例将打印唯一变量和站点的数量。

如果我们预测所有站点的所有变量，我们可以看到 39 个目标变量远小于（12 * 14）168。

Unique Variables: 12
Unique Sites: 14

让我们仔细看看目标变量的数据。

大块目标的时间结构

我们可以从每个块的目标结构和分布开始。

所有八天观察的前几个块的块大小为 1,3 和 5。

我们可以枚举所有目标列，并为每个列创建一个线图。这将为每个目标变量创建一个时间序列线图，以大致了解它如何随时间变化。

我们可以针对几个块重复这一点，以大致了解时间结构如何在块之间变化。

下面的函数名为 plot_chunk_targets（），以块格式和块 ID 列表绘制。它将创建一个包含 39 个线图的图形，每个目标变量一个，每个图块 n 行，每个块一个。

# plot all targets for one or more chunk ids
def plot_chunk_targets(chunks, c_ids):
	pyplot.figure()
	targets = range(56, 95)
	for i in range(len(targets)):
		ax = pyplot.subplot(len(targets), 1, i+1)
		ax.set_xticklabels([])
		ax.set_yticklabels([])
		column = targets[i]
		for chunk_id in c_ids:
			rows = chunks[chunk_id]
			pyplot.plot(rows[:,column])
	pyplot.show()

下面列出了完整的示例。

# plot targets for a chunk
from numpy import unique
from pandas import read_csv
from matplotlib import pyplot

# split the dataset by 'chunkID', return a dict of id to rows
def to_chunks(values, chunk_ix=1):
	chunks = dict()
	# get the unique chunk ids
	chunk_ids = unique(values[:, chunk_ix])
	# group rows by chunk id
	for chunk_id in chunk_ids:
		selection = values[:, chunk_ix] == chunk_id
		chunks[chunk_id] = values[selection, :]
	return chunks

# plot all targets for one or more chunk ids
def plot_chunk_targets(chunks, c_ids):
	pyplot.figure()
	targets = range(56, 95)
	for i in range(len(targets)):
		ax = pyplot.subplot(len(targets), 1, i+1)
		ax.set_xticklabels([])
		ax.set_yticklabels([])
		column = targets[i]
		for chunk_id in c_ids:
			rows = chunks[chunk_id]
			pyplot.plot(rows[:,column])
	pyplot.show()

# load data
dataset = read_csv('AirQualityPrediction/TrainingData.csv', header=0)
# group data by chunks
values = dataset.values
chunks = to_chunks(values)
# plot targets for some chunks
plot_chunk_targets(chunks, [1])

运行该示例将为块标识符“1”创建一个包含 39 个线图的单个图形。

这些图很小，但大致了解变量的时间结构。

我们可以看到，有多个变量没有这个块的数据。这些不能直接预测，也可能不是间接预测。

这表明除了没有所有站点的所有变量之外，甚至列标题中指定的变量也可能不存在于某些块中。

我们还可以在一些系列中看到缺失值的中断。这表明即使我们可能对块中的每个时间步进行观察，我们也可能没有块中所有变量的连续序列。

许多情节都有一个循环结构。大多数都有八个峰值，很可能对应于块内八天的观察。这种季节性结构可以直接建模，也可以在建模时从数据中删除，并添加回预测的时间间隔。

该系列似乎没有任何趋势。

Parallel Time Series Line Plots For All Target Variables for 1 Chunk

1 个块的所有目标变量的并行时间序列线图

我们可以重新运行该示例并使用完整数据绘制前三个块的目标变量。

# plot targets for some chunks
plot_chunk_targets(chunks, [1, 3 ,5])

运行该示例将创建一个图形，其中包含 39 个图形和每个图形的三个时间序列，一个用于每个块的目标。

绘图很忙，您可能希望增加绘图窗口的大小，以便更好地查看目标变量的块之间的比较。

对于许多具有循环日常结构的变量，我们可以看到在整个块中重复的结构。

这是令人鼓舞的，因为它表明为站点建模变量可能对块有所帮助。

此外，曲线 3 至 10 对应于七个不同位置的变量 11。这些图中时间结构的字符串相似性表明，对跨站点使用的每个变量的数据建模可能是有益的。

Parallel Time Series Line Plots For All Target Variables for 3 Chunks

3 个块的所有目标变量的并行时间序列线图

目标变量的箱线图分布

查看目标变量的分布也很有用。

我们可以首先通过为每个目标变量创建框和晶须图来查看每个目标变量的分布。

可以为每个目标并排创建单独的箱图，允许在相同比例下直接比较值的形状和范围。

# boxplot for target variables for a chuck
def plot_chunk_targets_boxplots(chunks, c_id):
	rows = chunks[c_id]
	pyplot.boxplot(rows[:,56:])
	pyplot.show()

下面列出了完整的示例。

# boxplots of targets for a chunk
from numpy import unique
from numpy import isnan
from numpy import count_nonzero
from pandas import read_csv
from matplotlib import pyplot

# split the dataset by 'chunkID', return a dict of id to rows
def to_chunks(values, chunk_ix=1):
	chunks = dict()
	# get the unique chunk ids
	chunk_ids = unique(values[:, chunk_ix])
	# group rows by chunk id
	for chunk_id in chunk_ids:
		selection = values[:, chunk_ix] == chunk_id
		chunks[chunk_id] = values[selection, :]
	return chunks

# boxplot for target variables for a chuck
def plot_chunk_targets_boxplots(chunks, c_id):
	rows = chunks[c_id]
	pyplot.boxplot(rows[:,56:])
	pyplot.show()

# load data
dataset = read_csv('AirQualityPrediction/TrainingData.csv', header=0)
# group data by chunks
values = dataset.values
chunks = to_chunks(values)
# boxplot for target variables
plot_chunk_targets_boxplots(chunks, 1)

运行该示例将创建一个包含 39 个箱图的图形，第一个块的 39 个目标变量中的每一个都有一个。

我们可以看到许多变量的中位数接近零或一;我们还可以看到大多数变量存在大的不对称差异，这表明变量可能与异常值存在偏差。

令人鼓舞的是，7 个站点的变量 11 的 4-10 的箱形图显示了类似的分布。这进一步证明了数据可以按变量分组并用于拟合可跨站点使用的模型。

Box and whisker plots of target variables for one chunk

一个块的目标变量的框和胡须图

我们可以使用跨所有块的数据重新创建此图，以查看数据集范围的模式。

下面列出了完整的示例。

# boxplots of targets for all chunks
from pandas import read_csv
from matplotlib import pyplot

# boxplot for all target variables
def plot_target_boxplots(values):
	pyplot.boxplot(values[:,56:])
	pyplot.show()

# load data
dataset = read_csv('AirQualityPrediction/TrainingData.csv', header=0)
# boxplot for target variables
values = dataset.values
plot_target_boxplots(values)

运行该示例将创建一个新图形，显示整个训练数据集的 39 个框和胡须图，而不管块是什么。

这有点乱，圆圈异常值掩盖了主要的数据分布。

我们可以看到异常值确实扩展到 5 到 10 个单位。这表明在建模时可能在标准化和/或重新缩放目标方面有一些用处。

也许最有用的发现是，有些目标没有任何（或非常多）数据而不管块数。这些列可能应该从数据集中排除。

Box and whisker plots of target variables for all training data

所有训练数据的目标变量的框和胡须图

显然空目标列

我们可以通过创建每列缺失数据比率的条形图来进一步研究明显的缺失数据，不包括开头的元数据列（例如前五列）。

下面的plot_col_percentage_missing()函数实现了这个功能。

# bar chart of the ratio of missing data per column
def plot_col_percentage_missing(values, ix_start=5):
	ratios = list()
	# skip early columns, with meta data or strings
	for col in range(ix_start, values.shape[1]):
		col_data = values[:, col].astype('float32')
		ratio = count_nonzero(isnan(col_data)) / len(col_data) * 100
		ratios.append(ratio)
		if ratio > 90.0:
			print(col, ratio)
	col_id = [x for x in range(ix_start, values.shape[1])]
	pyplot.bar(col_id, ratios)
	pyplot.show()

下面列出了完整的示例。

# summarize missing data per column
from numpy import isnan
from numpy import count_nonzero
from pandas import read_csv
from matplotlib import pyplot

# bar chart of the ratio of missing data per column
def plot_col_percentage_missing(values, ix_start=5):
	ratios = list()
	# skip early columns, with meta data or strings
	for col in range(ix_start, values.shape[1]):
		col_data = values[:, col].astype('float32')
		ratio = count_nonzero(isnan(col_data)) / len(col_data) * 100
		ratios.append(ratio)
		if ratio > 90.0:
			print(ratio)
	col_id = [x for x in range(ix_start, values.shape[1])]
	pyplot.bar(col_id, ratios)
	pyplot.show()

# load data
dataset = read_csv('AirQualityPrediction/TrainingData.csv', header=0)
# plot ratio of missing data per column
values = dataset.values
plot_col_percentage_missing(values)

如果比率高于 90％，则首先运行示例打印列 id（零偏移）和缺失数据的比率。

我们可以看到，实际上没有非 NaN 数据为零的列，但可能有二十多个（12）具有超过 90％的缺失数据。

有趣的是，其中七个是目标变量（指数 56 或更高）。

11 91.48885539779488
20 91.48885539779488
29 91.48885539779488
38 91.48885539779488
47 91.48885539779488
58 95.38880516115385
66 96.9805134713519
68 95.38880516115385
72 97.31630575606145
86 95.38880516115385
92 95.38880516115385
94 95.38880516115385

创建列索引号与缺失数据比率的条形图。

我们可以看到，对于缺失数据的比率可能存在一些分层，群集低于 10％，群集约为 70％，群集高于 90％。

我们还可以看到输入变量和目标变量之间的分离，前者非常规则，因为它们显示了在不同站点测量的相同变量类型。

某些目标变量的这些少量数据表明除了过去的观察之外还需要利用其他因素来做出预测。

Bar Chart of Percentage of Missing Data Per Column

每列缺失数据百分比的条形图

目标变量的直方图分布

目标变量的分布不整齐，至少可以是非高斯分布，或者最差的是高度多模态。

我们可以通过查看目标变量的直方图来查看单个块的数据。

matplotlib 中hist()函数的一个问题是它对 NaN 值不稳健。我们可以通过在绘制之前检查每列是否具有非 NaN 值并排除具有 NaN 值的行来克服此问题。

下面的函数执行此操作，并为一个或多个块的每个目标变量创建一个直方图。

# plot distribution of targets for one or more chunk ids
def plot_chunk_targets_hist(chunks, c_ids):
	pyplot.figure()
	targets = range(56, 95)
	for i in range(len(targets)):
		ax = pyplot.subplot(len(targets), 1, i+1)
		ax.set_xticklabels([])
		ax.set_yticklabels([])
		column = targets[i]
		for chunk_id in c_ids:
			rows = chunks[chunk_id]
			# extract column of interest
			col = rows[:,column].astype('float32')
			# check for some data to plot
			if count_nonzero(isnan(col)) < len(rows):
				# only plot non-nan values
				pyplot.hist(col[~isnan(col)], bins=100)
	pyplot.show()

下面列出了完整的示例。

# plot distribution of targets for a chunk
from numpy import unique
from numpy import isnan
from numpy import count_nonzero
from pandas import read_csv
from matplotlib import pyplot

# split the dataset by 'chunkID', return a dict of id to rows
def to_chunks(values, chunk_ix=1):
	chunks = dict()
	# get the unique chunk ids
	chunk_ids = unique(values[:, chunk_ix])
	# group rows by chunk id
	for chunk_id in chunk_ids:
		selection = values[:, chunk_ix] == chunk_id
		chunks[chunk_id] = values[selection, :]
	return chunks

# plot distribution of targets for one or more chunk ids
def plot_chunk_targets_hist(chunks, c_ids):
	pyplot.figure()
	targets = range(56, 95)
	for i in range(len(targets)):
		ax = pyplot.subplot(len(targets), 1, i+1)
		ax.set_xticklabels([])
		ax.set_yticklabels([])
		column = targets[i]
		for chunk_id in c_ids:
			rows = chunks[chunk_id]
			# extract column of interest
			col = rows[:,column].astype('float32')
			# check for some data to plot
			if count_nonzero(isnan(col)) < len(rows):
				# only plot non-nan values
				pyplot.hist(col[~isnan(col)], bins=100)
	pyplot.show()

# load data
dataset = read_csv('AirQualityPrediction/TrainingData.csv', header=0)
# group data by chunks
values = dataset.values
chunks = to_chunks(values)
# plot targets for some chunks
plot_chunk_targets_hist(chunks, [1])

运行该示例将创建一个包含 39 个直方图的图形，一个用于第一个块的每个目标变量。

情节难以阅读，但是大量的垃圾箱会显示变量的分布。

可以公平地说，也许没有一个目标变量具有明显的高斯分布。许多人的右侧尾部可能有偏斜的分布。

其他变量具有看似非常离散的分布，可能是所选测量设备或测量尺度的人工产物。

Histograms for each target variable for one chunk

一个块的每个目标变量的直方图

我们可以使用所有块的目标变量重新创建相同的图。

下面列出了完整的示例。

# plot distribution of all targets
from numpy import isnan
from numpy import count_nonzero
from pandas import read_csv
from matplotlib import pyplot

# plot histogram for each target variable
def plot_target_hist(values):
	pyplot.figure()
	targets = range(56, 95)
	for i in range(len(targets)):
		ax = pyplot.subplot(len(targets), 1, i+1)
		ax.set_xticklabels([])
		ax.set_yticklabels([])
		column = targets[i]
		# extract column of interest
		col = values[:,column].astype('float32')
		# check for some data to plot
		if count_nonzero(isnan(col)) < len(values):
			# only plot non-nan values
			pyplot.hist(col[~isnan(col)], bins=100)
	pyplot.show()

# load data
dataset = read_csv('AirQualityPrediction/TrainingData.csv', header=0)
# plot targets for all chunks
values = dataset.values
plot_target_hist(values)

运行该示例将创建一个包含 39 个直方图的图形，一个用于训练数据集中的每个目标变量。

我们可以看到更全面的分布，更具洞察力。

第一批情节可能显示出高度倾斜的分布，其核心可能是也可能不是高斯分布。

我们可以看到许多具有间隙的类高斯分布，这表明对高斯分布的连续变量施加了离散测量。

我们还可以看到一些显示指数分布的变量。

总之，这表明使用功率变换来探索将数据重塑为更高斯，和/或使用不依赖于变量的高斯分布的非参数建模方法。例如，可以预期经典线性方法具有困难时间。

Histograms for each target variable for the entire training dataset

整个训练数据集的每个目标变量的直方图

与目标变量的皱痕

比赛结束后，提供数据的人 David Chudzicki 总结了 12 个输出变量的真实含义。

这是在标题为“目标变量确实是”的表格帖子中提供的，部分复制如下：

Description						Target Variable
Carbon monoxide, 				8
Sulfur dioxide, 				4
SO2 max 5-min avg, 				3
Nitric oxide (NO),				10
Nitrogen dioxide (NO2),			14
Oxides of nitrogen (NOx), 		9
Ozone,							11
PM10 Total 0-10um STP, 			5
OC CSN Unadjusted PM2.5 LC TOT,	15
Total Nitrate PM2.5 LC, 		2
EC CSN PM2.5 LC TOT, 			1
Total Carbon PM2.5 LC TOT, 		7
Sulfate PM2.5 LC, 				8
PM2.5 Raw Data, 				4
PM2.5 AQI & Speciation Mass, 	3

这很有意思，因为我们可以看到目标变量本质上是气象的，并且与竞争的名称所暗示的空气质量有关。

问题是数据集中有 15 个变量，只有 12 种不同类型的目标变量。

此问题的原因是数据集中的相同目标变量可用于表示不同的目标变量。特别：

  • 目标 8 可以是'_ 一氧化碳 '或' 硫酸盐 PM2.5LC_ '的数据。
  • 目标 4 可以是'_ 二氧化硫 _'或' _PM2.5 原始数据 _'的数据。
  • 目标 3 可以是' _SO2 最大 5 分钟平均 _'或' _PM2.5 AQI＆amp;形态质量 _'。

根据变量的名称，将数据加倍到相同的目标变量，使用具有不同化学特征的变量，甚至可能使用例如测量的变量。它似乎是偶然的而不是战略性的。

目前尚不清楚，但很可能目标代表一个块中的一个变量，但可能代表块之间的不同变量。或者，变量可能在每个块内的站点之间不同。在前一种情况下，这意味着期望跨块的这些目标变量的一致性的模型，这是一个非常合理的假设，可能有困难。在后者中，模型可以将变量位点组合视为不同的变量。

可以通过比较这些变量在块之间的分布和比例来梳理差异。

这是令人失望的，并且取决于模型技能的重要性，可能需要从数据集中删除这些变量，这些变量是很多目标变量（39 个中的 20 个）。

关于建模的思考

在本节中，我们将利用我们发现的有关该问题的方法，并提出一些建模此问题的方法。

我喜欢这个数据集;它是凌乱的，现实的，抵制朴素的方法。

本节分为四个部分;他们是：

  • 取景。
  • 数据准备。
  • 造型。
  • 评价。

取景

该问题通常被视为多变量多步骤时间序列预测问题。

此外，需要跨多个站点预测多个变量，这是时间序列预测问题的常见结构细分，例如，预测不同物理位置（例如商店或车站）的可变物。

让我们来看看数据的一些可能的框架。

按变量和站点单变量

第一种方法可能是将每个站点的每个变量视为单变量时间序列预测问题。

对模型进行 8 天的每小时观察，并要求预测三天，从中获取并使用或评估预测提前期的特定子集。

在一些选择的情况下，这可能是可能的，这可以通过一些进一步的数据分析来确认。

然而，数据通常抵抗这种框架，因为并非所有块都对每个目标变量进行了 8 天的观察。此外，目标变量的时间序列可能非常不连续，即使不是大多数（90％到 95％）不完整。

我们可以放松对模型所需的先前数据的结构和数量的期望，设计模型以利用可用的任何东西。

这种方法每块需要 39 个模型，总共需要（208 * 39）或 8,112 个独立模型。这听起来有可能，但可能比我们从工程角度更喜欢的可扩展性。

可变站点组合可以跨块建模，仅需要 39 个模型。

单变量变量

目标变量可以跨站点聚合。

我们还可以放松使用滞后时间做出预测，并使用零填充或输入缺失值来显示可用的内容，或者甚至是滞后观察而忽略前置时间。

然后，我们可以根据给定变量的先前观察结果对问题进行框架，预测接下来的三天。

这些模型可能有更多可以使用，但会忽略基于站点的任何变量差异。这可能是也可能不是没有理由的，可以通过比较跨站点的变量分布来检查。

有 12 个独特的变量。

我们可以为每个块建模每个变量，给出（208 * 12）或 2,496 个模型。在块上建模 12 个变量可能更有意义，只需要 12 个模型。

多变量模型

也许一个或多个目标变量依赖于一个或多个气象变量，甚至依赖于其他目标变量。

这可以通过调查每个目标变量和每个输入变量之间的相关性以及其他目标变量来探索。

如果存在或可以假设这样的依赖性，则不仅可以预测具有更完整数据的变量，还可以预测具有 90％以上缺失数据的目标变量。

这些模型可以使用先前气象观测的一些子集和/或目标变量观测作为输入。数据的不连续性可能需要放宽输入变量的传统滞后时间结构，允许模型使用可用于特定预测的任何内容。

数据准备

根据模型的选择，输入和目标变量可能会受益于某些数据准备，例如：

  • 标准化。
  • 正常化。
  • 功率变换，高斯。
  • 季节性差异，存在季节性结构。

为了解决缺失数据，在某些情况下，可能需要通过简单的持久性或平均来进行估算。

在其他情况下，并且取决于模型的选择，可以直接从 NaN 值学习作为观察（例如 XGBoost 可以这样做）或填充 0 值并屏蔽输入（例如 LSTM 可以这样做）。

研究降尺度输入到 2,4 或 12 小时数据或类似数据以试图填补不连续数据中的空白，例如，可能是有趣的。每小时从 12 小时数据预测。

造型

建模可能需要一些原型来发现在方法和选择的输入观察方面有效的方法。

经典方法

可能存在具有完整数据的块的罕见示例，其中诸如 ETS 或 SARIMA 的经典方法可用于单变量预测。

通常，该问题抵制经典方法。

机器学习

一个很好的选择是使用非线性机器学习方法，这些方法与输入数据的时间结构无关，利用任何可用的东西。

此类模型可用于递归或直接策略以预测提前期。直接策略可能更有意义，每个所需的提前期有一个模型。

有 10 个交付周期和 39 个目标变量，在这种情况下，直接策略需要（39 * 10）或 390 个模型。

对问题建模的直接方法的缺点是模型无法利用预测区间中目标变量之间的任何依赖关系，特别是跨站点，跨变量或跨越提前期。如果这些依赖存在（并且肯定会存在），则可以在使用第二层集合模型时添加它们的味道。

特征选择可用于发现变量和/或滞后时间，这可以在预测每个目标变量和提前期时提供最大价值。

这种方法将提供很大的灵活性，正如竞争中所示，决策树的集合在很少调整的情况下表现良好。

深度学习

与机器学习方法一样，深度学习方法可能能够使用任何可用的多变量数据来做出预测。

对于这个问题，两类神经网络可能值得探索：

  • 卷积神经网络或 CNN。
  • 循环神经网络，特别是长短期记忆网络或 LSTM。

CNN 能够将多变量输入时间序列数据的长序列提取到小特征映射中，并且实质上从与预测最相关的序列中学习特征。它们在输入序列中处理噪声和特征不变性的能力可能是有用的。与其他神经网络一样，CNN 可以输出向量以预测预测的提前期。

LSTM 旨在处理序列数据，并可通过屏蔽直接支持缺失数据。它们也能够从长输入序列自动进行特征学习，单独或与 CNN 结合可以很好地解决这个问题。与编解码器架构一起，LSTM 网络可用于本地预测多个交付周期。

评估

一种反映竞争中使用的朴素方法可能最适合评估模型。

也就是说，将每个块分成训练和测试集，在这种情况下使用前五天的数据进行训练，其余三个用于测试。

通过在整个数据集上进行训练并将预测提交给 Kaggle 网站以评估所保持的测试数据集，最终确定模型可能是有趣的。

进一步阅读

如果您希望深入了解，本节将提供有关该主题的更多资源。

摘要

在本教程中，您发现并探索了空气质量预测数据集，该数据集代表了具有挑战性的多变量，多站点和多步骤时间序列预测问题。

具体来说，你学到了：

  • 如何加载和探索数据集的块结构。
  • 如何探索和可视化数据集的输入和目标变量。
  • 如何使用新的理解来概述一系列方法来构建问题，准备数据和建模数据集。

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

如何从智能手机数据建模人类活动

原文： machinelearningmastery.com/how-to-model-human-activity-from-smartphone-data/

人类活动识别是将由专用线束或智能电话记录的加速度计数据序列分类为已知的明确定义的运动的问题。

鉴于每秒产生大量观测结果，观测的时间性质以及缺乏将加速度计数据与已知运动联系起来的明确方法，这是一个具有挑战性的问题。

该问题的经典方法涉及基于固定大小的窗口和训练机器学习模型（例如决策树的集合）的时间序列数据中的手工制作特征。困难在于此功能工程需要该领域的深厚专业知识。

最近，诸如循环神经网络和一维卷积神经网络（CNN）之类的深度学习方法已经被证明在很少或没有数据特征工程的情况下提供具有挑战性的活动识别任务的最新结果。

在本教程中，您将发现用于时间序列分类的'_ 活动识别使用智能手机 _'数据集，以及如何加载和探索数据集以使其为预测性建模做好准备。

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

  • 如何下载数据集并将其加载到内存中。
  • 如何使用线图，直方图和箱线图来更好地理解运动数据的结构。
  • 如何建模问题，包括框架，数据准备，建模和评估。

让我们开始吧。

How to Model Human Activity From Smartphone Data

如何从智能手机数据模拟人类活动 照片由摄影师，保留一些权利。

教程概述

本教程分为 10 个部分;他们是：

  1. 人类活动识别
  2. 使用智能手机数据集进行活动识别
  3. 下载数据集
  4. 加载数据
  5. 活动类平衡
  6. 绘制一个主题的时间序列数据
  7. 绘制每个主题的直方图
  8. 绘制每个活动的直方图
  9. 绘制活动持续时间箱图
  10. 建模方法

1.人类活动识别

人类活动识别，或简称为 HAR，是基于使用传感器的移动痕迹来预测人正在做什么的问题。

运动通常是正常的室内活动，例如站立，坐着，跳跃和上楼梯。

传感器通常位于主体上，例如智能手机或背心，并且经常以三维（x，y，z）记录加速度计数据。

人类活动识别（HAR）旨在识别由一个人对他/她自己和周围环境进行一组观察所执行的动作。识别可以通过利用从各种来源（例如环境或身体佩戴的传感器）检索的信息来实现。

这个想法是，一旦主体的活动被识别和知道，智能计算机系统就可以提供帮助。

这是一个具有挑战性的问题，因为没有明确的分析方法将传感器数据与一般方式的特定动作联系起来。由于收集了大量的传感器数据（例如，每秒数十或数百次观察），并且在开发预测模型时从这些数据中经典使用手工制作的特征和启发式，因此在技术上具有挑战性。

最近，深度学习方法已经成功地证明了 HAR 问题，因为它们能够自动学习更高阶的特征。

基于传感器的活动识别从大量低水平传感器读数中寻找关于人类活动的深刻的高级知识。传统的模式识别方法在过去几年中取得了巨大的进步。然而，这些方法通常严重依赖于启发式手工特征提取，这可能会妨碍它们的泛化表现。 [...]最近，深度学习的最新进展使得可以执行自动高级特征提取，从而在许多领域实现了有希望的表现。

2.使用智能手机数据集识别活动

标准人类活动识别数据集是 2012 年提供的“_ 活动识别使用智能手机 _”数据集。

它由 Davide Anguita 等人准备并提供。来自意大利热那亚大学的 2013 年论文“使用智能手机进行人类活动识别的公共领域数据集”中对该数据集进行了全面描述。该数据集在他们的 2012 年论文中用机器学习算法建模，标题为“使用多类硬件友好支持向量机在智能手机上进行人类活动识别。“

数据集可用，可以从 UCI 机器学习库免费下载：

该数据来自 30 名年龄在 19 至 48 岁之间的受试者，他们进行 6 项标准活动中的一项，同时佩戴记录运动数据的腰部智能手机。记录执行活动的每个受试者的视频，并从这些视频手动标记移动数据。

以下是在记录其移动数据的同时执行活动的主体的示例视频。

<iframe allow="autoplay; encrypted-media" allowfullscreen="" frameborder="0" height="375" src="www.youtube.com/embed/XOEN9…" width="500"></iframe>

进行的六项活动如下：

  1. 步行
  2. 走上楼
  3. 走楼下
  4. 坐在
  5. 常设
  6. 铺设

记录的运动数据是来自智能手机的 x，y 和 z 加速度计数据（线性加速度）和陀螺仪数据（角速度），特别是三星 Galaxy S II 。

以 50Hz（即每秒 50 个数据点）记录观察结果。每个受试者进行两次活动，一次是左侧设备，另一次是右侧设备。

选择了一组 30 名志愿者，年龄从 19 岁到 48 岁不等。每个人都被要求遵循活动协议，同时佩戴腰部安装的三星 Galaxy S II 智能手机。六个选定的 ADL 是站立，坐下，放下，走路，走楼下和楼上。每个受试者执行两次协议：在第一次试验中，智能手机固定在皮带的左侧，而第二次试验由用户自己放置为首选

原始数据不可用。相反，可以使用预处理版本的数据集。

预处理步骤包括：

  • 使用噪声滤波器预处理加速度计和陀螺仪。
  • 将数据拆分为 2.56 秒（128 个数据点）的固定窗口，重叠率为 50％。
  • 将加速度计数据分割为重力（总）和身体运动分量。

使用中值滤波器和具有 20Hz 截止频率的 3 阶低通巴特沃斯滤波器对这些信号进行预处理以降低噪声。 [...]加速度信号，具有重力和身体运动成分，使用另一个巴特沃斯低通滤波器分离成身体加速度和重力。

特征工程应用于窗口数据，并且提供具有这些工程特征的数据的副本。

从每个窗口提取在人类活动识别领域中常用的许多时间和频率特征。结果是 561 元素的特征向量。

根据受试者的数据，将数据集分成训练（70％）和测试（30％）组。训练 21 个，测试 9 个。

这表明问题的框架，其中一系列运动活动被用作输入来预测当前正在进行的活动的部分（2.56 秒），其中使用已知受试者训练的模型来预测新受试者的运动活动。 。

使用旨在用于智能手机的支持向量机（例如定点算术）的早期实验结果导致测试数据集的预测准确度为 89％，实现与未修改的 SVM 实现类似的结果。

该方法适用于标准支持向量机（SVM）并利用定点算法来降低计算成本。与传统 SVM 的比较表明，计算成本方面有显着改善，同时保持了相似的精度[...]

现在我们已经熟悉了预测问题，我们将看看加载和探索这个数据集。

3.下载数据集

该数据集是免费提供的，可以从 UCI 机器学习库下载。

数据以单个 zip 文件的形式提供，大小约为 58 兆字节。此下载的直接链接如下：

下载数据集并将所有文件解压缩到当前工作目录中名为“HARDataset”的新目录中。

检查解压缩的内容，您会注意到以下几点：

  • 存在“train”和“test”文件夹，其包含用于建模的数据的分割部分（例如，70％/ 30％）。
  • 有一个“README.txt”文件，其中包含数据集的详细技术说明和解压缩文件的内容。
  • 有一个“features.txt”文件，其中包含工程特性的技术说明。

“train”和“test”文件夹的内容类似（例如文件夹和文件名），尽管它们包含的特定数据存在差异。

检查“train”文件夹显示了一些重要元素：

  • 包含预处理数据的“ Inertial Signals ”文件夹。
  • “X_train.txt”文件，包含用于拟合模型的工程特征。
  • “y_train.txt”包含每个观察的类标签（1-6）。
  • “subject_train.txt”文件，其中包含数据文件中每一行的映射及其主题标识符（1-30）。

每个文件中的行数匹配，表示每行是每个数据文件中的一个记录。

“ Inertial Signals ”目录包含 9 个文件。

  • _x，y 和 z 轴的重力加速度 _ 数据文件：total_acc_x_train.txt，total_acc_y_train.txt，total_acc_z_train.txt。
  • _x，y 和 z 轴的身体加速度 _ 数据文件：body_acc_x_train.txt，body_acc_y_train.txt，body_acc_z_train.txt。
  • _x，y 和 z 轴的体陀螺 _ 数据文件：body_gyro_x_train.txt，body_gyro_y_train.txt，body_gyro_z_train.txt。

该结构在“test”目录中进行镜像。

我们将把注意力集中在“_ 惯性信号 _”中的数据，因为这是开发可以学习合适表示的机器学习模型中最有趣的，而不是使用特定于域的特征工程。

检查数据文件显示列由空格分隔，值显示为缩放到-1,1。此缩放可以通过数据集随附的README.txt文件中的注释确认。

现在我们知道了我们拥有的数据，我们可以弄清楚如何将其加载到内存中。

4.加载数据

在本节中，我们将开发一些代码来将数据集加载到内存中。

首先，我们需要加载一个文件。

我们可以使用 read_csv（）Pandas 函数来加载单个数据文件，并指定该文件没有标题并使用空格分隔列。

dataframe = read_csv(filepath, header=None, delim_whitespace=True)

我们可以将它包装在名为load_file()的函数中。下面列出了此功能的完整示例。

# load dataset
from pandas import read_csv

# load a single file as a numpy array
def load_file(filepath):
	dataframe = read_csv(filepath, header=None, delim_whitespace=True)
	return dataframe.values

data = load_file('HARDataset/train/Inertial Signals/total_acc_y_train.txt')
print(data.shape)

运行该示例加载文件'total_acc_y_train.txt'，返回 NumPy 数组，并打印数组的形状。

我们可以看到训练数据由 7,352 行或数据窗口组成，其中每个窗口有 128 个观察值。

(7352, 128)

接下来，将一组文件（例如所有身体加速度数据文件）作为单个组加载将是有用的。

理想情况下，在处理多变量时间序列数据时，以下列格式构建数据非常有用：

[samples, timesteps, features]

这有助于分析，并且是对诸如卷积神经网络和循环神经网络的深度学习模型的期望。

我们可以通过多次调用上面的load_file()函数来实现这一点，对于组中的每个文件一次。

一旦我们将每个文件作为 NumPy 数组加载，我们就可以将所有三个数组组合或堆叠在一起。我们可以使用 dstack（）NumPy 函数来确保每个数组的堆叠方式使得要素在第三维中分离，正如我们所希望的那样。

函数load_group()对文件名列表实现此行为，如下所示。

# load a list of files, such as x, y, z data for a given variable
def load_group(filenames, prefix=''):
	loaded = list()
	for name in filenames:
		data = load_file(prefix + name)
		loaded.append(data)
	# stack group so that features are the 3rd dimension
	loaded = dstack(loaded)
	return loaded

我们可以通过加载所有总加速文件来演示此功能。

下面列出了完整的示例。

# load dataset
from numpy import dstack
from pandas import read_csv

# load a single file as a numpy array
def load_file(filepath):
	dataframe = read_csv(filepath, header=None, delim_whitespace=True)
	return dataframe.values

# load a list of files, such as x, y, z data for a given variable
def load_group(filenames, prefix=''):
	loaded = list()
	for name in filenames:
		data = load_file(prefix + name)
		loaded.append(data)
	# stack group so that features are the 3rd dimension
	loaded = dstack(loaded)
	return loaded

# load the total acc data
filenames = ['total_acc_x_train.txt', 'total_acc_y_train.txt', 'total_acc_z_train.txt']
total_acc = load_group(filenames, prefix='HARDataset/train/Inertial Signals/')
print(total_acc.shape)

运行该示例将打印返回的 NumPy 数组的形状，显示数据集的三个要素 x，y 和 z 的预期样本数和时间步长。

(7352, 128, 3)

最后，我们可以使用到目前为止开发的两个函数来加载训练和测试数据集的所有数据。

给定训练和测试文件夹中的并行结构，我们可以开发一个新函数来加载给定文件夹的所有输入和输出数据。该函数可以构建要加载的所有 9 个数据文件的列表，将它们作为一个具有 9 个特征的 NumPy 数组加载，然后加载包含输出类的数据文件。

下面的load_dataset()函数实现了这种行为。它可以被称为“train”组或“test”组，作为字符串参数传递。

# load a dataset group, such as train or test
def load_dataset(group, prefix=''):
	filepath = prefix + group + '/Inertial Signals/'
	# load all 9 files as a single array
	filenames = list()
	# total acceleration
	filenames += ['total_acc_x_'+group+'.txt', 'total_acc_y_'+group+'.txt', 'total_acc_z_'+group+'.txt']
	# body acceleration
	filenames += ['body_acc_x_'+group+'.txt', 'body_acc_y_'+group+'.txt', 'body_acc_z_'+group+'.txt']
	# body gyroscope
	filenames += ['body_gyro_x_'+group+'.txt', 'body_gyro_y_'+group+'.txt', 'body_gyro_z_'+group+'.txt']
	# load input data
	X = load_group(filenames, filepath)
	# load class output
	y = load_file(prefix + group + '/y_'+group+'.txt')
	return X, y

下面列出了完整的示例。

# load dataset
from numpy import dstack
from pandas import read_csv

# load a single file as a numpy array
def load_file(filepath):
	dataframe = read_csv(filepath, header=None, delim_whitespace=True)
	return dataframe.values

# load a list of files, such as x, y, z data for a given variable
def load_group(filenames, prefix=''):
	loaded = list()
	for name in filenames:
		data = load_file(prefix + name)
		loaded.append(data)
	# stack group so that features are the 3rd dimension
	loaded = dstack(loaded)
	return loaded

# load a dataset group, such as train or test
def load_dataset(group, prefix=''):
	filepath = prefix + group + '/Inertial Signals/'
	# load all 9 files as a single array
	filenames = list()
	# total acceleration
	filenames += ['total_acc_x_'+group+'.txt', 'total_acc_y_'+group+'.txt', 'total_acc_z_'+group+'.txt']
	# body acceleration
	filenames += ['body_acc_x_'+group+'.txt', 'body_acc_y_'+group+'.txt', 'body_acc_z_'+group+'.txt']
	# body gyroscope
	filenames += ['body_gyro_x_'+group+'.txt', 'body_gyro_y_'+group+'.txt', 'body_gyro_z_'+group+'.txt']
	# load input data
	X = load_group(filenames, filepath)
	# load class output
	y = load_file(prefix + group + '/y_'+group+'.txt')
	return X, y

# load all train
trainX, trainy = load_dataset('train', 'HARDataset/')
print(trainX.shape, trainy.shape)
# load all test
testX, testy = load_dataset('test', 'HARDataset/')
print(testX.shape, testy.shape)

运行该示例将加载训练和测试数据集。

我们可以看到测试数据集有 2,947 行窗口数据。正如预期的那样，我们可以看到训练和测试装置中的窗口大小匹配，并且每个训练和测试用例中的输出（y）大小与样本数量相匹配。

(7352, 128, 9) (7352, 1)
(2947, 128, 9) (2947, 1)

既然我们知道如何加载数据，我们就可以开始探索它了。

5.活动类平衡

对数据进行良好的首次检查是调查每项活动的平衡。

我们相信，30 个科目中的每一个都进行了六项活动。

确认此期望将检查数据是否确实平衡，使模型更容易，并确认我们正确加载和解释数据集。

我们可以开发一个函数来总结输出变量的细分，例如： y 变量。

下面的函数class_breakdown()实现了这种行为，首先将提供的 NumPy 数组包装在 DataFrame 中，按类值对行进行分组，并计算每个组的大小（行数）。然后总结结果，包括计数和百分比。

# summarize the balance of classes in an output variable column
def class_breakdown(data):
	# convert the numpy array into a dataframe
	df = DataFrame(data)
	# group data by the class value and calculate the number of rows
	counts = df.groupby(0).size()
	# retrieve raw rows
	counts = counts.values
	# summarize
	for i in range(len(counts)):
		percent = counts[i] / len(df) * 100
		print('Class=%d, total=%d, percentage=%.3f' % (i+1, counts[i], percent))

总结训练和测试数据集中类的细分以确保它们具有类似的细分可能是有用的，然后将结果与组合数据集的细分进行比较。

下面列出了完整的示例。

# summarize class balance
from numpy import array
from numpy import vstack
from pandas import read_csv
from pandas import DataFrame

# load a single file as a numpy array
def load_file(filepath):
	dataframe = read_csv(filepath, header=None, delim_whitespace=True)
	return dataframe.values

# summarize the balance of classes in an output variable column
def class_breakdown(data):
	# convert the numpy array into a dataframe
	df = DataFrame(data)
	# group data by the class value and calculate the number of rows
	counts = df.groupby(0).size()
	# retrieve raw rows
	counts = counts.values
	# summarize
	for i in range(len(counts)):
		percent = counts[i] / len(df) * 100
		print('Class=%d, total=%d, percentage=%.3f' % (i+1, counts[i], percent))

# load train file
trainy = load_file('HARDataset/train/y_train.txt')
# summarize class breakdown
print('Train Dataset')
class_breakdown(trainy)

# load test file
testy = load_file('HARDataset/test/y_test.txt')
# summarize class breakdown
print('Test Dataset')
class_breakdown(testy)

# summarize combined class breakdown
print('Both')
combined = vstack((trainy, testy))
class_breakdown(combined)

首先运行该示例总结了训练集的细分。我们可以看到每个类的非常相似的分布在数据集的 13％到 19％之间徘徊。

测试集和两个数据集上的结果看起来非常相似。

使用数据集可能是安全的，假设每个训练和测试集以及每个主题可以平衡类的分布。

Train Dataset
Class=1, total=1226, percentage=16.676
Class=2, total=1073, percentage=14.595
Class=3, total=986, percentage=13.411
Class=4, total=1286, percentage=17.492
Class=5, total=1374, percentage=18.689
Class=6, total=1407, percentage=19.138

Test Dataset
Class=1, total=496, percentage=16.831
Class=2, total=471, percentage=15.982
Class=3, total=420, percentage=14.252
Class=4, total=491, percentage=16.661
Class=5, total=532, percentage=18.052
Class=6, total=537, percentage=18.222

Both
Class=1, total=1722, percentage=16.720
Class=2, total=1544, percentage=14.992
Class=3, total=1406, percentage=13.652
Class=4, total=1777, percentage=17.254
Class=5, total=1906, percentage=18.507
Class=6, total=1944, percentage=18.876

6.绘制一个主题的时间序列数据

我们正在处理时间序列数据，因此导入检查是创建原始数据的线图。

原始数据由每个变量的时间序列数据窗口组成，窗口确实有 50％的重叠。这表明我们可能会在观察中看到一些重复作为线图，除非重叠被删除。

我们可以从使用上面开发的函数加载训练数据集开始。

# load data
trainX, trainy = load_dataset('train', 'HARDataset/')

接下来，我们可以在'train'目录中加载'subject_train.txt'，该目录提供行到它所属主题的映射。

我们可以使用load_file()函数加载这个文件。加载后，我们还可以使用unique()NumPy 函数来检索训练数据集中的唯一主题列表。

sub_map = load_file('HARDataset/train/subject_train.txt')
train_subjects = unique(sub_map)
print(train_subjects)

接下来，我们需要一种方法来检索单个主题的所有行，例如科目编号 1。

我们可以通过查找属于给定主题的所有行号并使用这些行号从训练数据集中加载的 X 和 y 数据中选择样本来完成此操作。

下面的data_for_subject()函数实现了这种行为。它将获取加载的训练数据，行号到主题的加载映射，以及我们感兴趣的主题的主题标识号，并将返回X和y仅针对该主题的数据。

# get all data for one subject
def data_for_subject(X, y, sub_map, sub_id):
	# get row indexes for the subject id
	ix = [i for i in range(len(sub_map)) if sub_map[i]==sub_id]
	# return the selected samples
	return X[ix, :, :], y[ix]

现在我们有一个主题的数据，我们可以绘制它。

数据由具有重叠的窗口组成。我们可以编写一个函数来消除这种重叠，并将窗口向下挤压给定变量，将其压缩成一个可以直接绘制为线图的长序列。

下面的to_series()函数对给定变量实现此行为，例如窗户数组。

# convert a series of windows to a 1D list
def to_series(windows):
	series = list()
	for window in windows:
		# remove the overlap from the window
		half = int(len(window) / 2) - 1
		for value in window[-half:]:
			series.append(value)
	return series

最后，我们有足够的情节来绘制数据。我们可以依次绘制主题的九个变量中的每一个，以及活动水平的最终图。

每个系列将具有相同数量的时间步长（x 轴的长度），因此，为每个变量创建一个子图并垂直对齐所有图可能很有用，这样我们就可以比较每个变量的运动。

下面的plot_subject()函数为单个主题的X和y数据实现此行为。该函数采用与load_dataset()函数中加载的变量（第 3 轴）相同的顺序。每个情节都会添加粗略的标题，因此我们不会轻易混淆我们正在看的内容。

# plot the data for one subject
def plot_subject(X, y):
	pyplot.figure()
	# determine the total number of plots
	n, off = X.shape[2] + 1, 0
	# plot total acc
	for i in range(3):
		pyplot.subplot(n, 1, off+1)
		pyplot.plot(to_series(X[:, :, off]))
		pyplot.title('total acc '+str(i), y=0, loc='left')
		off += 1
	# plot body acc
	for i in range(3):
		pyplot.subplot(n, 1, off+1)
		pyplot.plot(to_series(X[:, :, off]))
		pyplot.title('body acc '+str(i), y=0, loc='left')
		off += 1
	# plot body gyro
	for i in range(3):
		pyplot.subplot(n, 1, off+1)
		pyplot.plot(to_series(X[:, :, off]))
		pyplot.title('body gyro '+str(i), y=0, loc='left')
		off += 1
	# plot activities
	pyplot.subplot(n, 1, n)
	pyplot.plot(y)
	pyplot.title('activity', y=0, loc='left')
	pyplot.show()

下面列出了完整的示例。

# plot all vars for one subject
from numpy import array
from numpy import dstack
from numpy import unique
from pandas import read_csv
from matplotlib import pyplot

# load a single file as a numpy array
def load_file(filepath):
	dataframe = read_csv(filepath, header=None, delim_whitespace=True)
	return dataframe.values

# load a list of files, such as x, y, z data for a given variable
def load_group(filenames, prefix=''):
	loaded = list()
	for name in filenames:
		data = load_file(prefix + name)
		loaded.append(data)
	# stack group so that features are the 3rd dimension
	loaded = dstack(loaded)
	return loaded

# load a dataset group, such as train or test
def load_dataset(group, prefix=''):
	filepath = prefix + group + '/Inertial Signals/'
	# load all 9 files as a single array
	filenames = list()
	# total acceleration
	filenames += ['total_acc_x_'+group+'.txt', 'total_acc_y_'+group+'.txt', 'total_acc_z_'+group+'.txt']
	# body acceleration
	filenames += ['body_acc_x_'+group+'.txt', 'body_acc_y_'+group+'.txt', 'body_acc_z_'+group+'.txt']
	# body gyroscope
	filenames += ['body_gyro_x_'+group+'.txt', 'body_gyro_y_'+group+'.txt', 'body_gyro_z_'+group+'.txt']
	# load input data
	X = load_group(filenames, filepath)
	# load class output
	y = load_file(prefix + group + '/y_'+group+'.txt')
	return X, y

# get all data for one subject
def data_for_subject(X, y, sub_map, sub_id):
	# get row indexes for the subject id
	ix = [i for i in range(len(sub_map)) if sub_map[i]==sub_id]
	# return the selected samples
	return X[ix, :, :], y[ix]

# convert a series of windows to a 1D list
def to_series(windows):
	series = list()
	for window in windows:
		# remove the overlap from the window
		half = int(len(window) / 2) - 1
		for value in window[-half:]:
			series.append(value)
	return series

# plot the data for one subject
def plot_subject(X, y):
	pyplot.figure()
	# determine the total number of plots
	n, off = X.shape[2] + 1, 0
	# plot total acc
	for i in range(3):
		pyplot.subplot(n, 1, off+1)
		pyplot.plot(to_series(X[:, :, off]))
		pyplot.title('total acc '+str(i), y=0, loc='left')
		off += 1
	# plot body acc
	for i in range(3):
		pyplot.subplot(n, 1, off+1)
		pyplot.plot(to_series(X[:, :, off]))
		pyplot.title('body acc '+str(i), y=0, loc='left')
		off += 1
	# plot body gyro
	for i in range(3):
		pyplot.subplot(n, 1, off+1)
		pyplot.plot(to_series(X[:, :, off]))
		pyplot.title('body gyro '+str(i), y=0, loc='left')
		off += 1
	# plot activities
	pyplot.subplot(n, 1, n)
	pyplot.plot(y)
	pyplot.title('activity', y=0, loc='left')
	pyplot.show()

# load data
trainX, trainy = load_dataset('train', 'HARDataset/')
# load mapping of rows to subjects
sub_map = load_file('HARDataset/train/subject_train.txt')
train_subjects = unique(sub_map)
print(train_subjects)
# get the data for one subject
sub_id = train_subjects[0]
subX, suby = data_for_subject(trainX, trainy, sub_map, sub_id)
print(subX.shape, suby.shape)
# plot data for subject
plot_subject(subX, suby)

运行该示例将打印训练数据集中的唯一主题，第一个主题的数据样本，并创建一个包含 10 个图的图形，每个图形对应九个输入变量和输出类别。

[ 1 3 5 6 7 8 11 14 15 16 17 19 21 22 23 25 26 27 28 29 30]
(341, 128, 9) (341, 1)

在图中，我们可以看到与活动 1,2 和 3 相对应的大运动周期：步行活动。对于编号较高的活动，4,5 和 6（坐，站立和铺设），我们还可以看到更少的活动（即相对直线）。

这是我们正确加载解释原始数据集的确认。

我们可以看到这个主题已经执行了两次相同的一般活动顺序，并且一些活动执行了两次以上。这表明，对于特定主题，我们不应假设可能已执行的活动或其顺序。

我们还可以看到一些固定活动的相对较大的运动，例如铺设。这些可能是异常值或与活动转换相关。可以将这些观察结果平滑或移除为异常值。

最后，我们在九个变量中看到了很多共性。很可能只需要这些迹线的一部分来开发预测模型。

Line plot for all variables for a single subject

单个主题的所有变量的线图

我们可以通过做一个小的改动来重新运行另一个主题的例子，例如：选择训练数据集中第二个主题的标识符。

# get the data for one subject
sub_id = train_subjects[1]

第二个主题的情节显示出类似的行为，没有任何意外。

活动的双重序列确实比第一个主题更加规律。

Line plot for all variables for a second single subject

第二个单个主题的所有变量的线图

7.绘制每个受试者的直方图

由于问题是框架，我们有兴趣使用一些科目的运动数据来预测其他科目的运动活动。

这表明跨学科的运动数据必须有规律性。我们知道数据已经在-1 和 1 之间缩放，可能是每个受试者，这表明检测到的运动的幅度将是相似的。

我们还期望运动数据的分布在不同主题之间是相似的，因为它们执行相同的动作。

我们可以通过绘制和比较受试者的运动数据的直方图来检查这一点。一种有用的方法是为每个受试者创建一个图并绘制给定数据的所有三个轴（例如总加速度），然后对多个受试者重复此操作。可以修改绘图以使用相同的轴并水平对齐，以便可以比较跨主题的每个变量的分布。

下面的plot_subject_histograms()函数实现了这种行为。该函数采用加载的数据集和行到主题的映射以及要绘制的最大主题数，默认情况下固定为 10。

为每个主题创建一个绘图，并将一个数据类型的三个变量绘制为具有 100 个二进制位的直方图，以帮助使分布明显。每个图共享相同的轴，该轴固定在-1 和 1 的边界。

# plot histograms for multiple subjects
def plot_subject_histograms(X, y, sub_map, n=10):
	pyplot.figure()
	# get unique subjects
	subject_ids = unique(sub_map[:,0])
	# enumerate subjects
	xaxis = None
	for k in range(n):
		sub_id = subject_ids[k]
		# get data for one subject
		subX, _ = data_for_subject(X, y, sub_map, sub_id)
		# total acc
		for i in range(3):
			ax = pyplot.subplot(n, 1, k+1, sharex=xaxis)
			ax.set_xlim(-1,1)
			if k == 0:
				xaxis = ax
			pyplot.hist(to_series(subX[:,:,i]), bins=100)
	pyplot.show()

下面列出了完整的示例。

# plot histograms for multiple subjects
from numpy import array
from numpy import unique
from numpy import dstack
from pandas import read_csv
from matplotlib import pyplot

# load a single file as a numpy array
def load_file(filepath):
	dataframe = read_csv(filepath, header=None, delim_whitespace=True)
	return dataframe.values

# load a list of files, such as x, y, z data for a given variable
def load_group(filenames, prefix=''):
	loaded = list()
	for name in filenames:
		data = load_file(prefix + name)
		loaded.append(data)
	# stack group so that features are the 3rd dimension
	loaded = dstack(loaded)
	return loaded

# load a dataset group, such as train or test
def load_dataset(group, prefix=''):
	filepath = prefix + group + '/Inertial Signals/'
	# load all 9 files as a single array
	filenames = list()
	# total acceleration
	filenames += ['total_acc_x_'+group+'.txt', 'total_acc_y_'+group+'.txt', 'total_acc_z_'+group+'.txt']
	# body acceleration
	filenames += ['body_acc_x_'+group+'.txt', 'body_acc_y_'+group+'.txt', 'body_acc_z_'+group+'.txt']
	# body gyroscope
	filenames += ['body_gyro_x_'+group+'.txt', 'body_gyro_y_'+group+'.txt', 'body_gyro_z_'+group+'.txt']
	# load input data
	X = load_group(filenames, filepath)
	# load class output
	y = load_file(prefix + group + '/y_'+group+'.txt')
	return X, y

# get all data for one subject
def data_for_subject(X, y, sub_map, sub_id):
	# get row indexes for the subject id
	ix = [i for i in range(len(sub_map)) if sub_map[i]==sub_id]
	# return the selected samples
	return X[ix, :, :], y[ix]

# convert a series of windows to a 1D list
def to_series(windows):
	series = list()
	for window in windows:
		# remove the overlap from the window
		half = int(len(window) / 2) - 1
		for value in window[-half:]:
			series.append(value)
	return series

# plot histograms for multiple subjects
def plot_subject_histograms(X, y, sub_map, n=10):
	pyplot.figure()
	# get unique subjects
	subject_ids = unique(sub_map[:,0])
	# enumerate subjects
	xaxis = None
	for k in range(n):
		sub_id = subject_ids[k]
		# get data for one subject
		subX, _ = data_for_subject(X, y, sub_map, sub_id)
		# total acc
		for i in range(3):
			ax = pyplot.subplot(n, 1, k+1, sharex=xaxis)
			ax.set_xlim(-1,1)
			if k == 0:
				xaxis = ax
			pyplot.hist(to_series(subX[:,:,i]), bins=100)
	pyplot.show()

# load training dataset
X, y = load_dataset('train', 'HARDataset/')
# load mapping of rows to subjects
sub_map = load_file('HARDataset/train/subject_train.txt')
# plot histograms for subjects
plot_subject_histograms(X, y, sub_map)

运行该示例将创建一个包含 10 个绘图的单个图形，其中包含总加速度数据的三个轴的直方图。

给定图上的三个轴中的每一个具有不同的颜色，具体地，x，y 和 z 分别是蓝色，橙色和绿色。

我们可以看到给定轴的分布确实看起来是高斯分布的，具有大的独立数据组。

我们可以看到一些分布对齐（例如，中间的主要组大约为 0.0），这表明对象之间的运动数据可能存在一些连续性，至少对于这些数据而言。

Histograms of the total acceleration data for 10 subjects

10 个受试者的总加速度数据的直方图

我们可以更新plot_subject_histograms()函数，接下来绘制身体加速度的分布。更新的功能如下所示。

# plot histograms for multiple subjects
def plot_subject_histograms(X, y, sub_map, n=10):
	pyplot.figure()
	# get unique subjects
	subject_ids = unique(sub_map[:,0])
	# enumerate subjects
	xaxis = None
	for k in range(n):
		sub_id = subject_ids[k]
		# get data for one subject
		subX, _ = data_for_subject(X, y, sub_map, sub_id)
		# body acc
		for i in range(3):
			ax = pyplot.subplot(n, 1, k+1, sharex=xaxis)
			ax.set_xlim(-1,1)
			if k == 0:
				xaxis = ax
			pyplot.hist(to_series(subX[:,:,3+i]), bins=100)
	pyplot.show()

运行更新的示例会创建具有非常不同结果的相同图表。

在这里，我们可以看到所有数据聚集在一个主题内和主题之间的轴上。这表明数据可能是中心的（零均值）。跨受试者的这种强一致性可能有助于建模，并且可能表明总加速度数据中受试者之间的差异可能不那么有用。

Histograms of the body acceleration data for 10 subjects

10 名受试者的身体加速度数据的直方图

最后，我们可以为陀螺仪数据生成一个最终图。

更新的功能如下所示。

# plot histograms for multiple subjects
def plot_subject_histograms(X, y, sub_map, n=10):
	pyplot.figure()
	# get unique subjects
	subject_ids = unique(sub_map[:,0])
	# enumerate subjects
	xaxis = None
	for k in range(n):
		sub_id = subject_ids[k]
		# get data for one subject
		subX, _ = data_for_subject(X, y, sub_map, sub_id)
		# body acc
		for i in range(3):
			ax = pyplot.subplot(n, 1, k+1, sharex=xaxis)
			ax.set_xlim(-1,1)
			if k == 0:
				xaxis = ax
			pyplot.hist(to_series(subX[:,:,6+i]), bins=100)
	pyplot.show()

运行该示例显示与身体加速度数据非常相似的结果。

我们看到每个轴上的每个轴的高斯分布的可能性很高，以 0.0 为中心。分布更宽一些，显示更加丰富的尾巴，但这对于跨主题的运动数据建模是一个令人鼓舞的发现。

Histograms of the body gyroscope data for 10 subjects

10 名受试者的身体陀螺仪数据的直方图

8.绘制每个活动的直方图

我们有兴趣根据活动数据区分活动。

最简单的情况是区分单个主题的活动。调查此问题的一种方法是按活动审查主题的移动数据分布。我们希望看到单个主题的不同活动的运动数据之间的分布存在一些差异。

我们可以通过创建每个活动的直方图来检查这一点，每个图上给定数据类型的三个轴。同样，可以水平排列图以比较每个数据轴的活动分布。我们希望看到各地活动的分布存在差异。

首先，我们必须按活动对主题的跟踪进行分组。下面的data_by_activity()函数实现了这种行为。

# group data by activity
def data_by_activity(X, y, activities):
	# group windows by activity
	return {a:X[y[:,0]==a, :, :] for a in activities}

我们现在可以为给定主题的每个活动创建绘图。

下面的plot_activity_histograms()函数为给定主题的跟踪数据实现此功能。

首先，按活动对数据进行分组，然后为每个活动创建一个子图，并将数据类型的每个轴添加为直方图。该函数仅枚举数据的前三个特征，即总加速度变量。

# plot histograms for each activity for a subject
def plot_activity_histograms(X, y):
	# get a list of unique activities for the subject
	activity_ids = unique(y[:,0])
	# group windows by activity
	grouped = data_by_activity(X, y, activity_ids)
	# plot per activity, histograms for each axis
	pyplot.figure()
	xaxis = None
	for k in range(len(activity_ids)):
		act_id = activity_ids[k]
		# total acceleration
		for i in range(3):
			ax = pyplot.subplot(len(activity_ids), 1, k+1, sharex=xaxis)
			ax.set_xlim(-1,1)
			if k == 0:
				xaxis = ax
			pyplot.hist(to_series(grouped[act_id][:,:,i]), bins=100)
			pyplot.title('activity '+str(act_id), y=0, loc='left')
	pyplot.show()

下面列出了完整的示例。

# plot histograms per activity for a subject
from numpy import array
from numpy import dstack
from numpy import unique
from pandas import read_csv
from matplotlib import pyplot

# load a single file as a numpy array
def load_file(filepath):
	dataframe = read_csv(filepath, header=None, delim_whitespace=True)
	return dataframe.values

# load a list of files, such as x, y, z data for a given variable
def load_group(filenames, prefix=''):
	loaded = list()
	for name in filenames:
		data = load_file(prefix + name)
		loaded.append(data)
	# stack group so that features are the 3rd dimension
	loaded = dstack(loaded)
	return loaded

# load a dataset group, such as train or test
def load_dataset(group, prefix=''):
	filepath = prefix + group + '/Inertial Signals/'
	# load all 9 files as a single array
	filenames = list()
	# total acceleration
	filenames += ['total_acc_x_'+group+'.txt', 'total_acc_y_'+group+'.txt', 'total_acc_z_'+group+'.txt']
	# body acceleration
	filenames += ['body_acc_x_'+group+'.txt', 'body_acc_y_'+group+'.txt', 'body_acc_z_'+group+'.txt']
	# body gyroscope
	filenames += ['body_gyro_x_'+group+'.txt', 'body_gyro_y_'+group+'.txt', 'body_gyro_z_'+group+'.txt']
	# load input data
	X = load_group(filenames, filepath)
	# load class output
	y = load_file(prefix + group + '/y_'+group+'.txt')
	return X, y

# get all data for one subject
def data_for_subject(X, y, sub_map, sub_id):
	# get row indexes for the subject id
	ix = [i for i in range(len(sub_map)) if sub_map[i]==sub_id]
	# return the selected samples
	return X[ix, :, :], y[ix]

# convert a series of windows to a 1D list
def to_series(windows):
	series = list()
	for window in windows:
		# remove the overlap from the window
		half = int(len(window) / 2) - 1
		for value in window[-half:]:
			series.append(value)
	return series

# group data by activity
def data_by_activity(X, y, activities):
	# group windows by activity
	return {a:X[y[:,0]==a, :, :] for a in activities}

# plot histograms for each activity for a subject
def plot_activity_histograms(X, y):
	# get a list of unique activities for the subject
	activity_ids = unique(y[:,0])
	# group windows by activity
	grouped = data_by_activity(X, y, activity_ids)
	# plot per activity, histograms for each axis
	pyplot.figure()
	xaxis = None
	for k in range(len(activity_ids)):
		act_id = activity_ids[k]
		# total acceleration
		for i in range(3):
			ax = pyplot.subplot(len(activity_ids), 1, k+1, sharex=xaxis)
			ax.set_xlim(-1,1)
			if k == 0:
				xaxis = ax
			pyplot.hist(to_series(grouped[act_id][:,:,i]), bins=100)
			pyplot.title('activity '+str(act_id), y=0, loc='left')
	pyplot.show()

# load data
trainX, trainy = load_dataset('train', 'HARDataset/')
# load mapping of rows to subjects
sub_map = load_file('HARDataset/train/subject_train.txt')
train_subjects = unique(sub_map)
# get the data for one subject
sub_id = train_subjects[0]
subX, suby = data_for_subject(trainX, trainy, sub_map, sub_id)
# plot data for subject
plot_activity_histograms(subX, suby)

运行该示例将创建包含六个子图的图，每个子图用于训练数据集中第一个主题的每个活动。总加速度数据的 x，y 和 z 轴中的每一个分别具有蓝色，橙色和绿色直方图。

我们可以看到每个活动都有不同的数据分布，大运动（前三个活动）与固定活动（最后三个活动）之间存在显着差异。前三个活动的数据分布看起来是高斯的，可能有不同的均值和标准偏差。后一活动的分布看起来是多模态的（即多个峰值）。

Histograms of the total acceleration data by activity

按活动计算的总加速度数据的直方图

我们可以使用plot_activity_histograms()的更新版本重新运行相同的示例，该版本将绘制车身加速度数据。

更新的功能如下所示。

# plot histograms for each activity for a subject
def plot_activity_histograms(X, y):
	# get a list of unique activities for the subject
	activity_ids = unique(y[:,0])
	# group windows by activity
	grouped = data_by_activity(X, y, activity_ids)
	# plot per activity, histograms for each axis
	pyplot.figure()
	xaxis = None
	for k in range(len(activity_ids)):
		act_id = activity_ids[k]
		# total acceleration
		for i in range(3):
			ax = pyplot.subplot(len(activity_ids), 1, k+1, sharex=xaxis)
			ax.set_xlim(-1,1)
			if k == 0:
				xaxis = ax
			pyplot.hist(to_series(grouped[act_id][:,:,3+i]), bins=100)
			pyplot.title('activity '+str(act_id), y=0, loc='left')
	pyplot.show()

运行更新的示例会创建一个新的图。

在这里，我们可以看到在动作与静止活动之间的活动中有更多类似的分布。在动态活动的情况下数据看起来是双峰的，并且在静止活动的情况下可能是高斯的或指数的。

我们通过活动看到的总体对体加速度分布的模式反映了我们在上一节中使用相同数据类型看到的对象。也许总加速度数据是区分活动的关键。

Histograms of the body acceleration data by activity

通过活动的身体加速度数据的直方图

最后，我们可以再次更新示例以绘制陀螺仪数据的每个活动的直方图。

更新的功能如下所示。

# plot histograms for each activity for a subject
def plot_activity_histograms(X, y):
	# get a list of unique activities for the subject
	activity_ids = unique(y[:,0])
	# group windows by activity
	grouped = data_by_activity(X, y, activity_ids)
	# plot per activity, histograms for each axis
	pyplot.figure()
	xaxis = None
	for k in range(len(activity_ids)):
		act_id = activity_ids[k]
		# total acceleration
		for i in range(3):
			ax = pyplot.subplot(len(activity_ids), 1, k+1, sharex=xaxis)
			ax.set_xlim(-1,1)
			if k == 0:
				xaxis = ax
			pyplot.hist(to_series(grouped[act_id][:,:,6+i]), bins=100)
			pyplot.title('activity '+str(act_id), y=0, loc='left')
	pyplot.show()

运行该示例将创建具有与身体加速度数据类似的模式的绘图，尽管可能显示胖尾高斯分布而不是动态活动的双峰分布。

Histograms of the body gyroscope data by activity

通过活动的身体陀螺仪数据的直方图

所有这些图都是为第一个主题创建的，我们期望在其他主题的活动中看到类似的运动数据分布和关系。

9.绘制活动持续时间箱图

需要考虑的最后一个方面是受试者在每项活动上花费的时间。

这与班级的平衡密切相关。如果活动（类）在数据集中通常是平衡的，那么我们期望特定主题在其跟踪过程中的活动平衡也将相当平衡。

我们可以通过计算每个主题在每个活动上花费的时间（样本或行数）并查看每个活动的持续时间分布来确认这一点。

审查这些数据的一种方便方法是将分布总结为箱线图，显示中位数（线），中间 50％（方框），数据的一般范围（四分位数间距）和异常值（以点为单位） 。

下面的函数plot_activity_durations_by_subject()通过首先按主题分割数据集，然后按活动分割主题数据并计算在每个活动上花费的行，然后最终创建持续时间测量的每个活动的箱线图来实现此行为。

# plot activity durations by subject
def plot_activity_durations_by_subject(X, y, sub_map):
	# get unique subjects and activities
	subject_ids = unique(sub_map[:,0])
	activity_ids = unique(y[:,0])
	# enumerate subjects
	activity_windows = {a:list() for a in activity_ids}
	for sub_id in subject_ids:
		# get data for one subject
		_, subj_y = data_for_subject(X, y, sub_map, sub_id)
		# count windows by activity
		for a in activity_ids:
			activity_windows[a].append(len(subj_y[subj_y[:,0]==a]))
	# organize durations into a list of lists
	durations = [activity_windows[a] for a in activity_ids]
	pyplot.boxplot(durations, labels=activity_ids)
	pyplot.show()

下面列出了完整的示例。

# plot durations of each activity by subject
from numpy import array
from numpy import dstack
from numpy import unique
from pandas import read_csv
from matplotlib import pyplot

# load a single file as a numpy array
def load_file(filepath):
	dataframe = read_csv(filepath, header=None, delim_whitespace=True)
	return dataframe.values

# load a list of files, such as x, y, z data for a given variable
def load_group(filenames, prefix=''):
	loaded = list()
	for name in filenames:
		data = load_file(prefix + name)
		loaded.append(data)
	# stack group so that features are the 3rd dimension
	loaded = dstack(loaded)
	return loaded

# load a dataset group, such as train or test
def load_dataset(group, prefix=''):
	filepath = prefix + group + '/Inertial Signals/'
	# load all 9 files as a single array
	filenames = list()
	# total acceleration
	filenames += ['total_acc_x_'+group+'.txt', 'total_acc_y_'+group+'.txt', 'total_acc_z_'+group+'.txt']
	# body acceleration
	filenames += ['body_acc_x_'+group+'.txt', 'body_acc_y_'+group+'.txt', 'body_acc_z_'+group+'.txt']
	# body gyroscope
	filenames += ['body_gyro_x_'+group+'.txt', 'body_gyro_y_'+group+'.txt', 'body_gyro_z_'+group+'.txt']
	# load input data
	X = load_group(filenames, filepath)
	# load class output
	y = load_file(prefix + group + '/y_'+group+'.txt')
	return X, y

# get all data for one subject
def data_for_subject(X, y, sub_map, sub_id):
	# get row indexes for the subject id
	ix = [i for i in range(len(sub_map)) if sub_map[i]==sub_id]
	# return the selected samples
	return X[ix, :, :], y[ix]

# convert a series of windows to a 1D list
def to_series(windows):
	series = list()
	for window in windows:
		# remove the overlap from the window
		half = int(len(window) / 2) - 1
		for value in window[-half:]:
			series.append(value)
	return series

# group data by activity
def data_by_activity(X, y, activities):
	# group windows by activity
	return {a:X[y[:,0]==a, :, :] for a in activities}

# plot activity durations by subject
def plot_activity_durations_by_subject(X, y, sub_map):
	# get unique subjects and activities
	subject_ids = unique(sub_map[:,0])
	activity_ids = unique(y[:,0])
	# enumerate subjects
	activity_windows = {a:list() for a in activity_ids}
	for sub_id in subject_ids:
		# get data for one subject
		_, subj_y = data_for_subject(X, y, sub_map, sub_id)
		# count windows by activity
		for a in activity_ids:
			activity_windows[a].append(len(subj_y[subj_y[:,0]==a]))
	# organize durations into a list of lists
	durations = [activity_windows[a] for a in activity_ids]
	pyplot.boxplot(durations, labels=activity_ids)
	pyplot.show()

# load training dataset
X, y = load_dataset('train', 'HARDataset/')
# load mapping of rows to subjects
sub_map = load_file('HARDataset/train/subject_train.txt')
# plot durations
plot_activity_durations_by_subject(X, y, sub_map)

运行该示例将创建六个箱形图，每个活动一个。

每个箱图总结了训练数据集中每个活动花费的时间（行数或窗口数）。

我们可以看到，受试者在静止活动（4,5 和 6）上花费的时间更多，在动作活动中花费的时间更少（1,2 和 3），3 的分布最小，或者时间花费最少。

活动的分布并不大，这表明不需要削减较长时间的活动或过度采样活动。尽管如果运动活动的预测模型的技能通常更差，这些方法仍然可用。

Boxplot of activity durations per subject on train set

训练组上每个受试者的活动持续时间的箱线图

我们可以使用以下附加行为训练数据创建类似的箱线图。

# load test dataset
X, y = load_dataset('test', 'HARDataset/')
# load mapping of rows to subjects
sub_map = load_file('HARDataset/test/subject_test.txt')
# plot durations
plot_activity_durations_by_subject(X, y, sub_map)

运行更新的示例显示了活动之间的类似关系。

这是令人鼓舞的，这表明测试和训练数据集确实可以合理地代表整个数据集。

Boxplot of activity durations per subject on test set

测试集上每个受试者的活动持续时间的箱线图

现在我们已经探索了数据集，我们可以建议一些关于如何建模的想法。

10.建模方法

在本节中，我们总结了一些建模活动识别数据集的方法。

这些想法分为项目的主题。

问题框架

第一个重要的考虑因素是预测问题的框架。

原始作品中描述的问题的框架是基于已知主体的运动数据和活动，根据运动数据预测新主体的活动。

我们可以总结为：

  • 给出一个移动数据窗口预测活动。

这是一个合理而有用的问题框架。

将提供的数据构建为预测问题的其他一些可能方法包括：

  • 在给定移动数据的时间步长的情况下预测活动。
  • 给出多个移动数据窗口的预测活动。
  • 给定多个移动数据窗口的预测活动序列。
  • 给出预分段活动的一系列移动数据的预测活动。
  • 给定移动数据的时间步长，预测活动停止或转换。
  • 给定移动数据窗口预测静止或非静止活动

其中一些框架可能过于具有挑战性或太容易。

然而，这些框架提供了探索和理解数据集的其他方法。

数据准备

在使用原始数据训练模型之前可能需要一些数据准备。

数据似乎已经缩放到范围[-1,1]。

在建模之前可以执行的一些其他数据转换包括：

  • 跨学科规范化。
  • 每个科目的标准化。
  • 跨学科标准化。
  • 轴功能选择。
  • 数据类型功能选择。
  • 信号异常检测和清除。
  • 删除过多代表活动的窗口。
  • 对代表性不足的活动进行过采样。
  • 将信号数据下采样到 1 / 4,1 / 2,1,2 或其他部分。

预测性建模

通常，问题是时间序列多分类问题。

正如我们所看到的，它也可以被构造为二分类问题和多步时间序列分类问题。

最初的论文探讨了在数据集版本中使用经典机器学习算法，其中从每个数据窗口设计了特征。具体地说，是一种改进的支持向量机

数据集的特征工程版本上的 SVM 结果可以提供问题表现的基线。

从这一点扩展，在该版本的数据集上评估多个线性，非线性和集合机器学习算法可以提供改进的基准。

问题的焦点可能是数据集的未设计或原始版本。

在这里，可以探索模型复杂性的进展，以确定最合适的问题模型;一些候选模型包括：

  • 常见的线性，非线性和集成机器学习算法。
  • 多层感知机。
  • 卷积神经网络，特别是 1D CNN。
  • 循环神经网络，特别是 LSTM。
  • CNN 和 LSTMS 的杂交体，例如 CNN-LSTM 和 ConvLSTM。

模型评估

在原始论文中对模型的评估涉及使用 70％和 30％比率的受试者对数据进行训练/测试分割。

对这种预先定义的数据分割的探索表明，这两组都能够合理地代表整个数据集。

另一种替代方法可以是每个受试者使用留一法交叉验证或 LOOCV。除了为每个受试者提供用作保留测试集的机会之外，该方法还将提供可以平均和总结的 30 个分数的群体，这可以提供更稳健的结果。

使用分类精度和混淆矩阵表示模型表现，两者都适用于预测问题的多类性质。

具体而言，混淆矩阵将有助于确定某些类别比其他类别更容易或更具挑战性，例如静止活动与涉及运动的活动相比。

进一步阅读

如果您希望深入了解，本节将提供有关该主题的更多资源。

文件

API

用品

摘要

在本教程中，您发现了使用智能手机数据集进行时间序列分类的活动识别，以及如何加载和浏览数据集以使其为预测性建模做好准备。

具体来说，你学到了：

  • 如何下载数据集并将其加载到内存中。
  • 如何使用线图，直方图和箱线图来更好地理解运动数据的结构。
  • 如何建模问题包括框架，数据准备，建模和评估。

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

avatar
文章
阅读
粉丝