本文是基於《Python資料分析與挖掘實戰》的實戰部分的第10章的資料——《家用電器使用者行為分析與事件識別》做的分析。

旨在補充原文中的細節程式碼，並給出文中涉及到的內容的完整程式碼；另外，原文中的資料處理部分排版先後順序個人感覺較為凌亂，在此給出梳理。

在作者所給程式碼的基礎上增加的內容包括：  

1）在資料規約部分： 書中提到：規約掉熱水器"開關機狀態"=="關"且”水流量”==0的資料，說明熱水器不處於工作狀態，資料記錄可以規約掉。但由後文知，此條件不能進行規約 因為，"開關機狀態"=="關"且”水流量”==0可能是一次用水中的停頓部分，刪掉後則無法準確計算關於停頓的資料

 2）在一次完整用水事件的劃分模型中： 將時間間隔列資料離散並面元，探索了不同時間間隔中，用水事件的個數； 畫用水停頓時間間隔頻率分佈直方圖； 確定一次用水事件停頓閾值，然後劃分一次完整用水事件。 

3）用水事件閾值尋優模型： 通過頻率分佈直方圖-確定閾值的變化與劃分得到的事件個數關係 通過影象中斜率指標-確定閾值的變化與劃分得到的事件個數關係

4）屬性構造中： 原書中只給出了需要構造的屬性的定義，並未給出具體程式碼，本文給出了具體的程式碼；並給出了兩種方法求用水事件的時間間隔

 5）模型構造： 添加了顯示混淆矩陣視覺化預測結果，檢視訓練結果正確率

1 背景與目標分析

    根據熱水器廠商提供的資料進行分析，對使用者的用水事件進行分析，判斷用水是否是洗浴事件，識別不同使用者的用水習慣，以提供個性化的服務。

    實質：二分類問題

2 資料探索

data = pd.read_excel('original_data.xls',encoding='gbk')
 

data[u'發生時間'] = pd.to_datetime(data[u'發生時間'], format = '%Y%m%d%H%M%S')#將該特徵轉成日期時間格式（***）
data = data[data[u'水流量'] > 0] # 只要流量大於0的記錄
# print len(data) #7679

data[u'用水停頓時間間隔']= data[u'發生時間'].diff()/ np.timedelta64(1, 'm') #將datetime64[ns]轉成 以分鐘為單位（*****）
data= data.fillna(0) # 替換掉data[u'用水停頓時間間隔']的第一個空值
 

2.1 資料質量分析

#-----第*1*步-----資料探索，檢視各數值列的最大最小和空值情況
data_explore = data.describe().T
data_explore['null'] = len(data)-data_explore['count']
explore = data_explore[['min','max','null']]
explore.columns = [u'最小值',u'最大值',u'空值數']
explore

2.2 資料特徵分析

2.2.1 分佈分析

#----第*2*步-----離散化與面元劃分
# 將時間間隔列資料劃分為0~0.1，0.1~0.2，0.2~0.3....13以上，由資料描述可知，
# data[u'用水停頓時間間隔']的最大值約為2094，因此取上限2100
Ti = list(data[u'用水停頓時間間隔'])#將要面元化的資料轉成一維的列表
timegaplist = [0.0,0.1,0.2,0.3,0.5,1,2,3,4,5,6,7,8,9,10,11,12,13,2100]# 確定劃分區間

cats = pd.cut(Ti,timegaplist,right=False) # 包擴區間左端,類似"[0,0.1)",（預設為包含區間右端）
x = pd.value_counts(cats)
x.sort_index(inplace = True)
dx = DataFrame(x,columns=['num'])
dx['fn'] = dx['num']/sum(dx['num'])
dx['cumfn'] = dx['num'].cumsum()/sum(dx['num'])

f1 = lambda x :'%.2f%%' %  (x*100)
dx[['f']]= dx[['fn']].applymap(f1)
dx

#-----第*3*步-----畫用水停頓時間間隔頻率分佈直方圖
fig = plt.figure()
ax = fig.add_subplot(1,1,1)

dx['fn'].plot(kind='bar')
plt.ylabel(u'頻率/組距')
plt.xlabel(u'時間間隔（分鐘）')
p = 1.0*dx['fn'].cumsum()/dx['fn'].sum()# 數值等於 dx['cumfn']，但型別是列表
dx['cumfn'].plot(color = 'r', secondary_y = True, style = '-o',linewidth = 2)
plt.annotate(format((p[4]), '.4%'), xy = (7, p[4]), xytext=(7*0.9, p[4]*0.95), arrowprops=dict(arrowstyle="->", connectionstyle="arc3,rad=.2")) #添加註釋，即85%處的標記。這裡包括了指定箭頭樣式。
plt.ylabel(u'累計頻率')

plt.title(u'用水停頓時間間隔頻率分佈直方圖')
plt.grid(axis='y',linestyle='--')

# fig.autofmt_xdate() #自動根據標籤長度進行旋轉
for label in ax.xaxis.get_ticklabels():   #此語句完成功能同上,但是可以自定義旋轉角度
       label.set_rotation(60)

plt.savefig('Water-pause-times.jpg')
plt.show()

3 資料預處理

3.1 資料規約

# 規約掉"熱水器編號"、"有無水流"、"節能模式"三個屬性
# 注意：
#書中提到：規約掉熱水器"開關機狀態"=="關"且”水流量”==0的資料，說明熱水器不處於工作狀態，資料記錄可以規約掉。但由後文知，此條件不能進行規約

# 因為，"開關機狀態"=="關"且”水流量”==0可能是一次用水中的停頓部分，刪掉後則無法準確計算關於停頓的資料

or_data = pd.read_excel('original_data.xls',encoding='gbk')
or_data.head()

data = or_data.drop(or_data.columns[[0,5,9]],axis=1) # 刪掉不相關屬性

data.to_excel('data_guiyue.xlsx')
data.head()

3.2 資料變換

3.2.1 用水事件閾值尋優模型

根據資料探索部分，初步可以推測，用水事件停頓閾值為四分鐘較為合適。但是為了進一步確定是否為四分鐘，為此接下來探索的是——用水事件閾值尋優模型。

#第一步：確定閾值與事件數的關係
#****************************
#@1  目標：確定閾值的變化與劃分得到的事件個數關係
#    方法：通過頻率分佈直方圖
#****************************
timedeltalist = np.arange(2.25,8.25,0.25)
# 從2.25到8.25間，以間隔為0.25，確定閾值即，閾值範圍為[2.25,2.5,2.75,3,...,7.75,8]
counts = [] # 記錄不同閾值下的事件個數
for i in range(len(timedeltalist)):
    threshold = pd.Timedelta(minutes = timedeltalist[i])#閾值為四分鐘
    d = data[u'發生時間'].diff() > threshold #  # 相鄰時間做差分，比較是否大於閾值
    data[u'事件編號'] = d.cumsum() + 1 # 通過累積求和的方式為事件編號
    temp = data[u'事件編號'].max()
    counts.append(temp)
coun = pd.Series(counts, index=timedeltalist)

# 畫頻率分佈直方圖
#將閾值與對應的事件數繪製成頻率分佈直方圖，以確定最優閾值

plt.rcParams['font.sans-serif']=['SimHei']
plt.rcParams['axes.unicode_minus']= False
plt.rc('figure', figsize=(8,6))
np.set_printoptions(precision=4)
fig = plt.figure()
fig.set(alpha=0.2)#設定圖示透明度
ax = fig.add_subplot(1,1,1)
# coun.plot(linestyle='-.',color='r',marker='<')
coun.plot(style='-.r*')#同上
ax.locator_params('x',nbins = int(len(coun)/2)+1)  # (****)
ax.set_xlabel(u'用水事件間隔閾值(分鐘)')
ax.set_ylabel(u'事件數（個）')
ax.grid(axis='y',linestyle='--') # (****)
plt.savefig('threshold_numofCase.jpg')
plt.show()

    由上圖可知，影象趨勢平緩說明使用者的停頓習慣趨於穩定，所以取該段時間開始作為閾值，既不會將短的用水時間合併，也不會將長的用水時間拆開，因此，最後選取一次用水時間間隔閾值為4分鐘

利用閾值的斜率指標來作為某點的斜率指標

# 第二步：閾值優化
#****************************
#@2  目標：確定閾值的變化與劃分得到的事件個數關係
#    方法：通過影象中斜率指標
#****************************

# 當存在閾值的斜率指標 k<KS :
#     取閾值最小的點A（可能存在多個閾值的斜率指標小於1）的橫座標x作為用水事件劃分的閾值（該值是經過實驗資料驗證的專家閾值）
# 當不存在閾值的斜率指標 k<KS：
#     找所有閾值中“斜率指標最小”的閾值t1：
#     若：該閾值t1對應的斜率指標小於KS2：
#         則取該閾值作為用水事件劃分的閾值
#     若：該閾值t1對應的斜率指標不小於KS2
#         則閾值取預設值——4分鐘
# 備註：
# KS是評價斜率指標用的專家閾值1
# KS是評價斜率指標用的專家閾值2

data = pd.read_excel('dataExchange_divideEvent.xlsx')
n = 4 #使用以後四個點的平均斜率
KS = 1 # 專家閾值1
KS2 = 5 # 專家閾值2

def event_num(ts):
    d = data[u'發生時間'].diff() > ts # 相鄰時間做差分，比較是否大於閾值（*****）
    return d.sum()+1 # 直接返回事件數（*****）

dt = [pd.Timedelta(minutes = i) for i in np.arange(1,9,0.25)]#（***）
h = DataFrame(dt,columns = [u'閾值']) # 定義閾值列（**）
h[u'事件數'] = h[u'閾值'].apply(event_num) # 計算每個閾值對應的事件數（*****）
h[u'斜率'] = h[u'事件數'].diff()/0.25 # 計算每個相鄰點對應的斜率（****）
h[u'斜率指標偏移前'] = pd.rolling_mean(h[u'斜率'].abs(), n)# 採用當前指標和後n個指標斜率的絕對值的平均作為當前指標的斜率

h[u'斜率指標'] = np.nan
h[u'斜率指標'][:-4] = h[u'斜率指標偏移前'][4:]


mink = h[u'斜率指標'][h[u'斜率指標'] < KS]# 斜率指標小於1的值的集合
mink1 = h[u'斜率指標'][h[u'斜率指標'] < KS2]# 斜率指標小於5的值的集合

if list(mink): # 斜率指標值小於1不為空時，即，存在斜率指標值小於1時
    minky = [h[u'閾值'][i] for i in mink.index]# 取“閾值最小”的點A所對應的間隔時間作為ts
    ts = min(minky) #取最小時間為ts
elif list(mink1):# 當不存在斜率指標值小於1時,找所有閾值中“斜率指標最小”的閾值
    t1 = h[u'閾值'][h[u'斜率指標'].idxmin()] #“斜率指標最小”的閾值t1
    # ts = h[u'閾值'][h[u'斜率指標偏移前'].idxmin() - n] #等價於前一行作用（*****）
    # 備註：用idxmin返回最小值的Index，由於rolling_mean自動計算的是前n個斜率的絕對值的平均，所以結果要平移-n，得到偏移後的各個值的斜率指標，注意：最後四個值沒有斜率指標因為找不出在它以後的四個更長的值
    if h[u'斜率指標'].min()<5:
        ts = t1#當該閾值的斜率指標小於5，則取該閾值作為用水事件劃分的閾值
    else:
        ts = pd.Timedelta(minutes = 4)# 當該閾值的斜率指標不小於5，則閾值取預設值——4分鐘

tm = ts/np.timedelta64(1, 'm')

print "當前時間最優時間間隔為%s分鐘" % tm

threshold = pd.Timedelta(minutes=4)#閾值為四分鐘
d = data[u'發生時間'].diff() > threshold # 相鄰時間做差分，比較是否大於閾值（*****）
data[u'事件編號'] = d.cumsum() + 1 # 通過累積求和的方式為事件編號（*****）

data.to_excel('dataExchange_divideEvent.xlsx')

3.2.3 屬性構造

   由於屬性構造內容較多，因此，本小節單獨放在下一篇部落格中進行分享。

備註：本章節完整程式碼詳見