李海林，夏燕燕，鄒金串

（1.華僑大學 信息管理系，福建 泉州 362021；2.華僑大學 體育學院，福建 泉州 362021；3.華僑大學 應用統(tǒng)計與大數據研究中心，福建 廈門 361021）

0 概述

依據《大數據白皮書（2019 年）》的報告，全球數據量在2021 年可提升至41 ZB。數據產生大多帶有時間屬性，這在一定程度上推動了時間序列數據挖掘理論與應用的發(fā)展。隨著運動醫(yī)學領域時間序列數據的大量積累與存儲，以及理論研究的日臻完善，分析這些數據以獲取隱含的、具有參考價值的信息已成為多數機構和學者的需求。心肺運動試驗（Cardiopulmonary Exercise Testing，CPET）能將人體的呼吸系統(tǒng)、心血管系統(tǒng)等綜合為一體，其不僅能夠體現(xiàn)受試者的有氧運動能力，評估受試者的心肺耐力（Cardiorespiratory Fitness，CRF）［1-2］，而且能夠以整體整合醫(yī)學的視角來探究受試者對運動的應激反應［3］。對CPET 數據進行聚類分析，就受試個體而言，利于從整體整合醫(yī)學視角進一步探究生理指標狀況，為評估身體健康狀態(tài)提供參考。對受試群體進行聚類分析，則能夠發(fā)現(xiàn)“離群值”或是實現(xiàn)預分群，為運動員選拔工作提供指導，并有利于針對不同特征組群進行針對性訓練。

在時間序列聚類挖掘過程中，相似性度量與聚類算法是兩大關鍵內容。一方面，選擇一種合適的相似性度量方法，將有利于聚類主題更加清晰明確。研究成果表明，合理有效的相似性度量能夠提升算法性能［4］。相似性度量主要選擇動態(tài)時間彎曲（Dynamic Time Warping，DTW）方法［5］，其度量時序距離具有更大的靈活性，但往往時序數據維度大，計算的時間復雜度也相對較高。為此，李海林等［6］提出一種基于關鍵形態(tài)特征的降維方法，通過數據預處理來加速時間序列的相似性度量。針對DTW 的高時間復雜度問題，文獻［7］提出一種高效優(yōu)化的DTW 度量方法，文獻［8］從整體上提出一種基于序列波動特征的相似性度量方法，該方法針對序列的波動特性，通過設置閾值選取關鍵極值點，再通過時間精準匹配、近鄰匹配得出兩序列間的相似性。另一方面，在多種傳統(tǒng)聚類算法思想中，由于層次聚類事先需要給定的先驗知識較少，對結果的分析往往能夠得到更多啟發(fā)，因而被廣泛使用。文獻［9］提出基于DTW 改進的層次聚類算法，并驗證了該算法所得聚類結果的有效性。時間序列聚類算法已被廣泛應用于金融股票［10-11］、交通客流［12-14］、環(huán)境監(jiān)測［15-16］、隧道施工［17］等領域。鑒于時間因素的重要性［5，18］，研究人員將時間序列聚類技術擴展應用到金融期貨［19］、農場肉類數據檢測［20］、電力客戶群市場細分［21］等各領域中。

在運動醫(yī)學領域，時間序列數據挖掘技術的應用尚未成熟。國內中長跑運動員選拔大多沿用“邊訓邊選”的理念，選拔與訓練工作是交替進行的，并采用定性與定量相結合的方法加以評定。在傳統(tǒng)意義上，定量評定主要采用運動測試期間某一靜止時點數據（如最大耗氧量VO2max、最大攝氧量速度vVO2max）或是一段時間內某類指標的差值（如vVO2max持續(xù)時間tlim）加以判別，忽略了整個運動測試期間數據的變動信息。CPET作為一種科學、整體整合理念的心肺功能測試，其數據應用價值潛力巨大。本文提出一種基于時間序列形態(tài)特征的算法，用于對CPET數據進行凝聚層次聚類分析。該算法能夠將測試期間的趨勢、極值、平臺、周期等各項變動特征加以綜合考慮，將測試期內指標變動較為相似的受試者聚為一類。通過對聚類結果的評估與分析，以實現(xiàn)為中長跑運動員選拔等工作提供科學化指導。

1 動態(tài)時間彎曲

時間序列數據是指通過等時間間隔獲取到的序列數據，帶有時間屬性。數據排列可反映時間先后順序，該數據一般具有時間維與變量維兩個維度。變量維為1 的序列稱為一元時間序列（Univariate Time Series，UTS），反之稱之為多元時間序列（Multivariate Time Series，MTS）。

定義1（一元時間序列數據）時間維為n的一元時間序列X表示為{x1,x2,…,xn}。

DTW 度量方法可對非等長序列進行距離度量，在時間序列領域被廣泛使用。另外，DTW 度量方法允許序列間出現(xiàn)平移、振幅、收縮等狀態(tài)，能最大程度地度量序列之間的相似性。

令d(xi,yj)為原始距離矩陣元素，D(xi,yj)為DTW 累積距離矩陣元素，則最終時間序列X={(x1；t1),(x2；t2),…,(xn；tn)}與Y={(y1；t1),(y2；t2),…,(ym；tm)}的動態(tài)彎曲距離為D(xn,ym)。DTW 公式表示如下：

DTW 數據點匹配狀況相比于歐式距離更加靈活，更大程度上度量了數據間的相似性。為防止序列出現(xiàn)過度彎曲匹配的狀況，一般采用加窗改進的DTW 度量方法，在式（1）中，將i≠0且j≠0 條件轉變?yōu)閕≠0且j≠0；i-c≤j≤i+c，其中c為常數。

2 基于心肺時序數據的聚類方法

CPET 屬力竭性測試，不同受試者運動時長不一，數據在時間維一般較短而屬性維較長。遞增負荷運動測試所得心肺數據序列具有一定波動特征，如通氣震蕩，存在明顯的上升或是下降趨勢。在受試者極限運動負荷下，相關指標數據可能會出現(xiàn)一個穩(wěn)定的“平臺”，這是CPET 數據的關鍵特征。同時，由于是縱向分析角度，不同指標閾值特征點序列長短不一較為顯著，如心率（Heart Rate，HR）指標波動較小，幾乎隨時間逐漸遞增，閾值特征點序列較短；而生理死腔與潮氣量（Dead Volume/Tidal Volume，VD/VT）比值震蕩波動較大，閾值特征點序列較長。若不考慮時長因素，則直接對數據進行DTW 度量不太準確。時長指序列數據之間的時間間隔跨度，n個時間間隔則計數為n。在縱向相似性度量方面，需考量時長因素以實現(xiàn)通過閾值特征點對原始序列進行線性擬合。橫向角度由于通過對應同類指標進行度量而不需要考慮時長因素。針對以上CPET 數據特征，故需選擇基于閾值特征點提取的降維策略，并結合DTW 思想進行距離度量。

2.1 形態(tài)特征表示

時間序列形態(tài)特征主要是指序列數據波動特征、變動趨勢、周期規(guī)律等特點。通過提取數據的主要形態(tài)特征點，能夠使聚類任務更加明確可行，并基于提取閾值極值點與閾值平穩(wěn)點來實現(xiàn)。

2.1.1 閾值極值點

傳統(tǒng)極值點的定義將數據點xi與其相鄰兩點xi-1和xi+1進行比較。在此基礎上，給定誤差范圍ε，允許數據點之間的比較存在ε范圍內的波動。

定義3（閾值極大值點）若xi-xi-1＞ε且xixi+1＞ε，則可判定xi為閾值極大值點。

定義4（閾值極小值點）若xi-1-xi＞ε且xi+1-xi＞ε，則可判定xi為閾值極小值點。

2.1.2 閾值平穩(wěn)點

傳統(tǒng)平穩(wěn)點是指連續(xù)的數據值相等，在曲線形態(tài)上是一條平行于X軸的直線。閾值平穩(wěn)點是指數據之間允許有ε的浮動范圍。在ε范圍內，將數據看成平穩(wěn)波動，其值變換為該時間段內數據的均值。

定義5（閾值平穩(wěn)點）若xi既不是閾值極大值點且非閾值極小值點，當|xi-xi-1|≤ε時，xi-1、xi為閾值平穩(wěn)點；當連續(xù)多個數據點滿足以上條件時，均可視為閾值平穩(wěn)點。在尋找該點時，有可能會將以ε范圍波動的遞增、遞減數據序列作為平穩(wěn)波動序列處理。當出現(xiàn)以上情況時，會增加線性擬合的誤差。因而依據數據特點，閾值需要謹慎選擇。ε值可結合最終數據線性擬合誤差、壓縮率來綜合確定。

對于時間序列X={x1,x2,···,xn}，經過一次遍歷后得到特征點序列X′={x′1,x′2,···,x′i}。依據X′對 原始序列數據進行擬合，得到等長線性擬合序列數據X″={x″1,x″2,···,x″n}。線性擬合誤差表示如下：

通過比較可知，基于閾值特征點提取的序列數據比傳統(tǒng)特征點提取的序列數據更短，并且能夠提取序列數據的非極值點，即關鍵轉折點、重要點。針對微小頻繁波動的數據，閾值特征點線性擬合表現(xiàn)更加靈活。因此，基于閾值特征點降維能夠保留數據的大體波動、周期等信息，而忽略過于頻繁、細小的變動信息。

2.2 相似性度量

縱向維度相似性度量DTW_SimT 與橫向維度的DTW_Sim 主要區(qū)別在于：縱向維度相似性度量需將提取的閾值特征點結合時長進行線性擬合，得到和原序列長度相等的序列；橫向維度相似性度量不考慮閾值特征點的時長因素。本文主要敘述橫向維度相似性度量算法。

對于多元時間序列數據集XM，需遍歷每個樣本個體，并求出每個數據集屬性列的閾值特征點序列；若為閾值平穩(wěn)點則進行標記，并以該段連續(xù)被標記為閾值平穩(wěn)點集的均值代替原來的數值。兩個多元時間序列樣本之間的相似性距離，由使用加窗改進的DTW 方法對應計算各變量維序列距離所得的均值表示。

算法1相似性度量算法DTW_Sim

2.3 聚類方法

通過DTW_Sim（或DTW_SimT）相似性度量方法，結合基于離差平方和類間鄰近性評價方法，實現(xiàn)對數據進行凝聚層次聚類，即聚類算法DTW_MFAC。

算法2時間序列聚類算法DTW_MFAC

3 數值實驗

3.1 實驗數據集

為驗證所提算法的有效性，從分類數據網站上下載了6 類含有類別標簽的數據，如表1 所示。數據均為一元時間序列，并選擇實驗數據的訓練集進行聚類效果驗證。由于DodgerLoopGame 訓練集中除去缺失數據僅有8 個，因此選擇其對應測試集數據進行實驗。

表1 實驗數據集Table 1 Experimental datasets

在執(zhí)行算法前，需通過平均線性擬合誤差、平均壓縮率為各數據尋找最佳閾值ε。初始ε變化范圍為［0.001，1.0］，以0.05 為步長逐漸增加。為方便綜合考慮最佳閾值，將平均線性擬合誤差及平均壓縮率做均值方差標準化處理。

選取各數據集的最佳閾值如表2 所示，結果保留兩位小數。此時，數據線性擬合誤差較小，同時壓縮率較低。

表2 實驗數據集最佳閾值Table 2 Optimal threshold of experimental datasets

3.2 算法比較與評估

為最大限度地減少誤差平方和，在凝聚層次聚類過程中，選擇離差平方和類間鄰近性度量方法，并與傳統(tǒng)的歐式距離度量方法做比較，如表3 所示。

表3 實驗數據集聚類準確率比較Table 3 Comparison of clustering accuracy of experimental datasets

通過不同數據集比較實驗可知，相比傳統(tǒng)基于歐式距離的層次聚類算法，DTW_MFAC 算法在對數據基本形態(tài)特征提取的過程中，一定程度上保證了聚類效果，總體而言聚類效果良好。同時，對于數據形態(tài)特征較為顯著、波動幅度較大的數據集，聚類準確率更高。

4 CPET 聚類分析與運動員選拔

4.1 聚類評估

輪廓系數最早于1986 年由ROUSSEEUW 提出，是以樣本類間距離與類內距離（也稱為凝聚度與分離度）為出發(fā)點。當數據的類別未知時，可使用輪廓系數對聚類效果進行度量評估。對于單個樣本xi，xi與同一簇內所有樣本距離的平均值記為ai；樣本到其他各簇Ci(1 ≤i≤k)各樣本xi(xi?ci)之間距離平均值的最小值，記為bi，則樣本xi的輪廓系數為：

單個樣本的輪廓系數相加求平均值，即為樣本數據集的輪廓系數S，如式（4）所示：

其中：N為樣本量。輪廓系數取值在［-1，1］區(qū)間，有3 個極端分別為-1、0 和1。輪廓系數si越接近-1，表示聚類效果越差，xi應分配到其他簇中；si越接近0，表示聚類效果處于中間狀態(tài)，將xi分到該簇或離其最近的簇中均可；si越接近1，表示聚類效果越好，xi屬于該簇成員。

4.2 CPET 數據

選取嚴格按實驗室要求做完CPET 試驗的15 名大學體院男生（年齡為（20.8±1.42）歲；身高為（174.93±3.87）cm；體重為（64.73±5.57）kg；BMI 為（21.13±1.37））數據作為分析對象。與專業(yè)級中長跑運動員不同，15 名受試者均為業(yè)余愛好者。每位受試者進行跑臺力竭測試的時間大約在10 min 不等。在實驗中，每15 s 增加一個負荷并以該時間間隔獲取受試者的30 項生理指標，包括最大耗氧量（VO2max）、最大心率（HRmax）等指標。通氣無氧閾（Ventilatory Threshold，AT）則依據以下規(guī)則進行判別：1）VCO2、通氣當量（Ventilatory Equivalents，VE）非線性升高的拐點；2）VCO2與VO2的相交點；3）VE/VO2相對升高而VE/VCO2未見降低的拐 點；4）呼吸商（Respiratory Quotient，RQ）相對升高的拐點。受試者耐力評定關鍵生理指標值如表4 所示。表4 中前3 列指標將數據按降序排列，可用于初步評估各ID 運動員的心肺耐量狀況?？梢钥闯?，運動員耐力狀況與CPET 測試時長有著直接聯(lián)系。直觀上也可以看出，測試時間越長，說明運動員心肺耐量越好。

表4 受試者耐力評定關鍵生理指標值Table 4 Value of key physiological indexes for endurance assessment of subjects

4.3 閾值和生理指標確定

縱向分析是指對某一ID 數據集進行聚類分析?？v向聚類分析為橫向維度分析提供閾值參照以及聚類指標選取建議。由于VO2max是評估受試者心肺耐力、運動受限的關鍵指標，文中選取ID 為9 的志愿受試者數據進行縱向分析，其VO2max（單位：ml/min）為4 062，消除體重量綱之后，與該列VO2max（單位：ml/min/kg）最大值58.88 僅相差0.01。首先需預確定ID9 數據集的最佳閾值范圍，如圖1 所示。

圖1 CPET 數據最佳閾值范圍選取Fig.1 Selection of best threshold range of CPET data

ID9 數據最佳閾值范圍為（0，0.15］。此時，數據的平均擬合誤差以及壓縮率均較小。本文中選取0.05 作為最佳閾值。此時ID9 各屬性列數據的平均擬合誤差為0.53，平均壓縮率為0.62。同時，使用該閾值對CPET 數據集進行橫向聚類。

使用DTW_MFAC 算法對ID9 數據集聚類，首先需對數據進行均值方差標準化，再使用選取的閾值對序列形態(tài)特征進行提取，并選取動脈血氧分壓（Partial Pressure of Oxygen，PaO2）指標的擬合狀況進行查看。從圖2 可以看出，ID9 數據的PaO2擬合狀態(tài)較好，保留了指標數據的大體波動特征以及平穩(wěn)、趨勢等狀態(tài)。

圖2 ID9 PaO2指標線性擬合圖Fig.2 Linear fitting diagram of ID9 PaO2 index

對于同一受試者而言，CPET 測試始與末均分別對應受試者正常及力竭狀態(tài)，需設置序列的首末兩點對應匹配。依據輪廓系數評估結果，各簇類數目k的系數均超過0.5，如圖3 所示。ID9 各項指標數據最佳分類數目為2，其值約為0.72，此時聚類效果較好?？蓪祿譃? 個簇群，聚類結果如圖4 所示。

圖3 縱向維度聚類（ID9）輪廓系數評估結果Fig.3 Profile coefficient evaluation results of vertical dimension clustering（ID9）

圖4 縱向維度聚類可視化（ID9）結果Fig.4 Vertical dimension clustering visualization（ID9）results

ID9 聚類結果主要分為兩大簇，可以看出：左簇群大多與耗氧、代謝以及能量轉化相關，而氧氣可直接影響人體能量的轉化；右簇群主要與每搏輸出量（Stroke Volume，SV）、二氧化碳轉化指標以及呼吸時間相關。若選取k=3 作為聚類結果劃分，右簇群進一步被細分為兩類：一類主要與CO2氣體相關；另一類則與PaO2和呼吸時間相關。由此可大致判斷，心肺運動測試機理以人體呼吸過程為基礎，探究O2、CO2的轉換與利用狀態(tài)。呼吸牽涉血液循環(huán)、肺循環(huán)等多個系統(tǒng)，對數據進行縱向維度解讀，有利于深入探究受試者各指標間的“簇群分布”特征，對于研究不同人群例如不同專項運動員、不同等級運動員或是健康人群與患者之間呼吸生理指標的內在相似度提供參考。

4.4 選拔應用

通過縱向聚類分析，說明各指標間數據的ε誤差可選擇0.05。聚類指標選擇左簇群耗氧量（VO2）、二氧化碳（VCO2）、心率（HR）、分鐘通氣當量（VE）、代謝當量（Metabolic Equivalents，METS）、生理死腔與潮氣量比值（VD/VT）、呼吸商（npRQ）以及右簇群每搏輸出量（SV），共8 類指標。選取的指標能夠體現(xiàn)運動員攝取、利用氧的效率、肺循環(huán)以及心功能等綜合狀況。

輪廓系數評估結果如圖5 所示，15 名業(yè)余中長跑受試者個體差異較大，應當單獨自劃分為一簇。以13 和14 作為簇類數目k劃分時，輪廓系數值分別達到0.758 和0.88，相應的層次聚類結果如圖6 所示。說明ID6 與ID7 以及ID14 與ID15 運動員綜合耐力評估狀態(tài)較為相似。

圖5 CPET 數據橫向聚類輪廓系數評估結果Fig.5 Evaluation results of profile coefficient of lateral clustering of CPET data

圖6 CPET 數據橫向聚類可視化結果Fig.6 Visualization results of CPET data horizontal clustering

通過表4 可初步判斷，與其他業(yè)余中長跑運動員相比，ID6、ID7、ID14、ID15 心肺耐量較差，故針對中長跑等耐力項目選拔時以上4 位測試者可不予考慮。其余運動員可通過對應VO2max、HRmax以及AT 確立訓練目標并進行專項訓練。在一段時間后，可再依據CPET 各項指標結合聚類算法思想進行再選拔。閾值選取可依據平均擬合誤差、平均壓縮率或是領域實踐經驗加以確定。

5 結果分析

中長跑作為一項以體能為主的項目，已逐漸成為全民性健身運動，大眾參與度高，業(yè)余參賽選手比重大。此類運動需高度重視基礎耐力評定與提高［22］，且已有多種可量化指標供評估參照。其中，VO2max是最為基礎和綜合的指標。已有研究表明，VO2max受遺傳因素影響較高，達到90%以上，在青少年初步選材時期應用更加普遍。另一方面，通氣無氧閾（AT）也能夠用于評定有氧運動能力，在運動訓練領域起到了關鍵性的作用，且受遺傳因素影響較?。?3］。針對中長跑耐力運動員的評定與選拔，國內外通常以最大耗氧量速度（vVO2max）、最大耗氧量平臺（VO2maxPD）、跑步經濟性（Running Economy，RE）或是vVO2max持續(xù)時間、達到VO2max后持續(xù)時間等值作為參考。此外，優(yōu)選男子800 m 參賽選手時，可參考最大沖刺速度（Maximum Sprint Speed，MSS）及厭氧速度儲備（Anaerobic Speed Reserve，ASR）指標，即vVO2max與MSS 的差值速度，兩者值的提高被認為是具備技能競爭能力的重要條件［24］；而精英級長跑運動員在增量測試中具有較高的呼氣末二氧化碳濃度（Postapneic end-tidal Carbon Dioxide Pressure，PetCO2）以及更低的分鐘通氣量（VE）［25］。在以上評估評定過程中，通常會選擇多個指標綜合考慮，同時不可排除個體差異的影響。

目前，國內大多沿用“定量與定性結合”、“邊訓邊選”的模式對運動員進行選拔。在多數情況下也需教練多年的教學經驗作為輔助決策依據。但在這一過程中，往往會面臨很多難題，如：需要提前了解運動員的身體素質以及各項生理指標；在訓練過程中能夠合理確定運動強度和運動量；能夠預防運動損傷狀況；能夠根據運動員自身能力不足之處做科學且有針對性的訓練計劃［26］等。然而，心肺運動試驗（CPET）作為一項整體整合式的試驗，以上指標均可通過測試結果分析得到，并能有效解決上述面臨的主要問題。

CPET 通?？勺鳛榕懿饺耸康尼t(yī)療預檢項目。鑒于CPET 測試的潛在價值，本文提出一種基于此類時間序列聚類的算法DTW_MFAC。該算法將以DTW 為核心的相似性度量算法DTW_Sim 與凝聚層次聚類算法思想相結合，DTW_Sim 能夠識別出各指標序列的形態(tài)特征，如極值點、平穩(wěn)點、趨勢、周期等，在更大程度上度量了序列間的相似性。實驗的結果表明，DTW_MFAC針對形態(tài)特征較為顯著的MTS 數據集聚類效果較好。本文以15 名大學生業(yè)余組中長跑運動員的CPET 數據作為分析對象，選取了表征有氧能力與心肺耐量的8 類關鍵指標進行聚類分析，結果發(fā)現(xiàn)受試者個體差異較大。同時依據輪廓系數評估結果，可初步剔除心肺耐量較差的4 位測試者。此外，不同閾值以及選取不同指標用于聚類分析，結果存在一定差異，可根據實際需要靈活選擇。

6 結束語

本文提出一種針對CPET 時間序列數據進行聚類分析的算法DTW_MFAC，使用該算法能夠識別出CPET 序列的閾值極值點與平穩(wěn)點。通過CPET 分析，選取可反映運動員心肺耐量的8 類指標進行聚類分析，發(fā)現(xiàn)個體差異較大，沒有出現(xiàn)明顯的簇群分布結果。依據輪廓系數評估準則，心肺耐量較差的測試者可不予考慮。為使聚類效果更佳，選取受試個體進行閾值確定或結合整體數據對象進行考慮將是進一步的研究重點。另外，本文選擇凝聚層次聚類思想構建算法，未針對基于密度的聚類、基于圖的聚類等目前主流聚類算法進行比較分析，后續(xù)研究可據此對算法進行優(yōu)化和改進。