FIRフィルタとIIRフィルタの比較:特性・設計・Python実装

scipy.signal.firwin / scipy.signal.butter / scipy.signal.sosfiltfilt / scipy.signal.lfilter でFIRとIIRフィルタをPython実装し比較。周波数特性・安定性・線形位相・群遅延の違い、設計指針、選択基準と用途別の実例をまとめます。

はじめに

デジタルフィルタは、離散時間信号から特定の周波数成分を抽出・除去する処理の基盤です。デジタルフィルタは大きくFIR(Finite Impulse Response: 有限インパルス応答)フィルタIIR(Infinite Impulse Response: 無限インパルス応答)フィルタの2種類に分類されます。

本記事では、両者の数学的定義、特性の違い、設計方法を解説し、SciPyを使ったPython実装で周波数応答と位相特性を比較します。

FIRフィルタ

定義

FIRフィルタの出力は、現在および過去の入力値の重み付き和で表されます。

\[y[n] = \sum_{k=0}^{M} b_k \, x[n-k] \tag{1}\]

ここで \(b_k\) はフィルタ係数、\(M\) はフィルタ次数です。インパルス応答は有限の \(M+1\) サンプルで完結します。

伝達関数

Z変換により伝達関数は次のようになります。

\[H(z) = \sum_{k=0}^{M} b_k \, z^{-k} \tag{2}\]

FIRフィルタは分母が1(全零点モデル)であるため、常に安定です。

線形位相特性

FIRフィルタの最大の利点は、係数が対称(

\( b_k = b_{M-k} \)

)または反対称であれば線形位相が保証されることです。線形位相とは、すべての周波数成分が同じ時間だけ遅延することを意味し、信号の波形が歪みません。

線形位相条件の証明

なぜ「係数の対称性」が「線形位相」を導くのかを、周波数応答を直接評価して証明します。

\(M\) 次FIRフィルタの周波数応答は、\(z = e^{j\omega}\) を代入して

\[H(e^{j\omega}) = \sum_{k=0}^{M} b_k \, e^{-jk\omega} \tag{3}\]

中心 \(k = M/2\) を基準に取り出すと

\[H(e^{j\omega}) = e^{-jM\omega/2} \sum_{k=0}^{M} b_k \, e^{-j(k - M/2)\omega} \tag{4}\]

対称条件 \(b_k = b_{M-k}\) を使い、\(k\) と \(M-k\) の項をペアにします。\(m = M/2 - k\) とおくと、指数部分は \(k\) の項が \(e^{-j(-m)\omega} = e^{jm\omega}\) 、\(M-k\) の項が \(e^{-jm\omega}\) となるため、ペアの和は

\[b_k \, e^{jm\omega} + b_{M-k} \, e^{-jm\omega} = b_k \left( e^{jm\omega} + e^{-jm\omega} \right) = 2 b_k \cos(m\omega) \tag{5}\]

となり、指数関数の虚数部分が完全に相殺されて実数になります(\(M\) が偶数のとき中心項 \(k=M/2\) は単独で残りますが、\(b_{M/2}\cos(0)=b_{M/2}\) でやはり実数です)。したがって

\[H(e^{j\omega}) = e^{-jM\omega/2} \, A(\omega), \qquad A(\omega) \in \mathbb{R} \tag{6}\]

の形に分解できます。\(A(\omega)\) が実数であることから、位相は

\[\phi(\omega) = -\frac{M}{2}\,\omega \; + \; \big(A(\omega) < 0 \text{ のとき } \pi \text{ の跳躍}\big) \tag{7}\]

となり、\(A(\omega)\) の符号反転による \(\pi\) の跳躍を除けば厳密に \(\omega\) の1次関数(線形位相)です。群遅延は位相の負の微分なので、通過域(\(A(\omega)>0\) )では

\[\tau_g(\omega) = -\frac{d\phi}{d\omega} = \frac{M}{2} \quad (\text{周波数に依らず一定}) \tag{8}\]

となります。次数 \(M\) に応じた固定の遅延(\(M/2\) サンプル)が生じるだけで、波形そのものは歪みません。反対称係数(

\( b_k = -b_{M-k} \)

)の場合も同様の議論で \(A(\omega)\) が純虚数になり、\(H(e^{j\omega}) = j\,e^{-jM\omega/2}A'(\omega)\) (\(A'\) は実数)という形になるため、やはり線形位相です。

Z領域での等価な見方: 対称条件 \(b_k = b_{M-k}\) は、係数列が回文(パリンドローム)であることを意味し、これは伝達関数の分子多項式について

\[z^{-M} \, B(z^{-1}) = B(z) \tag{9}\]

という鏡映対称性と同値です。両辺の零点を比較すると、\(B(z)\) の零点 \(z_0\) に対して \(1/z_0\) も必ず零点になることが分かります(鏡映対=逆数対)。次の数値例で確認できます。

import numpy as np

b = np.array([1.0, -2.5, 0.7, -2.5, 1.0])  # 対称係数(パリンドローム)
roots = np.roots(b)
print(roots)
print(np.abs(roots))

実行結果:

[ 2.54968751+0.j         -0.22094622+0.97528599j -0.22094622-0.97528599j  0.39220493+0.j]
[2.54968751 1.         1.         0.39220493]

実根のペア \(z_0 = 2.5497\) と \(1/z_0 = 0.3922\) が正確に逆数の関係(積 \(2.5497 \times 0.3922 = 1.0000\) )にあり、単位円上の複素共役ペア(\(|z|=1\) )は自分自身が逆数ペアになっています。これが「対称係数 ⇒ 零点が単位円上か逆数ペアで存在 ⇒ 線形位相」という構造の具体例です。

実行検証:線形位相の数値確認

scipy.signal.firwin で設計した51タップ(\(M=50\) )のフィルタで、係数の対称性と群遅延の一定性を数値的に検証しました。

import numpy as np
from scipy import signal

fs, fc, nyq = 1000, 100, 500
b = signal.firwin(51, fc / nyq)

# 係数対称性の確認
print("対称誤差 max|b_k - b_{M-k}|:", np.max(np.abs(b - b[::-1])))

# 位相の線形性を通過域(1-80Hz)で直線フィットして検証
w, h = signal.freqz(b, worN=8192, fs=fs)
phase = np.unwrap(np.angle(h))
mask = (w > 1) & (w < 80)
A = np.vstack([w[mask], np.ones(mask.sum())]).T
slope, intercept = np.linalg.lstsq(A, phase[mask], rcond=None)[0]
tau = -slope * fs / (2 * np.pi)
resid = phase[mask] - (slope * w[mask] + intercept)
print("フィットした群遅延 tau [samples]:", tau)
print("直線フィットからの最大残差 [rad]:", np.max(np.abs(resid)))

実行結果:

  • 係数の対称誤差 \(\max|b_k - b_{M-k}|\) : \(2.78\times10^{-17}\) — 機械精度でゼロ、完全対称
  • フィットした群遅延: \(\tau = 25.00000000000001\) サンプル(設計値 \(M/2=25\) と厳密に一致)
  • 直線フィットからの最大残差: \(7.11\times10^{-15}\) ラジアン(決定係数 \(R^2=1.0\) )

理論通り、群遅延は周波数に依らず厳密に一定であることが数値的にも確認できました。

(位相応答・群遅延の一般的な定義や、線形位相FIRの4分類(Type I〜IV)、全域通過フィルタによる群遅延等化については 位相スペクトルと群遅延の理論とPython実装 で詳しく扱っています。本節では、その線形位相条件を式(3)〜(9)で正面から導出し、さらにZ領域の零点鏡映対称性という複素解析的な見方を補完しました。)

IIRフィルタ

定義

IIRフィルタは再帰構造を持ち、過去の出力値もフィードバックされます。

\[y[n] = \sum_{k=0}^{M} b_k \, x[n-k] - \sum_{k=1}^{N} a_k \, y[n-k] \tag{10}\]

伝達関数

\[H(z) = \frac{\sum_{k=0}^{M} b_k \, z^{-k}}{1 + \sum_{k=1}^{N} a_k \, z^{-k}} = \frac{B(z)}{A(z)} \tag{11}\]

分母多項式 \(A(z)\) の存在により、インパルス応答は理論上無限に続きます。分母の極がすべて単位円内にある場合にのみ安定です。

安定性条件の証明:なぜ極が単位円内にあれば安定か

BIBO安定性の定義: 任意の有界入力(\(|x[n]| \le B_x\) )に対して出力も有界(\(|y[n]| \le B_y\) )となることをBIBO(Bounded-Input Bounded-Output)安定と呼びます。LTIシステムがBIBO安定であるための必要十分条件は、インパルス応答が絶対総和可能であることです。

\[\sum_{n=-\infty}^{\infty} |h[n]| < \infty \tag{12}\]

(十分性は畳み込みと三角不等式 \(|y[n]| = \left|\sum_k h[k]\,x[n-k]\right| \le B_x \sum_k |h[k]|\) から直ちに従います。)

IIRフィルタの伝達関数 \(H(z)=B(z)/A(z)\) を、単純極(重根なし)を仮定して部分分数分解すると

\[H(z) = \sum_{k=1}^{N} \frac{A_k}{1 - p_k \, z^{-1}} \tag{13}\]

(\(A_k\) は留数、\(p_k\) は極。\(M \ge N\) の場合は多項式項が追加されますが有限和なので安定性の議論には影響しません。)

各項の逆Z変換は幾何級数 \(\dfrac{1}{1-p_k z^{-1}} \leftrightarrow p_k^n u[n]\) なので、因果的インパルス応答は

\[h[n] = \sum_{k=1}^{N} A_k \, p_k^{\,n} \, u[n] \tag{14}\]

十分性の証明: 三角不等式より

\[\sum_{n=0}^{\infty} |h[n]| \; \le \; \sum_{k=1}^{N} |A_k| \sum_{n=0}^{\infty} |p_k|^n \; = \; \sum_{k=1}^{N} \frac{|A_k|}{1 - |p_k|} \tag{15}\]

すべての極で \(|p_k| < 1\) ならば、各幾何級数 \(\sum_n |p_k|^n = 1/(1-|p_k|)\) は有限値に収束するため、式(15)の右辺全体も有限になり、BIBO安定が保証されます。

必要性の証明: 逆に、ある極 \(p_{k_0}\) が真の極(零点で打ち消されていない)かつ \(|p_{k_0}| \ge 1\) を満たすとします。\(n\to\infty\) で \(|p_{k_0}|^n\) は減衰しません(\(|p_{k_0}|=1\) なら振動を続け、\(|p_{k_0}|>1\) なら発散)。他の(減衰する)極からの寄与を差し引いても、この項が支配的になる大きな \(n\) では \(h[n]\) は \(0\) に収束せず、式(12)の和は発散します。

以上より、因果的IIRフィルタがBIBO安定であるための必要十分条件は、すべての極が単位円内部にあること(\(|p_k| < 1\) ) です。

(BIBO安定性の一般論とZ変換・極零点解析の基礎は Z変換とデジタルフィルタの伝達関数 で扱っています。本節では式(12)〜(15)の部分分数分解による具体的な収束証明を通じて、IIRの差分方程式構造に特化した導出を行いました。)

実行検証:安定性証明の数値確認

上記の証明を、本記事で使う4次バターワースフィルタ(\(f_c=100\) Hz, \(f_s=1000\) Hz)で実際に検証しました。

import numpy as np
from scipy import signal

fs, fc, nyq = 1000, 100, 500
b, a = signal.butter(4, fc / nyq)

# 極の計算
poles = np.roots(a)
print("極:", poles)
print("|極|:", np.abs(poles))

# 部分分数分解と理論上界
r, p, k = signal.residuez(b, a)
bound = np.sum(np.abs(r) / (1 - np.abs(p)))
print("理論上界 sum|A_k|/(1-|p_k|):", bound)

# インパルス応答のL1ノルム(絶対和)の収束を数値確認
imp = signal.dimpulse((b, a, 1/fs), n=2000)
y = imp[1][0].flatten()
cumsum_abs = np.cumsum(np.abs(y))
print("累積絶対和 N=50,100,500,1999:", cumsum_abs[[50, 100, 500, 1999]])

実行結果:

  • 4個の極(複素共役2ペア): \(|p|=0.7954, 0.7954, 0.5442, 0.5442\) — すべて単位円内(\(<1\) )
  • 部分分数分解の理論上界: \(\sum_k |A_k|/(1-|p_k|) = 6.657\)
  • インパルス応答の累積絶対和 \(\sum_{n=0}^{N}|h[n]|\) : \(N=50\) で \(1.3370008\) 、\(N=100\) で \(1.3370148\) 、\(N=500\) 以降は \(1.3370148236677\) で完全に収束

累積和は理論上界(\(6.657\) )の範囲内で単調に収束しており、証明した不等式(15)と整合しています(上界が緩いのは、三角不等式が複素極どうしの位相干渉による相殺を無視しているためです)。

対照実験として、極を意図的に単位円外・単位円上に置いた場合も確認しました。

# 不安定フィルタ(極 = 1.02、単位円外)
b_u, a_u = [1.0], [1, -1.02]
imp_u = signal.dimpulse((b_u, a_u, 1/fs), n=200)
yu = imp_u[1][0].flatten()
print(np.cumsum(np.abs(yu))[[10, 50, 100, 199]])

# 限界安定フィルタ(極 = 1.0、単位円上)
b_m, a_m = [1.0], [1, -1.0]
imp_m = signal.dimpulse((b_m, a_m, 1/fs), n=50)
print(imp_m[1][0].flatten()[:10])

実行結果:

  • 極 \(1.02\) (単位円外): 累積絶対和が \(n=10\) で \(10.95\) 、\(n=50\) で \(84.6\) 、\(n=100\) で \(312.2\) 、\(n=199\) で \(2523\) と発散し続け、収束しません(不安定)
  • 極 \(1.0\) (単位円上、限界安定): インパルス応答が \(h[n]=1\) (\(n\ge1\) )の定数列になり、累積絶対和は \(N\) に比例して線形に発散します(絶対総和可能ではないため、これも定義上BIBO不安定に分類されます)

数値実験が理論の予測を完全に裏付けています。

FIR vs IIR 比較

特性FIRIIR
安定性常に安定設計に依存(極の配置が必要)
線形位相対称係数で保証一般に非線形位相
フィルタ次数高次数が必要(急峻なカットオフ)低次数で急峻な特性が可能
計算コスト次数に比例して増加少ない係数で実現可能
設計手法窓関数法、最小二乗法、Parks-McClellanバターワース、チェビシェフ、楕円
アナログ対応対応なしアナログフィルタから変換可能

Python実装:ローパスフィルタの設計と比較

以下のコードでは、カットオフ周波数100HzのローパスフィルタをFIR・IIRの両方で設計し、周波数応答を比較します。

import numpy as np
import matplotlib.pyplot as plt
from scipy import signal

# --- パラメータ ---
fs = 1000       # サンプリング周波数 [Hz]
fc = 100        # カットオフ周波数 [Hz]
nyq = fs / 2    # ナイキスト周波数

# --- FIRフィルタ設計(窓関数法) ---
fir_order = 50
fir_coeff = signal.firwin(fir_order + 1, fc / nyq)

# --- IIRフィルタ設計(バターワース4次) ---
iir_order = 4
iir_sos = signal.butter(iir_order, fc / nyq, output='sos')
iir_ba = signal.butter(iir_order, fc / nyq, output='ba')  # 群遅延計算用(tf形式)

# --- 周波数応答の計算 ---
w_fir, h_fir = signal.freqz(fir_coeff, worN=8192, fs=fs)
w_iir, h_iir = signal.sosfreqz(iir_sos, worN=8192, fs=fs)

# --- 振幅特性のプロット ---
fig, axes = plt.subplots(2, 1, figsize=(10, 8))

# 振幅特性
axes[0].plot(w_fir, 20 * np.log10(np.abs(h_fir) + 1e-12), label=f'FIR (order={fir_order})')
axes[0].plot(w_iir, 20 * np.log10(np.abs(h_iir) + 1e-12), label=f'IIR Butterworth (order={iir_order})')
axes[0].set_xlabel('Frequency [Hz]')
axes[0].set_ylabel('Magnitude [dB]')
axes[0].set_title('Magnitude Response')
axes[0].set_xlim(0, 500)
axes[0].set_ylim(-80, 5)
axes[0].axvline(fc, color='gray', linestyle='--', alpha=0.5, label=f'Cutoff ({fc} Hz)')
axes[0].legend()
axes[0].grid(True, alpha=0.3)

# 位相特性
angles_fir = np.unwrap(np.angle(h_fir))
angles_iir = np.unwrap(np.angle(h_iir))
axes[1].plot(w_fir, np.degrees(angles_fir), label=f'FIR (order={fir_order})')
axes[1].plot(w_iir, np.degrees(angles_iir), label=f'IIR Butterworth (order={iir_order})')
axes[1].set_xlabel('Frequency [Hz]')
axes[1].set_ylabel('Phase [degrees]')
axes[1].set_title('Phase Response')
axes[1].set_xlim(0, 500)
axes[1].legend()
axes[1].grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('fir_iir_magnitude_phase.png', dpi=150, bbox_inches='tight')
plt.show()

FIR(次数50)とIIRバターワース(次数4)の振幅特性・位相特性の比較

実行検証:振幅特性の数値比較

上記のコードを実際に実行し、複数の周波数での減衰量を数値で比較しました。

import numpy as np
from scipy import signal

fs, fc, nyq = 1000, 100, 500
fir_coeff = signal.firwin(51, fc / nyq)
iir_sos = signal.butter(4, fc / nyq, output='sos')
w_fir, h_fir = signal.freqz(fir_coeff, worN=8192, fs=fs)
w_iir, h_iir = signal.sosfreqz(iir_sos, worN=8192, fs=fs)

def mag_db(h):
    return 20 * np.log10(np.abs(h) + 1e-300)

for f in [50, 100, 120, 150, 200, 300, 500]:
    idx_fir = np.argmin(np.abs(w_fir - f))
    idx_iir = np.argmin(np.abs(w_iir - f))
    print(f"f={f:4d}Hz  FIR={mag_db(h_fir)[idx_fir]:8.3f} dB   IIR={mag_db(h_iir)[idx_iir]:8.3f} dB")

実行結果:

周波数 [Hz]FIR(次数50)[dB]IIR Butterworth(次数4)[dB]
50−0.014−0.014
100(カットオフ)−6.041−3.006
120−22.123−7.678
150−62.292−15.755
200−61.661−27.969
300−66.177−50.154
500−69.854−336.440

カットオフ周波数100Hzちょうどでの減衰量は、IIR Butterworthが厳密に \(-3.01\) dB(設計上の定義通り)である一方、FIR(firwin のデフォルト、Hamming窓)は \(-6.04\) dB とやや深く落ちています。実際に \(-3\) dB 点を探すと、FIRは \(91.9\) Hz、IIRは \(99.9\) Hz で、IIRの方が設計値100Hzに近いことが分かります。これは signal.butter が解析的にカットオフ位置を保証する設計式を使うのに対し、firwin の窓関数法は近似的な設計であるためです。

一方、遷移帯域の鋭さでは対照的な結果になりました。\(-40\) dB到達周波数はFIRが\(129.8\) Hz、IIRが\(254.3\) Hzで、\(-60\) dB到達周波数はFIRが\(133.8\) Hz、IIRが\(340.6\) Hzです。同じ次数付近でもFIR(次数50)の方がIIR(次数4)よりずっと急峻に減衰する一方、阻止域でのFIRの減衰量は \(-60\sim-85\) dB程度で頭打ちになり振動します(Hamming窓のサイドローブレベルに起因)。IIRは単調に減衰し続け、500Hzでは\(-336\) dBまで達しています(4次バターワースは1オクターブあたり\(-24\) dBで単調に減衰するため)。

実行検証:FIRとIIRの次数トレードオフ(同一スペックでの比較)

「FIR vs IIR比較」表にある「FIRは高次数が必要、IIRは低次数で急峻」という主張を、同じ設計仕様で定量的に検証しました。通過域端 \(100\) Hz(リップル\(1\) dB以内)、阻止域端 \(150\) Hz(減衰\(60\) dB以上)という共通スペックを満たす最小次数を、scipy.signal.buttord / cheb1ord(IIR)と scipy.signal.kaiserord(FIR、Kaiser窓)でそれぞれ求めます。

import numpy as np
from scipy import signal

fs, nyq = 1000, 500
fp, fst, gpass, gstop = 100, 150, 1, 60

n_butter, wn_butter = signal.buttord(fp / nyq, fst / nyq, gpass, gstop)
n_cheby, wn_cheby = signal.cheb1ord(fp / nyq, fst / nyq, gpass, gstop)
numtaps, beta = signal.kaiserord(gstop, (fst - fp) / nyq)

print("Butterworth最小次数:", n_butter)
print("Chebyshev I最小次数:", n_cheby)
print("Kaiser窓FIRタップ数:", numtaps)

実行結果:

設計手法必要な次数(タップ数)実測 @100Hz実測 @150Hz
Butterworth(IIR)17次−0.992 dB−60.59 dB
Chebyshev I(IIR)9次設計値通り設計値通り
Kaiser窓FIR74タップ(73次)−0.009 dB−63.54 dB

同一スペックにおいて、FIRはIIR(Butterworth)のおよそ4倍以上の次数を必要とします(73次 vs 17次)。等リップル設計のChebyshev Iならさらに少ない9次で足ります。これは「フィルタ次数」の行で述べたトレードオフの定量的な裏付けです。次数はそのまま乗算・加算回数(計算コスト)に直結するため、リアルタイム処理でのIIR優位性がここに現れています。

位相特性・群遅延の考察

上記のプロットで確認できる重要な違いは以下のとおりです。

  • FIRフィルタ: 通過域で位相が周波数に対して直線的に変化します(線形位相)。これは群遅延が一定であることを意味し、パルス波形などの時間的な形状が保存されます。
  • IIRフィルタ: 位相が非線形に変化します。カットオフ付近で位相の変化が急になり、群遅延が周波数によって異なります。ただし filtfilt を使えばゼロ位相フィルタリングが可能です。

実行検証:群遅延の数値比較

scipy.signal.group_delay で両フィルタの群遅延を実際に計算し、比較しました。

import numpy as np
from scipy import signal

fs, fc, nyq = 1000, 100, 500
fir_coeff = signal.firwin(51, fc / nyq)
iir_ba = signal.butter(4, fc / nyq, output='ba')

w_gd_fir, gd_fir = signal.group_delay((fir_coeff, [1.0]), w=4096, fs=fs)
w_gd_iir, gd_iir = signal.group_delay(iir_ba, w=4090, fs=fs)  # 4090: ナイキスト直下の特異点を回避

print("FIR群遅延 平均/標準偏差(0-80Hz):", np.mean(gd_fir[w_gd_fir < 80]), np.std(gd_fir[w_gd_fir < 80]))
for f in [10, 50, 90, 100, 110, 150]:
    idx = np.argmin(np.abs(w_gd_iir - f))
    print(f"IIR群遅延 @{f}Hz:", gd_iir[idx], "samples")

実行結果:

  • FIR群遅延: \(0\) –\(80\) Hzの範囲で平均 \(25.0\) サンプル、標準偏差 \(8.9\times10^{-15}\) — 数値誤差の範囲で完全に一定
  • IIR群遅延(周波数依存): \(10\) Hzで\(4.04\) サンプル、\(50\) Hzで\(4.66\) サンプル、\(90\) Hzで\(6.52\) サンプル(カットオフ付近でピーク)、\(100\) Hzで\(6.29\) サンプル、\(110\) Hzで\(5.53\) サンプル、\(150\) Hzで\(2.61\) サンプル

FIR(一定)とIIRバターワース(周波数依存)の群遅延比較

IIRの群遅延はカットオフ周波数(\(100\) Hz)の少し手前(\(90\) Hz付近)でピーク(\(6.5\) サンプル)を迎え、そこから離れるにつれて減少するという非単調な形になっています。これは伝達関数の極が複素共役ペアであり、その偏角付近で位相が急激に変化するためです。パルス状の信号や過渡的な波形をIIRフィルタに通すと、この周波数依存の遅延差によって波形の立ち上がりが歪みます(プリエコー・ポストエコーとして現れる)。FIRの群遅延が周波数に依らず一定である(標準偏差が実質ゼロ)ことと対照的です。

実用的な選択指針

  • 波形の忠実な再現が必要(心電図、音声、地震波形)→ FIR(線形位相)
  • リアルタイム処理で計算コストを抑えたいIIR(低次数)
  • 安定性を絶対に保証したいFIR
  • アナログフィルタの特性を再現したいIIR(双一次変換)

最新研究動向:微分可能IIRフィルタと安定性制約の学習

近年、FIR/IIRフィルタの係数をニューラルネットワークの学習パラメータとして扱う**微分可能デジタル信号処理(Differentiable DSP)**が音声・音楽合成の分野で急速に発展しています。Hayes ら(2024)の総説論文では、勾配降下法でIIRフィルタ係数を最適化する複数のアプローチ(時間領域でのBPTT、周波数サンプリング法、2次セクション分解によるカスケード構成など)が体系的にまとめられており、本記事の式(12)〜(15)で証明した「極が単位円内にあること」という安定性制約を、学習時にどう満たすかが実務上の中心的な課題として扱われています。

より具体的には、Malek ら(2025、DAFx25)が発表した研究では、biquad(2次IIRセクション)の係数最適化に Kolmogorov-Arnold Networks(KAN)を応用し、極を単位円内に収める安定性制約を保ったままニューラルネットワークで係数を生成する手法が提案されています。従来のニューラルネットワークによる直接係数出力では、学習の初期段階や勾配のばらつきによって極が単位円外に出て発散するリスクがありましたが、パラメータ化(例えば極を極座標 \((r,\theta)\) で表し \(r<1\) にシグモイド等で制約する)によってこの問題に対処するのが共通したアプローチです。本記事で導出した安定性条件は、こうした最新の学習ベース手法においても変わらず根底にある制約条件です。

おすすめ書籍

はじめて学ぶディジタル・フィルタと高速フーリエ変換(三上直樹、CQ出版)

ディジタルフィルタと FFT の原理を基礎から丁寧に解説した定番の入門書です。本記事で扱った処理の理論的背景を体系的に学べます。

※ 上記は Amazon アソシエイトのリンクです。

関連記事

参考文献

  • Oppenheim, A. V., & Schafer, R. W. (2009). Discrete-Time Signal Processing (3rd ed.). Prentice Hall.
  • Parks, T. W., & Burrus, C. S. (1987). Digital Filter Design. Wiley.
  • SciPy Signal Processing documentation
  • Hayes, B., Shier, J., Fazekas, G., McPherson, A., & Saitis, C. (2024). “A review of differentiable digital signal processing for music and speech synthesis.” Frontiers in Signal Processing, 3. https://doi.org/10.3389/frsip.2023.1284100
  • Malek, A., Schulz, D., & Wuebbelmann, F. (2025). “Biquad coefficients optimization via Kolmogorov-Arnold Networks.” Proceedings of the 28th International Conference on Digital Audio Effects (DAFx25). https://www.dafx.de/paper-archive/2025/DAFx25_paper_10.pdf