[[178542]]

本文發(fā)表于新加坡政府開放數(shù)據(jù)門戶站的博客，經(jīng)授權(quán)由InfoQ中文站翻譯并分享，本文由InfoQ社區(qū)編輯劉志勇翻譯，感謝譯者的辛苦奉獻。

撰文: Daniel Sim | 分析: Lee Shangqian、Daniel Sim、Clarence Ng

編者按：大數(shù)據(jù)正在滲透各行各業(yè)，甚至能跟你考試能力測試、患上某種疾病的機率等非常生活化的場景應(yīng)用都發(fā)生緊密的聯(lián)系。今后大數(shù)據(jù)在我們的生活中就像是水和電一樣，讓社會整個信息質(zhì)量更好、讓信息利用效率更高效。

世界著名未來學家托夫勒曾說改變這個世界的力量有三種暴力、知識、金錢，而如今我們的世界正在被第四種力量改變，那就是大數(shù)據(jù)!大數(shù)據(jù)不管應(yīng)用在哪個行業(yè)它的核心都是通過技術(shù)來獲知事情發(fā)展的真相，最終利用這個“真相”來更加合理的配置資源。具體來說，要實現(xiàn)大數(shù)據(jù)的核心價值，還需要前兩個重要的步驟，***步是通過“眾包”的形式收集海量數(shù)據(jù)，第二步是通過大數(shù)據(jù)的技術(shù)途徑進行“全量數(shù)據(jù)挖掘”，***利用分析結(jié)果進行“資源優(yōu)化配置”。說白了，大數(shù)據(jù)最終的落地就是資源優(yōu)化配置。

本文揭示了新加坡政府是如何利用大數(shù)據(jù)技術(shù)來捕獲引發(fā)地鐵被中斷的反常列車，我們得以再一次見識大數(shù)據(jù)技術(shù)的神奇力量。 

最近幾個月，新加坡的地鐵環(huán)線(MRT Circle Line)遭到了一連串的神秘中斷，對數(shù)以千計的乘客造成了很大的混亂和痛苦。

同大多數(shù)同事一樣，我每天早晨搭乘環(huán)線地鐵到辦公室。因此，在11月5日，當我所在的團隊有調(diào)查原因的機會時，我就毫不猶豫地自告奮勇參加了。

根據(jù)新加坡地鐵公司(SMRT)和新加坡陸路交通管理局(Land Transport Authority，LTA)的先前調(diào)查，我們知道這些事件是由于某種形式的信號干擾造成的，導致了一些列車的信號丟失。信號丟失會觸發(fā)那些列車中的制動安全功能，并使它們沿著軌道隨機停止。

但是這起***次發(fā)生在八月份的事件——似乎是隨機發(fā)生的，使調(diào)查小組很難找到確切的原因。

我們獲得了由SMRT編譯的數(shù)據(jù)集，其中包含以下信息：

• 每個事件的日期和時間
• 事件的位置
• 涉及的列車編號
• 列車的方向

我們開始清理數(shù)據(jù)，在Jupyter Notebook中進行工作，這是一個流行的編寫和記錄Python代碼的工具。

像往常一樣，***步是導入一些有用的Python庫。

import math 
import xlrd 
import itertools as it 
import numpy as np 
import pandas as pd 
from datetime import datetime  

片段1

然后我們從原始數(shù)據(jù)中提取有用的部分。

dfincidents_0 = pd.read_excel('CCL EVAC E-brake occurrences hourly update_mod.xlsx', sheetname='Aug Sep') 
dfincidents_1 = pd.read_excel('CCL EVAC E-brake occurrences hourly update_mod.xlsx', sheetname='Nov') 
# Incident data for Nov had different format 
dfincidents_1['Time'] = dfincidents_1['Time'].str.strip('hrs').str.strip(' ') 
dfincidents_1['Time'] = pd.to_datetime(dfincidents_1['Time'], format='%H%M').dt.time 
dfincidents = pd.concat([dfincidents_0, dfincidents_1]) 
# Reset the index because they were concatenated from two data sources 
dfincidents.reset_index(inplace=True, drop=True)  

片段2

我們將日期和時間列合并為一個標準列，以便更容易地將數(shù)據(jù)可視化：

def datetime_from_date_and_time(row): 
    """ 
    Combines the date column and time column into a single column 
    """ 
    d = row['Date'] 
    t = row['Time'] 
    return datetime( 
        d.year, d.month, d.day, 
        t.hour, t.minute, t.second 
    ) 
# Add DateTime to the data for easier visualization 
dfincidents['DateTime'] = dfincidents.apply(datetime_from_date_and_time, axis=1)  

片段3

這就產(chǎn)生了如下表格： 

[[178544]] 

截圖1：初始處理的輸出

最初的可視化沒有明確的答案

我們在初步探索性的分析中找不到任何明顯的答案，如下圖所示：

1.事件發(fā)生在一天之中，整天的事件數(shù)量反映了高峰和非高峰旅行時間。 

圖1 出現(xiàn)反映高峰和非高峰旅行時間的次數(shù)

2.事故發(fā)生在環(huán)線上的各個地點，西側(cè)發(fā)生的事件略多一些。 

圖2 干擾的原因似乎與位置無關(guān)

3.只有一輛或兩輛列車的話，信號干擾并不會產(chǎn)生影響，但在這條環(huán)線上有許多列車。“PV”是“客運車輛”(Passenger Vehicle)的縮寫。 

圖3 60輛不同列車受到信號干擾

Marey圖表：顯示時間、位置和方向

我們的下一步是將多維度納入探索性分析。

我們的靈感來自Marey圖表，這是1983年Edward Tufte出版的經(jīng)典著作《定量信息的視覺顯示》(《The Visual Display of Quantitative Information》)最近，它被Mike Barry和Brian Card應(yīng)用在波士頓地鐵系統(tǒng)的大規(guī)?？梢暬椖浚?nbsp;

截圖2 摘自http://mbtaviz.github.io/

在該圖表中，垂直軸表示時間——按時間順序從上到下;而水平軸表示沿著列車線路的車站。對角線則表示列車運行。

在我們的Marey圖表中開始繪制軸： 

圖4 環(huán)線版本的一個空白Marey圖表

在正常情況下，在HarbourFront和Dhoby Ghaut之間運行的列車將在類似與此的線路上移動，每次單程行程只需要一個小時以上： 

圖5 環(huán)線上列車運行的程式化表示

我們的目的是在這個圖表上繪制事件——是點而不是線。

準備可視化數(shù)據(jù)

首先，我們將站名從三字母代碼轉(zhuǎn)換為數(shù)字：

• Marina Bay到Promenade之前：0到1.5;
• Dhoby Ghaut到HarbourFront：2到29。

如果事故發(fā)生在兩個站之間，則將其表示為0.5加上兩個站號中的較小數(shù)。例如，如果事件發(fā)生在HarbourFront(29)和Telok Blangah(28)之間，則位置將是“28.5”。這樣我們就很容易繪制沿水平軸的點。 

stations=("MRB,BFT,DBG,BBS,EPN,PMN,NCH,SDM,MBT,DKT,PYL,MPS,TSG,BLY,SER," 
 
          "LRC,BSH,MRM,CDT,BTN,FRR,HLV,BNV,ONH,KRG,HPV,PPJ,LBD,TLB,HBF").split(',')  
def loc_id(station1, station2 = None):  
    """  
    Translates a 3-letter station code to a number,  
    or a pair of 3-letter station codes to a number.  
    Single stations are represented as whole numbers.  
    Locations between stations are represented with a .5.  
    Example:  
    loc_id('MRB')         # 0 (Marina Bay)  
    loc_id('MRB', 'BFT')  # 0.5 (Between Marina Bay and Bayfront) 
    loc_id('DBG')         # 2 (Dhoby Ghaut)  
    loc_id('HBF')         # 29  
    loc_id('HBF', 'TLB')  # 28.5 (Between Harbourfront and Telok Blangah)  
    loc_id('HBF', 'DBG')  # throws and error, because these stations are not adjacent  
    """  
    if station2 == None or station2 == 'nan' or (type(station2) is float and math.isnan(station2)):  
        # Single stations  
        return stations.index(station1)  
    else: # Pairs of stations -- take the average to get the 0.5  
        stn1_index = stations.index(station1)  
        stn2_index = stations.index(station2)  
        # Handle the branch at Promenade  
        if (set(['PMN', 'EPN']) == set([station1, station2])):  
            return float(stations.index('EPN')) + 0.5  
        elif set(['PMN', 'BFT']) == set([station1, station2]):  
            return float(stations.index('BFT')) + 0.5  
        else:  
            # Require station pairs to be adjacent stations  
            assert(math.fabs(stn1_index - stn2_index) == 1)  
            return float(stn1_index + stn2_index) / 2 

片段 4

然后我們計算位置ID的數(shù)字……

def loc_id_from_stations(row): 
    try: 
        # This handles entries with both "Station from" and "Station to" 
        # and entries with only "Station from" 
        return loc_id(row['Station from'], str(row['Station to'])) 
    except ValueError: 
        # Some entries only have "Station to" but no 
        # "Station from" 
        return loc_id(row['Station to'])  

片段 5

并添加到數(shù)據(jù)集：

# Select only some columns that we are interested in 
sel_dfincidents = dfincidents[['DateTime', 'PV', 'Bound', 'Station from', 'Station to', 'Event', 'Remarks']] 
# Add the location ID into the dataset 
sel_dfincidents['LocID'] = sel_dfincidents.apply(loc_id_from_stations, axis=1)  

片段6

然后我們得到如下表格： 

截圖3 添加位置ID后的輸出表

通過數(shù)據(jù)處理，我們能夠創(chuàng)建所有緊急制動事件的散點圖。這里的每個點代表一個事件。我們還是無法發(fā)現(xiàn)任何明顯的事故模式。 

圖6 信號干擾事件表示為散點圖

接下來，我們通過將每個事件表示為指向左側(cè)或右側(cè)的三角形，而不是點，將列車方向添加到圖表中： 

圖7 方向由箭頭和顏色表示。

它看起來相當隨機，但當我們放大到如下的圖表，一個模式似乎似乎浮出了水面： 

圖8 上午6點到10點之間的事件

如果你仔細閱讀圖表，你會注意到，故障似乎按順序發(fā)生。當一趟列車受到干擾時，另一趟在同一方向行駛的列車很快就會受到波及。

信號如何互相干擾?

在這一點上依然不明，一趟列車是罪魁禍首。

我們已經(jīng)證實的是，似乎有一個隨時間和位置相關(guān)的模式：事件一個接一個地發(fā)生，與上一個事件的方向相反。似乎有一條“破壞的痕跡”，它會不會是導致數(shù)據(jù)集中那些事件的誘因?

事實上，連接事件的假想線看上去與截圖2的Marey圖表可疑地類似。干擾的原因會不會是對面軌道的列車? 

圖9 它會是一趟相反方向行駛的列車嗎?

我們決定檢驗這個“反常列車”的假說。

我們已經(jīng)知道，沿著環(huán)線的車站之間的行駛時間在兩到四分鐘之間。這意味著如果發(fā)生4分鐘的間隔，我們可以將所有緊急制動事件分組在一起。 

def same_cascade(i, j):  
    """  
    Given a pair of incidents (i,j), returns true if:  
        t <= d * 4 mins  
    where t is  the time difference between occurrences  
    and d is the distance (measured by difference in location ID).  
    Moreover, we consider the track direction, and only consider  
    incidents that are "moving backwards".  
    """  
    # If trains are not travelling in the same direction  
    # they cannot be due to the same "backward moving" interference  
    # source.  
    # (Note: This was the hypothesis when this code was written. 
     # It turned out that the rogue train could affect all  
    # trains in the vicinity, not just in the opposite track)  
    if i["Bound"] != j["Bound"] or \  
        i["Bound"] not in ['IT', 'OT']:  
        return False  
    # time difference in minutes  
    time_difference = (i["DateTime"] - j["DateTime"]) / np.timedelta64(1, 'm')  
    location_difference = i["LocID"] - j["LocID"]  
    if location_difference == 0:  
        return False  
    ratio = time_difference / location_difference  
    if i["Bound"] == 'OT':  
        return ratio > 0 and ratio < 4  
    elif i["Bound"] == 'IT':  
        return ratio < 0 and ratio > -4  

片段 7

我們發(fā)現(xiàn)了滿足這個條件的所有事件對：

incidents = sel_dfincidents.to_records()  
# (a, b, c, d, ...) --> ((a,b), (a,c), (a,d), ..., (b,c), (b,d), ..., (c,d), ...)  
incident_pairs = list(it.combinations(incidents, 2))  
related_pairs = [ip for ip in incident_pairs if same_cascade(*ip)]  
related_pairs = [(i[0], j[0]) for i,j in related_pairs]  

片段 8

然后，我們使用并查集的數(shù)據(jù)結(jié)構(gòu)將所有相關(guān)的事件對分組成更大的集合。這使我們可以將可能關(guān)聯(lián)到同一“反常列車”的事件分組。 

def pairs_to_clusters(pairs): 
    """ 
    A quick-and-dirty disjoint-set data structure. But this works fast enough for 200+ records 
    Could be better. 
    Example input: 
    (1,2), (2,3), (4,5) 
    Output: 
    1: {1,2,3} 
    2: {1,2,3} 
    3: {1,2,3} 
    4: {4,5} 
    5: {4,5} 
    """ 
    the_clusters = dict() 
    for i,j in pairs: 
        if i not in the_clusters: 
            if j in the_clusters: 
                the_clusters[j].add(i) 
                the_clusters[i] = the_clusters[j] 
            else: 
                the_clusters[i] = set(list([i, j])) 
                the_clusters[j] = the_clusters[i] 
        else: 
            if j in the_clusters: 
                if the_clusters[i] is not the_clusters[j]: # union the two sets 
                    for k in the_clusters[j]: 
                        the_clusters[i].add(k) 
                        the_clusters[k] = the_clusters[i] 
                else: # they are already in the same set 
                    pass 
            else: 
                the_clusters[i].add(j) 
                the_clusters[j] = the_clusters[i] 
    return the_clusters  

片段9

然后將我們的算法應(yīng)用于數(shù)據(jù)：

clusters = pairs_to_clusters(related_pairs) 
# Show each set only once 
clusters = [v for k,v in clusters.items() if min(v) == k] 
clusters[0:10]  

片段10

這些是我們確定的一些集群：

[{0, 1},

{2, 4},

{5, 6, 7},

{8, 9},

{18, 19, 20},

{21, 22, 24, 26, 27},

{28, 29, 30, 31, 32, 33, 34},

{42, 44, 45},

{47, 48},

{51, 52, 53, 56}]

接下來，我們計算了可以通過我們的聚類算法解釋的事件的百分比。 

# count % of incidents occurring in a cluster 
all_clustered_incidents = set() 
for i,clust in enumerate(clusters): 
    all_clustered_incidents |= clust 
(len(all_clustered_incidents), 
 len(incidents), 
 float(len(all_clustered_incidents)) / len(incidents))  

片段11

結(jié)果是： 

(189, 259, 0.7297297297297297) 

它表達的意思是：在數(shù)據(jù)集中的259個緊急制動事件中，189個案例(其中73個)可以通過“反常列車”假說來解釋。我們覺得我們真的走對了路。

我們基于聚類結(jié)果對事件圖進行了著色。具有相同顏色的三角形在同一個集群中。 

圖10 通過我們的算法聚類的事件

有多少反常列車?

如圖5所示，在環(huán)線上的每個端到端行程需要大約1小時。我們通過事件圖和與圖5密切匹配的線繪制***線。這強烈暗示只有一個“反常列車”。 

圖11 事件集群的時間強烈暗示干擾能夠關(guān)聯(lián)到單趟列車

我們還觀察到，不明的“反常列車”本身似乎沒有遇到任何信號問題，因為它沒有出現(xiàn)在我們的散點圖。

我們相信我們有一個很好的例子，決定進一步調(diào)查。

捕獲反常列車

日落之后，我們?nèi)チ薑im Chuan站以確定“反常列車”。我們不能檢查詳細的列車日志，因為SMRT需要更多的時間來提取數(shù)據(jù)。因此，我們決定用傳統(tǒng)方式通過查看在事件發(fā)生時到達和離開每個車站的列車的視頻記錄來識別列車。

上午3點，車隊出現(xiàn)了頭號嫌犯：PV46，一輛自2015年起投入使用的列車。

檢驗假說

11月6日(星期日)，LTA和SMRT測試，如果PV46是通過運行列車在非高峰時間問題的來源，那我們就對了——PV46確實導致了附近列車之間的通信丟失，并觸發(fā)那些列車上的緊急制動器。在PV46未投入使用的那天之前，并沒有這樣的事件發(fā)生。

在11月7日(星期一)，我們的團隊處理了PV46的歷史位置數(shù)據(jù)，并得出結(jié)論，從8月到11月的所有事件中，超過95%可以用我們的假說來解釋。其余的事件，可能是由于偶發(fā)在正常情況下的信號丟失。

在某些日子，如9月1日，這種模式特別明顯。你可以很容易地看到，當PV46投入運行時 ，在其運行期間或者前后，常常發(fā)生發(fā)生干擾事件。 

LTA和SMRT最終于11月11日發(fā)布了一份聯(lián)合新聞稿，與公眾分享調(diào)查結(jié)果。

小結(jié)

當我們開始時，我們曾經(jīng)希望能夠找到跨機構(gòu)調(diào)查組可能感興趣的模式，其中包括LTA、SMRT和DSTA的許多官員。SMRT和LTA提供的整潔的事件日志幫助我們開了一個好頭，因為在導入和分析數(shù)據(jù)之前需要盡可能少的清理。我們還對LTA和DSTA的有效跟蹤調(diào)查表示滿意，因為他們證實了PV46存在硬件問題。

從數(shù)據(jù)科學的角度來看，我們很幸運，事件發(fā)生得如此接近。這使我們能夠在如此短的時間內(nèi)識別問題和罪魁禍首。如果事件更加孤立，則Z字形模式將不那么明顯，并且它將讓我們耗費更多的時間和數(shù)據(jù)來解決這個謎。

當然，令人***興的是，現(xiàn)在所有人可以坐環(huán)線地鐵再次充滿信心地工作。

注意：本文提到的代碼是在2016年11月5日編寫的——我們處理SMRT數(shù)據(jù)以確定環(huán)線地鐵事件原因的實際日期。我們承認可能存在不足之處。您可以在這里下載我們的Jupyter Notebook的副本(源碼)。