メインコンテンツ

findchangepts

信号の急激な変化の検出

説明

ipt = findchangepts(x) は、x の平均値が最も大幅に変化するインデックスを返します。

  • xN 個の要素を持つベクトルである場合、findchangepts は、x を 2 つの領域 x(1:ipt-1)x(ipt:N) に分割します。これらは、各領域の局所的な平均値からの残差 (二乗) 誤差の和を最小化します。

  • xMN 列の行列である場合、findchangepts は、x を 2 つの領域 x(1:M,1:ipt-1)x(1:M,ipt:N) に分割し、各領域の局所的な M 次元の平均値からの残差誤差の和を最小化する列インデックスを返します。

ipt = findchangepts(x,Name,Value) は、名前と値の引数を使用して追加オプションを指定します。オプションには、レポートする変化点の数や、平均値の代わりに測定する統計量があります。詳細については、変化点の検出を参照してください。

[ipt,residual] = findchangepts(___) は、前述のいずれかの指定を組み込んでモデル化された変化に対する信号の残差誤差も返します。

出力引数なしの findchangepts(___) は、信号および検出された任意の変化点をプロットします。詳細については、Statisticを参照してください。

メモ

関数 findchangepts は、プロットの前に現在の Figure をクリアします (clf)。信号および検出された変化点をサブプロットにプロットするには、プロット関数を使用します。オーディオ ファイルの分割を参照してください。

すべて折りたたむ

8192 Hz でサンプリングされた列車の警笛の録音を含むファイルを読み込みます。信号の平方根平均二乗レベルが最も大きく変化する点を 10 個検出します。

load train

findchangepts(y,MaxNumChanges=10,Statistic="rms")

Figure contains an axes object. The axes object with title Number of changepoints = 10 Total log weighted dispersion = -34290.1456 contains 2 objects of type line.

信号の短時間パワー スペクトル密度を計算します。信号を 128 サンプルのセグメントに分割し、各セグメントにハミング ウィンドウを適用します。隣り合ったセグメント間のオーバーラップを 120 サンプル、DFT 点を 128 に指定します。パワー スペクトル密度が最も大きく変化する点を 10 個検出します。

[s,f,t,pxx] = spectrogram(y,128,120,128,Fs);

findchangepts(pow2db(pxx),MaxNumChanges=10)

Figure contains 2 axes objects. Axes object 1 contains 131 objects of type line. Axes object 2 with title Number of changepoints = 10 Total residual error = 2820745.8453 contains 11 objects of type image, line.

再現可能な結果が必要な場合は、乱数発生器をリセットします。次の条件でランダムな信号を生成します。

  • 平均値は 7 つの各領域内では一定だが、領域から領域にかけて急激に変化する。

  • 分散は 5 つの各領域内では一定だが、領域から領域にかけて急激に変化する。

rng("default")

lr = 20;

mns = [0 1 4 -5 2 0 1];
nm = length(mns);

vrs = [1 4 6 1 3];
nv = length(vrs);

v = randn(1,lr*nm*nv)/2;

f = reshape(repmat(mns,lr*nv,1),1,lr*nm*nv);
y = reshape(repmat(vrs,lr*nm,1),1,lr*nm*nv);

t = v.*y+f;

作成の手順を強調しながら信号をプロットします。

subplot(2,2,1)
plot(v)
title("Original")
xlim([0 700])

subplot(2,2,2)
plot([f;v+f]')
title("Means")
xlim([0 700])

subplot(2,2,3)
plot([y;v.*y]')
title("Variances")
xlim([0 700])

subplot(2,2,4)
plot(t)
title("Final")
xlim([0 700])

Figure contains 4 axes objects. Axes object 1 with title Original contains an object of type line. Axes object 2 with title Means contains 2 objects of type line. Axes object 3 with title Variances contains 2 objects of type line. Axes object 4 with title Final contains an object of type line.

信号の平均値が最も大きく変化する点を 5 つ検出します。

figure
findchangepts(t,MaxNumChanges=5)

Figure contains an axes object. The axes object with title Number of changepoints = 5 Total residual error = 1989.3814 contains 3 objects of type line.

信号の平方根平均二乗レベルが最も大きく変化する点を 5 つ検出します。

findchangepts(t,MaxNumChanges=5,Statistic="rms")

Figure contains an axes object. The axes object with title Number of changepoints = 5 Total log weighted dispersion = 871.1003 contains 2 objects of type line.

信号の平均値と標準偏差が最も大きく変化する点を検出します。

findchangepts(t,Statistic="std")

Figure contains an axes object. The axes object with title Number of changepoints = 1 Total log weighted dispersion = 1263.4625 contains 3 objects of type line.

Fs=7418Hz でサンプリングされた音声信号を読み込みます。ファイルには、"MATLAB®" という単語を発声している女性の録音音声が含まれています。

load mtlb

信号の分散が大きく変化している点を検出することによって、単語内の母音と子音を識別します。変化点の数を 5 つに制限します。

numc = 5;

[q,r] = findchangepts(mtlb,Statistic="rms",MaxNumChanges=numc);

変化点インデックスに基づいて音声信号の信号マスクを作成します。信号マスクの使用の詳細については、signalMaskを参照してください。

t = (0:length(mtlb)-1)/Fs;
roitable = ([[1;q] [q;length(mtlb)]]);

x = ["M" "A" "T" "L" "A" "B"]';
c = categorical(x,unique(x,"stable"));

msk = signalMask(table(t(roitable),c),SampleRate=Fs,RightShortening=1);
roimask(msk)
ans=6×2 table
           Var1            c
    ___________________    _

          0    0.017525    M
    0.01766     0.10461    A
    0.10475     0.22162    T
    0.22176     0.33675    L
    0.33688     0.46535    A
    0.46549     0.53909    B

音声信号および検出された変化点を関心領域とともに信号マスクからサブプロットにプロットします。

  • 上のサブプロットで、関数 plotsigroi を使用して信号マスク領域を可視化します。カラー バーが上に表示されるように設定を調整します。

  • 下のサブプロットで、元の音声信号をプロットし、検出された変化点を垂直線として追加します。

subplot(2,1,1)
plotsigroi(msk,mtlb)
colorbar("off")
nc = numel(c)-1;
colormap(gca,lines(nc));
colorbar(TickLabels=categories(c),Ticks=1/2/nc:1/nc:1, ...
    TickLength=0,Location="northoutside")
xlabel("")

subplot(2,1,2)
plot(t,mtlb)
hold on
xline(q/Fs)
hold off
xlim([0 t(end)])
xlabel("Seconds")

Figure contains 2 axes objects. Axes object 1 contains 6 objects of type line. Axes object 2 with xlabel Seconds contains 6 objects of type line, constantline.

各セグメントの後に一時停止しながら音声を再生するには、次の行のコメントを解除します。

% for k = 1:length(roitable)
%     intv = roitable(k,1):roitable(k,2);
%     soundsc(mtlb(intv).*hann(length(intv)),Fs)
%     pause(.5)
% end

可変振幅をもつ 2 つの正弦波と線形トレンドで構成される信号を作成します。

vc = sin(2*pi*(0:201)/17).*sin(2*pi*(0:201)/19).* ...
    [sqrt(0:0.01:1) (1:-0.01:0).^2]+(0:201)/401;

信号の平均値が最も大きく変化している点を検出します。この場合、'Statistic' の名前と値の引数はオプションです。残差誤差の最小改善値に 1 を指定します。

findchangepts(vc,'Statistic','mean','MinThreshold',1)

Figure contains an axes object. The axes object with title Number of changepoints = 2 Total residual error = 9.3939 contains 3 objects of type line.

信号の平方根平均二乗レベルが最も大きく変化する点を検出します。残差誤差の最小改善値に 6 を指定します。

findchangepts(vc,'Statistic','rms','MinThreshold',6)

Figure contains an axes object. The axes object with title Number of changepoints = 4 Total log weighted dispersion = -436.5368 contains 2 objects of type line.

信号の標準偏差が最も大きく変化する点を検出します。残差誤差の最小改善値に 10 を指定します。

findchangepts(vc,'Statistic','std','MinThreshold',10)

Figure contains an axes object. The axes object with title Number of changepoints = 26 Total log weighted dispersion = -1110.8065 contains 3 objects of type line.

信号の平均値と傾きが最も急激に変化している点を検出します。残差誤差の最小改善値に 0.6 を指定します。

findchangepts(vc,'Statistic','linear','MinThreshold',0.6)

Figure contains an axes object. The axes object with title Number of changepoints = 3 Total residual error = 7.9824 contains 3 objects of type line.

20 のランダムな制御点を使用して、1000 のサンプルの 2 次元ベジエ曲線を生成します。ベジエ曲線は次のように定義されます。

C(t)=k=0m(mk)tk(1-t)m-kPk,

ここで、Pkm 個の制御点の k 番目、t の範囲は 0 ~ 1、(mk) は二項係数です。曲線と制御点をプロットします。

m = 20;
P = randn(m,2);
t = linspace(0,1,1000)';

pol = t.^(0:m-1).*(1-t).^(m-1:-1:0);
bin = gamma(m)./gamma(1:m)./gamma(m:-1:1);
crv = bin.*pol*P;

plot(crv(:,1),crv(:,2),P(:,1),P(:,2),"o:")

Figure contains an axes object. The axes object contains 2 objects of type line.

各セグメントの点とセグメントの平均値からの距離が最小になるように、曲線を 3 つのセグメントに分割します。

findchangepts(crv',MaxNumChanges=3)

Figure contains 2 axes objects. Axes object 1 contains 5 objects of type line. Axes object 2 with title Number of changepoints = 2 Total residual error = 158.2579 contains 3 objects of type line.

曲線を、直線に当てはめた 20 のセグメントに分割します。

findchangepts(crv',Statistic="linear",MaxNumChanges=19)

Figure contains 2 axes objects. Axes object 1 contains 5 objects of type line. Axes object 2 with title Number of changepoints = 19 Total residual error = 0.090782 contains 21 objects of type line.

20 のランダムな制御点を使用して、3 次元ベジエ曲線を生成してプロットします。

P = rand(m,3);
crv = bin.*pol*P;

plot3(crv(:,1),crv(:,2),crv(:,3),P(:,1),P(:,2),P(:,3),"o:")
xlabel("x")
ylabel("y")

Figure contains an axes object. The axes object with xlabel x, ylabel y contains 2 objects of type line.

曲線を上から可視化します。

view([0 0 1])

Figure contains an axes object. The axes object with xlabel x, ylabel y contains 2 objects of type line.

各セグメントの点とセグメントの平均値からの距離が最小になるように、曲線を 3 つのセグメントに分割します。

findchangepts(crv',MaxNumChanges=3)

Figure contains 2 axes objects. Axes object 1 contains 7 objects of type line. Axes object 2 with title Number of changepoints = 2 Total residual error = 7.2855 contains 3 objects of type line.

曲線を、直線に当てはめた 20 のセグメントに分割します。

findchangepts(crv',Statistic="linear",MaxNumChanges=19)

Figure contains 2 axes objects. Axes object 1 contains 7 objects of type line. Axes object 2 with title Number of changepoints = 19 Total residual error = 0.0075362 contains 21 objects of type line.

入力引数

すべて折りたたむ

入力信号。実数値のベクトルまたは行列で指定します。

  • xN 個の要素を持つベクトルである場合、findchangepts は、x を 2 つの領域 x(1:ipt-1)x(ipt:N) に分割します。これらは、Statistic で指定された統計量の局所値からの各領域の残差 (二乗) 誤差の和を最小化します。

  • xMN 列の行列である場合、findchangepts は、x を 2 つの領域 x(1:M,1:ipt-1)x(1:M,ipt:N) に分割し、Statistic で指定された統計量の局所的な M 次元の平均値から各領域の残差誤差の和を最小化する列インデックスを返します。

例: reshape(randn(100,3)+[-3 0 3],1,300) は、平均値に 2 回の急激な変化を伴うランダムな信号です。

例: reshape(randn(100,3).*[1 20 5],1,300) は、平方根平均二乗レベルに 2 回の急激な変化を伴うランダムな信号です。

データ型: single | double

名前と値の引数

すべて折りたたむ

オプションの引数のペアを Name1=Value1,...,NameN=ValueN として指定します。ここで、Name は引数名で、Value は対応する値です。名前と値の引数は他の引数の後に指定しなければなりませんが、ペアの順序は重要ではありません。

例: MaxNumChanges=3,Statistic="rms",MinDistance=20 は、平方根平均二乗レベルの変化が最も大きく、かつ少なくとも 20 サンプルで隔てられている点を最大で 3 つ検出します。

R2021a より前では、コンマを使用して名前と値をそれぞれ区切り、Name を引用符で囲みます。

例: 'MaxNumChanges',3,'Statistic',"rms",'MinDistance',20 は、平方根平均二乗レベルの変化が最も大きく、かつ少なくとも 20 サンプルで隔てられている点を最大で 3 つ検出します。

返す大きな変化の最大数。整数スカラーとして指定します。最も大きな変化のある点を検出した後、findchangepts は、指定した最大値を超えない範囲でより多くの変化点を含めるために、徐々にその検索基準を緩和します。何らかの探索設定によってこの最大値を超えて返されている場合、それ以降、関数は何も返しません。MaxNumChanges が指定されていない場合、関数は最も大きな変化のある点を返します。MinThresholdMaxNumChanges を同時に指定することはできません。

例: findchangepts([0 1 0]) は 2 番目のサンプルのインデックスを返します。

例: findchangepts([0 1 0],MaxNumChanges=1) は空行列を返します。

例: findchangepts([0 1 0],MaxNumChanges=2) は 2 番目の点と 3 番目の点のインデックスを返します。

データ型: single | double

検出する変化のタイプ。次の値のいずれかとして指定します。

  • "mean" — 平均値の変化を検出。出力引数なしで findchangepts を呼び出すと、関数は、信号、変化点、連続変化点によって囲まれる各セグメントの平均値をプロットします。

  • "rms" — 平方根平均二乗レベルの変化を検出。出力引数なしで findchangepts を呼び出すと、関数は信号および変化点をプロットします。

  • "std" — ガウス対数尤度を使用して、標準偏差の変化を検出。出力引数なしで findchangepts を呼び出すと、関数は、信号、変化点、連続変化点によって囲まれる各セグメントの平均値をプロットします。

  • "linear" — 平均値と傾きの変化を検出。出力引数なしで findchangepts を呼び出すと、関数は、信号、変化点、および連続変化点によって囲まれる信号の各部分に最適なラインをプロットします。

例: findchangepts([0 1 2 1],Statistic="mean") は 2 番目のサンプルのインデックスを返します。

例: findchangepts([0 1 2 1],Statistic="rms") は 3 番目のサンプルのインデックスを返します。

変化点間のサンプルの最小数。整数スカラーとして指定します。この数値を指定しない場合、既定の設定は、平均値の変化については 1 に、その他の変化については 2 になります。

例: findchangepts(sin(2*pi*(0:10)/5),MaxNumChanges=5,MinDistance=1) は 5 インデックスを返します。

例: findchangepts(sin(2*pi*(0:10)/5),MaxNumChanges=5,MinDistance=3) は 2 インデックスを返します。

例: findchangepts(sin(2*pi*(0:10)/5),MaxNumChanges=5,MinDistance=5) はインデックスを返しません。

データ型: single | double

各変化点に対する残差誤差の合計の最小改善値。ペナルティを表す実数スカラーとして指定します。このオプションは、予想される各変化点に追加のペナルティを適用することによって、返される大きな変化の数を制限するように作用します。MinThresholdMaxNumChanges を同時に指定することはできません。

例: findchangepts([0 1 2],MinThreshold=0) は 2 インデックスを返します。

例: findchangepts([0 1 2],MinThreshold=1) は 1 インデックスを返します。

例: findchangepts([0 1 2],MinThreshold=2) はインデックスを返しません。

データ型: single | double

出力引数

すべて折りたたむ

変化点の位置。整数インデックスのベクトルとして返されます。

モデル化された変化に対する信号の残差誤差。ベクトルとして返されます。

詳細

すべて折りたたむ

参照

[1] Killick, Rebecca, Paul Fearnhead, and Idris A. Eckley. “Optimal detection of changepoints with a linear computational cost.” Journal of the American Statistical Association. Vol. 107, No. 500, 2012, pp. 1590–1598.

[2] Lavielle, Marc. “Using penalized contrasts for the change-point problem.” Signal Processing. Vol. 85, August 2005, pp. 1501–1510.

拡張機能

すべて展開する

バージョン履歴

R2016a で導入