一種面向多模態手術軌跡的快速無監督分割方法

邵振洲，趙紅發，渠 瀛，施智平，關 永，袁慧梅

1(首都師范大學 信息工程學院，北京 100048)

2(首都師范大學 輕型工業機器人與安全驗證北京市重點實驗室，北京 100048)

3(成像技術北京市高精尖創新中心，北京 100048)

4(田納西大學 諾克斯維爾分校工程學院，田納西 美國 37996)

1 引 言

醫療手術軌跡分割技術是機器人輔助微創手術(RMIS)中的重要研究方向，在示范學習[1]、技能評估[2-4]和復雜任務自動化[5]等諸多領域有廣泛的應用.以著名的“達芬奇”醫療手術機器人為例，在醫生手術過程中可以產生運動學和視頻記錄，每一個手術任務的過程軌跡(即運動學數據和視頻記錄)都可以分解為若干個可視為原始軌跡的子軌跡段(如圖1所示)，這些子軌跡段復雜度低、方差小、容易剔除異常值，可以通過相互組合來實現新的手術任務.如何更加快速準確地進行手術軌跡分割在機器人輔助微創手術領域具有十分重要的研究意義.

圖1 手術軌跡分割示意圖Fig.1 Diagram of surgical trajectory segmentation

目前，解決手術軌跡分割問題的方法主要包括基于機器學習的有監督和無監督方法.其中有監督的方法從人工標注中學習子軌跡段到預定義軌跡段的匹配關系.例如，Lin等人[6]采用基于線性判別分析(LDA)的簡單貝葉斯分類器學習手術軌跡和標注的分類匹配關系；受到隱馬爾科夫模型(HMM)在語音識別領域成功的啟發，文獻[7]提出了采用專家手術示范數據驅動的方法完成HMM的建模，對手術過程進行分割和識別；然而，這類方法沒有考慮細粒度軌跡的臨近關系，所以捕獲手術軌跡中的長程狀態轉換[8]有更好的效果.但是這些模型在訓練前需要預先完成人工標注，而人工標注是非常耗時的.于是為了擺脫對人工標注的依賴，一些無監督方法被提出，通常首先提取手術軌跡數據(包括運動學、視覺數據)的特征，然后通過聚類的方法，估計出一系列的時間分割點或者分割子軌跡段.例如，Zhou等人[9]提出了一種無監督分割方法，使用人為設計的特征(速度，方向，夾角)，采用基于LDA分類的分割方法，得到軌跡的轉折點，這種方法只考慮了軌跡數據中的部分信息.同樣是無監督的分割方法，Sang[10]采用主成分分析(PCA)提取特征，然后利用高斯混合模型(GMM)對特征進行聚類的分割方法達到了更好的效果.近年來，Krishnan等人[11]提出了一種對轉移狀態進行聚類的方法(TSC)，通過構造手術軌跡臨近時間序列的轉移表示，保留軌跡序列的局部時間相關性信息，提升了分割算法性能.最近，受到卷積神經網絡(CNN)在計算機視覺領域獲得巨大成功的啟發，文獻[12]將上述工作擴展到深度學習領域，提出一種采用VGG網絡提取手術視頻特征，然后結合運動學特征進行分割的方法，提升了分割質量.

然而，無監督方法仍存在一些不足.首先，在視頻特征提取階段，由于視頻本身數據量大并且所采用的神經網絡結構復雜，提取的特征數據量很大(1Gb以上)，導致在視頻特征的提取與編碼階段時間耗費較長.其次，動作細致是微創手術的特點，其運動學數據(即末端執行器的位置，方向，線速度，角速度，夾子的夾角)易受到噪聲干擾，由于醫生專業水平的差異性，難免存在一些手術偏差，這樣就導致一些噪聲的引入，即使同一個手術專家遠程操作重復同一個手術任務，手術軌跡數據也有很大不同，對手術軌跡進行平滑處理有助于分割.然而，已有的分割方法對運動學數據處理只是簡單的歸一化處理，沒有對其做進一步的去噪處理.

針對上述問題，本文提出了一種快速手術軌跡無監督分割方法(如圖2所示)，首先通過堆疊卷積自編碼器(SCAE)提取手術過程視頻特征，在將視頻數據映射到低維空間表示的同時保留視頻中的主要特征，其次采用小波對手術過程的運動學數據進行多尺度濾波處理，平滑短程軌跡，得到運動學特征.為了減少由于采用單一特征進行分割導致的信息缺失，將視覺特征及運動學特征進行融合，然后對融合特征進行時間軸加窗處理得到轉移表示，對轉移進行聚類找出轉移狀態.最后在視覺空間、運動學空間和時間上對轉移狀態進行層次化聚類得到聚類簇，并對其中的異常簇進行修剪，得到最終的分割結果分割子軌跡段.本文對涉及到的數學符號進行了統一匯總如表1所示.

圖2 快速手術軌跡無監督分割算法框圖Fig.2 Algorithm block diagram of fast surgical trajectory unsupervised segmentation

表1 數學符號意義Table 1 Meaning of mathematical symbols

2 基于SCAE的視覺特征提取

(1)

(2)

為了進一步降低編碼器輸出特征的維度，在編碼器的最后一層(式3)和解碼器的第一層(式4)采用卷積核為1×1大小的卷積層，在降低維度的同時保留了空間信息：

(3)

(4)

圖3 SCAE網絡結構示意圖Fig.3 Diagram of SCAE network structure

(5)

為了避免由于傳統轉置卷積帶來的棋盤效應，本文采用雙線性插值的方法進行上采樣，解碼器第l上采樣層可由式(6)計算.

(6)

網絡的訓練通過最小化均方誤差L式(7)進行.

(7)

視頻圖像特征提取網絡的偽代碼如下：

算法1.基于SCAE特征提取網絡訓練與預測

網絡訓練輸入：手術視頻圖像幀I

1. FORlin encoder layer

//編碼器卷積層，3×3卷積核

4. END FOR

//編碼器末層，1×1卷積核

//解碼器首層，1×1卷積核

7. FORlin decoder layer

//解碼器卷積層，3×3卷積核

10. END FOR

網絡預測輸出：特征圖H4

1. FORlin encoder layer

4. END FOR

3 基于視頻和運動學軌跡的層次化非參聚類

3.1 基于小波變換的運動學數據尺度平滑處理

小波變換可以將手術軌跡分解為低頻部分與高頻部分，高頻部分反映了軌跡在細節上的差異，而低頻部分可以看作是原始手術軌跡的平滑逼近[14].采用多尺度濾波處理手術軌跡可以在一定程度上消除噪聲，減少短程軌跡的影響.用雙通道多采樣率濾波器組表示的小波分解過程如圖4所示.

圖4 小波分解示意圖Fig.4 Diagram of wavelet decomposition

本文采用db10小波對運動學數據進行5級濾波處理，以去除由于專家遠程操作手術所引入的小尺度噪聲.小波濾波前后的右臂末端的位置隨時間的變化如圖5所示，縱坐標為末端所處的位置以米(m)為單位，橫坐標為時間以幀數表示(30幀/s)，x、y、z表示末端執行器的位置，xw、yw、zw為對應的小波濾波結果.

圖5 手術軌跡小波濾波前后對比圖Fig.5 Comparison of surgery trajectory before and after filtering by wavelet

3.2 層次化非參混合模型

本文采用狄利克雷高斯混合模型(DPGMM)進行聚類.高斯混合模型是生成模型，考慮一個數據集X=xi，i=1…n，GMM模型將X視為從k個高斯組分中生成的，即：X～GMM(μi，σi)，i=1…k.而如何確定參數k是比較困難的，通常通過構建狄里克雷混合模型進行解決，即認為GMM的參數是從狄里克雷過程產生：(μi，σi，k)～DP(α，H).DPGMM采用EM算法進行擬合.

DPGMM仍需要定義最大組分個數kmax，本文在非參聚類過程中，首先采用K-means對輸入特征進行聚類(k=kmax)，將聚類結果作為DPGMM的初始化，然后采用EM算法對DPGMM進行細致的迭代優化，能夠在一定程度上減少迭代次數.

結合視覺特征和運動學特征的快速手術軌跡無監督分割算法偽代碼如下：

算法2.快速手術軌跡無監督分割算法

輸入：手術視頻圖像幀I，手術軌跡運動學數據K

輸出：一組時間分割點Segs

//使用訓練好的SCAE提取視覺特征，見算法1

1.V←SCAEtest(I)

2.X←wavelet(K)//使用小波變換處理運動學數據

//K-means對轉移進行聚類，計算聚類中心cs

4.cs←Kmeans(D)

//DPGMM對轉移進行聚類，用cs初始化

5.CD←DPGMM(D，cs)

6.S←Dt；?t，Ct≠Ct+1//轉移狀態識別

//Kmeans初始化的DPGMM在視覺空間聚類

7.CSV←KDPGMM(VS)

//Kmeans初始化的DPGMM在運動學空間聚類

8.CSVX←KDPGMM(XC)；?c∈CSV

//Kmeans初始化的DPGMM在時間軸聚類

9.CSVXT←KDPGMM(T)；T=t(c)；?c∈CSVX

//計算聚類簇中心

10.Segs←mean(T(c))；?c∈CSVXT

4 實驗結果與分析

為了評估快速手術軌跡無監督分割算法的性能，本文進行了2組實驗.第一組實驗主要驗證本算法自身的通用性和探究使用視頻特征的必要性.在第二組實驗中，將本文算法與同樣使用了視頻數據和運動學數據進行分割的TSC-DL[12]，以及僅使用運動學數據進行分割的TSC[11]，GMM[10]，HMMGMM[15]，HSMM[16]進行了軌跡分割性能對比.

本文實驗采用約翰霍普金斯大學提供的JIGSAWS數據集[17].該數據集采集自達芬奇通用醫療機器人系統，共分為視頻數據和運動學數據兩部分.其中視頻數據為固定攝像機對手術區域錄制的640×480大小的視頻；運動學數據共76維，包括遠程端雙臂的末端執行器38維數據和病人端雙臂末端執行器38維數據；38維數據(每個臂19維)包括：雙臂末端的3維笛卡爾坐標、9維表示的旋轉矩陣、3維的速度、3維的角速度和1維的夾持器角度.視頻和運動學數據采樣頻率均為30Hz.該數據集包括不同熟練程度的8名外科醫生完成的三種手術任務：縫合，針傳遞和打結.

表2 實驗環境配置

本實驗選取JIGSAWS數據集的一個子集進行實驗，共包括2種手術任務(針傳遞和針縫合)，每種手術任務包括5個專家(E)，3個中級專家(I)，3個非專家(N)的數據.運動學數據只使用病人端的38維數據.對運動學數據和視覺數據進行3倍降采樣以降低計算量.在特征提取之前分別為針傳遞和針縫合訓練一個視頻特征提取網絡(SCAE)，SCAE的輸入為640×480的彩色圖像，輸出為一個70維的特征向量，采用Adam[18]進行訓練.本實驗的機器配置見表2.

4.1 評估方法

本文采用歸一化互信息(NMI)和整體運行時間對算法進行評估.

歸一化互信息：給定預測聚類結果A和真實標記B，歸一化互信息用于衡量A，B兩個聚類結果的相似度，可由式(8)計算：

(8)

其中H(A)，H(B)分別為A，B的信息熵，I(A，B)為A和B的互信息，NMI范圍在0-1之間，0表示兩個聚類結果沒有相關關系，而1表示完全相關.

整體運行時間：統計分割所有手術示范過程中特征提取及進行聚類分割的總運行時間，沒有使用視頻數據的只統計聚類分割的運行時間.

4.2 本文提出方法的性能評估

為了探究使用視頻數據對本文方法的影響，本文在同一數據集上對僅使用運動學數據、僅使用視頻數據和使用運動學數據及視頻數據三種軌跡分割方法進行了實驗，NMI見表3.使用的數據集是JIGSAWS數據集的子集，包括全專家E(5個專家示范)，專家及中級專家E+I(3個專家示范和3個中等專家示范)，所有級別專家E+I+N(3個專家示范、3個中等專家示范及3個非專家示范).

表3 采用不同特征進行分割的NMI對比Table 3 NMI of segmentation using different features

對于縫合任務，從表3中，我們可以看出，結合運動學和視頻數據進行分割的NMI比僅使用運動學數據提高了至少4.7%；對于僅使用運動學數據來說，NMI隨著非專家示范的比例增長呈現逐漸下降的趨勢，是因為專家示范通常比非專家完成的更加順利和迅速.然而在使用運動學數據和視頻數據進行分割時，卻沒有一致的趨勢.我們認為這與縫合任務的特殊性和新手的一些不一致的錯誤有關.

對于針傳遞任務，從表3中，我們可以看出，結合運動學和視頻數據進行分割的NMI比僅使用運動學數據提高了至少37.9%，這比縫合任務要顯著很多.與縫合任務相似，對于僅使用運動學數據的分割方法來說，NMI隨著非專家示范的比例增長呈現逐漸下降的趨勢.然而在使用運動學數據和視頻數據進行分割時，與縫合任務不同，針傳遞任務NMI表現出與僅使用運動學數據進行分割相同的趨勢.我們猜測這是因為針傳遞任務的軌跡相比于縫合任務的軌跡更加規則.

圖6對無監督聚類的分割結果進行了可視化，從圖中我們發現使用視覺和運動學數據進行分割(圖6(b)，圖6(d))的效果要好于僅使用運動學數據進行分割(圖6(a)，圖6(c))；僅使用運動學數的分割方法通常會導致過度分割，從而導致一些分割片段無法被正確分割；而結合視覺和運動學數據的分割方法結合了兩種互補的信息，能夠從整體進行分割，而減少過度分割的發生.

4.3 不同軌跡分割方法的對比

實驗中分別采用HSMM，HMMGMM，GMM，TSC，TSC-DL以及本文算法對同一手術示范數據集進行分割；手術示范數據集包括全專家E(5個專家示范)，專家及中級專家E+I(3個專家示范和3個中等專家示范)，所有級別專家E+I+N(3個專家示范、3個中等專家示范及3個非專家示范).不同方法的分割精度和運行時間如表4和表5所示.表5對不同分割方法的總運行時間進行了統計，對于僅使用運動學數據的分割方法(HSMM、HMMGMM、GMM、TSC)，運行時間為其聚類分割的運行時間.對于使用視覺數據和運動學數據的分割方法(TSC-DL、本文方法)，運行時間為視頻特征提取和聚類分割運行時間的總和.

表4 不同分割方法的NMI對比Table 4 NMI of segmentation using different approachs

從表4中我們可以看出，使用了視頻數據的軌跡分割方法(本文算法和TSC-DL)比沒有使用視頻數據的軌跡分割方法(HSMM，HMMGMM，GMM，TSC)的NMI有明顯的提升.在所有的軌跡分割任務中，本文的方法達到了最好的NMI.在針傳遞任務中，本文方法比同樣使用了視頻數據的TSC-DL算法的NMI提高了8.9%-18.7%；在針縫合任務中，本文方法比TSC-DL算法的NMI提高了2.0%-7.5%.因為本文采用小波對運動學數據進行多尺度平滑處理，濾除小尺度噪聲而降低了短程軌跡的干擾，從而間接提升了整體分割準確性；基于卷積自編碼的特征提取器能夠提取圖像中的主要特征，而不受細節干擾.

我們也注意到在針傳遞任務中，大多數算法(包括本文算法)的NMI隨著非專家手術示范比例的增加都有相同的變化趨勢，即NMI隨著非專家手術示范比例的增加而越來越低，其中E+I的比僅有E的要好一些的原因可能是選取的數據集有一定特殊性，而TSC-DL因為采用了基于專家等級進行加權的聚類簇裁剪機制，能夠在非專家手術示范比例的增加的情況下表現的越來越好.

表5 不同分割方法運行時間對比(單位：秒)Table 5 Comparison of running time using different segmentation methods(unit：s)

圖7 不同分割方法運行時間柱狀圖Fig.7 Bar graph of the running time of segmentation using different methods

從圖7中可以明顯看出，結合視頻數據和運動學數據的分割方法都會比僅使用運動學數據的分割方法慢10倍左右，因為結合視頻數據和運動學數據的軌跡分割方法往往需要從視頻中提取特征，然而這一步是比較耗時的.基于層次化聚類的分割方法(TSC-DL、TSC、本文方法)比基于單次聚類的分割方法(GMM、HMMGMM、HSMM)要慢.

在結合視頻數據和運動學數據的分割方法中，本文提出的方法幾乎是TSC-DL速度的10倍，甚至優于僅使用運動學數據進行分割的TSC算法，因為本文方法采用了一個更簡單的無監督深度模型對視頻數據進行特征提取，大大降低了視頻特征提取階段的時間消耗.盡管采用K-means對混合模型進行初始化有助于速度的提升，但是提升有限.

本實驗表明，相較于TSC-DL方法，本文方法能夠將結合視覺特征和運動學特征進行分割的速度提升10倍，并且分割結果得到較高的NMI.

5 總 結

本文針對手術軌跡分割中視頻空間的高維特征和分割效果容易受到短程軌跡干擾的問題，提出了一種快速手術軌跡無監督分割方法.該算法通過堆疊卷積自編碼器提取手術軌跡視頻特征，在對視頻特征進行降維的同時學習到視頻的主要特征，然后結合小波多尺度濾波處理后的運動學軌跡進行聚類分割，減少了短程軌跡對分割的干擾.在JIGSAWS數據集上實驗表明，該算法能夠在保證分割精度的情況下，能夠有效縮短分割時間.該分割算法可以大幅縮短外科醫生的技能評估用時、推進機器人手術自動化，在機器人輔助醫療領域具有廣泛的應用場景.