Main Content

最新のリリースでは、このページがまだ翻訳されていません。 このページの最新版は英語でご覧になれます。

深いネットワークを使用した ECG 信号の QRS 群と R ピークのラベル付け

この例では、信号ラベラーのカスタム自動ラベル付け関数を使用して、心電図 (ECG) 信号の QRS 群と R ピークにラベル付けする方法を示します。1 つのカスタム関数は、学習済みの再帰型深層学習ネットワークを使用し、QRS 群を特定して検索します。もう 1 つのカスタム関数は、単純なピークの検索を使用して R ピークを検索します。この例では、ネットワークは、ネットワークの学習とテストのプロセスから完全に独立した 2 つの信号の QRS 群にラベル付けします。

ECG 波形の 3 つの偏位で構成される QRS 群は、心臓の左右の心室の脱分極を反映しています。また、QRS は人間の心拍の最大振幅セグメントでもあります。QRS 群の調査は、人の心臓の全体的な健康状態と異常の有無を評価するために役立ちます [1]。特に、QRS 群内の R ピークを検出し、連続するピークの時間間隔を調べることで、診断医は患者の心拍変動を計算して、心不整脈を見つけることができます。

この例の深層学習ネットワークは深層学習を使用した波形セグメンテーションで導入されており、公開されている QT データベースの ECG 信号を使用して学習させています [2] [3]。データは、合計 105 人の患者からの約 15 分間の ECG の記録で構成され、250 Hz でサンプリングされています。各記録を取得するために、検査員は患者の胸部の異なる場所に 2 つの電極を配置して、2 チャネル信号にしています。データベースは、自動化されたエキスパート システムによって生成される信号領域ラベルを提供します [1]。追加されたラベルにより、深層学習ネットワークにデータを使用できるようになっています。

信号ラベラーへのデータの読み込み、リサンプリング、およびインポート

この例では、MIT-BIH Arrhythmia Database の信号をラベル付けします [4]。データベース内の各信号は、平均 360 Hz のレートで不規則にサンプリングされ、2 人の心臓専門医によって注釈が付けられ、結果の検証が可能になっています。

記録 200 および 203 に対応する 2 つの MIT データベース信号を読み込みます。信号をサンプル時間 1/250 秒の等間隔グリッドにリサンプリングします。これは、QT データベース データの定格サンプル レートに対応します。

load mit200
y200 = resample(ecgsig,tm,250);

load mit203
y203 = resample(ecgsig,tm,250);

信号ラベラー を開きます。[ラベラー] タブで、[インポート] ▼ をクリックし、[Members From Workspace] を選択します。ダイアログ ボックスで、信号 y200 および y203 を選択します。時間情報を追加します。ドロップダウン リストから Time を選択し、サンプル レートを 250 Hz に指定します。[インポート] をクリックし、ダイアログ ボックスを閉じます。[ラベル付き信号セット] ブラウザーに信号が表示されます。名前の隣にあるチェック ボックスをオンにして信号をプロットします。

ラベルの定義

ラベルを定義して、信号に付加します。

  1. QRS 群のカテゴリカル関心領域 (ROI) ラベルを定義します。[ラベラー] タブの [定義の追加] をクリックします。[ラベル名] として QRSregions を指定し、[ラベル タイプ]ROI を選択して、[データ型] として categorical と入力し、2 つの [カテゴリ]QRS および n/a を 1 行に 1 つずつ追加します。

  2. QRSregions のサブラベルを R ピークの数値の点ラベルとして定義します。[ラベルの定義] ブラウザーで QRSregions をクリックして選択します。[定義の追加] ▼ をクリックし、Add sublabel definition を選択します。[ラベル名]Rpeaks を指定して、[ラベル タイプ]Point を選択し、[データ型] として numeric と入力します。

カスタム自動ラベル付け関数の作成

2 つのカスタム ラベル付け関数を作成します。1 つは QRS 群を検索してラベル付けし、もう 1 つは各 QRS 群内の R ピークを検索してラベル付けします。(関数 findQRScomputeFSSTp2qrs、および findRpeaks のコードは、後ほどこの例で示します)。各関数を作成するために、[ラベラー] タブの [値の自動処理] ▼ をクリックして、[カスタム関数の追加] を選択します。信号ラベラーにより、関数の名前、説明、およびラベル タイプを求めるダイアログ ボックスが表示されます。

  1. QRS 群を検索する関数については、名前フィールドに findQRS と入力し、[ラベル タイプ] として ROI を選択します。説明フィールドは空のままにすることも、独自の説明を入力することもできます。

  2. R ピークを検索する関数については、名前フィールドに findRpeaks と入力し、[ラベル タイプ] として Point を選択します。説明フィールドは空のままにすることも、独自の説明を入力することもできます。

関数を書き込み済みでその関数が現在のフォルダーまたは MATLAB® パスにある場合、信号ラベラーはその関数をギャラリーに追加します。関数をまだ書き込んでいない場合、コードの入力または貼り付けができるように、信号ラベラーはエディターで空白のテンプレートを開きます。ファイルを保存します。ファイルを保存すると、関数がギャラリーに表示されます。

QRS 群と R ピークのラベル付け

入力信号の QRS 群を見つけてラベル付けします。

  1. [ラベル付き信号セット] ブラウザーで、y200 の横にあるチェック ボックスをオンにします。

  2. [ラベルの定義] ブラウザーで QRSregions を選択します。

  3. [値の自動処理] ギャラリーで、findQRS を選択します。

  4. [自動ラベル付け] をクリックし、Auto-Label Signals を選択します。表示されるダイアログ ボックスで [OK] をクリックします。

信号ラベラーは、すべての信号の QRS 群を検索してラベル付けしますが、チェック ボックスをオンにした信号のラベルのみを表示します。QRS 群が、プロットとラベル ビューアーの座標軸に影付きの領域として表示されます。[表示] タブの [パナー] をクリックしてパナーをアクティブにし、ラベル付き信号の領域を拡大します。

QRS 群に対応する R ピークを見つけてラベル付けします。

  1. [ラベルの定義] ブラウザーで Rpeaks を選択します。

  2. [ラベラー] タブに戻ります。[値の自動処理] ギャラリーで、findRpeaks を選択します。

  3. [自動ラベル付け] をクリックし、Auto-Label Signals を選択します。表示されるダイアログ ボックスで [OK] をクリックします。

ラベルとその数値が、プロットとラベル ビューアーの座標軸に表示されます。

ラベル付き信号のエクスポートと心拍変動の計算

ラベル付き信号をエクスポートして、各患者の心拍変動を比較します。[ラベラー] タブで、[エクスポート] ▼ をクリックし、Labeled Signal Set To File を選択します。表示されるダイアログ ボックスで、ラベル付き信号セットに HeartRates.mat という名前を付け、オプションで簡単な説明を追加します。[エクスポート] をクリックします。

MATLAB® コマンド ウィンドウに戻ります。ラベル付き信号セットを読み込みます。セット内の各信号について、連続する心拍間の時間差の標準偏差として心拍変動を計算します。差のヒストグラムをプロットし、心拍変動を表示します。

load HeartRates

nms = getMemberNames(heartrates);

for k = 1:heartrates.NumMembers
    
    v = getLabelValues(heartrates,k,{'QRSregions','Rpeaks'});
    
    hr = diff(cellfun(@(x)x.Location,v));
    
    subplot(2,1,k)
    histogram(hr,0.5:.025:1.5)
    legend(['hrv = ' num2str(std(hr))])
    ylabel(nms(k))
    ylim([0 6])

end

関数 findQRS: QRS 群の検索

関数 findQRS は、入力信号の QRS 群を見つけてラベル付けします。

この関数は、computeFSSTp2qrs の 2 つの補助関数を使用します。(両方の補助機能のコードは、後ほどこの例で示します)。同じディレクトリ内の別のファイルに関数を保存するか、最後の end ステートメントの前に関数を挿入して、関数 findQRS 内で関数を入れ子にすることができます。

computeFSSTp2qrs の呼び出しの間で、findQRS は関数 classify と学習済みの深層学習ネットワーク net を使用して、QRS 領域を特定します。classify を呼び出す前に、深層学習を使用した波形セグメンテーションで説明されているように、findQRS はデータを net で期待される形式に変換します。

  • 各信号は 250 Hz でサンプリングされ、2 行 N 列の cell 配列のスタックに分割されなければなりません。ここで、各行はチャネルに対応し、N は 5000 の倍数です。実際の分割とスタックは関数 computeFSST で行われます。

  • リサンプリングされた各 MIT 信号には 6945 個のサンプルがあり、5000 の倍数ではありません。各信号のすべてのデータを保持するには、信号を乱数でパディングします。プロセスの後半で、関数 p2qrs がこの乱数を QRS 群に属さないものとしてラベル付けし、破棄します。

function [labelVals,labelLocs] = findQRS(x,t,parentLabelVal,parentLabelLoc,varargin)
% This is a template for creating a custom function for automated labeling
%
%  x is a matrix where each column contains data corresponding to a
%  channel. If the channels have different lengths, then x is a cell array
%  of column vectors.
%
%  t is a matrix where each column contains time corresponding to a
%  channel. If the channels have different lengths, then t is a cell array
%  of column vectors.
%
%  parentLabelVal is the parent label value associated with the output
%  sublabel or empty when output is not a sublabel.
%  parentLabelLoc contains an empty vector when the parent label is an
%  attribute, a vector of ROI limits when parent label is an ROI or a point
%  location when parent label is a point.
%
%  labelVals must be a column vector with numeric, logical or string output
%  values.
%  labelLocs must be an empty vector when output labels are attributes, a
%  two column matrix of ROI limits when output labels are ROIs, or a column
%  vector of point locations when output labels are points.

labelVals = [];
labelLocs = [];

Fs = 250;

load('trainedQTSegmentationNetwork','net')

for kj = 1:size(x,2)

    sig = x(:,kj);
    
    % Create 10000-sample signal expected by the deep network
    
    sig = [sig;randn(10000-length(sig),1)/100]';
    
    % Resize input and compute synchrosqueezed Fourier transforms

    mitFSST = computeFSST(sig,Fs);
    
    % Use trained network to predict which points belong to QRS regions
    
    netPreds = classify(net,mitFSST,'MiniBatchSize',50);
    
    % Convert stack of cell arrays into a single vector
    
    Location = [1:length(netPreds{1}) length(netPreds{1})+(1:length(netPreds{2}))]';
    Value = [netPreds{1} netPreds{2}]';
    
    % Label QRS complexes as regions of interest and discard non-QRS data
    
    [Locs,Vals] = p2qrs(table(Location,Value));
    
    labelVals = [labelVals;Vals];
    labelLocs = [labelLocs;Locs/Fs];

end

% Insert computeFSST and p2qrs here if you want to nest them inside
% queryQRS instead of including them as separate functions in the folder.

end

関数 computeFSST: 入力のサイズ変更とフーリエ シンクロスクイーズド変換の計算

この関数は、入力データを net で期待される形式に整形し、関数fsstを使用して入力のフーリエ シンクロスクイーズド変換 (FSST) を計算します。深層学習を使用した波形セグメンテーションで、ネットワークのパフォーマンスは、各学習信号またはテスト信号の時間-周波数マップを入力として使用した場合に最も高くなります。FSST は変換が元の入力と同じ時間分解能を持つため、再帰型ネットワークに特に役立つ一連の機能を発揮します。関数は以下を実行します。

  • 適切な周波数分解能を与えるために、長さが 128 のカイザー ウィンドウを指定します。

  • 0.5 Hz ~ 40 Hz の周波数範囲でデータを抽出します。

  • 各信号の平均で減算し、標準偏差で除算します。

  • FSST の実数部と虚数部を個別の特徴として扱います。

function signalsFsst = computeFSST(xd,Fs)

targetLength = 5000;
signalsOut = {};

for sig_idx = 1:size(xd,1)

    current_sig = xd(sig_idx,:)';

    % Compute the number of targetLength-sample chunks in the signal
    
    numSigs = floor(length(current_sig)/targetLength);

    % Truncate to a multiple of targetLength
    
    current_sig = current_sig(1:numSigs*targetLength);

    % Create a matrix with as many columns as targetLength signals
    
    xM = reshape(current_sig,targetLength,numSigs);

    % Vertically concatenate into cell arrays
    
    signalsOut = [signalsOut; mat2cell(xM.',ones(numSigs,1))];

end

signalsFsst = cell(size(signalsOut));

for idx = 1:length(signalsOut)

   [s,f] = fsst(signalsOut{idx},Fs,kaiser(128));

   % Extract data over the frequency range from 0.5 Hz to 40 Hz
   
   f_indices = (f > 0.5) & (f < 40);
   signalsFsst{idx}= [real(s(f_indices,:)); imag(s(f_indices,:))];

   % Subtract the mean and divide by the standard deviation
   
   signalsFsst{idx} = (signalsFsst{idx}-mean(signalsFsst{idx},2)) ...
       ./std(signalsFsst{idx},[],2);

end

end

関数 p2qrs: 関心領域としての QRS 群のラベル付け

深層ネットワークは、入力信号のすべての点を、P 領域、QRS 群、T 領域に属するとして、またはいずれにも属さないとしてラベル付けする categorical 配列を出力します。この関数は、これらの点ラベルを QRS 関心領域ラベルに変換します。

  • 変換を実行するために、この関数は整数の数値をカテゴリに割り当て、関数findchangeptsを使用して数値配列の値が変化する点を見つけます。

  • これらの各変化点はカテゴリカル領域の左の端点であり、配列内でその前にある点は、前の領域の右の端点です。

  • このアルゴリズムは、1e-6 を右の端点に追加して、1 サンプル領域の持続時間がゼロになることを防ぎます。

  • df パラメーターは、持続時間が df サンプルよりも長い QRS 群のみを関心領域として選択します。

function [locs,vals] = p2qrs(k)

fc = 1e-6;
df = 20;

ctgs = categories(k.Value);
levs = 1:length(ctgs);
for jk = levs
   cat2num(k.Value == ctgs{jk}) = levs(jk);
end
chpt = findchangepts(cat2num,'MaxNumChanges',length(cat2num));
locs = [[1;chpt'] [chpt'-1;length(cat2num)]+fc];

vals = categorical(cat2num(locs(:,1))',levs,ctgs);
locs = locs+round(k.Location(1))-1;

qrs = find(vals=='QRS' & diff(locs,[],2)>df);

vals = categorical(string(vals(qrs)),["QRS" "n/a"]);

locs = locs(qrs,:);

end

関数 findRpeaks: R ピークの検索

この関数は、findQRS により見つかった QRS 関心領域の最も突出したピークを検索します。この関数は、MATLAB® 関数 islocalmax を、findQRS により見つかった区間の信号の絶対値に適用します。

function [labelVals,labelLocs] = findRpeaks(x,t,parentLabelVal,parentLabelLoc,varargin)

labelVals = zeros(size(parentLabelLoc,1),1);
labelLocs = zeros(size(parentLabelLoc,1),1);

for kj = 1:size(parentLabelLoc,1)
    tvals = t>=parentLabelLoc(kj,1) & t<=parentLabelLoc(kj,2);
    ti = t(tvals);
    xi = x(tvals);
    lc = islocalmax(abs(xi),'MaxNumExtrema',1);
    labelVals(kj) = xi(lc);
    labelLocs(kj) = ti(lc);
end

end

参考文献

[1] Laguna, Pablo, Raimon Jané, and Pere Caminal."Automatic Detection of Wave Boundaries in Multilead ECG Signals: Validation with the CSE Database." Computers and Biomedical Research.Vol. 27, No. 1, 1994, pp. 45–60.

[2] Goldberger, Ary L., Luis A. N. Amaral, Leon Glass, Jeffery M. Hausdorff, Plamen Ch. Ivanov, Roger G. Mark, Joseph E. Mietus, George B. Moody, Chung-Kang Peng, and H. Eugene Stanley. "PhysioBank, PhysioToolkit, and PhysioNet: Components of a New Research Resource for Complex Physiologic Signals." Circulation. Vol. 101, No. 23, 2000, pp. e215–e220. [Circulation Electronic Pages: http://circ.ahajournals.org/content/101/23/e215.full].

[3] Laguna, Pablo, Roger G. Mark, Ary L. Goldberger, and George B. Moody."A Database for Evaluation of Algorithms for Measurement of QT and Other Waveform Intervals in the ECG."Computers in Cardiology.Vol. 24, 1997, pp. 673–676.

[4] Moody, George B., and Roger G. Mark."The impact of the MIT-BIH Arrhythmia Database."IEEE Engineering in Medicine and Biology Magazine.Vol. 20, No. 3, May–June 2001, pp. 45–50.

参考

アプリ

関数

関連する例

詳細