株価チャートを見たとき、日々の値動きのギザギザではなく「25日移動平均線」に注目した経験はないでしょうか。センサーで温度や加速度を計測すると、どんなに精密な機器でも細かなノイズが混入します。このような「ギザギザ」を滑らかに除去したいとき、まず最初に手を伸ばすべきツールが移動平均フィルタです。
移動平均フィルタは、デジタル信号処理で最もシンプルなFIR（有限インパルス応答）ローパスフィルタです。設計パラメータは「窓長 $M$」のたった1つ。それでもなお、このフィルタは驚くほど多くの場面で活躍します。
移動平均フィルタが登場する具体的な場面を2つ挙げます。
- IoTセンサーのノイズ平滑化: 温度センサー・IMU（慣性計測ユニット）・圧力センサーの生データには高周波ノイズが含まれます。移動平均は計算量がほぼゼロで、マイコンのリアルタイム処理に組み込めます。
- 金融時系列の平滑化: 株価・為替・出来高の移動平均線は、短期的な変動を除いてトレンドを把握するために世界中で使われています。
本記事では以下の内容を解説します。
- 数学的定義 $y[n] = \frac{1}{M}\sum_{k=0}^{M-1} x[n-k]$ の意味
- インパルス応答が「矩形窓」になる理由
- 周波数応答 $H(e^{j\omega})$ の厳密な導出（sinc 形状）
- 窓長 $M$ と遮断周波数・平滑度のトレードオフ
- 線形位相特性と群遅延 $(M-1)/2$
- 端効果（エッジ効果）とその対処法
- 累積和による O(1) 更新の高速化
- 単純移動平均（SMA）と指数移動平均（EMA）の違い
- Python での実装（スクラッチ・NumPy・SciPy）
前提知識
以下の記事を先に読んでおくと理解が深まります。
移動平均フィルタとは何か
上の図が移動平均フィルタの本質を表しています。左のパネルのギザギザした信号（グレー）に対して、移動平均を適用すると右のパネルのように滑らかな曲線（赤）が得られます。青の破線が「真の信号」ですが、移動平均出力はそれに非常に近い形を保ちながらノイズを取り除いています。
直感的なイメージ
移動平均の直感は「スライドする窓」です。長さ $M$ の窓を信号上を1サンプルずつスライドさせながら、窓の中に入っている $M$ 個のサンプルの平均値を出力します。
たとえば $M=5$ の場合、時刻 $n=10$ の出力は $x[10], x[9], x[8], x[7], x[6]$ の5つの平均です。時刻 $n=11$ の出力は $x[11], x[10], x[9], x[8], x[7]$ の5つの平均です。窓がひとつズレるたびに、一番古いサンプルが抜け、新しいサンプルが入ります。
高周波のノイズはランダムに正と負に振れます。それを平均すると、正と負が打ち消し合って小さくなります。一方、低周波のゆっくりしたトレンド成分は、隣り合うサンプルで大きく変化しないため、平均しても大きな値が残ります。これが「移動平均がローパスフィルタとして機能する」直感的な理由です。
数学的定義
移動平均フィルタの出力 $y[n]$ は、入力 $x[n]$ を用いて次のように定義されます。
$$ y[n] = \frac{1}{M} \sum_{k=0}^{M-1} x[n-k] $$
ここで $M$ は窓長（ウィンドウ長）と呼ばれる正の整数です。$k=0$ は現在のサンプル、$k=1$ は1サンプル前、$k=M-1$ は $M-1$ サンプル前を表します。$1/M$ で割るのは、出力のスケールを入力と揃えるための正規化です。
この式を展開すると次のようになります。
$$ y[n] = \frac{1}{M}\bigl(x[n] + x[n-1] + x[n-2] + \cdots + x[n-M+1]\bigr) $$
線形時不変（LTI）システムの観点から書き直せば、これは入力信号 $x[n]$ とフィルタ係数列 $h[n]$ の畳み込み（たたみ込み）です。
$$ y[n] = \sum_{k=-\infty}^{\infty} h[k]\, x[n-k] $$
ここで、フィルタ係数（インパルス応答）は次のように与えられます。
$$ h[k] = \begin{cases} \frac{1}{M} & (0 \le k \le M-1) \\ 0 & (\text{それ以外}) \end{cases} $$
この定義が何を意味するか、次のセクションで視覚的に確認しましょう。
インパルス応答：矩形窓の形をしたFIRフィルタ
図を見ると、フィルタ係数 $h[n]$ が $n=0$ から $n=M-1$ の $M$ 個の区間でのみ $1/M$ の値をとり、それ以外では $0$ であることがわかります。この形は「矩形」または「箱型」に見えることから、矩形窓（rectangular window）と呼ばれます。
なぜこれが FIR フィルタなのか
FIR（Finite Impulse Response）フィルタとは、「有限な長さのインパルス応答を持つフィルタ」という意味です。移動平均フィルタのインパルス応答は $M$ 個の係数だけからなり、それ以降は必ず0になります。これは有限長ですから、移動平均フィルタは定義からして FIR フィルタです。
IIR（Infinite Impulse Response）フィルタは過去の出力を再帰的に使うため、理論上は無限の長さのインパルス応答を持ちます。移動平均フィルタには過去出力へのフィードバックがないため、安定性の問題が生じません。
FIR であることの利点
FIR フィルタには理論的に重要な利点があります。過去の出力を使わないため、フィルタが無条件に安定です。また、後に詳しく説明しますが、係数が対称（この場合はすべて等しい）であることから、線形位相特性が保証されます。線形位相は「フィルタを通しても信号の波形の形が変わらない」という性質で、音声処理や計測において特に重要です。
次は、このフィルタが周波数の観点からどのように動作するかを掘り下げます。
周波数応答の導出：sinc 形状のフィルタ特性
フィルタが「ローパス」として機能することを数式で示しましょう。フィルタの周波数応答（周波数特性）$H(e^{j\omega})$ は、インパルス応答 $h[n]$ の離散時間フーリエ変換（DTFT）として定義されます。
$$ H(e^{j\omega}) = \sum_{n=-\infty}^{\infty} h[n]\, e^{-j\omega n} $$
移動平均フィルタでは $h[n] = 1/M$（$0 \le n \le M-1$ のとき）ですから、非ゼロ項だけを取り出して次のように書けます。
$$ H(e^{j\omega}) = \frac{1}{M} \sum_{n=0}^{M-1} e^{-j\omega n} $$
これは公比 $e^{-j\omega}$ の等比数列の和です。$e^{-j\omega} \neq 1$（すなわち $\omega \neq 0$）のとき、等比数列の和の公式 $\sum_{n=0}^{M-1} r^n = \frac{1-r^M}{1-r}$ を適用します。
$$ H(e^{j\omega}) = \frac{1}{M} \cdot \frac{1 – e^{-j\omega M}}{1 – e^{-j\omega}} $$
ここで分子・分母をそれぞれ整理します。指数の中を $-j\omega M/2$ で括り出すことがポイントです。
分子について、$e^{-j\omega M/2}$ を括り出すと、
$$ 1 – e^{-j\omega M} = e^{-j\omega M/2}\bigl(e^{j\omega M/2} – e^{-j\omega M/2}\bigr) = e^{-j\omega M/2} \cdot 2j\sin\!\left(\frac{M\omega}{2}\right) $$
同様に分母について、$e^{-j\omega/2}$ を括り出すと、
$$ 1 – e^{-j\omega} = e^{-j\omega/2} \cdot 2j\sin\!\left(\frac{\omega}{2}\right) $$
これらを代入すると、
$$ H(e^{j\omega}) = \frac{1}{M} \cdot \frac{e^{-j\omega M/2} \cdot 2j\sin\!\bigl(\frac{M\omega}{2}\bigr)}{e^{-j\omega/2} \cdot 2j\sin\!\bigl(\frac{\omega}{2}\bigr)} $$
$2j$ が約分され、指数部をまとめると、
$$ \boxed{H(e^{j\omega}) = \frac{1}{M} \cdot \frac{\sin\!\bigl(\frac{M\omega}{2}\bigr)}{\sin\!\bigl(\frac{\omega}{2}\bigr)} \cdot e^{-j\omega(M-1)/2}} $$
これが移動平均フィルタの周波数応答の厳密な表式です。この式を分解すると、次の3つの要素からなることがわかります。
- 振幅部分 $\dfrac{1}{M} \cdot \dfrac{|\sin(M\omega/2)|}{|\sin(\omega/2)|}$：sinc 関数に似た形状（デジタルsinc）
- 位相部分 $e^{-j\omega(M-1)/2}$：線形位相（純粋な時間遅れ）
この周波数応答の形を視覚で確認しましょう。
左のパネル（線形スケール）から、$\omega = 0$（直流）で振幅が $1.0$（入力がそのまま通過）になり、周波数が高くなるにつれて振幅が下がっていることが確認できます。ただし減衰は単調ではなく、ゼロになってから再び上がる「サイドローブ」が存在します。右のパネル（dB スケール）では、サイドローブが $-13$ dB 程度と、理想的なローパスフィルタに比べてかなり大きいことがわかります。
ここで自然な疑問が生まれます——窓長 $M$ を変えると、この特性はどのように変化するのでしょうか？
窓長 M と遮断周波数のトレードオフ
$M = 2, 4, 8, 16$ の4つの場合を比較したのが上の図です。図から明確なパターンが読み取れます。
$M$ が大きいほど遮断周波数が低くなります。つまり、より低い周波数から信号を削り始め、より強い平滑化（ノイズ除去）効果が得られます。
遮断周波数 $\omega_c$（振幅が $1/\sqrt{2}$（$-3\,\mathrm{dB}$）になる角周波数）の近似値は、
$$ \omega_c \approx \frac{2\pi}{M} \cdot \frac{1}{2} \cdot C = \frac{0.886 \cdot 2\pi}{M} $$
実用上の近似として、$\omega_c \approx \frac{2\pi}{M}$（ラジアン/サンプル）が使われます。サンプリング周波数 $f_s$ [Hz] の場合、Hz 単位の遮断周波数は次のように表されます。
$$ f_c \approx \frac{f_s}{M} $$
たとえば $f_s = 1000\,\mathrm{Hz}$、$M = 10$ なら、$f_c \approx 100\,\mathrm{Hz}$ です。
主要なトレードオフ
| パラメータ | $M$ を大きくすると | $M$ を小さくすると |
|---|---|---|
| 遮断周波数 | 低くなる（強い平滑化） | 高くなる（弱い平滑化） |
| ノイズ除去 | 強くなる | 弱くなる |
| 群遅延 | 大きくなる（$(M-1)/2$ サンプル） | 小さくなる |
| 立ち上がり応答 | 遅くなる | 速くなる |
| サイドローブ | 変化なし（常に $-13$ dB 程度） | 変化なし |
左のパネルでは同じノイズ信号に $M=3, 8, 20$ の移動平均を適用しています。$M=3$ ではギザギザが残り、$M=20$ ではほぼ完全に滑らかになりますが、信号の立ち上がりや谷の部分がなまります。右のパネルのトレードオフ曲線では、$M$ が増えるにつれてノイズ残留量（RMSE）が下がる一方、遅延が直線的に増加することが視覚的に確認できます。「どこで妥協するか」はアプリケーションの要件次第です。
移動平均フィルタの欠点として、サイドローブが比較的高く（約 $-13\,\mathrm{dB}$）、阻止域のノイズ除去能力が完全ではないことが挙げられます。より急峻な遮断特性が必要な場合は、ハミング窓やカイザー窓を使った設計（Parks-McClellan 法など）が必要です。しかし、サンプリング周波数に対してノイズが十分に高周波な場合は、移動平均で実用上十分なことが多いです。
移動平均フィルタの重要な性質として、次に「信号の形が変わらない」という線形位相特性を確認しましょう。
線形位相と群遅延：信号の形を保つ
周波数応答の式 $H(e^{j\omega}) = \frac{1}{M}\frac{\sin(M\omega/2)}{\sin(\omega/2)} \cdot e^{-j\omega(M-1)/2}$ において、位相成分は
$$ \angle H(e^{j\omega}) = -\frac{M-1}{2}\,\omega \quad (\text{通過帯域内}) $$
となります。これは $\omega$ の線形関数です。この「線形位相」という性質が、フィルタ後の波形が元の形を保つために非常に重要です。
群遅延とは
群遅延（group delay）$\tau_g(\omega)$ は位相の負の傾きとして定義されます。
$$ \tau_g(\omega) = -\frac{d\angle H(e^{j\omega})}{d\omega} $$
移動平均フィルタでは、
$$ \tau_g(\omega) = -\frac{d}{d\omega}\left(-\frac{M-1}{2}\omega\right) = \frac{M-1}{2} \quad [\text{サンプル}] $$
群遅延がすべての周波数で一定（$\omega$ によらず $(M-1)/2$）であることがわかります。
群遅延一定の意味を直感的に説明すると、信号の各周波数成分がすべて同じ時間だけ遅延します。異なる周波数成分が異なる時間だけ遅れると、合成した波形は元の形と異なってしまいます（位相歪み）。移動平均フィルタでは群遅延が一定なので、通過した信号の形が保たれ、ただし $(M-1)/2$ サンプル分だけ遅れて出てくるわけです。
左のパネルでは位相特性（青）が線形近似（赤破線）とほぼ一致しています。右のパネルでは、通過帯域内の群遅延が $(M-1)/2 = 3.5$ サンプルで一定になっていることが確認できます。これはインパルス応答の係数が対称（すべて $1/M$ で等しい）であることから生じる対称 FIR フィルタの性質です。
線形位相の根拠：係数の対称性
FIR フィルタのインパルス応答が
$$ h[n] = h[M-1-n] $$
という対称性を満たすとき（Type I / Type II FIR フィルタ）、周波数応答は常に線形位相になることが証明できます。移動平均フィルタでは $h[n] = 1/M$（$0 \le n \le M-1$）で、確かに $h[n] = h[M-1-n]$ が成り立ちます。
この性質は実用的に非常に大切です。音声・映像・計測など、信号の時間的な「形」を正確に保ちたい場面では、線形位相 FIR フィルタが強く好まれます。
さて、線形位相の性質はわかりましたが、フィルタを実際に使うと「端」で問題が生じます。次のセクションで確認しましょう。
端効果：信号の始端・終端での歪み
移動平均を信号の始端・終端近くに適用すると、問題が生じます。たとえば $M=9$ のとき、$n=0$ の出力を求めようとすると $x[0], x[-1], x[-2], \ldots, x[-8]$ の9個が必要です。しかし $x[-1], x[-2], \ldots$ は存在しません。
この問題を「端効果」または「境界効果」と呼びます。コンピュータ上での実装では、境界でどう処理するかによって結果が変わります。
左のパネルの赤い領域（始端・終端各 $\lfloor M/2 \rfloor$ サンプル）が端効果の影響を受ける範囲です。この部分では、本来のフィルタ係数が揃わないため、出力が歪みます。
対処法：パディング
端効果への主な対処法は「パディング（padding）」です。信号の前後に仮想的なサンプルを補充することで、窓が常に $M$ 個のサンプルにかかるようにします。
1. ゼロパディング（zero padding）
前後に0を詰めます。単純で実装が容易ですが、信号が突然0になる仮定のため、始端・終端で「引っ張り」効果が出ます（右パネルの赤破線）。NumPy の np.convolve(x, h, mode='same') のデフォルト動作がこれに相当します。
2. 反射パディング（reflect padding） 信号を鏡像で折り返して補充します。$x[-1] = x[1]$、$x[-2] = x[2]$ のように使います。境界近くでの連続性が高く、ゼロパディングより歪みが小さくなります（右パネルの緑の線）。
3. 定数パディング（constant padding） 端のサンプル値で補充します（例：$x[-k] = x[0]$）。段差のない信号に向いています。
SciPy では scipy.signal.lfilter や scipy.ndimage.uniform_filter1d（mode オプション付き）で各種パディングを指定できます。
import numpy as np
from scipy.ndimage import uniform_filter1d
# 反射パディングによる移動平均
x = np.array([1.0, 2.0, 3.0, 2.0, 1.0, 0.0, -1.0])
M = 3
# reflect モード（端での歪みが最小）
y_reflect = uniform_filter1d(x, size=M, mode='reflect')
print("反射パディング:", y_reflect)
# -> [1.333 2.000 2.333 2.000 1.000 0.000 -0.667]
# nearest モード（端のサンプルで埋める）
y_nearest = uniform_filter1d(x, size=M, mode='nearest')
print("最近傍パディング:", y_nearest)
均一フィルタ uniform_filter1d は反射パディングを使っており、端での歪みを抑えた移動平均を簡単に計算できます。mode='reflect' を指定すると境界の連続性が最も自然になります。
端効果の問題を理解したところで、次は「大量データに移動平均を高速に適用する方法」を見ていきましょう。
累積和による高速化：O(1) の更新
移動平均の定義どおりに実装すると、1サンプルの出力を更新するたびに $M$ 回の加算と1回の除算が必要で、$N$ サンプルの信号全体では $O(NM)$ の計算量になります。$M$ が大きいほどコストが高くつきます。
スライディングウィンドウの工夫
$y[n]$ から $y[n+1]$ への更新を考えます。
$$ y[n] = \frac{1}{M}\bigl(x[n] + x[n-1] + \cdots + x[n-M+1]\bigr) $$
$$ y[n+1] = \frac{1}{M}\bigl(x[n+1] + x[n] + \cdots + x[n-M+2]\bigr) $$
差を取ると、
$$ y[n+1] = y[n] + \frac{1}{M}\bigl(x[n+1] – x[n-M+1]\bigr) $$
「新しく入ってきたサンプルを足して、窓から外れた古いサンプルを引く」だけで更新できます。これは加算1回・減算1回・除算1回の $O(1)$ 操作です。全体で $O(N)$ となり、$M$ が大きくても計算量が変わりません。
このアイデアは累積和（cumulative sum, prefix sum）を使って簡潔に実装できます。
import numpy as np
def moving_average_fast(x: np.ndarray, M: int) -> np.ndarray:
"""累積和を使った O(N) 移動平均フィルタ。
Args:
x: 入力信号 (1D array)
M: 窓長（正の整数）
Returns:
移動平均出力（valid 部分のみ、長さ N-M+1）
"""
# 累積和を計算（先頭に 0 を追加して差分を取りやすくする）
cumsum = np.cumsum(np.concatenate([[0], x]))
# cumsum[n+M] - cumsum[n] = x[n] + x[n+1] + ... + x[n+M-1]
return (cumsum[M:] - cumsum[:-M]) / M
# 動作確認
x = np.array([1, 3, 5, 7, 9, 7, 5, 3, 1], dtype=float)
y = moving_average_fast(x, M=3)
print("入力:", x)
print("出力 (M=3):", y)
# -> [3.000 5.000 7.000 7.667 7.000 5.000 3.000]
累積和 cumsum を使うと、どんなに大きな $M$ でも1回の np.cumsum と1回の引き算でまとめて計算できます。$M=100$ でも $M=1000$ でも、速度はほぼ変わりません。マイコンでリアルタイム処理するときは、整数の累積和レジスタを1つ持ってスライドさせるだけで済みます。
Python での標準的な実装
実際のコードでは np.convolve または scipy.signal を使うのが一般的です。
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import lfilter
# サンプル信号の生成
rng = np.random.default_rng(42)
N = 200
t = np.linspace(0, 10, N)
clean = np.sin(2 * np.pi * 0.5 * t)
noisy = clean + rng.normal(0, 0.4, N)
# 方法1: np.convolve（フィルタ係数を明示）
M = 10
h = np.ones(M) / M # 矩形窓（フィルタ係数）
y_conv = np.convolve(noisy, h, mode='same') # 'same' で入力と同じ長さ
# 方法2: scipy.ndimage.uniform_filter1d（反射パディング）
from scipy.ndimage import uniform_filter1d
y_uniform = uniform_filter1d(noisy, size=M, mode='reflect')
# 方法3: pandas（金融データ向け）
import pandas as pd
s = pd.Series(noisy)
y_pandas = s.rolling(window=M).mean().values
print(f"RMSE（フィルタなし）: {np.sqrt(np.mean((noisy - clean)**2)):.4f}")
print(f"RMSE（np.convolve）: {np.sqrt(np.mean((y_conv - clean)**2)):.4f}")
print(f"RMSE（uniform_filter1d）: {np.sqrt(np.mean((y_uniform - clean)**2)):.4f}")
np.convolve は最も素直で、フィルタ係数を直接渡すだけです。ただし mode='same' でもゼロパディングが行われるため端には注意が必要です。uniform_filter1d は反射パディングがデフォルトで、端の歪みが小さく計測データに向いています。pandas.rolling().mean() は金融データの分析で最も広く使われています。
3つの方法の RMSE（真の信号との二乗平均平方根誤差）は $M=10$ でほぼ同等になります。方法によって端の処理が異なるため、端を除いた中央部分では完全に一致します。
周波数応答の数値計算
フィルタの周波数特性を数値的に求めるには scipy.signal.freqz が便利です。
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import freqz
fig, axes = plt.subplots(1, 2, figsize=(12, 4))
for M, color in zip([4, 8, 16], ['#2196F3', '#E91E63', '#4CAF50']):
h = np.ones(M) / M # 移動平均フィルタ係数
b = h # 分子多項式係数（FIRなので分母は [1]）
a = [1] # 分母多項式係数
w, H = freqz(b, a, worN=2048) # w は 0〜π のラジアン/サンプル
axes[0].plot(w / np.pi, np.abs(H), label=f'M={M}', color=color, linewidth=2)
axes[1].plot(w / np.pi, 20 * np.log10(np.abs(H) + 1e-10),
label=f'M={M}', color=color, linewidth=2)
for ax in axes:
ax.legend()
ax.grid(True, alpha=0.3)
ax.set_xlabel('正規化角周波数 ω/π')
axes[0].set_ylabel('|H(e^jω)|')
axes[0].set_title('振幅特性（線形）')
axes[1].set_ylabel('振幅 [dB]')
axes[1].set_title('振幅特性（dB）')
axes[1].axhline(-3, color='gray', linestyle='--', linewidth=1)
plt.tight_layout()
plt.show()
freqz(b, a) は FIR フィルタの係数 b（分子）と a（分母、FIR なら $[1]$）から周波数応答 $H(e^{j\omega})$ を計算します。返ってくる w は $[0, \pi]$ の角周波数、H は複素数の周波数応答です。振幅は np.abs(H)、位相は np.angle(H)、dB スケールは 20 * np.log10(np.abs(H)) で得られます。
$M$ が大きいほど遮断周波数が低くなる（曲線が左に移動する）ことが確認できます。また dB スケールでのサイドローブが約 $-13\,\mathrm{dB}$ でほぼ変わらないことも確認できます。これが移動平均フィルタの限界であり、より厳密な阻止域が必要な場合は Parks-McClellan 法などを使う必要があります。
単純移動平均（SMA）と指数移動平均（EMA）の比較
移動平均には2種類の主要な形式があります。ここまで解説してきた「単純移動平均（SMA）」と、金融分析で広く使われる「指数移動平均（EMA）」です。
EMA の定義
EMA は次の漸化式で定義されます。
$$ y[n] = \alpha\, x[n] + (1-\alpha)\, y[n-1], \quad 0 < \alpha \le 1 $$
$\alpha$ は「平滑化係数」と呼ばれます。$\alpha$ が大きいほど最新の入力 $x[n]$ の重みが大きく（素早く反応）、$\alpha$ が小さいほど過去の値の重みが大きい（ゆっくり追従）です。
この漸化式を繰り返し展開すると、
$$ y[n] = \sum_{k=0}^{\infty} \alpha(1-\alpha)^k\, x[n-k] $$
となり、過去のサンプルに対して 指数的に減衰する重み が付いていることがわかります。直近のサンプルが最も重く、遠い過去ほど軽くなります。
SMA と EMA の比較
| 性質 | 単純移動平均（SMA） | 指数移動平均（EMA） |
|---|---|---|
| インパルス応答 | 有限長（矩形窓） | 無限長（指数減衰） |
| フィルタ種別 | FIR | IIR |
| 線形位相 | あり（群遅延一定） | なし（位相歪みあり） |
| 計算コスト | O(M) または O(1)（累積和） | O(1)（1サンプルごと） |
| メモリ | M 個のバッファ | 前回出力 1 個のみ |
| 段差追従 | M サンプルの遅延 | 指数的に追従（残留偏差あり） |
| マイコン実装 | やや複雑（バッファ管理） | 非常に簡単 |
上の図の左パネル（段差信号の追従）を見ると、SMA は段差後 $M$ サンプルかけて完全に追従するのに対し、EMA は指数的に追従しますが理論上は完全には到達しません（残留偏差）。実用上は $\alpha=0.2$ なら数サンプルで実用的な近さに到達します。
右パネル（インパルス応答の形状）では、SMA が矩形（有限長）なのに対し、EMA が指数減衰（無限長）になっています。EMA はメモリ1個で実装できる超軽量フィルタです。
等価な窓長の関係
SMA の窓長 $M$ と EMA の $\alpha$ を「おおよそ等価」にする変換式として、
$$ \alpha = \frac{2}{M+1} $$
がよく使われます（ただし厳密に等価ではありません）。$M=10$ なら $\alpha = 2/11 \approx 0.182$ となります。
import numpy as np
def ema(x: np.ndarray, alpha: float) -> np.ndarray:
"""指数移動平均（EMA）の計算。
Args:
x: 入力信号
alpha: 平滑化係数 (0 < alpha <= 1)
Returns:
EMA の出力
"""
y = np.zeros_like(x, dtype=float)
y[0] = x[0]
for i in range(1, len(x)):
y[i] = alpha * x[i] + (1 - alpha) * y[i - 1]
return y
# 動作確認
x = np.array([0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 0], dtype=float)
alpha = 0.3
y = ema(x, alpha)
print(f"EMA (α={alpha}):", np.round(y, 3))
# -> [0. 0. 0. 0.3 0.51 0.657 0.76 0.832 0.582 0.408 0.285 0.2 ]
EMA は段差応答の追従が指数関数的で、立ち上がりは速いが完全には到達しません。一方、SMA は窓の外に段差が出るまでの $M$ サンプル間、直線的に上昇して完全に追従します。どちらが良いかは「素早い反応が必要か（EMA 向き）」か「信号の形を完全に保ちたいか（SMA 向き）」によります。
EMA と SMA の違いを直感的に把握したところで、実際のノイズ除去効果を確認してみましょう。
ノイズ除去デモ：時間領域とスペクトル
時間領域での比較
この図では、$y = \sin(2\pi \cdot 0.5 t) + 0.5\sin(2\pi \cdot 1.0 t)$ という2成分の信号に標準偏差 0.5 のガウスノイズを加えた信号（①上段）に対し、$M=12$ の移動平均を適用した結果（②中段）を示しています。
中段のパネルでは、ノイズが大幅に除去され、移動平均出力（赤）が真の信号（青破線）に非常に近くなっています。RMSE が改善し、明らかに品質が向上しています。③下段は除去されたノイズ成分（入力−出力の差分）で、ランダムなノイズのみが残っており、信号成分がほぼ含まれていないことが確認できます。
スペクトル領域での比較
左のパネルでは、フィルタ前後のパワースペクトル密度（PSD）を比較しています。入力（グレー）では高周波全体に広がるノイズフロアが存在しますが、移動平均後は高周波成分が大幅に減少しています。1 Hz と 5 Hz の信号成分（黒点線）は両方とも残っています。
右のパネルは移動平均フィルタの理論的な周波数応答（dB スケール）で、$M=8$ ではカットオフが約 $f_s/M$ Hz、$M=20$ ではさらに低いカットオフになっています。実測のスペクトルと理論特性が一致していることが確認できます。
import numpy as np
import matplotlib.pyplot as plt
# パワースペクトル密度の計算と比較
rng = np.random.default_rng(42)
N = 1024
fs = 100.0 # サンプリング周波数 [Hz]
t = np.arange(N) / fs
# 信号生成（2Hz と 10Hz の正弦波 + ノイズ）
clean = np.sin(2 * np.pi * 2.0 * t) + 0.5 * np.sin(2 * np.pi * 10.0 * t)
noisy = clean + rng.normal(0, 0.5, N)
# 移動平均フィルタ適用
M = 8
h = np.ones(M) / M
filtered = np.convolve(noisy, h, mode='same')
# FFT によるスペクトル計算
freqs = np.fft.rfftfreq(N, d=1/fs)
S_noisy = 20 * np.log10(np.abs(np.fft.rfft(noisy)) / N + 1e-12)
S_filt = 20 * np.log10(np.abs(np.fft.rfft(filtered)) / N + 1e-12)
print(f"RMSE（フィルタなし）: {np.sqrt(np.mean((noisy - clean)**2)):.4f}")
print(f"RMSE（M={M}） : {np.sqrt(np.mean((filtered - clean)**2)):.4f}")
print(f"SNR 改善（概算） : {20*np.log10(np.sqrt(np.mean((noisy-clean)**2)) / (np.sqrt(np.mean((filtered-clean)**2)) + 1e-10)):.1f} dB")
このコードで RMSE の改善量と SNR 改善量が定量的に確認できます。ノイズ成分が高周波に広がっているほど、移動平均による SNR 改善効果が大きくなります。計算された SNR 改善量は正の値になるはずで、これがローパスフィルタとしての効果を端的に表しています。
移動平均フィルタの実装まとめ
これまで学んだ内容を一つの関数にまとめます。
import numpy as np
from scipy.ndimage import uniform_filter1d
def moving_average(x: np.ndarray, M: int,
mode: str = 'reflect') -> np.ndarray:
"""移動平均フィルタ（scikit-style API）。
Args:
x : 入力信号（1D）
M : 窓長（奇数推奨 → 整数遅延）
mode : パディング方式 ('reflect', 'nearest', 'wrap', 'constant')
Returns:
移動平均出力（入力と同じ長さ）
例::
>>> y = moving_average(np.array([1,3,5,7,9], dtype=float), M=3)
>>> print(y)
[2. 3. 5. 7. 8.]
"""
return uniform_filter1d(x.astype(float), size=M, mode=mode)
def moving_average_fast(x: np.ndarray, M: int) -> np.ndarray:
"""累積和 O(N) 移動平均フィルタ（valid 部分のみ）。
Returns:
長さ len(x) - M + 1 の出力配列
"""
cs = np.cumsum(np.concatenate([[0.0], x]))
return (cs[M:] - cs[:-M]) / M
def ema_filter(x: np.ndarray, alpha: float) -> np.ndarray:
"""指数移動平均（EMA）フィルタ。
Args:
alpha: 平滑化係数 (0 < alpha <= 1)。小さいほど強い平滑化
"""
y = np.empty_like(x, dtype=float)
y[0] = x[0]
for i in range(1, len(x)):
y[i] = alpha * x[i] + (1 - alpha) * y[i - 1]
return y
# --- 動作確認 ---
rng = np.random.default_rng(0)
x = np.sin(np.linspace(0, 4 * np.pi, 100)) + rng.normal(0, 0.3, 100)
y_sma = moving_average(x, M=9)
y_fast = moving_average_fast(x, M=9)
y_ema = ema_filter(x, alpha=2/(9+1)) # M=9 相当
# 中央部分で完全一致を確認（端効果を除く）
center = slice(5, 90)
print(f"SMA vs fast（中央部）最大差: {np.max(np.abs(y_sma[center] - y_fast[:len(y_fast[center])])):.6f}")
# -> 0.000000 程度（パディングが異なるため端は一致しない）
print(f"SMA RMSE: {np.sqrt(np.mean((y_sma - np.sin(np.linspace(0, 4*np.pi, 100)))**2)):.4f}")
print(f"EMA RMSE: {np.sqrt(np.mean((y_ema - np.sin(np.linspace(0, 4*np.pi, 100)))**2)):.4f}")
このコードでは3つの実装（uniform_filter1d・累積和・EMA）を比較しています。中央部分では SMA の2実装がほぼ一致します（パディングの違いで端は異なります）。SMA と EMA の RMSE は同程度になりますが、EMA は位相遅れが生じるため波形の一致度は異なります。
まとめ
本記事では移動平均フィルタについて、定義から周波数特性の厳密な導出、実装の高速化まで体系的に解説しました。
- 定義: $y[n] = \frac{1}{M}\sum_{k=0}^{M-1} x[n-k]$。直近 $M$ サンプルの相加平均
- インパルス応答: 矩形窓（有限長 FIR フィルタ）
- 周波数応答: $H(e^{j\omega}) = \frac{\sin(M\omega/2)}{M\sin(\omega/2)}\,e^{-j\omega(M-1)/2}$（sinc 形状）
- 遮断周波数: $\omega_c \approx 2\pi/M$（ラジアン/サンプル）、$f_c \approx f_s/M$（Hz）
- 線形位相: 群遅延 $(M-1)/2$ が全周波数で一定 → 波形形状を保存
- 端効果: 反射パディングで歪みを軽減
- 高速化: 累積和で $O(1)$ 更新（$O(N)$ 全体）
- EMA との違い: EMA は IIR で位相歪みあり・メモリ1個で超軽量
移動平均フィルタはシンプルさゆえに「最初に試すフィルタ」として広く使われますが、サイドローブ（$-13\,\mathrm{dB}$ 程度）と急峻でない遮断特性が限界です。より高品質な遮断特性が必要な場合は、窓関数（ハミング・カイザー等）を使った FIR フィルタ設計や、Parks-McClellan 法（等リプルフィルタ）を検討してください。
