信号の特徴抽出と機械学習を使用した行動認識

この例では次を使用します。

この例では、スマートフォンの加速度計信号から特徴を抽出し、機械学習アルゴリズムを使用して人の行動を分類する方法を説明します。データの特徴抽出は、signalTimeFeatureExtractor オブジェクトおよび signalFrequencyFeatureExtractor オブジェクトを使用して行います。この特徴を、サポートベクターマシン (SVM) モデルの学習に使用します。

データセット

Sensor HAR (行動認識) アプリ [1] を使用して、[2] の加速度計の生の信号を収集しました。5 種類の身体活動を行う間、被験者はスマートフォンを身に着けていました。データセットをバッファーし、特定の活動に対応する 44 サンプル長の信号を取得しました。データセットの "ダンス" 活動および y 方向と z 方向の加速度計信号は除外し、BufferedHumanActivity データセットを作成して、この例で使用する BufferedHumanActivity.mat ファイルに保存しました。

BufferedHumanActivity データセットを読み込みます。

load BufferedHumanActivity.mat

このデータセットには 7776 個の x 方向の加速度計信号が格納されています。各信号の持続時間は 44 サンプル分あり、4 種類の身体活動 ("座る"、"立つ"、"歩く"、"走る") のいずれかに対応しています。このデータセットには、以下の変数が含まれています。

atx — x 方向の固定長バッファー済み加速度計センサーデータ (44 行 7776 列の行列)
actid — 応答ベクトル。整数で表される活動 ID を含みます。整数 1、2、3、4 は、それぞれ "座る"、"立つ"、"歩く"、"走る" を意味します。
actnames — 各活動 ID に対応する活動名のリスト
fs — 加速度計センサーデータのサンプルレート

特徴抽出

加速度計信号には、2 つの主成分が含まれていると考えられます。その 1 つが、身体のダイナミクス (被験者の身体の動き) によって発生する "速い" 経時変動です。もう 1 つが、垂直な重力場に対する身体の位置によって発生する "遅い" 経時変動です。

遅い信号変動から速い信号変動を分離するため、元の信号に対してハイパスフィルターを適用します。signalTimeFeatureExtractor オブジェクトと signalFrequencyFeatureExtractor オブジェクトを使用して、フィルター処理された信号とフィルター処理されていない信号からさまざまな特徴を抽出します。これらのオブジェクトにより、1 回の関数呼び出しで時間領域と周波数領域の複数の特徴を効率よく計算できます。

% Filter the signals with a highpass filter
atxFiltered = highpass(atx,0.7,fs);

時間の特徴の場合は、2 つの signalTimeFeatureExtractor オブジェクトが構成されます。1 つを使用して、フィルター処理されていない信号の平均を抽出し (meanFE)、もう 1 つを使用して、フィルター処理された信号の平方根平均二乗、形状係数、ピーク値、クレストファクター、クリアランス係数、インパルス係数を抽出します (timeFE)。

meanFE = signalTimeFeatureExtractor("Mean",true,"SampleRate",fs);
timeFE = signalTimeFeatureExtractor("RMS",true,...
    "ShapeFactor",true,...
    "PeakValue",true,...
    "CrestFactor",true,...
    "ClearanceFactor",true,...    
    "ImpulseFactor",true,...
    "SampleRate",fs);

周波数の特徴の場合は、signalFrequencyFeatureExtractor を使用して、フィルター処理された信号の周波数平均、帯域パワー、パワーが半分の帯域幅、ピーク振幅、およびピーク位置を抽出します。

freqFE = signalFrequencyFeatureExtractor("PeakAmplitude",true,...
    "PeakLocation",true,...
    "MeanFrequency",true,...
    "BandPower",true,...
    "PowerBandwidth",true,...
    "SampleRate",fs);

さらにパラメーターを設定することで、スペクトルピークの計算を調整できます。たとえば、ピークの最大数が 6、各スペクトルピーク間の最小距離が 0.25 Hz に設定されているとします。さらに、256 の FFT 長および 44 サンプル (すなわち、信号長) の箱型ウィンドウを選択し、スペクトル推定を計算します。

fftLength = 256;
window = rectwin(size(atx,1));
setExtractorParameters(freqFE,"WelchPSD","FFTLength",fftLength,"Window",window);
mindist_xunits = 0.25;
minpkdist = floor(mindist_xunits/(fs/fftLength));
setExtractorParameters(freqFE,"PeakAmplitude","MaxNumExtrema",6,"MinSeparation",minpkdist);
setExtractorParameters(freqFE,"PeakLocation","MaxNumExtrema",6,"MinSeparation",minpkdist);

すべての信号の特徴の計算は、変換後の配列データストアを使用して並列化できます。特徴抽出器オブジェクトの関数 extract を使用して、データストアが各行列の列を読み取り、特徴を計算します。

meanFeatureDs = arrayDatastore(atx,"IterationDimension",2);
meanFeatureDs = transform(meanFeatureDs,@(x)meanFE.extract(x{:}));
timeFeatureDs = arrayDatastore(atxFiltered,"IterationDimension",2);
timeFeatureDs = transform(timeFeatureDs,@(x)timeFE.extract(x{:}));
freqFeatureDs = arrayDatastore(atxFiltered,"IterationDimension",2);
freqFeatureDs = transform(freqFeatureDs,@(x)freqFE.extract(x{:}));

Parallel Computing Toolbox がインストールされている場合は、"UseParallel" オプションを true に設定して変換後のデータストアの readall メソッドを呼び出し、ワーカーのプール間で計算を分散させます。結果として得られる計算された特徴は結合され、7776 個の信号観測値のそれぞれについて 22 個の特徴が得られることになります。

meanFeatures = readall(meanFeatureDs,"UseParallel",true);

Starting parallel pool (parpool) using the 'Processes' profile ...
Connected to parallel pool with 8 workers.

timeFeatures = readall(timeFeatureDs,"UseParallel",true);
freqFeatures = readall(freqFeatureDs,"UseParallel",true);

features = [meanFeatures timeFeatures freqFeatures];

抽出された特徴を使用した SVM 分類器の学習

特徴と活動ラベルを分類学習器アプリにインポートして、SVM 分類器を学習させることができます。あるいは、特徴 (予測子) と活動ラベル (応答) を含む、以下のような特徴テーブルを使用して、SVM テンプレートと分類器を作成することもできます。

まず、予測子と応答を含むテーブルを作成します。

featureTable = array2table(features);
actioncats = categorical(actnames)';
featureTable.ActivityID = actioncats(actid);
head(featureTable)

    features1    features2    features3    features4    features5    features6    features7    features8    features9    features10    features11    features12    features13    features14    features15    features16    features17    features18    features19    features20    features21    features22    ActivityID
    _________    _________    _________    _________    _________    _________    _________    _________    _________    __________    __________    __________    __________    __________    __________    __________    __________    __________    __________    __________    __________    __________    __________

     -73.408      0.10678      1.2695       0.24501      2.2946       3.5282        2.913       3.1208       0.011402     0.22658      0.0037348     0.0043388      0.0049913      0.014314    0.0032949     0.0042457      0.74219        1.6797        3.2031        3.8281        4.2188        4.5703       Sitting  
      -73.43      0.06735      1.2521       0.13304      1.9753       2.9553       2.4733       1.8959       0.004536     0.18083      0.0078615      0.001071      0.0046553     0.0023938    0.0017709      0.002039      0.74219          1.25        1.5625        2.3438        3.5938        3.9453       Sitting  
      -73.41       0.0626       1.303       0.15569       2.487       3.9603       3.2407       2.4191      0.0039188     0.18783      0.0036916      0.001265     0.00086816    0.00098286    0.0029621     0.0044119      0.74219        1.4062        2.2266        2.7734        3.0859        4.6094       Sitting  
     -73.393     0.072056      1.3457       0.20023      2.7788       4.6384       3.7394       2.9361       0.005192     0.21444      0.0028194     0.0016623      0.0028484     0.0018433     0.003666     0.0026144      0.89844        1.6797        2.3047        3.2422        4.0234        4.6484       Sitting  
     -73.409     0.080133      1.3786       0.21548       2.689        4.602       3.7069       3.2548      0.0064214      0.2053      0.0035392     0.0015361      0.0061205     0.0010848    0.0072086     0.0055945       1.5625        2.3828        3.0469        3.5156        3.8672        4.6484       Sitting  
      -73.43     0.071148      1.1902       0.13832      1.9441       2.6268       2.3139       3.0519      0.0050621     0.25175      0.0022982     0.0027692      0.0040954     0.0045089    0.0016846      0.003589      0.82031        2.3047        3.1641        3.9062        4.2188        4.5312       Sitting  
     -73.441     0.091667       1.169       0.19139      2.0879       2.7515       2.4408       2.8127      0.0084028     0.25907      0.0021497     0.0029254      0.0035706     0.0018514     0.015439     0.0030516      0.89844        2.1094        2.3828        2.6562        3.0859        4.5703       Sitting  
     -73.419      0.10858      1.1976       0.20506      1.8886       2.5625       2.2619       2.3954       0.011789     0.17288       0.010823     0.0088772      0.0078451     0.0071845    0.0066219     0.0024052      0.74219        1.5625        2.2656        3.0469        3.8281        4.5312       Sitting

信号の 75% を学習用、25% をテスト用に割り当て、データセットを分割します。関数 cvpartition を使用して、各区画に活動ラベルが同じような比率で含まれるようにします。

% Extract predictors and response
predictors = featureTable(:, 1:end-1);
response = featureTable.ActivityID;

% For reproducible results
rng default

% Partition the data and extract training predictors and response data
cvp = cvpartition(response,'Holdout',0.25);
trainingPredictors = predictors(cvp.training, :);
trainingResponse = response(cvp.training, :);

% Train the classifier
template = templateSVM(...
    'KernelFunction', 'polynomial', ...
    'PolynomialOrder', 2, ...
    'KernelScale', 'auto', ...
    'BoxConstraint', 1, ...
    'Standardize', true);
classificationSVM = fitcecoc(...
    trainingPredictors, ...
    trainingResponse, ...
    'Learners', template, ...
    'Coding', 'onevsone', ...
    'ClassNames',actioncats);

テスト区画の分類器をテストし、その分類精度を解析します。

% Extract test predictors and response data
testPredictors = predictors(cvp.test, :);
testResponse = response(cvp.test, :);

% Predict activity on the test data
testPredictions = predict(classificationSVM,testPredictors);

% Plot the confusion matrix to analyze performance of the classifier
figure
cm = confusionchart(testResponse, testPredictions, ...
    ColumnSummary="column-normalized", RowSummary="row-normalized");

accuracy = trace(cm.NormalizedValues)/sum(cm.NormalizedValues, "all");
fprintf("The classification accuracy on the test partition is %2.1f%%", accuracy*100)

The classification accuracy of the classifier on the test partition is 95.0%

誤差の大半は、走るを歩くに、立つを座るに誤分類した場合に発生します。

まとめ

この例では、signalTimeFeatureExtractor および signalFrequencyFeatureExtractor を使用して、スマートフォンのセンサー信号に基づいて人の行動の特徴を抽出する方法を確認しました。抽出された特徴を使用して SVM モデルの学習を行ったところ、約 95% の精度が得られました。別の方法として、featureInput 層を使用して深層学習分類器の学習を確認することもできます。

参考文献

[1] El Helou, A. Sensor HAR recognition App. MathWorks File Exchange https://www.mathworks.com/matlabcentral/fileexchange/54138-sensor-har-recognition-app

[2] El Helou, A. Sensor Data Analytics. MathWorks File Exchange https://www.mathworks.com/matlabcentral/fileexchange/54139-sensor-data-analytics-french-webinar-code