【作者主页】Francek Chen
【专栏介绍】 ⌈ ⌈ ⌈机器学习与数据挖掘实战 ⌋ ⌋ ⌋ 机器学习是人工智能的一个分支，专注于让计算机系统通过数据学习和改进。它利用统计和计算方法，使模型能够从数据中自动提取特征并做出预测或决策。数据挖掘则是从大型数据集中发现模式、关联和异常的过程，旨在提取有价值的信息和知识。机器学习为数据挖掘提供了强大的分析工具，而数据挖掘则是机器学习应用的重要领域，两者相辅相成，共同推动数据科学的发展。本专栏介绍机器学习与数据挖掘的相关实战案例。
【GitCode】专栏资源保存在我的GitCode仓库：https://gitcode.com/Morse_Chen/ML-DM_cases。
文章目录
- 一、背景与挖掘目标
- （一）背景
- （二）目标分析
- 二、数据探索分析
- 三、数据预处理
- （一）属性归约
- （二）划分用水事件
- （三）确定单次用水事件时长阈值
- （四）属性构造
- （五）筛选候选洗浴事件
- 四、分析与建模
- 五、模型评估与检验
- 小结
一、背景与挖掘目标
（一）背景
随着国内大家电品牌的进入和国外品牌的涌入，电热水器相关技术在过去20年间得到了快速发展，屡屡创新。首次提出封闭式电热水器的概念到水电分离技术的研发。漏电保护技术的应用及出水断电技术和防电墙技术专利的申请突破。如今高效能技术颠覆了业内对电热水器“高能耗”的认知。
当下的热水器行业也并非一片太平盛世，行业内正在上演一幕弱肉强食的“丛林法则”戏码，市场份额逐步向龙头企业集中，尤其是那些在资金、渠道和品牌影响力等方面拥有实力的综合家电品类巨头，正在不断蚕食鲸吞市场蛋糕。要想在该行业立足，只能走产品差异化的路线，提升技术实力和产品质量，在功能卖点、外观等方面做出自身特色。
在热水器用户行为分析过程中，用水事件识别是最为关键的环节。根据该热水器生产厂商提供的数据热水器用户用水事件划分与识别项目的整体目标如下。
- 根据热水器采集到的数据，划分一次完整用水事件。
- 在划分好的一次完整用水事件中，识别出洗浴事件。
（二）目标分析
热水器用户用水事件划分与识别主要包括以下5个步骤。
- 对热水用户的历史用水数据进行选择性抽取，构建专家样本。
- 对步骤（1）形成的数据集进行数据探索分析与预处理，包括探索水流量的分布情况，删除冗余属性，识别用水数据的缺失值，并对缺失值作处理，根据建模的需要进行属性构造等。根据以上处理，对用水样本数据建立用水事件时间间隔识别模型和划分一次完整的用水事件模型，再在一次完整用水事件划分结果的基础上，剔除短暂用水事件缩小识别范围等。
- 在步骤2得到的建模样本数据基础上，建立洗浴事件识别模型，对洗浴事件识别模型进行模型分析评价。
- 对步骤3形成的模型结果应用并对洗浴事件划分进行优化。
- 调用洗浴事件识别模型，对实时监控的热水器流水数据进行洗浴事件自动识别。
热水器用户用水事件划分与识别的总体流程如图所示。
二、数据探索分析
在热水器的使用过程中，热水器的状态会经常发生改变，比如开机和关机、由加热转到保温、由无水流到有水流、水温由50℃变为49℃等。而智能热水器在状态发生改变或者水流量非零时，每两秒会采集一条状态数据。由于数据的采集频率较高，并且数据来自大量用户，数据总量非常大。
对原始数据采用无放回随机抽样法抽取200家热水器用户从2014年1月1日至2014年12月31日的用水记录作为原始建模数据。由于用户不仅使用热水器来洗浴，而且包括了洗手、洗脸、刷牙、洗菜、做饭等用水行为，所以热水器采集到的数据来自各种不同的用水事件。
热水器采集的用水数据包含12个属性：热水器编码，发生时间，开关机状态，加热中，保温中，有无水流，实际温度，热水量，水流量，节能模式，加热剩余时间和当前设置温度。其解释说明下表所示。
| 属性名称 | 说明 |
|---|---|
| 热水器编码 | 热水器出厂编号 |
| 发生时间 | 记录热水器处于某状态的时刻 |
| 开关机状态 | 热水器是否开机 |
| 加热中 | 热水器处于对水进行加热的状态 |
| 保温中 | 热水器处于对水进行保温的状态 |
| 有无水流 | 热水水流量大于等于 10L/min 为有水，否则为无 |
| 实际温度 | 热水器中热水的实际温度 |
| 热水量 | 热水器热水的含量 |
| 水流量 | 热水器热水的水流速度，单位：L/min |
| 节能模式 | 热水器的一种节能工作模式 |
| 加热剩余时间 | 加热到设定温度还需多长时间 |
| 当前设置温度 | 热水器加热时热水能够到达的最大温度 |
探索分析热水器的水流量状况，其中有无水流和水流量属性最能直观体现热水器的水流量情况，对这两个属性进行探索分析，得到不同水流状态的记录的条形图，无水流状态的记录明显比有水流状态的记录要多。
import pandas as pd
import matplotlib.pyplot as pltinputfile = '../data/original_data.xls' # 输入的数据文件
data = pd.read_excel(inputfile) # 读取数据# 查看有无水流的分布
# 数据提取
lv_non = pd.value_counts(data['有无水流'])['无']
lv_move = pd.value_counts(data['有无水流'])['有']# 绘制条形图
fig = plt.figure(figsize = (6 ,5)) # 设置画布大小
plt.rcParams['font.sans-serif'] = 'SimHei' # 设置中文显示
plt.rcParams['axes.unicode_minus'] = False
plt.bar(x=range(2), height=[lv_non,lv_move], width=0.4, alpha=0.8, color='skyblue')
plt.xticks([index for index in range(2)], ['无','有'])
plt.xlabel('水流状态')
plt.ylabel('记录数')
plt.title('不同水流状态记录数')
plt.show()
plt.close()
不同水流状态的记录的箱线图，箱体贴近0，说明无水流量的记录较多，水流量的分布与水流状态的分布一致。
# 查看水流量分布
water = data['水流量']
# 绘制水流量分布箱型图
fig = plt.figure(figsize = (5 ,8))
plt.boxplot(water, patch_artist=True,labels = ['水流量'], # 设置x轴标题boxprops = {'facecolor':'lightblue'}) # 设置填充颜色
plt.title('水流量分布箱线图')
# 显示y坐标轴的底线
plt.grid(axis='y')
plt.show()
用水停顿时间间隔定义为一条水流量不为0的流水记录同下一条水流量不为0的流水记录之间的时间间隔。根据现场实验统计，两次用水的过程的用水停顿的间隔时长一般在不大于4分钟。为了探究用户真实用水停顿时间间隔的分布情况，统计用水停顿的时间间隔并做频率分布表。通过频率分布表分析用户用水停顿时间间隔的规律性，具体的数据如下表所示。
用水停顿时间间隔频数分布表：
| 间隔时长 | 0~0.1 | 0.1~0.2 | 0.2~0.3 | 0.3~0.5 | 0.5~1 | 1~2 | 2~3 | 3~4 | 4~5 |
|---|---|---|---|---|---|---|---|---|---|
| 停顿频率 | 78.71% | 9.55% | 2.52% | 1.49% | 1.46% | 1.29% | 0.74% | 0.48% | 0.26% |
| 间隔时长 | 5~6 | 6~7 | 7~8 | 8~9 | 9~10 | 10~11 | 11~12 | 12~13 | 13以上 |
|---|---|---|---|---|---|---|---|---|---|
| 停顿频率 | 0.27% | 0.19% | 0.17% | 0.12% | 0.09% | 0.09% | 0.10% | 0.11% | 2.36% |
停顿时间间隔为0~0.3分钟的频率很高，根据日常用水经验可以判断其为一次用水时间中的停顿；停顿时间间隔为6~13分钟的频率较低，分析其为两次用水事件之间的停顿间隔。两次用水事件的停顿时间间隔分布在3~7分钟。根据现场实验统计用水停顿的时间间隔近似。
三、数据预处理
（一）属性归约
由于热水器采集的用水数据属性较多，做以下处理，因分析的主要对象为用户，分析的主要目标为用户的洗浴行为的一般规律，所以“热水器编号”属性可以去除；因热水器采集的数据中，“有无水流”属性可以通过“水流量”属性反应，“节能模式”属性取值相同均为“关”，对分析无作用，可以去除。
import pandas as pd
import numpy as np
data = pd.read_excel('../data/original_data.xls')
print('初始状态的数据形状为：', data.shape)
# 删除热水器编号、有无水流、节能模式属性
data.drop(labels=["热水器编号","有无水流","节能模式"],axis=1,inplace=True)
print('删除冗余属性后的数据形状为：', data.shape)
data.to_csv('../tmp/water_heart.csv',index=False)
删除冗余属性“热水器编号”“有无水流”“节能模式”，删除冗余属性后得到用来建模的属性（../tmp/water_heart.csv）如下表所示。
| 发生时间 | 开关机状态 | 加热中 | 保温中 | 实际温度 | 热水量 | 水流量 | 加热剩余时间 | 当前设置温度 |
|---|---|---|---|---|---|---|---|---|
| 20141019161042 | 开 | 开 | 关 | 48°C | 25% | 0 | 1分钟 | 50°C |
| 20141019161106 | 开 | 开 | 关 | 49°C | 25% | 0 | 1分钟 | 50°C |
| 20141019161147 | 开 | 开 | 关 | 49°C | 25% | 0 | 0分钟 | 50°C |
| 20141019161149 | 开 | 关 | 开 | 50°C | 100% | 0 | 0分钟 | 50°C |
| 20141019172319 | 开 | 关 | 开 | 50°C | 50% | 0 | 0分钟 | 50°C |
| 20141019172321 | 关 | 关 | 关 | 50°C | 50% | 62 | 0分钟 | 50°C |
| 20141019172323 | 关 | 关 | 关 | 50°C | 50% | 63 | 0分钟 | 50°C |
（二）划分用水事件
用户的用水数据存储在数据库中，记录了各种各样的用水事件，包括洗浴，洗手，刷牙，洗脸，洗衣，洗菜等，而且一次用水事件由数条甚至数千条的状态记录组成。首先需要在大量的状态记录中划分出哪些连续的数据是一次完整的用水事件。
用水状态记录中，水流量不为0表明用户正在使用热水；而水流量为0时用户用热水发生停顿或者用热水结束。对于任一个用水记录，如果它的向前时差超过阈值 T T T，则将它记为事件的开始编号；如果向后时差超过阈值 T T T，则将其记为事件的结束编号。划分模型的符号说明如下表所示。
| 符号 | 释义 |
|---|---|
| t 1 t_1 t1 | 所有水流量不为0的用水行为的发生时间 |
| T T T | 时间间隔阈值 |
一次完整用水事件的划分步骤如下。
- 读取数据记录，识别到所有水流量不为0的状态记录，将它们的发生时间记为序列 t 1 t_1 t1。
- 对序列 t 1 t_1 t1构建其向前时差列和向后时差列，并分别与阈值进行比较。向前时差超过阈值 T T T，则将它记为新的用水事件的开始编号；如果向后时差超过阈值 T T T，则将其记为用水事件的结束编号。
循环执行步骤（2）直到向前时差列和向后时差列与均值比较完毕，结束事件划分。
用水事件划分主要分为两个步骤：
- 确定单次用水时长间隔
- 计算两条相邻记录的时间
# 读取数据
data = pd.read_csv('../tmp/water_heart.csv')
# 划分用水事件
threshold = pd.Timedelta('4 min') # 阈值为4分钟
data['发生时间'] = pd.to_datetime(data['发生时间'], format = '%Y%m%d%H%M%S') # 转换时间格式
data = data[data['水流量'] > 0] # 只要流量大于0的记录
sjKs = data['发生时间'].diff() > threshold # 相邻时间向前差分，比较是否大于阈值
sjKs.iloc[0] = True # 令第一个时间为第一个用水事件的开始事件
sjJs = sjKs.iloc[1:] # 向后差分的结果
sjJs = pd.concat([sjJs,pd.Series(True)]) # 令最后一个时间作为最后一个用水事件的结束时间
# 创建数据框，并定义用水事件序列
sj = pd.DataFrame(np.arange(1,sum(sjKs)+1),columns = ["事件序号"])
sj["事件起始编号"] = data.index[sjKs == 1]+1 # 定义用水事件的起始编号
sj["事件终止编号"] = data.index[sjJs == 1]+1 # 定义用水事件的终止编号
print('当阈值为4分钟的时候事件数目为：',sj.shape[0])
sj.to_csv('../tmp/sj.csv',index = False)
对用户的用水数据进行划分结果（../tmp/sj.csv）如下表所示。
| 发生时间 | …… | 事件编号 | |
|---|---|---|---|
| 2 | 2014-10-19 07:01:56 | …… | 1 |
| 56 | 2014-10-19 07:38:16 | …… | 2 |
| 381 | 2014-10-19 09:46:38 | …… | 3 |
| 382 | 2014-10-19 09:46:40 | …… | 3 |
| 384 | 2014-10-19 09:47:15 | …… | 3 |
| 404 | 2014-10-19 11:50:17 | …… | 4 |
| …… | …… | …… | …… |
（三）确定单次用水事件时长阈值
对某热水器用户的数据根据不同的阈值划分用水事件，得到了相应的事件个数，阈值变化与划分得到事件个数如下表所示。
| 阈值（分钟） | 2.25 | 2.5 | 2.75 | 3 | 3.25 | 3.5 | 3.75 | 4 | 4.25 | 4.5 | 4.75 | 5 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 事件个数 | 650 | 644 | 626 | 602 | 588 | 565 | 533 | 530 | 530 | 530 | 522 | 520 |
| 阈值（分钟） | 5.25 | 5.5 | 5.75 | 6 | 6.25 | 6.5 | 6.75 | 7 | 7.25 | 7.5 | 7.75 | 8 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 事件个数 | 510 | 506 | 503 | 500 | 480 | 472 | 466 | 462 | 460 | 460 | 460 | 460 |
阈值与划分事件个数关系如下图所示。
上图为阈值与划分事件个数的散点图，图中某段阈值范围内，下降趋势明显，说明在该段阈值范围内，用户的停顿习惯比较集中。如果趋势比较平缓，则说明用户的停顿热水的习惯趋于稳定，所以取该段时间开始的时间点作为阈值，既不会将短的用水事件合并，又不会将长的用水事件拆开。
用户停顿热水的习惯在方框的位置趋于稳定，说明热水器用户的用水的停顿习惯用方框开始的时间点作为划分阈值会有一个好的效果。
曲线在上图（散点图）中方框趋于稳定时，其方框开始的点的斜率趋于一个较小的值。为了用程序来识别这一特征，将这一特征提取为规则。根据用水停顿时间间隔频数分布表说明如何识别上图中的方框中起始的时间。
每个阈值对应一个点，给每个阈值计算得到一个斜率指标，如上图所示，其中， 点是要计算的斜率指标点。为了直观的展示，如下表所示。
| 符号名称 | 符号说明 |
|---|---|
| k A i k_{Ai} kAi | 与 i i i点的斜率的绝对值， i ∈ { B , C , D , E } i ∈ \{B, C, D, E\} i∈{B,C,D,E} |
| k k k | 任意两点 ( x 1 , y 1 ) (x_1, y_1) (x1,y1), ( x 2 , y 2 ) (x_2, y_2) (x2,y2)的斜率的绝对值 |
| k k k | 五个点的斜率之和的平均值 |
| ( x i , y i ) (x_i, y_i) (xi,yi) | i i i点的坐标， i ∈ { B , C , D , E } i ∈ \{B, C, D, E\} i∈{B,C,D,E} |
首先计算出 k A B 、 k A C 、 k A D 、 k A E k_{AB}、k_{AC}、k_{AD}、k_{AE} kAB、kAC、kAD、kAE四个斜率。然后可以计算出四个斜率之和的平均值 K K K。
K = y 1 − y 2 x 1 − x 2 K = \frac{y_1 - y_2}{x_1 - x_2} K=x1−x2y1−y2 K = k A B + k A C + k A D + k A E 4 K = \frac{k_{AB} + k_{AC} + k_{AD} + k_{AE}}{4} K=4kAB+kAC+kAD+kAE
将 K K K作为 A A A点的斜率指标。特别指出，横坐标上的最后四个点没有斜率指标，因为找不出在它以后的4个更长的阈值。但这不影响对最优阈值的寻找，因为可以提高阈值的上限，以使最后的4个阈值不在考虑范围内。
阈值优化的结果如下：
- 当存在一个阈值的斜率指标 K < 1 K<1 K<1时，取阈值最小的点 A A A（可能存在多个阈值的斜率指标小于1）的横坐标 x A x_A xA作为用水事件划分的阈值。其中， K < 1 K<1 K<1中的“1”是经过实际数据验证的一个专家阈值。
- 当不存在 K < 1 K<1 K<1时，则找所有阈值中斜率指标最小的阈值；如果该阈值的斜率指标小于5，则取该阈值作为用水事件划分的阈值；如果该阈值的斜率指标不小于5，则阈值取默认值4分钟。其中，斜率指标小于5中的“5”是经过实际数据验证的一个专家阈值。
# 确定单次用水事件时长阈值
n = 4 # 使用以后四个点的平均斜率
threshold = pd.Timedelta(minutes=5) # 专家阈值
data['发生时间'] = pd.to_datetime(data['发生时间'], format='%Y%m%d%H%M%S')
data = data[data['水流量'] > 0] # 只要流量大于0的记录
# 自定义函数：输入划分时间的时间阈值，得到划分的事件数
def event_num(ts):d = data['发生时间'].diff() > ts # 相邻时间作差分，比较是否大于阈值return d.sum() + 1 # 这样直接返回事件数
dt = [pd.Timedelta(minutes=i) for i in np.arange(1, 9, 0.25)]
h = pd.DataFrame(dt, columns=['阈值']) # 转换数据框，定义阈值列
h['事件数'] = h['阈值'].apply(event_num) # 计算每个阈值对应的事件数
h['斜率'] = h['事件数'].diff()/0.25 # 计算每两个相邻点对应的斜率
h['斜率指标']= h['斜率'].abs().rolling(4).mean() # 往前取n个斜率绝对值平均作为斜率指标
ts = h['阈值'][h['斜率指标'].idxmin() - n]
# 用idxmin返回最小值的Index，由于rolling_mean()计算的是前n个斜率的绝对值平均
# 所以结果要进行平移（-n）
if ts > threshold:ts = pd.Timedelta(minutes=4)
print('计算出的单次用水时长的阈值为：',ts)
（四）属性构造
1. 构建用水时长与频率属性
不同用水事件的用水时长是基础特征之一。根据用水时长这一特征可以构建下表所示的事件开始时间、事件结束时间、洗浴时间点、用水时长、总用水时长和用水时长/总用水时长这6个特征。
| 特征 | 构建方法 | 说明 |
|---|---|---|
| 事件开始时间 | 事件开始时间 = (起始数据的时间 - 发送阈值) / 2 | 热水事件开始发生的时间。 |
| 事件结束时间 | 事件结束时间 = (结束数据的时间 + 发送阈值) / 2 | 热水事件结束发生的时间。 |
| 洗浴时间点 | 洗浴时间点 = 事件开始时间的小时点，如时间为“20:00:10”，则洗浴时间点为“20” | 开始用水的时间点。 |
| 用水时长 | 用水时长 = 每条用水数据时长的和 = 和下条数据的间隔时间 / 2 + 和下条数据的间隔时间 / 2。 | 一次用水过程中有热水流出的时长。 |
| 总用水时长 | 从划分出的用水事件，起始数据的时间到终止数据的时间间隔 + 发送阈值。 | 记录整个用水阶段的时长。 |
| 用水时长 / 总用水时长 | 用水时长 / 总用水时长 | 判断用水时长占总用水时长的比重。 |
构建用水开始时间或结束的时间两个特征时分别减去或加上了发送阈值（发送阈值是指热水器传输数据的频率的大小）。其原因如下，取平均值会导致很大的偏差。综合分析构建“用水开始时间”为起始数据的时间减去“发送阈值”的一半。
用水时长相关的特征只能够区分出一部分用水事件，不同用水事件的用水停顿和频率也不同。主要用水频率类属性构建说明。
| 特征 | 构建方法 | 说明 |
|---|---|---|
| 停顿时长 | 一次完整用水事件中，对水流量为 0 的数据做计算 | 停顿时长 = 每条用水停顿数据时长的和 = （和下条数据的间隔时间 / 2 + 和上条数据的间隔时间 / 2）的和。标记一次完整用水事件中的每次用水停顿的时长。 |
| 总停顿时长 | 一次完整用水事件中的所有停顿时长之和 | 标记一次完整用水事件中的总停顿时长。 |
| 平均停顿时长 | 一次完整用水事件中的所有停顿时长的平均值 | 标记一次完整用水事件中的停顿的平均时长。 |
| 停顿次数 | 一次完整用水事件中断用水的次数之和 | 帮助识别洗浴及连续洗浴事件。 |
data = pd.read_excel('../data/water_hearter.xlsx',encoding='gbk') # 读取热水器使用数据记录
sj = pd.read_csv('../tmp/sj.csv') # 读取用水事件记录
# 转换时间格式
data["发生时间"] = pd.to_datetime(data["发生时间"],format="%Y%m%d%H%M%S")# 构造特征：总用水时长
timeDel = pd.Timedelta("0.5 sec")
sj["事件开始时间"] = data.iloc[sj["事件起始编号"]-1,0].values- timeDel
sj["事件结束时间"] = data.iloc[sj["事件终止编号"]-1,0].values + timeDel
sj['洗浴时间点'] = [i.hour for i in sj["事件开始时间"]]
sj["总用水时长"] = np.int64(sj["事件结束时间"] - sj["事件开始时间"])/1000000000 + 1# 构造用水停顿事件
# 构造特征“停顿开始时间”、“停顿结束时间”
# 停顿开始时间指从有水流到无水流，停顿结束时间指从无水流到有水流
for i in range(len(data)-1):if (data.loc[i,"水流量"] != 0) & (data.loc[i + 1,"水流量"] == 0) :data.loc[i + 1,"停顿开始时间"] = data.loc[i +1, "发生时间"] - timeDelif (data.loc[i,"水流量"] == 0) & (data.loc[i + 1,"水流量"] != 0) :data.loc[i,"停顿结束时间"] = data.loc[i , "发生时间"] + timeDel# 提取停顿开始时间与结束时间所对应行号，放在数据框Stop中
indStopStart = data.index[data["停顿开始时间"].notnull()]+1
indStopEnd = data.index[data["停顿结束时间"].notnull()]+1
Stop = pd.DataFrame(data={"停顿开始编号":indStopStart[:-1],"停顿结束编号":indStopEnd[1:]})
# 计算停顿时长，并放在数据框stop中，停顿时长=停顿结束时间-停顿结束时间
Stop["停顿时长"] = np.int64(data.loc[indStopEnd[1:]-1,"停顿结束时间"].values-data.loc[indStopStart[:-1]-1,"停顿开始时间"].values)/1000000000
# 将每次停顿与事件匹配,停顿的开始时间要大于事件的开始时间，
# 且停顿的结束时间要小于事件的结束时间
for i in range(len(sj)):Stop.loc[(Stop["停顿开始编号"] > sj.loc[i,"事件起始编号"]) & (Stop["停顿结束编号"] < sj.loc[i,"事件终止编号"]),"停顿归属事件"]=i+1# 删除停顿次数为0的事件
Stop = Stop[Stop["停顿归属事件"].notnull()]# 构造特征 用水事件停顿总时长、停顿次数、停顿平均时长、
# 用水时长，用水/总时长
stopAgg = Stop.groupby("停顿归属事件").agg({"停顿时长":sum,"停顿开始编号":len})
sj.loc[stopAgg.index - 1,"总停顿时长"] = stopAgg.loc[:,"停顿时长"].values
sj.loc[stopAgg.index-1,"停顿次数"] = stopAgg.loc[:,"停顿开始编号"].values
sj.fillna(0,inplace=True) # 对缺失值用0插补
stopNo0 = sj["停顿次数"] != 0 # 判断用水事件是否存在停顿
sj.loc[stopNo0,"平均停顿时长"] = sj.loc[stopNo0,"总停顿时长"]/sj.loc[stopNo0,"停顿次数"]
sj.fillna(0,inplace=True) # 对缺失值用0插补
sj["用水时长"] = sj["总用水时长"] - sj["总停顿时长"] # 定义特征用水时长
sj["用水/总时长"] = sj["用水时长"] / sj["总用水时长"] # 定义特征 用水/总时长
print('用水事件用水时长与频率特征构造完成后数据的特征为：\n',sj.columns)
print('用水事件用水时长与频率特征构造完成后数据的前5行5列特征为：\n', sj.iloc[:5,:5])
2. 构建用水量与波动属性
除了用水时长，停顿和频率外，用水量也是识别该事件是否为洗浴事件的重要特征。可以构建出下表所示的两个用水量特征。
| 特征 | 构建方法 | 说明 |
|---|---|---|
| 总用水量 | 总用水量 = 每条有水流数据的用水量 = 持续时间 × 水流大小 | 一次用水过程中使用的总的水量，单位为L。 |
| 平均水流量 | 平均水流量 = 总用水量 / 有水流时间 | 一次用水过程中，开花洒时平均水流量大小（为热水），单位为L/min。 |
同时用水波动也是区分不同用水事件的关键。根据不同用水事件的这一特征可以构建下表所示的水流量波动和停顿时长波动两个特征。
| 特征 | 构建方法 | 说明 |
|---|---|---|
| 水流量波动 | 水流量波动 = ∑ ( 单次水流的值 − 平均水流值 ) 2 × 持续时间 总的有水流量的时间 \sum\frac{(单次水流的值-平均水流值)^2\times 持续时间}{总的有水流量的时间} ∑总的有水流量的时间(单次水流的值−平均水流值)2×持续时间 | 一次用水过程中，开花洒时水流量的波动大小 |
| 停顿时长波动 | 停顿时长波动 = ∑ ( 单次水流的值 − 平均水流值 ) 2 × 持续时间 总停顿时长 \sum\frac{(单次水流的值-平均水流值)^2\times 持续时间}{总停顿时长} ∑总停顿时长(单次水流的值−平均水流值)2×持续时间 | 一次用水过程中，用水停顿时长的波动情况 |
data["水流量"] = data["水流量"] / 60 # 原单位L/min，现转换为L/sec
sj["总用水量"] = 0 # 给总用水量赋一个初始值0
for i in range(len(sj)):Start = sj.loc[i,"事件起始编号"]-1End = sj.loc[i,"事件终止编号"]-1if Start != End:for j in range(Start,End):if data.loc[j,"水流量"] != 0:sj.loc[i,"总用水量"] = (data.loc[j + 1,"发生时间"] - data.loc[j,"发生时间"]).seconds* \data.loc[j,"水流量"] + sj.loc[i,"总用水量"]sj.loc[i,"总用水量"] = sj.loc[i,"总用水量"] + data.loc[End,"水流量"] * 2else:sj.loc[i,"总用水量"] = data.loc[Start,"水流量"] * 2sj["平均水流量"] = sj["总用水量"] / sj["用水时长"] # 定义特征 平均水流量
# 构造特征：水流量波动
# 水流量波动=∑(((单次水流的值-平均水流量)^2)*持续时间)/用水时长
sj["水流量波动"] = 0 # 给水流量波动赋一个初始值0
for i in range(len(sj)):Start = sj.loc[i,"事件起始编号"] - 1End = sj.loc[i,"事件终止编号"] - 1for j in range(Start,End + 1):if data.loc[j,"水流量"] != 0:slbd = (data.loc[j,"水流量"] - sj.loc[i,"平均水流量"])**2slsj = (data.loc[j + 1,"发生时间"] - data.loc[j,"发生时间"]).secondssj.loc[i,"水流量波动"] = slbd * slsj + sj.loc[i,"水流量波动"]sj.loc[i,"水流量波动"] = sj.loc[i,"水流量波动"] / sj.loc[i,"用水时长"] # 构造特征：停顿时长波动
# 停顿时长波动=∑(((单次停顿时长-平均停顿时长)^2)*持续时间)/总停顿时长
sj["停顿时长波动"] = 0 # 给停顿时长波动赋一个初始值0
for i in range(len(sj)):if sj.loc[i,"停顿次数"] > 1: # 当停顿次数为0或1时，停顿时长波动值为0，故排除for j in Stop.loc[Stop["停顿归属事件"] == (i+1),"停顿时长"].values:sj.loc[i,"停顿时长波动"] = ((j - sj.loc[i,"平均停顿时长"])**2) * j + sj.loc[i,"停顿时长波动"]sj.loc[i,"停顿时长波动"] = sj.loc[i,"停顿时长波动"] / sj.loc[i,"总停顿时长"]print('用水量和波动特征构造完成后数据的特征为：\n',sj.columns)
print('用水量和波动特征构造完成后数据的前5行5列特征为：\n',sj.iloc[:5,:5])
（五）筛选候选洗浴事件
洗浴事件的识别是建立在一次用水事件识别的基础上，也就是从已经划分好的一次用水事件中识别出哪些一次用水事件是洗浴事件。
可以使用3个比较宽松的条件筛选掉那些非常短暂的用水事件，确定不可能为洗浴事件的数据去除掉，剩余的事件称为“候选洗浴事件”。这三个条件是“或”的关系，也就是说，只要一次完整的用水事件满足任意一个条件，就被判定为短暂用水事件，即会被筛选掉。3个筛选条件如下：
- 一次用水事件中总用水量小于5升。
- 用水时长小于100秒。
- 总用水时长小于120秒。
sj_bool = (sj['用水时长'] >100) & (sj['总用水时长'] > 120) & (sj['总用水量'] > 5)
sj_final = sj.loc[sj_bool,:]
sj_final.to_excel('../tmp/sj_final.xlsx',index=False)
print('筛选出候选洗浴事件前的数据形状为：',sj.shape)
print('筛选出候选洗浴事件后的数据形状为：',sj_final.shape)
筛选前用水事件数目总共172个，经过筛选后，余下75个用水事件。结合日志，最终用于建模的属性的总数为11个，其基本状况如下：
| 特征名称 | 均值 | 中位数 | 标准差 |
|---|---|---|---|
| 洗浴时间点 | 19.000000 | 20.000000 | 3.263227 |
| 总用用水时长 | 529.506667 | 503.000000 | 261.902621 |
| 总停顿时长 | 57.893333 | 4.000000 | 95.050566 |
| 停顿次数 | 1.213333 | 1.000000 | 1.544767 |
| 平均停顿时长 | 34.167302 | 2.000000 | 51.083390 |
| 用水时长 | 471.613333 | 461.000000 | 206.411416 |
| 用水/总时长 | 0.921799 | 0.989899 | 0.116112 |
| 总用水量 | 241.015556 | 235.116667 | 127.539757 |
| 平均水流量 | 0.497794 | 0.498853 | 0.118436 |
| 水流量波动 | 0.155609 | 0.019534 | 0.728971 |
| 停顿时长波动 | 619.675823 | 0.000000 | 1999.449248 |
四、分析与建模
根据建模样本数据建立BP神经网络模型识别洗浴事件。由于洗浴事件与普通用水事件在特征上存在不同，而且这些不同的特征在特征上被体现出来。于是，根据用户提供的用水日志，将其中洗浴事件的数据状态记录作为训练样本训练BP神经网络。然后根据训练好的网络来检验新采集到的数据，具体过程如下图所示。
在训练神经网络的时候，选取了“候选洗浴事件”的11个属性作为网络的输入，分别为：洗浴时间点，总用水时长，总停顿时长，平均停顿时长，停顿次数，用水时长，用水时长/总用水时长，总用水量，平均水流量，水流量波动和停顿时长波动。
训练BP网络时给定的输出（教师信号）为1与0，其中1代表该次事件为洗浴事件，0表示该次事件不是洗浴事件。是否为洗浴事件的标签是根据热水器的用水记录日志得到。
构建神经网络模型需要注意数据本身属性之间的存在量级差异，因此需要进行标准化，消除量级差异。另外，为了便于后续应用模型，可以用joblib.dump函数保存模型。
import pandas as pd
from sklearn.preprocessing import StandardScaler
from sklearn.neural_network import MLPClassifier
import joblib# 读取数据
Xtrain = pd.read_excel('../tmp/sj_final.xlsx')
ytrain = pd.read_excel('../data/water_heater_log.xlsx')
test = pd.read_excel('../data/test_data.xlsx')
# 训练集测试集区分。
x_train, x_test, y_train, y_test = Xtrain.iloc[:,5:],test.iloc[:,4:-1],ytrain.iloc[:,-1],test.iloc[:,-1]
# 标准化
stdScaler = StandardScaler().fit(x_train)
x_stdtrain = stdScaler.transform(x_train)
x_stdtest = stdScaler.transform(x_test)
# 建立模型
bpnn = MLPClassifier(hidden_layer_sizes = (17,10), max_iter = 200, solver = 'lbfgs',random_state=50)
bpnn.fit(x_stdtrain, y_train)
# 保存模型
joblib.dump(bpnn,'../tmp/water_heater_nnet.m')
print('构建的模型为：\n',bpnn)
在训练BP神经网络时，对神经网络的参数进行了寻优，发现含2个隐层的神经网络训练效果较好，其中2个隐层的隐节点数分别为17和10时训练的效果较好。
根据样本，得到训练好的神经网络后，就可以用来识别对应的用户家的洗浴事件，其中待检测的样本的11个属性作为输入，输出层输出一个值在[-1,1]范围内，如果该值小于0，则该事件不是洗浴事件，如果该值大于0，则该事件是洗浴事件。某热水器用户记录了两周的热水器用水日志，将前一周的数据作为训练数据，后一周的数据作为测试数据，代入上述模型进行测试。
五、模型评估与检验
结合模型评价相关的知识，使用精确率（precision）、召回率（recall）和f1值来做模型评价的效果先顾地较为客观、准确。同时结合ROC曲线，可以进一步更加直观地评价模型的效果，得到模型的ROC曲线如下图所示， ROC曲线覆盖的面积较大，说明模型的识别效果较好。
from sklearn.metrics import classification_report
from sklearn.metrics import roc_curve
import joblib
import matplotlib.pyplot as pltbpnn = joblib.load('../tmp/water_heater_nnet.m') # 加载模型
y_pred = bpnn.predict(x_stdtest) # 返回预测结果
print('神经网络预测结果评价报告：\n',classification_report(y_test,y_pred))
# 绘制roc曲线图
plt.rcParams['font.sans-serif'] = 'SimHei' # 显示中文
plt.rcParams['axes.unicode_minus'] = False # 显示负号
fpr, tpr, thresholds = roc_curve(y_pred,y_test) # 求出TPR和FPR
plt.figure(figsize=(6,4)) # 创建画布
plt.plot(fpr,tpr) # 绘制曲线
plt.title('用户用水事件识别ROC曲线') # 标题
plt.xlabel('FPR') # x轴标签
plt.ylabel('TPR') # y轴标签
plt.savefig('../tmp/用户用水事件识别ROC曲线.png') # 保存图片
plt.show() # 显示图形
根据该热水器用户提供的用水日志判断事件是否为洗浴与多层神经网络模型识别结果报告，如下表所示。
| precision | recall | f1-score | support | |
|---|---|---|---|---|
| 0 | 0.47 | 0.75 | 0.58 | 12 |
| 1 | 0.90 | 0.73 | 0.81 | 37 |
| avg/total | 0.80 | 0.73 | 0.75 | 49 |
根据模型评估报告表可以看出，在洗浴事件的识别上精确率（precision）非常高，达到了90%，同时召回率（recall）也达到了70%以上。综合上述结果，可以确定此次创建的模型是有效并且效果良好的能够用于实际的洗浴事件的识别中。
小结
- 本案例针对电热水器用户行为分析，提出通过数据挖掘技术识别洗浴事件的方法。
- 首先，基于热水器状态数据，采用时间间隔阈值法划分用水事件，并通过斜率指标优化确定4分钟为最佳停顿阈值。在特征工程阶段，构造了用水时长、停顿时长波动、总用水量等11个特征，筛选出候选洗浴事件。
- 随后，利用BP神经网络模型（隐层结构17-10）进行训练，输入标准化后的特征数据，输出洗浴事件分类。模型评估显示，洗浴事件识别精确率达90%，召回率73%，F1值0.81，ROC曲线AUC表现良好，验证了模型有效性。
- 该方案系统性地结合了时序分析、特征工程与机器学习，为家电智能化提供了数据驱动的解决方案，未来可进一步优化用户个性化差异的适应性。
附：以上文中的数据集及相关资源下载地址：
链接：https://pan.quark.cn/s/5eaba1ae893e
提取码：nFec
