このページの内容は最新ではありません。最新版の英語を参照するには、ここをクリックします。
欠損サンプルのある信号の処理
うるう年の 2012 年に (ポンド単位で) 記録された特定の人の体重について考えます。この人は自分の体重を毎日は記録していませんでした。信号の周期性を調べたいとします。ただし、それを行う前に、欠損データの処理をしなければなりません。
データを読み込み、測定値をキログラムに変換します。欠損している読み取りが NaN
に設定されます。いくつの点が欠損しているかを判定します。
load('weight2012.dat') wgt = weight2012(:,2)/2.20462; daynum = 1:length(wgt); missing = isnan(wgt); fprintf('Missing %d samples of %d\n',sum(missing),max(daynum))
Missing 27 samples of 366
resample
を使用して欠損点に値を代入します。既定の設定では、resample
は線形内挿を使用して推定値を生成します。元の読み取りと内挿された読み取りをプロットします。欠損点のほぼ半数が含まれている第 200 ~ 250 日を拡大します。
wgt_orig = wgt; wgt = resample(wgt,daynum); plot(daynum,wgt_orig,'.-',daynum,wgt,'o') xlabel('Day') ylabel('Weight (kg)') axis([200 250 73 77]) legend('Original','Interpolated') grid
信号が周期的かどうかを周波数領域で解析して判定します。サイクルの期間を週数で測定し、求めます。変動に集中できるよう、平均値を減算します。
Fs = 7;
[p,f] = pwelch(wgt-mean(wgt),[],[],[],Fs);
plot(f,p)
xlabel('Frequency (week^{-1})')
grid
該当者の体重が週次でどのように変動したかに注目します。週単位で注目すべきパターンがあるでしょうか。年の最後の 2 日間は 52 週にするため消去します。曜日別に測定値を並べ替えます。
wgd = reshape(wgt(1:7*52),[7 52]); plot(wgd') xlabel('Week') ylabel('Weight (kg)') dweek = datetime([repmat([2012 1],7,1) (1:7)'],'Format','eeee'); legend(string(dweek),'Location','northwest') grid
データのサブセットを低次多項式に近似するフィルターを使用して、変動を滑らかにします。特に、7 日間のセットを 3 次多項式に近似するように設定します。
wgs = sgolayfilt(wgd',3,7); plot(wgs) xlabel('Week') ylabel('Weight (kg)') title('Smoothed Weight Fluctuations') legend(string(dweek),'Location','southeast'); grid
該当者は、週末に他の曜日よりも多く食べるため体重も増加する傾向があります。日次平均値を計算して確認します。
for jk = 1:7 fprintf('%s mean: %5.1f kg\n',dweek(jk),mean(wgd(jk,:))) end
Sunday mean: 76.2 kg Monday mean: 75.7 kg Tuesday mean: 75.2 kg Wednesday mean: 74.9 kg Thursday mean: 75.1 kg Friday mean: 75.3 kg Saturday mean: 75.8 kg