【实验】基于随机森林的气温预测(上)

简介: 【实验】基于随机森林的气温预测

一、天气最高温度预测任务


实验数据以及源代码 https://github.com/wyn-365/temp-predict.git


1.1 任务


我们要完成三项任务: 使用随机森林算法完成基本建模任务

基本任务需要我们处理数据,观察特征,完成建模并进行可视化展示分析

观察数据量与特征个数对结果影响

在保证算法一致的前提下,加大数据个数,观察结果变换。重新考虑特征工程,引入新特征后观察结果走势。

对随机森林算法进行调参,找到最合适的参数

掌握机器学习中两种经典调参方法,对当前模型进行调节

# 数据读取
import pandas as pd
features = pd.read_csv('data/temps.csv')
features.head(5)

数据表中

year,moth,day,week分别表示的具体的时间

temp_2:前天的最高温度值

temp_1:昨天的最高温度值

average:在历史中,每年这一天的平均最高温度值

actual:这就是我们的标签值了,当天的真实最高温度

friend:这一列可能是凑热闹的,你的朋友猜测的可能值,咱们不管它就好了


1.2 数据大小

print('The shape of our features is:', features.shape)
# 统计指标
features.describe()
# 对于时间数据,我们也可以进行一些转换,目的就是有些工具包在绘图或者计算的过程中,需要标准的时间格式:
# 处理时间数据
import datetime
# 分别得到年,月,日
years = features['year']
months = features['month']
days = features['day']
# datetime格式
dates = [str(int(year)) + '-' + str(int(month)) + '-' + str(int(day)) for year, month, day in zip(years, months, days)]
dates = [datetime.datetime.strptime(date, '%Y-%m-%d') for date in dates]
dates[:5]

结果显示:The shape of our features is: (348, 9),表示我们的数据一共有348条记录,每个样本有9个特征。如果你想观察一下各个指标的统计特性,还可以用.describe()来直接展示一下:


1.3 数据展示

# 准备画图
import matplotlib.pyplot as plt
%matplotlib inline
# 指定默认风格
plt.style.use('fivethirtyeight')
# 接着我们设计画图的布局,这里我们需要展示4项指标,分别为最高气温的标签值,前天,昨天,朋友预测的气温最高值。既然是4个图,那不妨就2*2的规模来画吧,这样会更清晰一点,对每个图再指定好其名字和坐标轴含义就可以了:
# 设置布局
fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(nrows=2, ncols=2, figsize = (10,10))
fig.autofmt_xdate(rotation = 45)
# 标签值
ax1.plot(dates, features['actual'])
ax1.set_xlabel(''); ax1.set_ylabel('Temperature'); ax1.set_title('Max Temp')
# 昨天
ax2.plot(dates, features['temp_1'])
ax2.set_xlabel(''); ax2.set_ylabel('Temperature'); ax2.set_title('Previous Max Temp')
# 前天
ax3.plot(dates, features['temp_2'])
ax3.set_xlabel('Date'); ax3.set_ylabel('Temperature'); ax3.set_title('Two Days Prior Max Temp')
# 我的逗逼朋友
ax4.plot(dates, features['friend'])
ax4.set_xlabel('Date'); ax4.set_ylabel('Temperature'); ax4.set_title('Friend Estimate')
plt.tight_layout(pad=2)
# 各项指标看起来都还算正常,由于是国外的天气数据所以跟咱们的统计标准有些区别。接下来就要考虑数据预处理问题了,原始数据中在week列中并不是一些数值特征,而是表示周几的字符串,这些计算机可不认识,需要我们来转换一下:

1.5 数据预处理

# 独热编码
features = pd.get_dummies(features)
features.head(5)
# 这样就完成了数据集中属性值的预处理工作,默认会把所有属性值都转换成独热编码的格式,并且还帮我们自动添加了后缀看起来更清晰了,这里我们其实也可以按照自己的方式来设置编码特征的名字的,如果大家遇到了一个不太熟悉的函数,想看一下其中的细节,有一个更直接的方法就是在notebook当中直接调help工具来看一下它的API文档,下面返回的就是其细节介绍,不光有各个参数说明,还有一些小例子,建议大家在使用的过程中一定要养成多练多查的习惯,查找解决问题的方法也是一个很重要的技能:
print (help(pd.get_dummies))
print('Shape of features after one-hot encoding:', features.shape)
# 数据与标签
import numpy as np
# 标签
labels = np.array(features['actual'])
# 在特征中去掉标签
features= features.drop('actual', axis = 1)
# 名字单独保存一下,以备后患
feature_list = list(features.columns)
# 转换成合适的格式
features = np.array(features)

1.6 训练集与测试集

# 数据集切分
from sklearn.model_selection import train_test_split
train_features, test_features, train_labels, test_labels = train_test_split(features, labels, test_size = 0.25,
                                                                           random_state = 42)
print('训练集特征:', train_features.shape)
print('训练集标签:', train_labels.shape)
print('测试集特征:', test_features.shape)
print('测试集标签:', test_labels.shape)

1.7 建立一个基础的随机森林模型


万事俱备,我们可以来建立随机森林模型啦,首先导入工具包,先建立1000个树试试吧,其他参数先用默认值,之后我们会再深入到调参任务中

# 导入算法
from sklearn.ensemble import RandomForestRegressor
# 建模
rf = RandomForestRegressor(n_estimators= 1000, random_state=42)
# 训练
rf.fit(train_features, train_labels)
# 由于数据样本量还是非常小的,所以很快就可以得到结果了,这里我们先用MAPE指标来进行评估,也就是平均绝对百分误差,其实对于回归任务,评估方法还是比较多,给大家列出来几种,很简单就可以实现出来,也可以选择其他指标来进行评估:

1.8 测试

# 预测结果
predictions = rf.predict(test_features)
# 计算误差
errors = abs(predictions - test_labels)
# mean absolute percentage error (MAPE)
mape = 100 * (errors / test_labels)
print ('MAPE:',np.mean(mape))

1.9 特征重要性

# 得到特征重要性
importances = list(rf.feature_importances_)
# 转换格式
feature_importances = [(feature, round(importance, 2)) for feature, importance in zip(feature_list, importances)]
# 排序
feature_importances = sorted(feature_importances, key = lambda x: x[1], reverse = True)
# 对应进行打印
[print('Variable: {:20} Importance: {}'.format(*pair)) for pair in feature_importances]

1.10 用最重要的特征再来试试

# 选择最重要的那两个特征来试一试
rf_most_important = RandomForestRegressor(n_estimators= 1000, random_state=42)
# 拿到这俩特征
important_indices = [feature_list.index('temp_1'), feature_list.index('average')]
train_important = train_features[:, important_indices]
test_important = test_features[:, important_indices]
# 重新训练模型
rf_most_important.fit(train_important, train_labels)
# 预测结果
predictions = rf_most_important.predict(test_important)
errors = abs(predictions - test_labels)
# 评估结果
mape = np.mean(100 * (errors / test_labels))
print('mape:', mape)
# 转换成list格式
x_values = list(range(len(importances)))
# 绘图
plt.bar(x_values, importances, orientation = 'vertical')
# x轴名字
plt.xticks(x_values, feature_list, rotation='vertical')
# 图名
plt.ylabel('Importance'); plt.xlabel('Variable'); plt.title('Variable Importances'); 

1.11 预测值与真实值之间的差异


# 日期数据
months = features[:, feature_list.index('month')]
days = features[:, feature_list.index('day')]
years = features[:, feature_list.index('year')]
# 转换日期格式
dates = [str(int(year)) + '-' + str(int(month)) + '-' + str(int(day)) for year, month, day in zip(years, months, days)]
dates = [datetime.datetime.strptime(date, '%Y-%m-%d') for date in dates]
# 创建一个表格来存日期和其对应的标签数值
true_data = pd.DataFrame(data = {'date': dates, 'actual': labels})
# 同理,再创建一个来存日期和其对应的模型预测值
months = test_features[:, feature_list.index('month')]
days = test_features[:, feature_list.index('day')]
years = test_features[:, feature_list.index('year')]
test_dates = [str(int(year)) + '-' + str(int(month)) + '-' + str(int(day)) for year, month, day in zip(years, months, days)]
test_dates = [datetime.datetime.strptime(date, '%Y-%m-%d') for date in test_dates]
predictions_data = pd.DataFrame(data = {'date': test_dates, 'prediction': predictions}) 
# 真实值
plt.plot(true_data['date'], true_data['actual'], 'b-', label = 'actual')
# 预测值
plt.plot(predictions_data['date'], predictions_data['prediction'], 'ro', label = 'prediction')
plt.xticks(rotation = '60'); 
plt.legend()
# 图名
plt.xlabel('Date'); plt.ylabel('Maximum Temperature (F)'); plt.title('Actual and Predicted Values');
plt.show()


看起来还可以,这个走势我们的模型已经基本能够掌握了,接下来我们要再深入到数据中了,考虑几个问题: 1.如果可以利用的数据量增大,会对结果产生什么影响呢? 2.加入新的特征会改进模型效果吗?此时的时间效率又会怎样?


二、更多的数据效果会不会更好呢?


2.1 数据导入

# 导入工具包
import pandas as pd
# 读取数据
features = pd.read_csv('data/temps_extended.csv')
features.head(5)
print('数据规模',features.shape)
# 统计指标
round(features.describe(), 2)

数据规模 (2191, 12)

新的数据中,数据规模发生了变化,数据量扩充到了2191条并且加入了新的天气指标:

ws_1:前一天的风速

prcp_1: 前一天的降水

snwd_1:前一天的积雪深度

既然有了新的特征,先来看看他们长什么样吧,同样的方式绘制就可以了:


2.2 数据时间处理

# 转换成标准格式
import datetime
# 得到各种日期数据
years = features['year']
months = features['month']
days = features['day']
# 格式转换
dates = [str(int(year)) + '-' + str(int(month)) + '-' + str(int(day)) for year, month, day in zip(years, months, days)]
dates = [datetime.datetime.strptime(date, '%Y-%m-%d') for date in dates]
# 绘图
import matplotlib.pyplot as plt
%matplotlib inline
# 风格设置
plt.style.use('fivethirtyeight')
# Set up the plotting layout
fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(nrows=2, ncols=2, figsize = (15,10))
fig.autofmt_xdate(rotation = 45)
# Actual max temperature measurement
ax1.plot(dates, features['actual'])
ax1.set_xlabel(''); ax1.set_ylabel('Temperature (F)'); ax1.set_title('Max Temp')
# Temperature from 1 day ago
ax2.plot(dates, features['temp_1'])
ax2.set_xlabel(''); ax2.set_ylabel('Temperature (F)'); ax2.set_title('Prior Max Temp')
# Temperature from 2 days ago
ax3.plot(dates, features['temp_2'])
ax3.set_xlabel('Date'); ax3.set_ylabel('Temperature (F)'); ax3.set_title('Two Days Prior Max Temp')
# Friend Estimate
ax4.plot(dates, features['friend'])
ax4.set_xlabel('Date'); ax4.set_ylabel('Temperature (F)'); ax4.set_title('Friend Estimate')
plt.tight_layout(pad=2)
# 设置整体布局
fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(nrows=2, ncols=2, figsize = (15,10))
fig.autofmt_xdate(rotation = 45)
# 平均最高气温
ax1.plot(dates, features['average'])
ax1.set_xlabel(''); ax1.set_ylabel('Temperature (F)'); ax1.set_title('Historical Avg Max Temp')
# 风速
ax2.plot(dates, features['ws_1'], 'r-')
ax2.set_xlabel(''); ax2.set_ylabel('Wind Speed (mph)'); ax2.set_title('Prior Wind Speed')
# 降水
ax3.plot(dates, features['prcp_1'], 'r-')
ax3.set_xlabel('Date'); ax3.set_ylabel('Precipitation (in)'); ax3.set_title('Prior Precipitation')
# 积雪
ax4.plot(dates, features['snwd_1'], 'ro')
ax4.set_xlabel('Date'); ax4.set_ylabel('Snow Depth (in)'); ax4.set_title('Prior Snow Depth')
plt.tight_layout(pad=2)

在数据分析和特征提取的过程中,我们的出发点都是尽可能多的选择有价值的特征,因为其实阶段我们能得到的信息越多,之后建模可以利用的信息也是越多的,比如在这份数据中,我们有完整日期数据,但是显示天气的变换肯定是跟季节因素有关的,但是在原始数据集中并没有体现出季节的指标,我们可以自己创建一个季节变量当做新的特征,无论是对之后建模还是分析都会起到帮助的: 有了季节特征之后,假如我想观察一下不同季节的时候上述各项指标的变换情况该怎么做呢?这里给大家推荐一个非常实用的绘图函数pairplot,需要我们先安装seaborn这个工具包,它相当于是在Matplotlib的基础上进行封装,说白了就是用起来更简单规范了:


2.3 画图

# 创建一个季节变量
seasons = []
for month in features['month']:
    if month in [1, 2, 12]:
        seasons.append('winter')
    elif month in [3, 4, 5]:
        seasons.append('spring')
    elif month in [6, 7, 8]:
        seasons.append('summer')
    elif month in [9, 10, 11]:
        seasons.append('fall')
# 有了季节我们就可以分析更多东西了
reduced_features = features[['temp_1', 'prcp_1', 'average', 'actual']]
reduced_features['season'] = seasons
# 导入seaborn工具包
import seaborn as sns
sns.set(style="ticks", color_codes=True);
# 选择你喜欢的颜色模板
palette = sns.xkcd_palette(['dark blue', 'dark green', 'gold', 'orange'])
# 绘制pairplot
sns.pairplot(reduced_features, hue = 'season', diag_kind = 'kde', palette= palette, plot_kws=dict(alpha = 0.7),
                   diag_kws=dict(shade=True)); 

可以看到,x轴和y轴都是我们这4项指标,不同颜色的点表示不同的季节,在主对角线上x轴和y轴都是相同特征表示其在不同季节时的数值分布情况,其他位置用散点图来表示两个特征之间的关系,例如在左下角temp_1和actual就呈现出了很强的相关性。


2.4 数据预处理

# 独热编码
features = pd.get_dummies(features)
# 提取特征和标签
labels = features['actual']
features = features.drop('actual', axis = 1)
# 特征名字留着备用
feature_list = list(features.columns)
# 转换成所需格式
import numpy as np
features = np.array(features)
labels = np.array(labels)
# 数据集切分
from sklearn.model_selection import train_test_split
train_features, test_features, train_labels, test_labels = train_test_split(features, labels, 
                                                                            test_size = 0.25, random_state = 0)
print('Training Features Shape:', train_features.shape)
print('Training Labels Shape:', train_labels.shape)
print('Testing Features Shape:', test_features.shape)
print('Testing Labels Shape:', test_labels.shape)
# 新的训练集由1643个样本组成,测试集有548个样本。先来看看老数据集的情况,这里由于我们重新打开了一个新的notebook,所有代码中重新读取了老温度数据集,并进行了相同的预处理操作:

2.5 先看看老数据的结果

# 工具包导入
import pandas as pd
# 为了剔除特征个数对结果的影响,这里特征统一只有老数据集中特征
original_feature_indices = [feature_list.index(feature) for feature in
                                      feature_list if feature not in
                                      ['ws_1', 'prcp_1', 'snwd_1']]
# 读取老数据集
original_features = pd.read_csv('data/temps.csv')
original_features = pd.get_dummies(original_features)
import numpy as np
# 数据和标签转换
original_labels = np.array(original_features['actual'])
original_features= original_features.drop('actual', axis = 1)
original_feature_list = list(original_features.columns)
original_features = np.array(original_features)
# 数据集切分
from sklearn.model_selection import train_test_split
original_train_features, original_test_features, original_train_labels, original_test_labels = train_test_split(original_features, original_labels, test_size = 0.25, random_state = 42)
# 同样的树模型进行建模
from sklearn.ensemble import RandomForestRegressor
# 同样的参数与随机种子
rf = RandomForestRegressor(n_estimators= 100, random_state=0)
# 这里的训练集使用的是老数据集的
rf.fit(original_train_features, original_train_labels);
# 为了测试效果能够公平,统一使用一致的测试集,这里选择了刚刚我切分过的新数据集的测试集
predictions = rf.predict(test_features[:,original_feature_indices])
# 先计算温度平均误差
errors = abs(predictions - test_labels)
print('平均温度误差:', round(np.mean(errors), 2), 'degrees.')
# MAPE
mape = 100 * (errors / test_labels)
# 这里的Accuracy为了方便观察,我们就用100减去误差了,希望这个值能够越大越好
accuracy = 100 - np.mean(mape)
print('Accuracy:', round(accuracy, 2), '%.')

可以看到,当我们把数据量增大之后,效果发生了一些提升,这也符合实际情况,在机器学习任务中,我们都是希望数据量能够越大越好,这样可利用的信息就更多了。 下面我们要再对比一下特征数量对结果的影响,之前这两次比较还没有加入新的特征,这回我们把降水,风速,积雪3特征加入训练集中,看看效果又会怎样:

目录
相关文章
|
6月前
|
机器学习/深度学习 图计算
R语言广义线性模型(GLM)、全子集回归模型选择、检验分析全国风向气候数据(2)
R语言广义线性模型(GLM)、全子集回归模型选择、检验分析全国风向气候数据(2)
|
5月前
|
机器学习/深度学习 算法 数据可视化
使用Logistic算法从疝气病症预测病马的死亡率
使用Logistic算法从疝气病症预测病马的死亡率
65 1
|
6月前
|
机器学习/深度学习
R语言广义线性模型(GLM)、全子集回归模型选择、检验分析全国风向气候数据(1)
R语言广义线性模型(GLM)、全子集回归模型选择、检验分析全国风向气候数据
|
6月前
|
机器学习/深度学习 算法 数据挖掘
R语言气象模型集成预报:神经网络、回归、svm、决策树用环流因子预测降雨降水数据
R语言气象模型集成预报:神经网络、回归、svm、决策树用环流因子预测降雨降水数据
|
6月前
|
机器学习/深度学习 算法 Serverless
数据分享|R语言武汉流动人口趋势预测:灰色模型GM(1,1)、ARIMA时间序列、logistic逻辑回归模型
数据分享|R语言武汉流动人口趋势预测:灰色模型GM(1,1)、ARIMA时间序列、logistic逻辑回归模型
|
6月前
|
机器学习/深度学习 图计算
R语言广义线性模型(GLM)、全子集回归模型选择、检验分析全国风向气候数据
R语言广义线性模型(GLM)、全子集回归模型选择、检验分析全国风向气候数据
|
6月前
|
机器学习/深度学习
R语言用线性回归模型预测空气质量臭氧数据
R语言用线性回归模型预测空气质量臭氧数据
|
6月前
|
数据可视化 计算机视觉
用回归和主成分分析PCA 回归交叉验证分析预测城市犯罪率数据
用回归和主成分分析PCA 回归交叉验证分析预测城市犯罪率数据
|
6月前
|
存储 数据库
R语言用多项式回归和ARIMA模型预测电力负荷时间序列数据
R语言用多项式回归和ARIMA模型预测电力负荷时间序列数据
|
6月前
|
存储 数据可视化 数据库
R语言泊松Poisson回归模型预测人口死亡率和期望寿命
R语言泊松Poisson回归模型预测人口死亡率和期望寿命