機械学習のための
信号処理

MathWorks
アプリケーションエンジニアリング部
竹本佳充© 2016 The MathWorks, Inc.
1

Scania Emergency Braking System
2

Scania’s Analytics-Driven
Emergency Braking System• Sensor fusion algorithms
developed using machinelearning trained with
80TB of driving data
• Algorithms deployed on a
production embedded
system using ModelBased Design“Sensor fusion brings the sensor data together and, algorithmically,
processes the combined sensor data to understand the real-timesituation and apply the right control strategy to ensure braking occurs
when necessary.”

3

Celestion社
楽器用ラウドスピーカーの分類に機械学習技術を適用

課題
楽器用ラウドスピーカーの評価における、マニュアルでの官能評価プロセスを自動化

ソリューション
ラウドスピーカーの高精度かつ高確度な分類のための
機械学習アルゴリズム構築にMATLABを使用

結果▪ プロセスを自動化することにより、マニュアルで工数のかかる作業を数秒で実現▪ 多種多様な機械学習および信号処理技術を評価し適用▪ エラー率を0.1%程度に低減

Celestion社のラウドスピーカーテスト用リスニングルーム“MATLABと Statistics and Machine
Learning Toolboxの使用により, 一度
データを用意すれば, 多種多様な分類手法を数分で評価することができました. もしそれらをマニュアルでコーディングしていたら、数年はかかっていたと思います. ”

Link to article

Andrew Harper
Celestion
4

2016 with
MATLAB

1977
Original
RecordingRemastering The Beatles’ Only Live Recording
5

“The Beatles Live at the Hollywood Bowl”
Remastered and Remarkable•

Used MATLAB, Signal Processing Toolbox, andrelated examples already written in MATLAB•

Developed tool that uses spectrum analysis,
optimization, and musical modelling“It doesn’t exist as a software program that is easy to use. It’s a lot of Matlab, more like a
research tool. There’s no graphical front end where you can just load a piece of audio up, paint a
track, and extract the audio. I write manual scripts, which I then put into the engine to process.”- James Clarke, Abbey Road Studios
6

CellScope’s Oto consumer otoscope
7

CellScope’s Oto Diagnoses
Ear Infections with a Smartphone• Visual analysis of theeardrum, by comparing
with a database of 10,000
expertly classified images• Developed using a full
spectrum of MATLAB tools:– Computer vision– Machine learning– Bioinformatics– Parallel computing– Code generation“What is great about having such an attachment is the power toimmediately share the image with the patient and better explain to them
what you are seeing...”8

Agenda 背景
 センサー信号の場合
 音声信号の場合
 画像の場合
 関連情報 まとめ

9Agenda 背景
 センサー信号の場合
 音声信号の場合
 画像の場合
 関連情報 まとめ

10

MATLABによるアナリティクスフロー
アクション

意思決定 & 設計• レポート & Apps
• 拡張性のある展開• 設計最適化

知識

理解
•
•
•
•

分析
周波数 & 時間領域
予測分析
推定

情報
Displacement Acceleration

フィルタリング
信号解析
データ処理
プロット

Weight

•
•
•
•

20
20
10
400
200

4000
2000

Horsepower

整理

MPG

40

200
150
100
50
20
MPG

40

10

20

Acceleration

200

400 2000

Displacement

4000

Weight50 100150200
Horsepow er

データ

観測

計測器

•
•
•
•

検出
収集
状態
データ取得
11

データアナリティクス・ワークフロー: 「機械学習」と「信号処理」
アクション

意思決定 & 設計• レポート & Apps
• 拡張性のある展開• 設計最適化

種類

センサー

サービス

防災
ヘルスケア
生産ライン
各種故障予知

知識

理解
•
•
•
•

分析
周波数 & 時間領域
予測分析
推定

音
モバイルアプリ公共サイン音
異音検査
定量評価

画像
セキュリティ
自動運転
医療診断

機械学習の適用範囲

新たな知見 固有の周期性
異常値
特異なパターン

音声認識
構文解析
音質

画像認識
オブジェクト分類
トラッキング

特徴量

平均値
標準偏差
AR係数
周波数特性

インパルス応答
フォルマント
LPC
MFCC

HOG
SURF
FREAK
BRISK

データ収集

各種センサーデバイス
測定器、データ収集I/O
テキスト、バイナリ

マイク
WAVE
FLAC
MPEG-4

カメラ
JPEG
MPEG
WMV

情報

整理
•
•
•
•

フィルタリング
信号解析
データ処理
プロット

データ

観測

計測器

•
•
•
•

検出
収集
状態
データ取得

信号処理の適用範囲

12

何故、センサー、音、画像？

メディカル分野の例（音＋画像）… the AirSonea® technology detects
wheeze patterns in an image createdfrom recorded breath sounds...
センサー

iSonea

レーダー
レーザー
ソナー

➢
➢
➢

外界環境
➢
➢
➢

白線, 標識天候
障害物 / 歩行者

自動運転システム（センサー＋画像）
車両ダイナミクス

カメラ

複合的な要素技術の応用が必要とされる時代
 MATLABは多分野の専門ツールを提供

物体検出 / 認識➢➢

Computer vision
Machine learning

判断➢

距離 / 速度 推定➢➢

Phased Array
RF

行動の決定

制御

センサー
フュージョン

ブレーキ /
アクセル

ハンドル

13

Agenda 背景
 センサー信号の場合
 音声信号の場合
 画像の場合
 関連情報 まとめ

14センサー信号の取得
無償サポートパッケージ
Data Acquisition Toolbox™
▪ データ収集ボードへのアクセス▪ アナログI/O▪ デジタルI/O▪ サポートベンダー▪ National Instruments など

Instrument Control Toolbox™
▪ 測定器へのアクセス
オシロスコープ▪ 信号発生器▪ スペアナ：▪

▪ シリアル通信(RS-232など）
15例: 人の活動解析とクラシフィケーション

FeatureExtraction

Classification

Dataset courtesy of:Davide Anguita, Alessandro Ghio, Luca Oneto, Xavier Parra and Jorge L. Reyes-Ortiz.
Human Activity Recognition on Smartphones using a Multiclass Hardware-Friendly Support Vector Machine.
International Workshop of Ambient Assisted Living (IWAAL 2012). Vitoria-Gasteiz, Spain. Dec 2012
http://archive.ics.uci.edu/ml/datasets/Human+Activity+Recognition+Using+Smartphones
16クラシフィケーションのための信号解析
加速度センサー情報による人の活動認識

直立の
状態

寝転んだ
状態

歩いている
状態

階段を上っている状態

17

Demo: クラシフィケーションのための信号解析
MATLAB環境

ワークスペース（定義された変数）

実行セル

編集用エディタ

plotコマンドによるセンサー波形（垂直方向の加速度）
18

Demo: クラシフィケーションのための信号解析
可視化（視覚的情報からの予想）

動的状態
➢ 歩行?
➢ 階段を上る?
➢ 階段を下る？

直立？
表示を
カスタマイズ

座？

寝？

静的状態

19

Demo: クラシフィケーションのための信号解析
ベーシックな特徴量：ヒストグラムの平均値、RMS等

歩く vs. 寝る↓
平均値で識別歩く vs. 立つ↓
標準偏差やRMSで識別

歩く vs. 階段を上る↓
何で識別する？

20

Demo: クラシフィケーションのための信号解析
前処理：重力の影響を除去（フィルタ設計環境）>>fdatool21

Demo: クラシフィケーションのための信号解析
フィルタ設計GUI

プログラム内でフィルタを適用する方法①エクスポート：
フィルタ係数を出力②MATLABコードを生成：
フィルタ処理関数の自動生成

フィルタの仕様
フィルタの種類22

Demo: クラシフィケーションのための信号解析
フィルタ設計環境(cont’d)
前頁「②MATLABコードを生成」を選択した場合自動生成された自作フィルタ関数

フィルタの仕様e.g.)% フィルタオブジェクト定義
H = myfilter;% フィルタ処理
y = filter(H,data);23

Demo: クラシフィケーションのための信号解析
体の動きによる成分の取得
動作変化時の過渡応答

オリジナル（ハイパスフィルタ処理前）

重力によるオフセット除去

体の動きによる寄与分（ハイパスフィルタ処理後）24

Demo: クラシフィケーションのための信号解析
周波数ドメインでの解析
ベーシックなfft

自作スクリプト

ウェルチ法

25

Demo: クラシフィケーションのための信号解析
周波数ドメインでの解析: 各動作の比較、特徴量の獲得「階段を上る」の方が全体的に左より→一周期の時間が長い→動作が遅い「階段を上る」の方が高調波の減衰が大きい→ピュアな正弦波に近い→動作がスムーズ

他のデータでも同傾向であることを、
3Dプロットで確認

26

Demo: クラシフィケーションのための信号解析
周波数ドメインでの解析: 特徴量の自動取得
制約なしの場合>>findpeaks(p);ピークの数、ピーク間の距離、
プロミネンス等の制約を与えた場合>>findpeaks(p,'minpeakprominence',5,...
'minpeakdistance',0.15);
27

Demo: クラシフィケーションのための信号解析
自己相関関数:基本周波数の取得等

基本周波数の取得「歩く」の自己相関関数

階段を「上る」と「下る」の自己相関関数

ピーク値に4サンプルのずれ
28

Demo: クラシフィケーションのための信号解析
3軸加速度センサからの、特徴量一括取得
特徴量一括取得用関数（総数66）
特徴量2
（RMS）
前処理
（フィルタリング）

特徴量1
（平均値）

特徴量3
（自己相関）
特徴量4（ハーモニクス）
特徴量5（エネルギー）29

Demo: クラシフィケーションのための信号解析
classficationLeanerアプリによる機械学習
学習アルゴリズム

混合行列表示App内で試行した学習モデルのエクスポート

試行済み
学習履歴>>classificationLearner30

Demo: クラシフィケーションのための信号解析
ニューラルネットワークによるクラシフィケーション>>net = patternnet(n);>>net = train(net, X, targets);>>plotconfusion(targets, outputs)
混合行列

学習の様子
31

補足：センサーデータ前処理に使える新機能

32

ギャップフィリング
ARモデリングによる欠落したデータの復元: fillgaps▪

Resample()で大規模欠落
データの補間▪

複数のギャップのサポート▪

赤池情報量規準を元に自動でARモデルの次数を決定▪

非定常モデルをオプションで選択可能» load clockex» x(29:42) = NaN;
» fillgaps(x)33局所統計量の変化点の検出
信号の統計的性質の変化を検出: findchangepts▪

局所的な統計量の変化に応じた信号分割を実施– 統計量として平均、分散、平均と分散、線形回帰をサポート▪

cusum は累積和を使って局所平均値の変化点を検出» load('engineRPM.mat','x')» findchangepts(x,'Statistic','linear', …
» 'MinThreshold',var(x))» load('meanshift','x')» cusum(x)34ダイナミックタイムワーピング
信号の伸長、アライン、比較: dtw▪
▪
▪

2つの信号の軌跡の類似性を測定
タイムアラインメントをオプションで指定可能
時系列データの分類に使われる» dtw(x,y)

35

実信号間編集距離
実信号間の調整: edr▪
▪
▪

動的タイム ワーピング(DTW) の拡張機能
局所的な外れ値があった場合に正しく信号をDTW
整列
ミスアライン
追加サポート– 複素信号– Ｎ次元信号» x = chirp(0:999,0,1000,1/100);» y = cos(2*pi*5*(0:199)/200);» x(400:410) = 7;
» y(100:115) = 7;
» edr(x,y,.1)

EDR

適切なアライン36

類似性サーチでデータ内の信号を見つける
信号の最適にマッチする箇所を探す: findsignal▪

任意の次元の実数/複素数ベクト
ルの比較

▪
▪

複数のマッチング
有名なメトリックス(Euclidean, Squared
Euclidean, Symmetric
Kullback-Liebler)» data = exp(»»▪

((1:300)/100).^2).*cos(2*pi*(1:300)/100);
signal = sin(2*pi*(1:100)/100);
findsignal(data,signal)

様々な正規化(Zero-mean, Unit norm,
average power)

37

Agenda 背景
 センサー信号の場合
 音声信号の場合
 画像の場合
 関連情報 まとめ38

Celestion社事例

http://jp.mathworks.com/company/newsletters/articles/applying-machine-learning-techniques-to-classify-musicalinstrument-loudspeakers.html?s_tid=srchtitle
39

Celestion社事例抜粋： 特徴量の選択、得られた教訓等評価した特徴量 I developed and examined more than 50 features over several
months. These included energy features, such as the energy center;
spectral features, such as high-frequency roll-off and brightness;and psychoacoustic features, such as sharpness.

心理音響値の定量化

 To compute psychoacoustic values, I used MATLAB, Signal Processing
Toolbox, and DSP System Toolbox™…. I applied functions forsmoothing, convolution, peak finding, and computing Welch’s
power spectral density (PSD) estimates, among others.

三つの教訓 The first is that there is value in applying simpler algorithms such as
linear logistic regression, even if they end up performing poorly
compared to SVM, neural networks, or ensemble methods.
 The second lesson is that with a long and open-ended process likefeature extraction, having access to readily available signal
processing functions is virtually a necessity
 The third lesson is that no single machine learning technique works
well with all datasets.
40

Celestion社事例抜粋： 使用ツール
Using MATLAB®, Statistics and Machine
Learning Toolbox™, and Signal
Processing Toolbox™, I developed an
application that uses machine learning
To compute psychoacoustic values, I usedalgorithms to classify one loudspeaker
MATLAB, Signal Processing Toolbox, andmodel as either good or bad.
DSP System Toolbox™ to implement
established formulas for psychoacoustics
from a range of textbooks.

I used Optimization Toolbox™ to findan optimal value forlambda andminimize the error on both sets.
I used Audio System Toolbox™ tostream real-time audio directly from
MATLAB, created the main UI for playing
and selection…

41

Audio System Toolbox™
▪

▪

オーディオシステムの開発とシミュレーション

オーディオシステムの開発– オーディオエフェクト用ライブラリ– ラウドネス測定– ASIOドライバサポートによる低レイテンシリアルタイムシミュレーション– マルチチャネルI/O– MIDIコントローラによるパラメータチューニング

audioTestBenchによる可視化と官能評価

パラメータ調整と官能評価– audioTestBenchで波形やスペクトルの可視化とパラメータチューニング
– MATLABコードからVST Plugin生成・取り込み

VSTプラグイン生成とDAW上での評価
42

例: AudioSystem Toolboxによるデータ取得%% オブジェクトの定義
v1 = audioDeviceReader%% サウンド入力パラメータ設定
v1.SampleRate= Fs;
v1.NumChannels = NumC;
v1.Driver = 'ASIO';
v1.Device = 'device name';%% 可視化オブジェクト
hts = dsp.TimeScope('SampleRate', Fs);
hsa = dsp.SpectrumAnalyzer('SampleRate',Fs);%% スタート
for ii = 1:500
data = step(v1);step(hsa, data)step(hts, data)end% 取得データの変数への代入

43

メル周波数ケプストラム係数：
MFCC(Mel-Frequency Cepstral Coefficients)

ケプストラム（英: Cepstrum）：（信号の）ケプストラムとは、（信号の）フーリエ変換の対数（位相アンラッピングを施したもの）をフーリエ変換したもの。スペクトルのスペクトルとも呼ばれる。
応用：… 人間の声や音楽の信号を表す特徴ベクトルとしても使われている。
この場合、スペクトルをまずメル尺度を使って変換する。その結果はメ
ル周波数ケプストラム係数（MFCC）と呼ばれる。これは話者認識
やピッチ抽出アルゴリズムなどに応用されている。…

メル尺度（実周波数と知覚音高の関係）（出展：フリー百科事典「ウィキペディア(Wikipedia)」, “ケプストラム”より抜粋）

44ソース・フィルターモデル
ソース：声帯励振（音源）
合成音

𝑌 𝑧 =𝐻 𝑧 𝑋 𝑧フィルタ：声道（MFCC）

肺からの空気が声帯を震わせその信号が声道を通り発音する

𝑌 𝑧

𝐻 𝑧
𝑋 𝑧
45

demo：MFCCとPitchを用いた話者の推定

46

音声”Two”の、有声音部と無声音部

子音/T/(無声音部)

母音/UW/(有声音部)

ゼロクロスの頻度大：無声音
ゼロクロスの頻度小：有声音
47

有声音部のピッチ推定: 相互相関法% Center-clip signalx(abs(x)<clipLevel) = 0;

1.音声信号にセンタークリップ処理
声道の応答に起因する微小な振動成分を除去
2.信号分解能確保のためアップサンプリング
信号を補間する
3.自己相関を計算
サンプルのラグを判定% Upsample the signal for
x = step(FIRInterpolatorFour,x(:,1));

4.ラグ値からピッチ周期を換算% Apply pitch detection to the word “two”
pD = audiopluginexample.SpeechPitchDetector;[~,pitch] = process(pD,audioIn);>> edit audiopluginexample.SpeechPitchDetector% Perform autocorrelationy = xcorr(x);48

データセット（Census Database: AN4）
http://www.speech.cs.cmu.edu/databases/an4/49

MFCC演算ステップ
メルフィルタバンク

50

音声全データ
fash

テスト用
学習用

fbbh
女性

fclc
fejs
ffmm
mblb
mblw

男性 mbmg
mcel

mcen
51

データセットから取得した特徴量（抜粋）Filename

Pitch

ラベル

MFCC1～13

ピッチ
MFCC1

MFCC2

MFCC3

MFCC4

MFCC5

MFCC6

MFCC7

MFCC8

MFCC9

MFCC10 MFCC11 MFCC12 MFCC13

Label

1 'an251-fash-b.raw'

210.8731 -4.0064-1.05022.21301.59440.16360.5905-0.1960

0.4518-0.7614-0.0169-0.0412-0.25250.0019 'fash'

2 'an251-fash-b.raw'

214.4054 -3.5244-1.7537

1.7664

1.87630.11160.3083-0.2133

0.3668-0.5598

0.2236-0.1668-0.3488-0.1594 'fash'

3 'an251-fash-b.raw'

213.3333 -2.9804-1.51771.45081.51710.04590.2228-0.15710.3317-0.5655

0.0566-0.2429-0.4489-0.2143 'fash'

4 'an251-fash-b.raw'

213.3333 -2.4698-1.4726

1.3716

1.52390.05300.1686-0.2833

0.1835-0.6909-0.0505-0.2472-0.3179-0.1572 'fash'

5 'an251-fash-b.raw'

208.8091 -2.3291-1.1027

1.31911.6805-0.35290.1277-0.5211-0.0254-0.6645-0.0556-0.3002-0.2029-0.2032 'fash'

6 'an251-fash-b.raw'

212.9784 -2.2898-0.8004

0.91951.6846-0.7346

0.1145-0.62510.1587-0.5014-0.0298-0.3634-0.2166-0.2122 'fash'

7 'an251-fash-b.raw'

207.4554 -2.3984-0.72930.52581.5897-0.7562-0.0039-0.69400.2292-0.5705

0.0363-0.3367-0.0192-0.0899 'fash'

8 'an251-fash-b.raw'

202.5316 -2.5082-0.7673

0.28551.3773-0.6114-0.3232-0.55420.3293-0.6243

0.1565-0.16160.2889-0.0642 'fash'

ピッチとMFCCの値のスケールが異なるため正規化が必要52

データセットから取得した特徴量（抜粋: 正規化後）
ピッチ
Filename

ラベル

MFCC1～13

Pitch

MFCC1

MFCC2

MFCC3

MFCC4

MFCC5

MFCC6

MFCC7

MFCC8

MFCC9

MFCC10 MFCC11 MFCC12 MFCC13

Label

1 'an251-fash-b.raw'

0.9923-1.7968-1.7525

2.42870.38830.90170.8702-0.2387

1.1208-2.1363-0.0567-0.2265-0.6421-0.0378 'fash'

2 'an251-fash-b.raw'

1.0660-1.5036-2.3170

1.95900.65190.82340.2628-0.2732

0.9087-1.5951

0.5719-0.6201-0.9597-0.6153 'fash'

3 'an251-fash-b.raw'

1.0437-1.1728-2.1276

1.62690.31610.72460.0788-0.16110.8213-1.61050.1352-0.8587-1.2898-0.8116 'fash'

4 'an251-fash-b.raw'

1.0437-0.8623-2.0915

1.5436

0.32240.7352-0.0378-0.41270.4515-1.9471-0.1446-0.8719-0.8579-0.6075 'fash'

5 'an251-fash-b.raw'

0.9492-0.7767-1.79461.4884

0.46880.1239-0.1259-0.8867-0.0698-1.8761-0.1581-1.0382-0.4783-0.7719 'fash'

6 'an251-fash-b.raw'

1.0363-0.7528-1.55211.06800.4727-0.4508-0.1544-1.09410.3896-1.4382-0.0906-1.2364-0.5235-0.8043 'fash'

7 'an251-fash-b.raw'

0.9210-0.8188-1.4950

0.65390.3840-0.4834-0.4092-1.2315

0.5656-1.6238

0.0824-1.15250.1274-0.3666 'fash'

8 'an251-fash-b.raw'

0.8182-0.8856-1.52550.40120.1854-0.2653-1.0963-0.9527

0.8153-1.7682

0.3964-0.6039

1.1439-0.2743 'fash'>> m = mean(features{:,2:15});>> s = std(features{:,2:15});>> features{:,2:15} = (features{:,2:15}-m)./s;
53

KNN 分類器（K最近傍分類器）
原理：近くにある K個の点を探して、多数決を採る

「Bさん」

×

「Aさん」

分類したい新しい特徴量

54KNN分類器（K最近傍分類器）
分類器： d = fitcknn(X, T, ‘NumNeighbors’, 4)

特徴量

ラベル

予測： C = predict(d, Y)

予測されるラベル

新しい特徴量55

混合行列% Train a classifierclassificationKNN = fitcknn(...
% Perform cross-validation
partitionedModel = crossval(...
% Compute validation accuracy
validationAccuracy = 1 - kfoldLoss(...

話者（正解）

HelperTrainKNNClassifier.m（主要関数抜粋）% Compute confusion matrix
confMatrix = confusionmat(...

話者（予測結果）

56

予測結果

合致
Filename

ActualSpeaker

Predicted
Speaker

Confidence
Percentage

HelperTestKNNClassifier.m（主要部抜粋）1'cen1-fash-b.raw'

fash

fash

73.18442'cen1-fbbh-b.raw'

fbbh

fbbh

56.0000

for idx = 1:length(Filenames)
T = features_test(Fidx==idx,2:end); % Rows that correspond to one file3

'cen1-fclc-b.raw'

fclc

fclc

77.06094'cen1-fejs-b.raw'

fejs

fejs

87.0536

predictedLabels = string(predict(trainedClassifier,T(:,1:15))); % Predict
totalVals = size(predictedLabels,1);

5'cen1-ffmm-b.raw'

ffmm

ffmm

51.34106'cen1-mblb-b.raw'

mblb

mblb

76.4368[predictedLabel, freq] = mode(predictedLabels); % Find most frequently predicted label
match = freq/totalVals*100;

7'cen1-mblw-b.raw'

mblw

mblw

36.71888'cen1-mbmg-b.raw' mbmg

mbmg

77.6596

result_file.Filename = Filenames(idx);
result_file.ActualSpeaker = T{1,15}{:};
result_file.PredictedSpeaker = char(predictedLabel);result_file.ConfidencePercentage = match;9'cen1-mcel-b.raw'

mcel

mcel

63.524610'cen1-mcen-b.raw'

mcen

mcen

82.967011'cen2-fash-b.raw'

fash

fash

35.5140

12'cen2-fbbh-b.raw'

fbbh

fbbh

64.7059

13'cen2-fclc-b.raw'

fclc

fclc

78.3784

14'cen2-fejs-b.raw'

fejs

fejs

86.9565

15'cen2-ffmm-b.raw'

ffmm

ffmm

48.011416'cen2-mblb-b.raw'

mblb

mblb

78.082217'cen2-mblw-b.raw'

mblw

mblw

44.2478

18'cen2-mbmg-b.raw' mbmg

mbmg

78.0645

19'cen2-mcel-b.raw'

mcel

mcel

50.349720'cen2-mcen-b.raw'

mcen

mcen

54.3624

result = [result; struct2table(result_file)]; %#ok
end

■分割された音声データについて10クラスに予測■各予測結果の頻度スコアのトップを予測結果とする■頻度のパーセント換算値を”信頼パーセンテージ”とする■全テストファイルが実際の話者と合致

57

Agenda 背景
 センサー信号の場合
 音声信号の場合
 画像の場合
 関連情報 まとめ58

イメージの取得
Image Acquisition Toolbox™
▪

業界標準のHWからの動画像取込み機能を提供– フレームグラバ (画像入力ボード)▪

Analog 入力▪

Camera Link® 入力– DCAM IIDC ® 1394 (FireWire)– GigE Vision ®互換カメラ– 一般的なUSB Webカメラ, IPカメラ(MATLABサポートパッケージ)– USB3 Vision互換カメラ(GenICam GenTL API経由)– Microsoft ® Kinect ® for Windows ®, Kinect ® for Windows ® v2▪
▪

カスタムアダプター開発キット付属
Simulinkカメラブロックセット
認識した人

深度画像

Kinect v2では手の開閉状態も検出

RGB画像＋骨格座標59

Demo: HSV vs HOG
瞬きの検出に適した特徴量は？

■HSV特徴量（明度のヒストグラム）>>hsv = rgb2hsv(img);>>f2 = imhist(hsv(:,:,3))';■HOG特徴量>>f1 = extractHOGFeatures(img);

60

画像処理で用いられる特徴量

画像における局所特徴：
 エッジ、コーナー、ブロブ等の特徴的なパターン
 材質や色等が周辺と異なる場合に顕著に現れる



良い特徴とは？
 検出が繰り返し可能
 特徴のパターンに特殊性がある 位置推定が可能

Computer Vision System Toolbox
以下、ドキュメント抜粋
http://jp.mathworks.com/help/vision/ug/local-feature-detection-and-extraction.html
61画像処理で用いられる特徴量（cont’d）
検出器と記述子の選択のための考慮事項
基準
特徴のタイプ目的:
 キー ポイントのマッチング
 分類
歪み（スケールの変化、回転）のタイプ

パフォーマンス要件:
 リアルタイム パフォーマンスが必要 正確度と処理速度

推奨
データに適した検出器を使用
 例1) バクテリア細胞のイメージ: ブロブ検出器
 例2) 都市部の航空写真: コーナー検出器
 HOG 、 SURF：分類用途
 BRISK 、FREAK：イメージ間の対応、レジストレーション
 スケールの変化が生じない場合：コーナー検出器
 スケールや回転などの大きな歪みがデータに含まれる場合： SURF 特徴検出器
および記述子（ただし計算量は多い）
バイナリ記述子は、勾配ベースの記述子と比べて処理は高速。正確度の面で劣る。
正確度を高めるには、同時に複数の検出器と記述子を使用

62

画像処理で用いられる特徴量（cont’d）
特徴のタイプに基づく検出関数の選択
検出器

特徴のタイプ
コーナー
コーナー
コーナー
ブロブ
コーナー
強度が一様な領域

FAST [1]最小固有値アルゴリズム [4]Harris [3]SURF [11]BRISK [6]MSER [8]関数スケール依存×××○○○

detectFASTFeaturesdetectMinEigenFeaturesdetectHarrisFeaturesdetectSURFFeaturesdetectBRISKFeaturesdetectMSERFeatures

記述子メソッドの選択
記述子

バイナリ

関数とメソッド

HOG
LBP
SURF

×
×
×

FREAK

不変要素

一般的な使用例

extractHOGFeatures(I, ...)
extractLBPFeatures(I, ...)
extractFeatures(I, 'Method','SURF')

スケール
×
×
○

回転
×
○
○

画像間の対応
×
なし
○

分類
○
○
○

○

extractFeatures(I, 'Method','FREAK')

○

○

○

×

BRISK

○

extractFeatures(I, 'Method','BRISK')

○

○

○

×

 ブロック キー ポイントの周りのピクセル近傍×extractFeatures(I, 'Method','Block')

×

×

○

○

[1] Rosten, E., and T. Drummond. “Machine Learning for High-Speed Corner Detection.” 9th European Conference on Computer Vision. Vol. 1, 2006, pp. 430–443.
[3] Harris, C., and M. J. Stephens. “A Combined Corner and Edge Detector.” Proceedings of the 4th Alvey Vision Conference. August 1988, pp. 147–152.
[4] Shi, J., and C. Tomasi. “Good Features to Track.” Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition. June 1994, pp. 593–600.
[6] Leutenegger, S., M. Chli, and R. Siegwart. “BRISK: Binary Robust Invariant Scalable Keypoints.” Proceedings of the IEEE International Conference. ICCV, 2011.
[8] Matas, J., O. Chum, M. Urba, and T. Pajdla. "Robust wide-baseline stereo from maximally stable extremal regions."Proceedings of British Machine Vision Conference. 2002, pp. 384–396.
[11] Bay, H., A. Ess, T. Tuytelaars, and L. Van Gool. “SURF: Speeded Up Robust Features.” Computer Vision and Image Understanding (CVIU). Vol. 110, No. 3, 2008, pp. 346–359.

63

画像処理で用いられる特徴量：例題（分類） 分類器の作成：trainCascadeDetector：
 Haar, LBP, HOGに対応
 HOG:
 オブジェクトの全体的な形状のキャプチャ向き 人、自動車の検出等
 Haar, LBP:
 スケールの細かいテクスチャ向き
 材質の識別、顔の検出等

レンガ

レンガ（回転）LBP特徴量によるテクスチャの判別
HOG特徴量カーペット

 黄：レンガ vsレンガ（回転） 青：レンガ vs カーペット

64

特徴点ベースのワークフロー
 特徴点のマッチングを用いた、物体認識
対象物体

SURF/Harris/FAST/BRISK等の特徴点検出SURF/FREAK/BRISK/HOG/LBP等の特徴量抽出RANSACを用いた対応点推定

物体の検出 特徴点ベースのレジストレーション（位置合せ）65

特徴量抽出の例題：（アドオンの追加）①”アドオンの入手”をクリック③”追加”をクリック②キーワードを入力

https://jp.mathworks.com/matlabcentral/fileexchange/55107-brain-mri-tumor-detection-and-classification

66

Brain MRI GUI
腫瘍部データの離散二次元Waveletから特徴量を抽出し、腫瘍の良性・悪性の判定に応用した例（ Wavelet Toolboxを使用）

https://jp.mathworks.com/matlabcentral/fileexchange/55107-brain-mri-tumor-detection-and-classification

67

学習用画像の拡張

左右反転

任意の
反転

68

学習用画像の拡張
物体の色味の変更

物体の外形の変形

画像操作のさまざまな関数で、容易に画像の拡張を行っていただけます69補足：機械学習アルゴリズムの実装
ユーザーインターフェース

ユーザーインターフェース
設計環境 (GUIDE)>> guide.fig ファイルが生成される用いた画像のソース： http://www.cdc.gov/dpdx/index.html(Centers for Disease Control and Prevention)70

補足：機械学習アルゴリズムの実装

MATLAB Compiler>> deploytool
もしくは アプリケーションタブから起動

MATLAB
Toolboxes

１：最上位のMATLABファイルを指定（スクリプトもしくは関数)

２：ビルド
ボタンを押すMATLAB Compiler
MATLAB Compiler SDK

.exe

.dll

ほぼすべての関数をサポート ・ GUIやFigureも対応

ロイヤリティーフリー

制限：http://www.mathworks.co.jp/help/compiler/limitations-about-what-may-be-compiled.html
対応Toolbox関数：http://www.mathworks.co.jp/products/compiler/supported/compiler_support.html

71

補足：機械学習アルゴリズムの実装
MATLAB
myapp.m

MATLAB Compiler

Compilemyapp.exe

Copy & Installmyapp.exe

myapp.exe

myapp.exe72

補足：機械学習アルゴリズムの実装

MATLAB
Production Server>> productionServerCompiler
もしくは アプリケーションタブから起動

MATLAB
Toolboxes

１：最上位のMATLABファイルを指定（関数)

２：パッケージ化
ボタンを押すMATLAB Compiler
MATLAB Compiler SDK.ctf

ほぼすべての関数をサポート

ロイヤリティーフリー

制限：http://www.mathworks.co.jp/help/compiler/limitations-about-what-may-be-compiled.html
対応Toolbox関数：http://www.mathworks.co.jp/products/compiler/supported/compiler_support.html

73

補足：機械学習アルゴリズムの実装

MATLAB
myfun.m

MATLAB
Production Server

Compile !

myfun_deployed.ctf

Copy !
・MATLAB関数をサーバ上で実行
・HTTP / HTTPS によるアクセス・サーバ上での一元管理
・メンテナンス性の大幅な向上

・ブラウザや軽量なクライアントアプリ（お客様開発）を使ったアクセル
・Java, C#, VB, C/C++, Python 等が利用可能

クライアントPC 群
MATLAB
Production Server

Access !

ClientApplication

ClientApplication

TCP / IP

ClientApplication/auto_deploy/myfun_deployed.ctf
74Agenda 背景
 センサー信号の場合
 音声信号の場合
 画像の場合
 関連情報 まとめ

75

エッジデバイス実装: Simulinkハードウエア連携機能(cont’d)

Arduino® Nano, Uno, Mega 2560

Raspberry PiTM
Raspberry PiTM2Apple® iOS Devices

LEGO® MINDSTORMS® EV3

KINECTTM for Windows®
LEGO MINDSTORMS NXT

BeagleBoard /PandaBoard

Digilent® Analog DiscoveryTM

RTL-SDR Radio

Samsung® GALAXY AndroidTM

TurtleBot™

KINECTTM for Windows® V2

サポートされているハードウェア一覧
http://jp.mathworks.com/hardware-support76

簡易IoTシステム概観

!

クラウド

遠隔監視/
警告

監視対象

IoTエッジ

ホストPC

77

Raspberry Piの概観、インターフェイス
音声出力

ホストPC
(Simulink)

モニター（HDMI出力）

シリアルバス
 I2C
 SPI

webカメラ本体

音声入力(USBマイク)

GPIO スイッチ
 LED
 センサ
78

Raspberry Piの歴史
年
RasPi トピック
2008  エベン・アプトン氏Raspberry Pi財団設立
2011  試作品完成
2012  初代Raspberry Pi リリース 初日二時間で1万台を売りきる
 初日だけで10万台の発注
2013  累計200万台
2014  累計300万台
 Compute moduleリリース
2015  Raspberry Pi 2, Zeroリリース
 累計500万台(Feb)2016  Raspberry Pi 3 リリース 愛知県稲沢市で製造開始
 累計1,100万台(Nov)
2017  Compute module 3, Zero Wリリース
 累計1,500万台(Jul)2018  累計1,900万台(Mar)

Simulink Support Package

 リリース当初は、コンピュータ
サイエンスの授業利用を想定

 リリース（R2013a）
 ThingSpeak Write (14a)






Pi 2対応 (15a)
MAC OS対応 (15a)Pi 3対応 (16a)Linux OS対応 (16b) 学生をモチベートする基本方針➢ プログラミング可能➢ 楽しい➢ 小型で丈夫➢ 安価 現状教育用と産業用は半々(産業用途が増加傾向）

 I2C,SPI,Serial（17a) Audio read, logging (17b) Pi Zero W対応 (18a)https://en.wikipedia.org/wiki/Raspberry_Pi

79

Raspberry Piによるトラフィックモニタ、分析

80

Raspberry Piによるトラフィックモニタ、分析（cont’d）
PCとWebカメラで動
作チェック

台数計測
雑音除去

フォアグランド
ディテクター
Blob解析

必要なデータをThingSpeakでアップロードRaspberry Pi 用
カメラペリフェラル
ブロック
81

Raspberry Piによるトラフィックモニタ、分析（cont’d）
トラフィック量のヒストグラム一週間の交通量（24時間で見た推移）

一週間の交通量（日ごとの推移）
桜通の1週間の交通量（西行きvs東行き）

日曜日の交通量（交通量の多い時間帯は？）82

データ収集&解析を手軽に！ ThingSpeak●無料ではじめられます！！
●準備は簡単３ステップ

データ収集＆解析Webサービス-アカウント/チャネル作成-関数群のインストール-専用コマンド/ブロックで簡単利用

thingSpeakWrite

収集

コマンド 1行
専用ブロック

可視化＆解析

-WEB上でデータの確認が可能-公開範囲Private/Public選択-MATLAB関数で解析＆可視化

アクション-ｃｓｖ形式でデータ保存-条件を満たしたらツイート等
様々なアクション機能●低価格HWとの相性が◎IoTプラットフォームThingSpeak利用のメリット●インターネットでどこからでもデータ確認&解析●Arduino/Raspberry Piとの活用で、低価格な実験系を構築●プログラムの完了/データ異常をスマホでチェック
83

ThingSpeakはいろいろな場面で活用できます
実験での活用測定器の自動制御にthingSpeakWrite
コマンド１行追加で無信号状態をすぐに把握。「数時間後気づいたらテストが進んでいなかった」を避け、離れた実験室との無駄な行き来もなくなります。

気象・農水産分野での活用
人やPCでやっていた温度/湿度等のデータ収集や植物の成長などの単純かつ長期間の計測。Arduino/Raspberry pi + ThingSpeakにより、数千円で自動観測システムを構築できます。

植物成長記録をWEB上で観測

情報・統計解析分野での活用
カメラ画像から車を検知する
コンピュータビジョンのアルゴリズムをRaspberry Piに実装し、道路の交通量を観測。ThingSpeakでは取得した統計データを「時間ごと」「曜日ごと」などのソートやMATLAB関数を用いて解析することができます。

医療・スポーツ科学分野での活用
生体信号モニタリングで必要な長時間のデータ収集もThingSpeak上で行うことができます。その他、加速度、温度、位置情報のデータ等関連パラメータの同時収集/解析に利用可能。異常値に対してアクション機能も活用できます。

84

例：温度センサーのモニタリング、アラートの配信～Raspberry Pi とセンサデバイス間の配線～
RasPiと温度センサ間の配線図
SDA
V+(3.3V)

GND

SCL
V+

V+

V+

ADD0

SDA

TMP102 ×2

RasPi

TMP102ADD0

TMP102

GND
RasPi

SCL
tmp102アドレスピン設定（同時に4台まで接続可能）
アドレス

ADD0ピン

0x48

GND

0x49

V+

0x4a

SDA

0x4b

SCL

85

例：温度センサーのモニタリング、アラートの配信～Simulinkモデル概観～
異常判定ロジック➢ 閾値以下:０➢ 閾値異常:1

I2C読み込み
専用ブロック

可視化およびクラウド送信部

86

例：温度センサーのモニタリング、アラートの配信～専用ブロックのパラメータ設定～

ThingSpeak Readブロック

I2C Master Read
ブロック

スレーブの
アドレス

書き込み用API
送信するデータの数87例：温度センサーのモニタリング、アラートの配信～Tweetの条件設定～

監視対象

条件

アクション
（Tweetを送信）

88例：温度センサーのモニタリング、アラートの配信～実行結果～

tmp102A

tmp102B
27°を超えた時
アラート
エクスターナルモードでの実行結果

クラウド（Thingspeak）での可視化とアラートの配信

89Agenda 背景
 センサー信号の場合
 音声信号の場合
 画像の場合
 関連情報 まとめ

90

使用ツール
種類

センサー

音

画像

配布、
データ・
アグリゲートMATLAB Compiler
MATLAB Production ServerThingspeak

機械学習

Statistics and Machine Learning Toolbox
Neural Network Toolbox

前処理、
特徴量抽出
データ収集

Signal Processing Toolbox
DSP System Toolbox
Data Acquisition Toolbox
Instrument Control Toolbox

Audio System
ToolboxImage Processing Toolbox
Computer Vision System Toolbox
Image Acquisition Toolbox

91

まとめ
 背景

 事例の紹介 センサー信号の場合

 信号処理関数で特徴量の試行錯誤を効率化
 音声信号の場合

 ユーザ事例、Audio System Toolboxのご紹介
 画像の場合

 学習用データの拡張

機械学習用途には信号処理ツールをお忘れなく 関連情報 エッジデバイス実装、Thingspeakのご紹介

92