2015-10-15 61 views
3

我有3種事件類型的數據,我想估計轉移概率Pij(1)。這些表明事件i後面有事件j發生的概率,假設事件發生了(所以我需要條件概率)。我也想知道Pij(2)和Pij(3),這是事件i之後的第二個(第三個)事件是事件j的條件概率。估計轉移概率(大熊貓)

看一看一些實物模型數據:

import pandas as pd 
import numpy as np 
np.random.seed(5) 
strings=list('ABC') 
events=[strings[i] for i in np.random.randint(0,3,20)] 
groups=[1,1,1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2] 
index=pd.date_range('2/2/2012',periods=20,freq='T') 
dfm=pd.DataFrame(data={'event':events,'group':groups},index=index) 
dfm.head() 

        event group 
2012-02-02 00:00:00 C 1 
2012-02-02 00:01:00 B 1 
2012-02-02 00:02:00 C 1 
2012-02-02 00:03:00 C 1 
2012-02-02 00:04:00 A 1 

到目前爲止,我還跟着一個非常不雅的和幼稚的戰略,用來shift,看看哪些事件發生在未來的時期:

#Create new columns containing the shifted values 
for i in range(1,4): 
    dfm['event_t%i'%i]=dfm.event.groupby(dfm.group).shift(-i) 
#Combine the columns with current and shifted values into one 
for i in range(1,4): 
    dfm['NEWevent_t%i'%i]=dfm['event']+' '+dfm['event_t%i'%i] 
    dfm = dfm.drop('event_t%i'%i, 1) 

#Count the number of times each combination occurs 
A=dfm['NEWevent_t1'].groupby(dfm.group).value_counts() 
B=dfm['NEWevent_t2'].groupby(dfm.group).value_counts() 
C=dfm['NEWevent_t3'].groupby(dfm.group).value_counts() 

merged=pd.concat([A, B, C], axis=1) 

這確實給出了每個組發生特定事件組合(例如AA,AB,..)的次數。繼續這樣做,我可以使用組變量和兩個字母對中的第一個字母作爲分組變量來做groupby。這種蠻力解決方案可能看起來像:

merged=merged.reset_index() 
merged['first']=merged['level_1'].apply(lambda x: x[0]) 
merged.columns=['group','i j','t1','t2','t3','first'] 
merged.groupby(['group','first'])['t1','t2','t3'].sum() 
sums=merged.groupby(['group','first'])['t1','t2','t3'].sum() 
merged=pd.merge(merged,sums,left_on=['group','first'],right_index=True) 
merged['Pij(1)']=merged.t1_x/merged.t1_y 
merged['Pij(2)']=merged.t2_x/merged.t2_y 
merged['Pij(3)']=merged.t3_x/merged.t3_y 
merged[['group','i j','Pij(1)','Pij(2)','Pij(3)']] 
merged.head() 

    group i j Pij(1) Pij(2)  Pij(3) 
0 1 A A 0.25 0.666667 0.666667 
1 1 A B 0.25 NaN   NaN 
2 1 A C 0.50 0.333333 0.333333 
3 1 B A 0.50 0.500000 0.500000 
4 1 B C 0.50 0.500000 0.500000 

我相信必須有一個更簡單的方法來實現這個?有關如何提高效率的任何建議?

注意:我的實際數據集包含500萬行,10個事件類型和100個組。

回答

4

提出轉換概率的最佳方法是在轉移矩陣中,其中T(i,j)是Ti前往Tj的概率。讓我們從您的數據開始:

import pandas as pd 
import numpy as np 

np.random.seed(5) 
strings=list('ABC') 
events=[strings[i] for i in np.random.randint(0,3,20)] 
groups=[1,1,1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2] 
index=pd.date_range('2/2/2012',periods=20,freq='T') 
dfm=pd.DataFrame(data={'event':events,'group':groups},index=index) 
for i in range(1,4): 
    dfm['event_t%i'%i]=dfm.event.groupby(dfm.group).shift(-i) 

我認爲您的換檔命令是好的,但那只是我。無論如何,從這裏限制到'group' == 1並填充轉換矩陣。最後,除以列以獲得轉換概率。

trans = pd.DataFrame(columns=strings, index=strings) 
g_dfm = dfm[dfm['group']==1] 

for s1 in strings: 
    for s2 in strings: 
     events = g_dfm[(g_dfm['event']==s1) & (g_dfm['event_t1']==s2)] 
     trans.ix[s1, s2] = len(events) 

trans = trans.astype(float).div(trans.sum(axis=1), axis=0) 
trans = trans.fillna(0) 

從那裏,你可以做一個熱圖:

import matplotlib.pyplot as plt 

fig, ax = plt.subplots(figsize=(3,3)) 
ax.pcolormesh(trans.values, cmap=plt.get_cmap('Blues'), vmin=0, vmax=1) 
ax.invert_yaxis() 
ax.set_yticks(np.arange(0, len(trans.index))+0.5) 
ax.set_xticks(np.arange(0, len(trans.columns))+0.5) 
ax.set_yticklabels(trans.index, fontsize=16, color='k') 
ax.set_xticklabels(trans.columns, fontsize=16, color='k') 
ax.tick_params(direction='out', pad=10) 
ax.set_frame_on(True) 

for tk1, tk2 in zip(ax.xaxis.get_major_ticks(), ax.yaxis.get_major_ticks()): 
    tk1.tick1On, tk2.tick1On, tk1.tick2On, tk2.tick2On = [False]*4 

plt.show() 

enter image description here

沖洗和重複所有的組和第二和第三過渡。

+1

謝謝。但是如果我是正確的,那麼行trans = trans/trans.sum()'有問題。 (1)這似乎給出了錯誤的結果,並且(2)不能通過零處理除法。解決方案將是trans = trans.astype(float).div(trans.sum(axis = 1),axis = 0) trans = trans.fillna(0)'。 – Pilik

+0

謝謝,我會解決它。它爲我工作,或者我不會發布它,但你的可能會更普遍。 – thefourtheye