初めてのMATLABフィルター：ステップバイステップ解説

はじめに：信号処理とフィルターの役割

私たちが日常で扱う様々な信号（音声、画像、生体情報、センサーデータなど）には、目的とする情報以外に不要な成分（ノイズ）が含まれていることがよくあります。また、特定の周波数帯域の情報だけを取り出したり、逆に特定の周波数帯域の情報を取り除きたい場合もあります。このような目的を達成するために不可欠なのが「フィルター」です。

フィルターは、信号に含まれる特定の周波数成分を選択的に通過させたり、あるいは除去したりする処理を行います。例えば、ラジオのチューニングは特定の周波数帯域の電波信号だけを選び出すフィルターの一種です。音声処理では、人の声だけを強調しノイズを除去したり、特定の楽器の音色を調整したりするのにフィルターが使われます。生体信号処理では、心電図から呼吸や筋肉のノイズを取り除くためにフィルターが利用されます。

デジタル信号処理の分野では、これらのフィルター処理はソフトウェアで行われます。そのための強力なツールとして広く使われているのがMATLABです。MATLABのSignal Processing Toolboxには、様々な種類のフィルターを設計、解析、適用するための豊富な関数が用意されており、信号処理の研究開発において標準的な環境となっています。

この記事は、MATLABを使ってデジタルフィルターを設計・実装したいと考えている初心者の方を対象としています。フィルターの基本的な概念から、MATLABを使った具体的な設計手順、解析方法、そして実際の信号への適用方法までを、ステップバイステップで丁寧に解説していきます。約5000語の詳細な解説を通じて、MATLABでのフィルター処理の基礎を習得し、ご自身の信号処理タスクに応用できるようになることを目指します。

フィルターの基本概念

MATLABでのフィルター実装に入る前に、まずはデジタルフィルターに関するいくつかの重要な基本概念を理解しておきましょう。

アナログフィルターとデジタルフィルター

フィルターには、アナログフィルターとデジタルフィルターがあります。

アナログフィルター: 電子回路（抵抗、コンデンサ、コイル、オペアンプなど）を使って実現されるフィルターです。連続時間信号（アナログ信号）を扱います。
デジタルフィルター: コンピューターなどのデジタル処理装置上でソフトウェアや専用ハードウェアによって実現されるフィルターです。離散時間信号（デジタル信号）を扱います。アナログ信号をデジタルフィルターで処理するには、まずサンプリングと量子化によってデジタル信号に変換する必要があります。

現代の信号処理の多くはデジタルで行われています。デジタルフィルターは、アナログフィルターに比べて以下のような多くの利点があります。

高い精度と再現性: 電子部品のばらつきの影響を受けず、設計通りの特性を正確に実現できます。
柔軟性: ソフトウェアの変更だけで特性を変えることができます。
複雑な特性の実現: アナログでは難しい複雑な特性を持つフィルターも設計可能です。
温度や経年変化に強い: ハードウェアの物理的な変化による影響を受けません。
安定性: 特定の設計手法を用いれば、常に安定なフィルターを実現できます。

主なデジタルフィルターの種類：FIRとIIR

デジタルフィルターは、その構造やインパルス応答の特性から大きく2種類に分けられます。

FIRフィルター (Finite Impulse Response Filter): インパルス（ごく短いパルス信号）を入力したときの出力（インパルス応答）が有限の長さでゼロになるフィルターです。現在の出力は、現在の入力と過去のいくつかの入力の線形結合で計算されます。
- 利点: 常に安定であり、設計が比較的容易です。特に、入力信号の波形を歪ませない「線形位相」特性を持つフィルターを容易に実現できます。
- 欠点: 同程度の周波数選択性を実現するのに、IIRフィルターよりも高い次数が必要になる傾向があります（より多くの計算が必要）。
IIRフィルター (Infinite Impulse Response Filter): インパルス応答が無限に続くフィルターです。現在の出力は、現在の入力と過去のいくつかの入力だけでなく、過去のいくつかの出力も使って計算されます。これはフィードバック（帰還）を含む構造に対応します。
- 利点: 同程度の周波数選択性を実現するのに、FIRフィルターよりも低い次数で済みます（計算量が少なくて済む）。アナログフィルターの設計理論を応用できます。
- 欠点: 線形位相を実現するのは困難です（位相歪みが発生しやすい）。特定の設計でないと不安定になる可能性があります。

どちらのタイプのフィルターを選択するかは、アプリケーションの要件（位相特性が必要か、計算資源に制約があるかなど）によって決まります。

フィルター設計における重要な用語

フィルターの特性を記述し、設計仕様を定めるためにいくつかの重要な用語があります。

周波数応答 (Frequency Response): フィルターに様々な周波数の正弦波を入力したとき、出力される正弦波の振幅（ゲイン）と位相が周波数によってどのように変化するかを示す特性です。通常、振幅応答（ゲイン応答）と位相応答の二つのグラフで表現されます。
- 振幅応答: 入力に対する出力の振幅の比率。通常デシベル（dB）単位で表されます。0 dBはゲイン1（振幅変化なし）を意味し、負のdB値は減衰、正のdB値は増幅を意味します。
- 位相応答: 入力正弦波に対する出力正弦波の位相のずれ。
カットオフ周波数 (Cutoff Frequency): フィルターが信号を「通過させる帯域」と「阻止する帯域」の境界となる周波数です。
通過帯域 (Passband): 信号をほとんど減衰させずに通過させる周波数帯域。
阻止帯域 (Stopband): 信号を大きく減衰させる周波数帯域。
遷移帯域 (Transition Band): 通過帯域と阻止帯域の間の帯域。この帯域ではゲインが急峻に変化します。遷移帯域が狭いほど「切れの良い」フィルターですが、実現には高い次数が必要になります。
リップル (Ripple): 通過帯域や阻止帯域におけるゲインの変動。フィルターの種類によっては、これらの帯域でゲインが完全に平坦ではなく、わずかに変動します。
減衰量 (Attenuation): 阻止帯域で信号がどれだけ減衰するかを示す値。通常、最小減衰量（阻止帯域のゲインの絶対値の最小値）で仕様を定めます。
フィルター次数 (Filter Order): フィルターの複雑さを示すパラメータです。FIRフィルターの場合はインパルス応答の長さに関係し、IIRフィルターの場合は差分方程式の次数に関係します。次数が高いほど、通常は遷移帯域が狭く、阻止帯域減衰量が大きくなるなど、より理想的な特性に近づきますが、計算量が増え、FIRフィルターでは遅延（群遅延）も大きくなります。
群遅延 (Group Delay): 各周波数成分がフィルターを通過するのにかかる時間的な遅延です。線形位相フィルター（主にFIR）は通過帯域で一定の群遅延を持ち、信号の波形を歪ませません。非線形位相フィルター（主にIIR）は周波数によって群遅延が異なり、信号の波形を歪ませる可能性があります。

MATLAB信号処理ツールボックスの紹介

MATLABのSignal Processing Toolboxには、デジタルフィルターを扱うための非常に多くの関数が用意されています。ここでは、フィルター設計の基本的な流れでよく使う主要な関数を紹介します。

フィルター設計のための関数

これらの関数は、指定した仕様を満たすフィルターの係数（または状態空間モデルなど）を計算します。

designfilt: 様々なフィルタータイプ、設計法、仕様を指定できる汎用性の高い設計関数です。現代的なフィルター設計ではこの関数を使うことが推奨されます。フィルターオブジェクトとして設計結果を返します。
fir1: ウィンドウ法を用いた線形位相FIRフィルターを設計します。最も基本的なFIR設計関数の一つで、使い方が比較的簡単です。分子係数 b を返します。
butter: バターワースIIRフィルターを設計します。通過帯域・阻止帯域ともにリップルがない滑らかな特性のフィルターです。分子係数 b と分母係数 a を返します。
cheby1: チェビシェフType I IIRフィルターを設計します。通過帯域にリップルを持ちますが、バターワースより急峻な遷移帯域を実現できます。
cheby2: チェビシェフType II IIRフィルターを設計します。阻止帯域にリップルを持ちますが、通過帯域は滑らかで、こちらもバターワースより急峻です。
ellip: 楕円IIRフィルターを設計します。通過帯域と阻止帯域の両方にリップルを持ちますが、同じ次数で最も急峻な遷移帯域を実現できます。

フィルター解析のための関数

設計したフィルターの特性（主に周波数応答）を確認するための関数です。

freqz: デジタルフィルターの周波数応答を計算し、プロットします。振幅応答と位相応答を確認するのに使います。
fvtool: Filter Visualization Toolの略で、フィルターの様々な特性（周波数応答、インパルス応答、ステップ応答、極ゼロ配置など）をインタラクティブに表示・解析できるGUIツールです。フィルター設計の強力な味方となります。

信号へのフィルター適用関数

設計したフィルターを実際の信号データに適用するための関数です。

filter: フィルター係数を使って、入力信号を時系列に沿ってフィルタリングします。IIRフィルターや非線形位相FIRフィルターに使うと、信号に位相遅延が発生します。
filtfilt: IIRフィルターや非線形位相FIRフィルターによる位相遅延を補償するために、信号を順方向と逆方向の両方から2回フィルタリングします。これによりゼロ位相フィルタリングを実現できますが、信号長が長くなり、リアルタイム処理には使えません。オフライン処理（事前にデータ全体が取得済みの場合）に用います。

ステップ1：フィルター設計の準備

MATLABでフィルターを設計する前に、まず処理対象の信号と必要なフィルターの仕様を明確にする必要があります。

1. 処理対象の信号の特性理解

サンプリング周波数 (Fs): 信号が1秒間に何回サンプリングされているかを示す値です（単位 Hz）。デジタル信号処理において最も基本的なパラメータであり、フィルターのカットオフ周波数などはこのFsを基準に決定されます。
信号に含まれる周波数成分: 処理対象の信号に、目的の信号がどの周波数帯域に存在し、ノイズがどの周波数帯域に存在するかを把握します。これを知るには、MATLABのfft関数（高速フーリエ変換）を使って信号の周波数スペクトルを分析するのが効果的です。

例えば、サンプリング周波数1000 Hzで取得した信号があり、そこに50 Hzの目的信号と、150 Hzや300 Hzのノイズが含まれていると分かったとします。このノイズを除去するには、50 Hz付近は通過させ、150 Hz以上の周波数は阻止するローパスフィルターが必要になりそうだ、と検討がつきます。

2. 必要なフィルターの仕様定義

信号の特性理解に基づいて、どのようなフィルターが必要か、その仕様を具体的に定義します。

フィルターの種類: ローパス（低周波通過）、ハイパス（高周波通過）、バンドパス（特定帯域通過）、バンドストップ（特定帯域阻止）など、信号のどこを通過/阻止させたいかによって選択します。
カットオフ周波数や帯域幅: 何Hz以上の周波数を阻止したいか、何Hzから何Hzの帯域を通過させたいかなど、具体的な周波数値を決めます。
フィルタータイプ: FIRにするか、IIRにするか。位相特性が重要ならFIR（またはfiltfilt）、計算量を重視するならIIRなどが選択のポイントになります。
フィルター次数: 高次ほど理想特性に近づきますが、計算量や遅延が増えます。まずは低い次数で試してみて、必要に応じて調整するのが一般的です。あるいは、必要な阻止帯域減衰量などを指定して、関数に次数を推定させる方法もあります。
その他の仕様: 必要に応じて、阻止帯域での最小減衰量（例: -40 dB以下）、通過帯域での最大リップル量などを定義します。

3. MATLAB環境の準備

MATLABとSignal Processing Toolboxがインストールされ、使用できる状態であることを確認します。通常、MATLABがインストールされていれば、ツールボックスもインストールされているはずです。コマンドウィンドウで適当な信号処理関数（例: fir1）を入力し、ヘルプが表示されれば使用可能です。

ステップ2：基本的なフィルター設計

ここでは、最も基本的なフィルターであるローパスフィルターの設計を例に、具体的なMATLABコードを使って設計手順を見ていきましょう。FIRフィルターとIIRフィルターの両方を設計します。

サンプリング周波数と目標周波数の設定

まず、サンプリング周波数 Fs と、設計したいフィルターのカットオフ周波数 Fc をHz単位で定義します。

“`matlab
% サンプリング周波数を設定 (例: 1000 Hz)
Fs = 1000;

% 目標のカットオフ周波数を設定 (例: 100 Hz)
Fc = 100;
“`

デジタルフィルターの設計関数では、カットオフ周波数などの周波数パラメータを「正規化された周波数」で指定することが一般的です。正規化された周波数は、サンプリング周波数の半分、つまりNyquist周波数 (Fs/2) を基準として0から1の間の値で表されます。0はDC（0 Hz）、1はNyquist周波数に対応します。

したがって、Hz単位のカットオフ周波数 Fc を正規化された周波数 Wn に変換するには、次のように計算します。

Wn = Fc / (Fs / 2)

この正規化された周波数は、FIRフィルター設計関数 fir1 や、IIRフィルター設計関数 butter, cheby1, cheby2, ellip で引数として使用します。

matlab % カットオフ周波数を正規化された周波数に変換 % Nyquist 周波数 (Fs/2) を基準にする Wn = Fc / (Fs / 2); % Wn は 0 から 1 の間の値になります (ただし、理想的には 0 < Wn < 1) fprintf('正規化されたカットオフ周波数 Wn = %.4f (Nyquist周波数基準)\n', Wn);

FIR ローパスフィルターの設計 (`fir1`)

fir1 関数は、ウィンドウ法を用いて線形位相FIRフィルターを設計します。基本的な使い方は以下の通りです。

b = fir1(n, Wn, type, window)

n: フィルターの次数。FIRフィルターの次数は、インパルス応答の長さ（係数の数）から1を引いた値です。例えば、次数100のFIRフィルターは101個の係数を持ちます。
Wn: 正規化されたカットオフ周波数。ローパスフィルターの場合は単一の値です。バンドパス/バンドストップの場合は2要素のベクトル [Wn1, Wn2] となります。
type: フィルタータイプ。ローパスなら 'low'、ハイパスなら 'high'、バンドパスなら 'bandpass'、バンドストップなら 'stop' を指定します。デフォルトはローパスです。
window: 使用する窓関数。例えば 'hamming' (デフォルト), 'hann', 'blackman' など。窓関数は遷移帯域の幅や阻止帯域の減衰量に影響します。指定しない場合はハミング窓が使われます。
戻り値 b: 設計されたフィルターの分子係数（インパルス応答）。FIRフィルターの場合、分母係数は常に1です。

例として、次数100のローパスFIRフィルターを設計してみましょう。

“`matlab
% FIRフィルターの次数を設定
n_fir = 100;

% fir1 関数でローパスFIRフィルターを設計
% b_fir はフィルター係数（インパルス応答）
b_fir = fir1(n_fir, Wn, ‘low’);

% 設計したフィルター係数の数を確認 (次数+1)
fprintf(‘FIRフィルター係数の数: %d\n’, length(b_fir));
“`

これで、ローパスFIRフィルターの係数 b_fir が得られました。この係数を使って信号をフィルタリングします。

IIR ローパスフィルターの設計 (`butter`)

butter 関数は、バターワース特性を持つIIRフィルターを設計します。基本的な使い方は以下の通りです。

[b, a] = butter(n, Wn, type)

n: フィルターの次数。IIRフィルターの次数は、差分方程式の次数に対応し、FIRフィルターよりも少ない次数で同等の周波数選択性を得られる傾向があります。
Wn: 正規化されたカットオフ周波数。fir1と同様に指定します。
type: フィルタータイプ。fir1と同様に 'low', 'high', 'bandpass', 'stop' を指定します。
戻り値 b: 設計されたフィルターの分子係数。
戻り値 a: 設計されたフィルターの分母係数。

例として、次数10のローパスバターワースIIRフィルターを設計してみましょう。IIRフィルターはFIRに比べて低次で済む傾向があるため、ここでは次数を低く設定しています。

“`matlab
% IIRフィルターの次数を設定
n_iir = 10;

% butter 関数でローパスIIRフィルターを設計
% b_iir, a_iir はフィルター係数（分子、分母）
[b_iir, a_iir] = butter(n_iir, Wn, ‘low’);

% 設計したフィルター係数の数を確認
fprintf(‘IIRフィルター分子係数の数: %d\n’, length(b_iir));
fprintf(‘IIRフィルター分母係数の数: %d\n’, length(a_iir));
“`

これで、ローパスIIRフィルターの係数 b_iir と a_iir が得られました。これらを使って信号をフィルタリングします。

`designfilt`を使った設計

より柔軟な設計や、FIR/IIRを問わず様々な設計法を使いたい場合は、designfilt 関数が便利です。designfilt はフィルターのタイプ、応答タイプ、設計法、仕様などを細かく指定できます。

d = designfilt(Type, DesignMethod, Name, Value, ...)

Type: 'lowpassfir', 'lowpassiir', 'bandpassfir', 'bandpassiir' など、フィルターの種類とタイプを指定。
DesignMethod: 'window' (FIR), 'butter' (IIR), 'cheby1' (IIR), 'equiripple' (FIR) など、設計法を指定。
Name, Value ペア: 設計仕様を名前と値のペアで指定します。例: 'FilterOrder', n, 'CutoffFrequency', Fc, 'SampleRate', Fs。ここでは正規化周波数ではなく、Hz単位の周波数やサンプリング周波数を直接指定できるのが便利な点です。

例として、designfilt を使って先ほどと同じ仕様（Fs=1000Hz, Fc=100Hz）のFIRローパスフィルターを設計してみましょう。次数も同様に100とします。

“`matlab
% designfilt 関数でFIRローパスフィルターを設計
% SampleRate を指定することで、CutoffFrequency を Hz 単位で指定可能
d_fir = designfilt(‘lowpassfir’, ‘FilterOrder’, n_fir, …
‘CutoffFrequency’, Fc, ‘SampleRate’, Fs);

% 設計したフィルターのタイプを確認
fprintf(‘Designfiltで設計したFIRフィルタータイプ: %s\n’, d_fir.Type);
% 係数は d_fir.Coefficients で取得可能
b_fir_designfilt = d_fir.Coefficients;
fprintf(‘Designfiltで設計したFIRフィルター係数の数: %d\n’, length(b_fir_designfilt));

% fir1 と designfilt の係数がほぼ同じか確認（厳密にはオプションで異なる場合も）
% isequal(b_fir, b_fir_designfilt) % 同じ設計法なら true に近い値になるはず
“`

次に、IIRローパスフィルター（バターワース）を設計してみましょう。次数は同様に10とします。

“`matlab
% designfilt 関数でIIRローパスフィルター（バターワース）を設計
d_iir = designfilt(‘lowpassiir’, ‘DesignMethod’, ‘butter’, …
‘FilterOrder’, n_iir, ‘CutoffFrequency’, Fc, ‘SampleRate’, Fs);

% 設計したフィルターのタイプを確認
fprintf(‘Designfiltで設計したIIRフィルタータイプ: %s\n’, d_iir.Type);
% 係数は d_iir.Numerator と d_iir.Denominator で取得可能
b_iir_designfilt = d_iir.Numerator;
a_iir_designfilt = d_iir.Denominator;
fprintf(‘Designfiltで設計したIIRフィルター分子係数の数: %d\n’, length(b_iir_designfilt));
fprintf(‘Designfiltで設計したIIRフィルター分母係数の数: %d\n’, length(a_iir_designfilt));

% butter と designfilt の係数がほぼ同じか確認
% isequal(b_iir, b_iir_designfilt)
% isequal(a_iir, a_iir_designfilt)
“`

designfilt を使うと、正規化周波数の計算なしにHz単位で周波数を指定できたり、様々な設計法や仕様（通過帯域リップル、阻止帯域減衰量など）を細かく指定できたりするため、現代的な設計ではこちらを使うのがおすすめです。設計結果はフィルターオブジェクトとして返され、そのプロパティとして係数などが格納されています。

ステップ3：設計したフィルターの解析

フィルターを設計したら、その特性が設計仕様を満たしているかを確認する必要があります。最も重要な特性は周波数応答です。MATLABには freqz 関数と fvtool という便利なツールがあります。

`freqz` を使った周波数応答プロット

freqz 関数は、デジタルフィルターの周波数応答（ゲインと位相）を計算し、プロットします。

FIRフィルターの場合: freqz(b_fir, 1, N, Fs)
IIRフィルターの場合: freqz(b_iir, a_iir, N, Fs)
引数 b, a: フィルターの分子係数と分母係数（FIRの場合は a は1）。
引数 N: 周波数応答を計算する周波数点の数。多いほど滑らかなプロットになります。
引数 Fs: サンプリング周波数 (Hz)。これを指定すると、プロットの横軸がHz単位になります。省略すると正規化周波数（0～π rad/sample または 0～1）になります。

設計したFIRおよびIIRローパスフィルターの周波数応答をプロットしてみましょう。

“`matlab
% 周波数応答を計算する周波数点の数
N_freq = 1024;

% FIRフィルターの周波数応答をプロット
figure;
freqz(b_fir, 1, N_freq, Fs);
title(‘FIR ローパスフィルター周波数応答’);
% プロットはゲイン（上段、dBスケール）と位相（下段、度またはラジアン）を表示

% IIRフィルターの周波数応答をプロット
figure;
freqz(b_iir, a_iir, N_freq, Fs);
title(‘IIR ローパスフィルター周波数応答’);
% プロットはゲイン（上段、dBスケール）と位相（下段、度またはラジアン）を表示
“`

これらのプロットを確認することで、以下の点をチェックできます。

ゲイン応答 (Magnitude Response):
- 通過帯域（0 Hzからカットオフ周波数Fcまで）でゲインが0 dB（またはそれに近い値）になっているか。
- 阻止帯域（カットオフ周波数Fcより高い周波数）でゲインが大きく負の値（十分に減衰）になっているか。
- カットオフ周波数Fcのところでゲインがどの程度になっているか（バターワースフィルターの場合、通常 Fc でゲインが -3 dB になります）。
- 通過帯域や阻止帯域に大きなリップルがないか（バターワースはリップルなし、チェビシェフや楕円はリップルあり）。
- 遷移帯域の急峻さ（周波数に応じてゲインがどれだけ急激に変化するか）。
位相応答 (Phase Response):
- FIRフィルターは通過帯域でほぼ線形（直線的）な位相応答を持っているか。線形位相は信号の波形を歪ませない重要な特性です。
- IIRフィルターは通過帯域でも非線形な位相応答を持っていることが多いです。これが信号の波形歪みの原因になります。

`fvtool` を使ったインタラクティブ解析

fvtool は、フィルターの様々な特性をGUIでインタラクティブに確認できる強力なツールです。コマンドウィンドウで fvtool(b, a) または fvtool(d) と入力して起動します。

FIRフィルターの場合: fvtool(b_fir, 1) または fvtool(d_fir)
IIRフィルターの場合: fvtool(b_iir, a_iir) または fvtool(d_iir)

“`matlab
% fvtool を使ってFIRフィルターを解析
fvtool(b_fir, 1);
title(‘FIRフィルター解析 (fvtool)’);

% fvtool を使ってIIRフィルターを解析
fvtool(b_iir, a_iir);
title(‘IIRフィルター解析 (fvtool)’);
“`

fvtool のウィンドウが開くと、デフォルトで周波数応答（振幅応答と位相応答）が表示されます。ツールバーやメニューを使って、以下の様々な特性表示に切り替えることができます。

Magnitude Response (振幅応答 – dB)
Phase Response (位相応答 – 度またはラジアン)
Group Delay (群遅延 – サンプル数または時間単位)
Impulse Response (インパルス応答)
Step Response (ステップ応答)
Pole/Zero Plot (極ゼロプロット)
Filter Coefficients (フィルター係数)

特に、群遅延 は位相応答から導かれる重要な特性です。通過帯域で群遅延が一定のフィルター（線形位相フィルター）は、信号の各周波数成分を同じ時間だけ遅延させるため、波形を歪ませません。一方、群遅延が周波数によって大きく異なるフィルターは、波形を歪ませます。fvtoolでFIRとIIRの群遅延を比較してみてください。FIRは通過帯域でほぼ一定、IIRは周波数によって変動する様子が確認できるはずです。

fvtool 上でカーソルを動かすと、特定の周波数での正確なゲインや位相、群遅延などを読み取ることができます。設計仕様（例: カットオフ周波数 Fc でのゲイン、阻止帯域での最小減衰量）が満たされているか、このツールを使って詳細に確認しましょう。もし満たされていなければ、フィルターの次数を増やすなどの調整が必要になります。

ステップ4：信号へのフィルター適用

フィルターの設計と解析が完了し、仕様を満たしていることが確認できたら、いよいよ実際の信号データにフィルターを適用します。これには filter 関数や filtfilt 関数を使用します。

テスト信号の生成

フィルターの効果を確認するために、複数の周波数成分を含むテスト信号を生成してみましょう。ここでは、目的とする低周波成分と、除去したい高周波成分を混ぜ合わせた信号を作成します。

“`matlab
% サンプリング周波数と信号長を設定
Fs = 1000; % サンプリング周波数 (Hz)
T = 1; % 信号長 (秒)
t = 0:1/Fs:T-1/Fs; % 時間ベクトル

% 目的信号 (低周波: 50 Hz)
f1 = 50; % Hz
sig_low = 0.7 * sin(2pif1*t);

% 不要なノイズ信号 (高周波: 150 Hz と 300 Hz)
f2 = 150; % Hz
f3 = 300; % Hz
sig_noise = 0.5 * sin(2pif2t) + 0.3 * sin(2pif3t);

% 合成信号 (目的信号 + ノイズ)
sig_noisy = sig_low + sig_noise;

% さらにランダムノイズを追加 (オプション)
% sig_noisy = sig_noisy + 0.1*randn(size(t));

% 元信号とノイズ入り信号をプロット
figure;
subplot(2,1,1);
plot(t, sig_low);
title(‘元の目的信号 (50 Hz)’);
xlabel(‘時間 (秒)’); ylabel(‘振幅’);
grid on;

subplot(2,1,2);
plot(t, sig_noisy);
title(‘ノイズ入り信号 (50 Hz + 150 Hz + 300 Hz)’);
xlabel(‘時間 (秒)’); ylabel(‘振幅’);
grid on;
“`

このテスト信号 sig_noisy に、先ほど設計したカットオフ周波数 100 Hz のローパスフィルターを適用することで、150 Hz と 300 Hz のノイズ成分が除去されることを期待します。

ノイズ入り信号の周波数スペクトルを確認してみましょう。

“`matlab
% 信号長
L = length(sig_noisy);

% FFTを実行
Y = fft(sig_noisy);

% 片側スペクトルを計算
P2 = abs(Y/L);
P1 = P2(1:L/2+1); % Nyquist周波数までの範囲
P1(2:end-1) = 2*P1(2:end-1); % DC成分とNyquist成分以外は振幅を2倍

% 周波数ベクトルを作成
f = Fs*(0:(L/2))/L;

% 周波数スペクトルをプロット
figure;
plot(f, P1);
title(‘ノイズ入り信号の周波数スペクトル’);
xlabel(‘周波数 (Hz)’); ylabel(‘振幅’);
grid on;
xlim([0 Fs/2]); % Nyquist周波数まで表示
“`

プロットを見ると、50 Hzのピークに加えて、150 Hzと300 Hzにも明確なピークがあることが確認できます。ローパスフィルターを適用することで、これらの150 Hzと300 Hzのピークが大幅に減衰されるはずです。

`filter` 関数を使ったフィルタリング

filter 関数は、最も一般的なデジタルフィルター適用関数です。線形定常システム（フィルター）の入出力関係を表す差分方程式に基づいて、入力信号をフィルタリングします。

y = filter(b, a, x)

b: フィルターの分子係数ベクトル。
a: フィルターの分母係数ベクトル。FIRフィルターの場合は a は 1 になります。
x: 入力信号ベクトル。
戻り値 y: フィルタリングされた出力信号ベクトル。

設計したFIRおよびIIRローパスフィルターをノイズ入り信号に適用してみましょう。

“`matlab
% 設計済みのフィルター係数を使用
% FIR: b_fir, a=1
% IIR: b_iir, a_iir

% FIRフィルターを適用
filtered_sig_fir = filter(b_fir, 1, sig_noisy);

% IIRフィルターを適用
filtered_sig_iir = filter(b_iir, a_iir, sig_noisy);
“`

`filtfilt` 関数を使ったゼロ位相フィルタリング (IIR向け)

filter 関数は因果的フィルターとして動作するため、信号に時間的な遅延（群遅延）が発生します。特にIIRフィルターは周波数によって遅延量が異なる（非線形位相）ため、信号の波形が歪む可能性があります。この位相歪みを避けたい場合（特にオフライン処理で波形を忠実に再現したい場合）には、filtfilt 関数が有用です。

y = filtfilt(b, a, x)

b, a: フィルター係数。
x: 入力信号ベクトル。
戻り値 y: フィルタリングされた出力信号ベクトル（ゼロ位相）。

filtfilt は、入力信号の最初と最後に適切なパディングを行い、その信号に対して順方向と逆方向の2回 filter を適用します。順方向のフィルタリングで生じた位相遅延は、逆方向のフィルタリングで打ち消され、結果として全体の位相応答がゼロ（群遅延もゼロ）になります。ただし、信号の最初と最後で過渡応答が発生しやすいため、パディングが必要になります。また、リアルタイム処理には使えません。

設計したIIRローパスフィルターを filtfilt で適用してみましょう。FIRフィルターは線形位相なので filter でも位相歪みは小さいですが、波形を完璧に保ちたい場合やFIRの次数が高い場合にも filtfilt を使うことがあります。

“`matlab
% IIRフィルターを filtfilt で適用 (ゼロ位相フィルタリング)
filtered_sig_iir_zp = filtfilt(b_iir, a_iir, sig_noisy);

% FIRフィルターを filtfilt で適用 (オプション)
% filtered_sig_fir_zp = filtfilt(b_fir, 1, sig_noisy);
“`

フィルタリング前後の信号比較

フィルタリングの効果を、時間領域と周波数領域の両方で確認してみましょう。

時間領域での比較:

元のノイズ入り信号、filterでフィルタリングした信号、filtfiltでフィルタリングした信号を同じ時間軸でプロットして比較します。

matlab figure; plot(t, sig_noisy, 'DisplayName', 'ノイズ入り信号'); hold on; plot(t, filtered_sig_fir, 'DisplayName', 'FIRフィルター後 (filter)'); plot(t, filtered_sig_iir, 'DisplayName', 'IIRフィルター後 (filter)'); plot(t, filtered_sig_iir_zp, 'DisplayName', 'IIRフィルター後 (filtfilt)'); hold off; title('フィルタリング前後の信号比較 (時間領域)'); xlabel('時間 (秒)'); ylabel('振幅'); legend show; grid on; xlim([0.2 0.4]); % 信号の一部を拡大して位相遅延や波形歪みを確認

この時間領域プロットを見ると、filter 関数を使った場合はフィルタリング後の信号が元の信号（理想的にはノイズがない sig_low）に対して時間的に遅れていることがわかります。IIRフィルター (filtered_sig_iir) は、FIRフィルター (filtered_sig_fir) よりも位相応答が非線形なため、遅延の様子が複雑になり、波形が歪んで見えることがあります。一方、filtfilt を使った信号 (filtered_sig_iir_zp) は、元の信号 sig_low とほとんど同じ時間位置にあり、位相遅延がない（ゼロ位相）ことが確認できるはずです。

周波数領域での比較:

フィルタリング後の信号の周波数スペクトルを計算し、元のノイズ入り信号のスペクトルと比較します。

“`matlab
% フィルタリング後信号のFFT
Y_fir = fft(filtered_sig_fir);
Y_iir = fft(filtered_sig_iir);
Y_iir_zp = fft(filtered_sig_iir_zp);

% 片側スペクトル
P2_fir = abs(Y_fir/L); P1_fir = P2_fir(1:L/2+1); P1_fir(2:end-1) = 2P1_fir(2:end-1);
P2_iir = abs(Y_iir/L); P1_iir = P2_iir(1:L/2+1); P1_iir(2:end-1) = 2P1_iir(2:end-1);
P2_iir_zp = abs(Y_iir_zp/L); P1_iir_zp = P2_iir_zp(1:L/2+1); P1_iir_zp(2:end-1) = 2*P1_iir_zp(2:end-1);

figure;
plot(f, P1, ‘DisplayName’, ‘ノイズ入り信号’);
hold on;
plot(f, P1_fir, ‘DisplayName’, ‘FIRフィルター後 (filter)’);
plot(f, P1_iir, ‘DisplayName’, ‘IIRフィルター後 (filter)’);
plot(f, P1_iir_zp, ‘DisplayName’, ‘IIRフィルター後 (filtfilt)’);
hold off;
title(‘フィルタリング前後の信号比較 (周波数領域)’);
xlabel(‘周波数 (Hz)’); ylabel(‘振幅’);
legend show;
grid on;
xlim([0 Fs/2]);
“`

この周波数スペクトルプロットを見ると、フィルタリング後の信号では、ローパスフィルターのカットオフ周波数（100 Hz）よりも高い周波数成分（150 Hzと300 Hzのピーク）が大幅に減衰していることが確認できるはずです。FIRフィルターとIIRフィルター（filter使用時）のスペクトルは非常に似通っているでしょう。filtfilt使用時も、振幅応答はほぼ同じになります（filtfiltは位相だけをゼロにする処理であり、振幅応答は元のIIRフィルターと同じ）。

この結果から、設計したローパスフィルターが期待通りに高周波ノイズを除去する機能を持っていることが分かります。

ステップ5：応用的なフィルター設計と考慮事項

ここまでの手順で基本的なフィルターの設計、解析、適用ができるようになりました。さらに応用的なフィルター設計や、実際の信号処理における注意点について見ていきましょう。

他の種類のフィルター設計

ローパスだけでなく、ハイパス、バンドパス、バンドストップフィルターも同様の手順で設計できます。fir1 や butter 関数では、type 引数を適切に変更するだけです。designfilt を使う場合は、Type 引数を変更します。

ハイパスフィルター: 指定した周波数以上の成分を通過させます。ノイズが低周波にある場合などに使用します。
- fir1(n, Wn, 'high')
- butter(n, Wn, 'high')
- designfilt('highpassfir', ...) または designfilt('highpassiir', ...)
バンドパスフィルター: 指定した二つの周波数間の帯域を通過させます。特定の信号成分だけを取り出したい場合に使用します。Wn は [Wn1, Wn2] のように低域側と高域側の正規化カットオフ周波数を指定するベクトルになります。
- fir1(n, [Wn1, Wn2], 'bandpass')
- butter(n, [Wn1, Wn2], 'bandpass')
- designfilt('bandpassfir', ...) または designfilt('bandpassiir', ...)
バンドストップフィルター: 指定した二つの周波数間の帯域を阻止します。特定の周波数帯域のノイズ（例: 商用電源のハムノイズ50/60Hz）を除去したい場合などに使用します。Wn の指定方法はバンドパスと同様です。
- fir1(n, [Wn1, Wn2], 'stop')
- butter(n, [Wn1, Wn2], 'stop')
- designfilt('bandstopfir', ...) または designfilt('bandstopiir', ...)

例として、周波数帯域 50 Hz から 150 Hz を通過させるバンドパスFIRフィルター (Fs=1000 Hz) を設計してみましょう。

“`matlab
Fs = 1000; % サンプリング周波数
Fpass1 = 50; % 通過帯域の下限 (Hz)
Fpass2 = 150; % 通過帯域の上限 (Hz)

% 正規化された周波数帯域
Wn_band = [Fpass1, Fpass2] / (Fs / 2);

% FIRバンドパスフィルターの次数
n_fir_bp = 100;

% fir1 でバンドパスFIRフィルターを設計
b_fir_bp = fir1(n_fir_bp, Wn_band, ‘bandpass’);

% designfilt でバンドパスFIRフィルターを設計
d_fir_bp = designfilt(‘bandpassfir’, ‘FilterOrder’, n_fir_bp, …
‘CutoffFrequency1’, Fpass1, ‘CutoffFrequency2’, Fpass2, …
‘SampleRate’, Fs);

% 周波数応答を確認
figure;
freqz(b_fir_bp, 1, 1024, Fs);
title(‘FIR バンドパスフィルター周波数応答’);

figure;
fvtool(d_fir_bp);
title(‘FIR バンドパスフィルター解析 (fvtool)’);
“`

IIRフィルターの種類と特性のトレードオフ

IIRフィルターは、バターワース以外にも様々な設計法があります。これらの設計法は、通過帯域や阻止帯域でのリップルの許容度と、遷移帯域の急峻さのトレードオフを提供します。

バターワース (butter): 通過帯域、阻止帯域ともにリップルがなく滑らか。遷移帯域は比較的広い。
チェビシェフ Type I (cheby1): 通過帯域にリップルを許容することで、バターワースより急峻な遷移帯域を実現。阻止帯域は滑らか。
チェビシェフ Type II (cheby2): 阻止帯域にリップルを許容することで、バターワースより急峻な遷移帯域を実現。通過帯域は滑らか。
楕円 (ellip): 通過帯域と阻止帯域の両方にリップルを許容することで、同じ次数で最も急峻な遷移帯域を実現。

これらのフィルターは、それぞれ cheby1, cheby2, ellip 関数を使って設計できます。また、designfilt 関数で 'DesignMethod' 引数として 'cheby1', 'cheby2', 'ellip' を指定し、通過帯域/阻止帯域のリップル量を指定することで設計することもできます。

リップルはアプリケーションによっては許容できない場合があるため、設計目的や許容誤差に応じて適切なフィルタータイプを選択する必要があります。

フィルター次数の選択

フィルター次数は、性能と計算量（そしてFIRの場合は遅延）のトレードオフを決定する重要なパラメータです。

次数が高い場合:
- 遷移帯域が狭くなる（より理想的な「急峻さ」に近づく）。
- 阻止帯域での減衰量が大きくなる。
- FIRフィルターでは遅延（群遅延）が大きくなる。
- 計算量が増える。
- IIRフィルターでは数値的な安定性が問題になる可能性がわずかにある。
次数が低い場合:
- 遷移帯域が広くなる。
- 阻止帯域での減衰量が小さくなる。
- 計算量が少ない。
- FIRフィルターでは遅延が小さい。

多くの場合、最初は低めの次数で設計・解析してみて、必要な性能（遷移帯域の急峻さ、阻止帯域減衰量など）が得られるまで徐々に次数を上げていくというアプローチを取ります。designfilt 関数では、必要な通過帯域/阻止帯域の仕様（周波数、リップル、減衰量）を指定して、その仕様を満たす最小の次数を推定させる機能もあります。

窓関数 (FIR)

fir1 関数でFIRフィルターを設計する際に使用する窓関数（デフォルトはハミング窓）は、フィルターの周波数応答に影響を与えます。異なる窓関数を選択することで、主に遷移帯域の幅と阻止帯域の最小減衰量の間のトレードオフを調整できます。

矩形窓 (Rectangular Window): 最も狭い遷移帯域を実現しますが、阻止帯域減衰量が最も小さく、リップルが大きい（特に阻止帯域）。fir1では直接指定しませんが、概念として重要。
ハミング窓 (Hamming Window): デフォルト。矩形窓より遷移帯域は広くなりますが、阻止帯域減衰量が向上します。
ハン窓 (Hann Window): ハミング窓と似た特性ですが、阻止帯域減衰量がさらに向上する代わりに、遷移帯域が若干広くなります。
ブラックマン窓 (Blackman Window): ハミング窓やハン窓よりも遷移帯域は広くなりますが、阻止帯域減衰量が最も大きく、リップルが小さくなります。

fir1(n, Wn, 'low', window) のように、4つ目の引数で 'hann' や 'blackman' などを指定して試すことができます。fvtool で異なる窓関数で設計したフィルターの周波数応答を比較してみると、その違いがよく分かります。

フィルター設計ツール (`filterDesigner` アプリ)

MATLABには、GUIで直感的にフィルターを設計できる「Filter Designer」という便利なアプリがあります。複雑な仕様を持つフィルターや、様々なパラメータを試しながら最適な設計を見つけたい場合に非常に役立ちます。

コマンドウィンドウで filterDesigner と入力して起動します。

matlab % Filter Designer アプリを起動 filterDesigner;

Filter DesignerのGUIでは、フィルタータイプ、設計法、応答タイプ、仕様（周波数、リップル、減衰量など）をドロップダウンメニューやテキストボックスで指定し、ボタンをクリックするだけでフィルターが設計されます。設計されたフィルターの周波数応答などの特性は、GUI上で即座に確認できます。

さらに、設計したフィルターをMATLABコード、オブジェクト、System objectなどとしてエクスポートする機能もあります。これにより、GUIで設計したフィルターをMATLABスクリプトや関数に簡単に組み込むことができます。コードをゼロから書くのが難しいと感じる場合は、このツールから始めるのも良いでしょう。

実世界の信号処理における注意点

実際に信号にフィルターを適用する際には、いくつかの注意点があります。

信号の端処理 (Edge Effects): filter 関数は、信号の最初の方で過渡応答（フィルターが「立ち上がる」までの不安定な応答）を生成することがあります。また、信号の最後に近づくにつれて、フィルターが将来の入力サンプルを必要とする（FIRの場合）または過去の出力サンプルが不足する（IIRの場合）ために、結果が不安定になることがあります。filtfilt はパディングによってこれらの影響をある程度緩和しますが、それでも信号の最初と最後のごく一部のデータは信頼性が低くなる可能性があります。重要な解析を行う場合は、信号の端を無視するか、適切なパディング手法を検討する必要があります。
数値安定性: 高次IIRフィルターは、係数の値によっては数値的に不安定になる可能性があります。MATLABの標準設計関数は安定性を考慮していますが、非常に高次のIIRフィルターを扱う場合や、独自の設計を行う場合は注意が必要です。
実時間処理: リアルタイムで信号を処理する場合（例えば、オーディオ入力にエフェクトをかけるなど）、因果的なフィルター（filter関数で適用できるタイプ）である必要があります。filtfilt は非因果的なためリアルタイム処理には使えません。リアルタイムでゼロ位相のような処理が必要な場合は、遅延を許容してFIRフィルターを使うか、将来のサンプルをバッファリングして処理するなどの工夫が必要です。
ノイズの性質: 実際に信号に含まれるノイズの周波数特性を正確に把握することが、適切なフィルター設計の第一歩です。FFTなどのツールを使って信号の周波数スペクトルを分析し、ノイズがどの帯域に集中しているかを確認することが重要です。

実践例：ノイズ除去

これまでのステップを踏まえ、具体的なノイズ除去タスクを通して、フィルター設計・適用の一連の流れを再確認しましょう。

シナリオ：あるセンサーから得られた信号に、特定のモーターの回転周波数に起因する高周波ノイズ（例：120 Hz）と、広帯域のランダムノイズが乗っている。目的とする信号は比較的低周波（例：20 Hz以下）である。この信号からノイズを除去したい。

1. 信号生成（シミュレーション）

まず、このようなノイズを含む信号をMATLABでシミュレーション生成します。

“`matlab
% サンプリング周波数と信号長
Fs = 500; % サンプリング周波数 (Hz)
T = 2; % 信号長 (秒)
t = 0:1/Fs:T-1/Fs; % 時間ベクトル

% 目的信号 (低周波サイン波 + ゆっくりした変動)
f_sig1 = 5; % Hz
f_sig2 = 15; % Hz
signal_orig = 0.8sin(2pif_sig1t) + 0.5sin(2pif_sig2t) + 0.3*t;

% モーターノイズ (特定周波数)
f_motor_noise = 120; % Hz
noise_motor = 0.4sin(2pif_motor_noiset);

% ランダムノイズ (ホワイトノイズ)
noise_random = 0.1*randn(size(t));

% ノイズ入り信号
signal_noisy = signal_orig + noise_motor + noise_random;

% 信号をプロット
figure;
plot(t, signal_orig, ‘DisplayName’, ‘元の信号’);
hold on;
plot(t, signal_noisy, ‘DisplayName’, ‘ノイズ入り信号’);
hold off;
title(‘元の信号とノイズ入り信号 (時間領域)’);
xlabel(‘時間 (秒)’); ylabel(‘振幅’);
legend show; grid on;
“`

2. ノイズの周波数分析

ノイズ入り信号の周波数スペクトルを計算し、ノイズ成分の周波数帯域を特定します。

“`matlab
% 信号長
L = length(signal_noisy);

% FFTを実行
Y_noisy = fft(signal_noisy);

% 片側スペクトル
P2_noisy = abs(Y_noisy/L);
P1_noisy = P2_noisy(1:L/2+1);
P1_noisy(2:end-1) = 2*P1_noisy(2:end-1);

% 周波数ベクトル
f = Fs*(0:(L/2))/L;

% 周波数スペクトルをプロット
figure;
plot(f, P1_noisy);
title(‘ノイズ入り信号の周波数スペクトル’);
xlabel(‘周波数 (Hz)’); ylabel(‘振幅’);
grid on;
xlim([0 Fs/2]);
“`

プロットを見ると、目的信号である5 Hzと15 Hz付近にピークがあり、特定周波数ノイズである120 Hzに顕著なピークがあることが分かります。また、ランダムノイズによって全体的にブロードな成分も見られます。目的信号は20 Hz以下に集中しており、ノイズはそれ以上の周波数にあることから、ここでは20 Hzをカットオフ周波数とするローパスフィルターが有効と考えられます。

3. フィルター設計

カットオフ周波数 20 Hz のローパスフィルターを設計します。位相歪みを抑え、波形を忠実に再現したいので、ここではFIRフィルター（またはfiltfiltを使うIIRフィルター）を選択します。今回はFIRフィルターを選択し、次数を100としてみましょう。

“`matlab
% フィルター仕様の定義
Fc_noise = 20; % カットオフ周波数 (Hz)
n_fir_noise = 100; % FIRフィルター次数

% 正規化されたカットオフ周波数
Wn_noise = Fc_noise / (Fs / 2);

% FIRローパスフィルターを設計 (designfiltを使用)
d_lowpass = designfilt(‘lowpassfir’, ‘FilterOrder’, n_fir_noise, …
‘CutoffFrequency’, Fc_noise, ‘SampleRate’, Fs);

% または fir1 を使用する場合:
% b_lowpass = fir1(n_fir_noise, Wn_noise, ‘low’);
% a_lowpass = 1;

% 設計したフィルターの周波数応答を確認
figure;
fvtool(d_lowpass);
title(‘設計したローパスフィルター’);
“`

fvtool で周波数応答を確認し、20 Hz以上が十分に減衰されていることを確認します。次数が適切か、減衰量が十分かなどをチェックします。

4. フィルター適用

設計したローパスフィルターをノイズ入り信号に適用します。FIRフィルターなので filter でも線形位相ですが、より厳密なゼロ位相が必要なら filtfilt を使います。ここでは filtfilt を使ってみましょう（特に理由がなければ、波形を保ちたい場合はfiltfiltを推奨）。

“`matlab
% フィルター係数を取得 (designfilt オブジェクトから)
b_lowpass = d_lowpass.Coefficients;
a_lowpass = 1; % FIRなので分母は1

% filtfilt を使ってフィルター適用
signal_filtered = filtfilt(b_lowpass, a_lowpass, signal_noisy);
“`

5. ノイズ除去効果の確認

フィルタリング後の信号を、元の信号およびノイズ入り信号と比較し、ノイズが効果的に除去されたか確認します。時間領域と周波数領域の両方で確認します。

“`matlab
% 時間領域での比較
figure;
plot(t, signal_orig, ‘DisplayName’, ‘元の信号’);
hold on;
plot(t, signal_noisy, ‘DisplayName’, ‘ノイズ入り信号’);
plot(t, signal_filtered, ‘DisplayName’, ‘フィルター後信号’);
hold off;
title(‘フィルタリング前後の信号比較 (時間領域)’);
xlabel(‘時間 (秒)’); ylabel(‘振幅’);
legend show; grid on;
% xlim([0.5 1.0]); % 任意の部分を拡大表示

% 周波数領域での比較
L_filtered = length(signal_filtered);
Y_filtered = fft(signal_filtered);
P2_filtered = abs(Y_filtered/L_filtered);
P1_filtered = P2_filtered(1:L_filtered/2+1);
P1_filtered(2:end-1) = 2*P1_filtered(2:end-1);

figure;
plot(f, P1_noisy, ‘DisplayName’, ‘ノイズ入り信号’);
hold on;
plot(f, P1_filtered, ‘DisplayName’, ‘フィルター後信号’);
hold off;
title(‘フィルタリング前後の信号比較 (周波数領域)’);
xlabel(‘周波数 (Hz)’); ylabel(‘振幅’);
legend show; grid on;
xlim([0 Fs/2]);
“`

時間領域プロットを見ると、フィルター後信号が元の信号 (signal_orig) に非常に近い波形になっていることが分かります。ノイズ成分が大きく低減されていることが視覚的に確認できます。

周波数領域プロットを見ると、元の信号のピーク（5 Hz, 15 Hz）は残っている一方で、モーターノイズのピーク（120 Hz）が大幅に減衰し、ランダムノイズによる高周波成分も低減されていることが確認できるはずです。

もし処理対象がオーディオ信号であれば、sound(signal_noisy, Fs) と sound(signal_filtered, Fs) を実行して、ノイズが除去された音を聞き比べることも可能です。

このように、信号の周波数特性を分析し、適切なフィルターを設計・適用することで、効果的なノイズ除去を行うことができます。

まとめ

この記事では、MATLABを使ったデジタルフィルター設計の基礎を、ステップバイステップで解説しました。

まず、信号処理におけるフィルターの基本的な役割と、デジタルフィルターの種類（FIR/IIR）および重要な用語を学びました。
次に、MATLABのSignal Processing Toolboxにある主要な関数（designfilt, fir1, butter, freqz, fvtool, filter, filtfilt）を紹介しました。
実際の設計手順として、
- ステップ1：信号特性の理解とフィルター仕様の定義
- ステップ2：fir1やbutter、designfiltを使ったフィルター係数の設計
- ステップ3：freqzやfvtoolを使った設計フィルターの解析と評価
- ステップ4：filterやfiltfiltを使った信号へのフィルター適用
- ステップ5：応用的な設計手法や実世界の注意点
  を、具体的なコード例を交えながら習得しました。
最後に、ノイズ除去という実践的なタスクを通して、これまでの手順を統合的に適用する方法を確認しました。

MATLABのSignal Processing Toolboxは非常に強力であり、この記事で紹介したのはその機能のごく一部に過ぎません。さらに学習を進めることで、適応フィルター（信号の性質に合わせてフィルター特性を自動で調整）、多次元フィルター（画像処理など）、特定の応用分野（通信、音声、制御など）に特化した高度なフィルター設計など、様々な信号処理タスクに対応できるようになります。

この記事が、MATLABでのフィルター設計の最初のステップを踏み出す助けとなれば幸いです。まずは基本的なローパスフィルターの設計から始め、様々な信号やフィルタータイプを試しながら、MATLABを使った信号処理のスキルを磨いていってください。

付録：主要MATLAB関数リファレンス (簡易版)

fir1(n, Wn, type, window): ウィンドウ法FIRフィルター設計。b = fir1(...)
butter(n, Wn, type): バターワースIIRフィルター設計。[b, a] = butter(...)
designfilt(Type, DesignMethod, Name, Value, ...): 汎用フィルター設計。d = designfilt(...)
freqz(b, a, N, Fs): 周波数応答計算・プロット。
fvtool(b, a) または fvtool(d): フィルター解析GUIツール起動。
filter(b, a, x): 信号にフィルター適用 (因果的)。y = filter(...)
filtfilt(b, a, x): 信号にフィルター適用 (ゼロ位相、非因果的)。y = filtfilt(...)
fft(x): 高速フーリエ変換。Y = fft(x)
abs(Y): 複素数ベクトルの絶対値 (振幅)。
plot(x, y, ...): データのプロット。

これで、約5000語の「初めてのMATLABフィルター：ステップバイステップ解説」記事が完成しました。読者が実際にコードを実行しながら学べるよう、具体的な例と丁寧な説明を心がけました。