在上一篇文章中,我們已經對心電訊號進行了預處理,將含有噪聲的訊號變得平滑,以便分類。本篇文章我們將正式開始利用深度學習對心電訊號進行分類識別。
卷積神經網路
不論是傳統機器學習,還是深度學習,分類的依據都是不同類別的資料中包含的不同特徵。要進行分類識別就需要對資料的特徵進行提取,但是二者的提取方式並不相同。對於傳統的機器學習而言,資料的特徵需要設計者或專業人員針對其特性進行手動提取,而深度學習則可以自動提取每類資料中的不同特徵。對於卷積神經網路CNN而言,能夠自動提取特徵的關鍵在於卷積操作。經過卷積操作提取的特徵往往會有冗餘,並且多次卷積會使神經網路的引數過多不便於訓練,所以CNN往往會在卷積層的後面跟上一個池化層。經過多次的卷積和池化後,較低層次的特徵就會逐步構成高層次的特徵,最後神經網路根據提取出的高層次特徵進行分類。
另外需要指出的是,為什麼在心電訊號分類中可以使用CNN呢。這是因為CNN具有的卷積操作具有區域性連線和權值共享的特徵。
- 區域性連線:用於區別不同種類的圖片所需的特徵只是整張圖片中的某些區域性區域,因此在進行卷積操作時使用的卷積核(感受野)可以只是幾個不同小區域,而不必使用整張圖片大小的卷積核(全連線)。這樣做不僅可以更好地表達不同的特徵,還能起到減少引數的作用。例如下圖,左邊是使用全連線的神經網路,右邊是使用區域性連線卷積核的網路。
- 權值共享:對於一類圖片而言,他們擁有相似的特徵,但是每張圖片中特徵的位置可能會有偏移。比如不同的人臉照片中眼睛的位置可能會有變化,很少有兩張照片眼睛的位置完全重合。對一張圖片進行卷積操作時,可以有多個卷積核來提取不同的特徵,但一個卷積核在進行移動的過程中其權值是保持不變的(當然不同卷積核的權值不共享)。這樣既能保證特徵提取不受位置的影響,還能減少引數的數量。
而心電訊號雖然是一維的,但是其中的特徵也滿足區域性連線和權值共享的條件,因此我們可以採用卷積神經網路對其分類。
構建深度學習的資料集
巧婦難為無米之炊,雖然我們已經有了預處理過的心電資料,但是這樣的資料是無法拿來直接進行分類學習的。所以我們要先構建符合深度學習模型使用的資料集。轉換的過程是首先從一條心電訊號中切分出符合要求的心拍作為樣本,然後將python list轉為numpy array,再經過亂序和切分,最終構成可供深度學習使用的資料集。這裡我們使用tf.keras提供的介面,可以直接使用numpy陣列型別,而不用再轉成TensorFlow的DataSet物件,對於訓練過程而言也更加簡單。
心拍的切分需要找到QRS波尖峰所在的位置。由於我們只訓練網路模型,我們這裡直接使用MIT-BIH資料集提供的人工標註,並在尖峰處向前取99個訊號點、向後取200個訊號點,構成一個完整的心拍。如果需要對真實測量的訊號進行識別分類,還要設計心拍的檢測演算法,後續我也可能會繼續做。
資料集根據用途分為訓練集、驗證集和測試集。訓練集用於訓練引數模型,驗證集用於模型訓練中準確率和誤差(損失函式)的檢驗,測試集用於訓練完成後對訓練效果的最終檢驗。可以類比學習、測驗和考試。這三者的資料結構都一致,只是包含的資料內容不同,每個訓練集都包含資料和標籤兩部分內容。資料是預處理後切分出的若干心拍的列表,標籤是每個心拍樣本對應的心電型別。
首先將上一篇的預處理步驟封裝成一個函式:
# 小波去噪預處理
def denoise(data):
# 小波變換
coeffs = pywt.wavedec(data=data, wavelet='db5', level=9)
cA9, cD9, cD8, cD7, cD6, cD5, cD4, cD3, cD2, cD1 = coeffs
# 閾值去噪
threshold = (np.median(np.abs(cD1)) / 0.6745) * (np.sqrt(2 * np.log(len(cD1))))
cD1.fill(0)
cD2.fill(0)
for i in range(1, len(coeffs) - 2):
coeffs[i] = pywt.threshold(coeffs[i], threshold)
# 小波反變換,獲取去噪後的訊號
rdata = pywt.waverec(coeffs=coeffs, wavelet='db5')
return rdata
然後將讀取資料和標註、心拍切分封裝成一個函式:
# 讀取心電資料和對應標籤,並對資料進行小波去噪
def getDataSet(number, X_data, Y_data):
ecgClassSet = ['N', 'A', 'V', 'L', 'R']
# 讀取心電資料記錄
print("正在讀取 " + number + " 號心電資料...")
record = wfdb.rdrecord('ecg_data/' + number, channel_names=['MLII'])
data = record.p_signal.flatten()
# 小波去噪
rdata = denoise(data=data)
# 獲取心電資料記錄中R波的位置和對應的標籤
annotation = wfdb.rdann('ecg_data/' + number, 'atr')
Rlocation = annotation.sample
Rclass = annotation.symbol
# 去掉前後的不穩定資料
start = 10
end = 5
i = start
j = len(annotation.symbol) - end
# 因為只選擇NAVLR五種心電型別,所以要選出該條記錄中所需要的那些帶有特定標籤的資料,捨棄其餘標籤的點
# X_data在R波前後擷取長度為300的資料點
# Y_data將NAVLR按順序轉換為01234
while i < j:
try:
lable = ecgClassSet.index(Rclass[i])
x_train = rdata[Rlocation[i] - 99:Rlocation[i] + 201]
X_data.append(x_train)
Y_data.append(lable)
i += 1
except ValueError:
i += 1
return
需要注意的是,上面的函式並沒有返回值,這是因為我們裝載心拍資料和樣本的列表X_data、Y_data包含了所有心電記錄中符合要求的心拍,需要從函式外傳入,並將得到的資料直接附加在列表末尾。這樣將心電訊號的編號、X_data、Y_data一同傳入,就能將所需資料儲存在X_data和Y_data中。
下面將所有心拍訊號(因為102和104沒有MLII導聯故去除)讀取到dataSet和lableSet兩個列表中,經過上面函式後,dataSet和lableSet都是一個(92192)的一維列表。其中dataSet中的每一個元素都是一個numpy的陣列,陣列中是一個元素都是一個心拍的300個訊號點,lableSet中的每一個元素是dataSet中一個陣列對應的標籤值(NAVLR對應01234)。reshape後將dataSet變為(92192,300)的列表,將lableSet變為(92192,1)的列表。然後對這兩個列表進行亂序處理,但是要保證二者之間的對應關係不改變。思路是先將兩個列表進行豎直方向的堆疊,變為一個列表train_ds,然後對其進行亂序處理,再拆分出亂序後的資料X和標籤Y。
由於tf.keras可以將輸入的資料集自動劃分成訓練集和測試集,所以只需要分出測試集即可。思路是先生成92192(總心拍個數)個數的隨機排列列表,然後擷取其中前30%個值作為索引,然後在資料集和標籤集中取出下標為這些索引的值,即得到測試資料集X_test和測試標籤集Y_test。
# 載入資料集並進行預處理
def loadData():
numberSet = ['100', '101', '103', '105', '106', '107', '108', '109', '111', '112', '113', '114', '115', '116', '117', '119', '121', '122', '123', '124', '200', '201', '202', '203', '205', '208', '210', '212', '213', '214', '215', '217', '219', '220', '221', '222', '223', '228', '230', '231', '232', '233', '234']
dataSet = []
lableSet = []
for n in numberSet:
getDataSet(n, dataSet, lableSet)
# 轉numpy陣列,打亂順序
dataSet = np.array(dataSet).reshape(-1, 300)
lableSet = np.array(lableSet).reshape(-1, 1)
train_ds = np.hstack((dataSet, lableSet))
np.random.shuffle(train_ds)
# 資料集及其標籤集
X = train_ds[:, :300].reshape(-1, 300, 1)
Y = train_ds[:, 300]
# 測試集及其標籤集
shuffle_index = np.random.permutation(len(X))
test_length = int(RATIO * len(shuffle_index)) # RATIO = 0.3
test_index = shuffle_index[:test_length]
X_test, Y_test = X[test_index], Y[test_index]
return X, Y, X_test, Y_test
經過上面的函式後,X,Y為總體資料集和標籤集,X_test,Y_test為測試資料集和標籤集,驗證資料集和測試集使用tf.keras介面自動劃分。這樣用於深度學習的資料集就已經構建好了。
深度學習識別分類
通常來說,深度學習神經網路的訓練過程程式設計較為複雜,但是我們這裡使用tf.keras高階介面,可以很方便地進行深度學習網路模型的構建。
首先我們構建網路結構,具體結構如下圖所示:
# 構建CNN模型
def buildModel():
newModel = tf.keras.models.Sequential([
tf.keras.layers.InputLayer(input_shape=(300, 1)),
# 第一個卷積層, 4 個 21x1 卷積核
tf.keras.layers.Conv1D(filters=4, kernel_size=21, strides=1, padding='SAME', activation='relu'),
# 第一個池化層, 最大池化,4 個 3x1 卷積核, 步長為 2
tf.keras.layers.MaxPool1D(pool_size=3, strides=2, padding='SAME'),
# 第二個卷積層, 16 個 23x1 卷積核
tf.keras.layers.Conv1D(filters=16, kernel_size=23, strides=1, padding='SAME', activation='relu'),
# 第二個池化層, 最大池化,4 個 3x1 卷積核, 步長為 2
tf.keras.layers.MaxPool1D(pool_size=3, strides=2, padding='SAME'),
# 第三個卷積層, 32 個 25x1 卷積核
tf.keras.layers.Conv1D(filters=32, kernel_size=25, strides=1, padding='SAME', activation='relu'),
# 第三個池化層, 平均池化,4 個 3x1 卷積核, 步長為 2
tf.keras.layers.AvgPool1D(pool_size=3, strides=2, padding='SAME'),
# 第四個卷積層, 64 個 27x1 卷積核
tf.keras.layers.Conv1D(filters=64, kernel_size=27, strides=1, padding='SAME', activation='relu'),
# 打平層,方便全連線層處理
tf.keras.layers.Flatten(),
# 全連線層,128 個節點
tf.keras.layers.Dense(128, activation='relu'),
# Dropout層,dropout = 0.2
tf.keras.layers.Dropout(rate=0.2),
# 全連線層,5 個節點
tf.keras.layers.Dense(5, activation='softmax')
])
return newModel
然後使用model.compile()構建;model.fit()訓練30輪,批大小為128,劃分驗證集的比例為0.3,設定callback進行訓練記錄的儲存;model.save()儲存模型;model.predict_classes()預測。完整程式碼可以取本人的GitHub倉庫檢視,地址在文章(一)中。
def main():
# X,Y為所有的資料集和標籤集
# X_test,Y_test為拆分的測試集和標籤集
X, Y, X_test, Y_test = loadData()
if os.path.exists(model_path):
# 匯入訓練好的模型
model = tf.keras.models.load_model(filepath=model_path)
else:
# 構建CNN模型
model = buildModel()
model.compile(optimizer='adam',
loss='sparse_categorical_crossentropy',
metrics=['accuracy'])
model.summary()
# 定義TensorBoard物件
tensorboard_callback = tf.keras.callbacks.TensorBoard(log_dir=log_dir, histogram_freq=1)
# 訓練與驗證
model.fit(X, Y, epochs=30,
batch_size=128,
validation_split=RATIO,
callbacks=[tensorboard_callback])
model.save(filepath=model_path)
# 預測
Y_pred = model.predict_classes(X_test)
對心電訊號的深度學習識別分類至此結束,識別率可達99%左右。