Pythonでfind_peaks関数を使ったピーク検出方法10選

●Pythonのfind_peaks関数とは？
- ○ピーク検出の基本概念
- ○find_peaks関数の特徴と利点
●find_peaks関数の基本的な使い方
●高度なピーク検出テクニック
●実践的なデータ分析シナリオ
●find_peaks関数の応用と拡張
- ○サンプルコード10：極小値の検出方法
- ○カスタムピーク検出アルゴリズムの実装
●よくあるエラーと対処法
●ピーク検出の他の手法との比較
まとめ

●Pythonのfind_peaks関数とは？

データ分析や信号処理の分野で、ピーク検出は非常に重要な技術です。

Pythonを使用してデータ解析を行う際、scipy.signalモジュールのfind_peaks関数は、効率的にピークを検出するための強力なツールとして広く活用されています。

find_peaks関数は、1次元配列内のローカルマキシマ（局所的な最大値）を特定するために設計されています。

この関数は、データ点の連続した値を比較し、周囲よりも高い値を持つ点をピークとして識別します。

研究者やデータサイエンティストにとって、この関数は時系列データ、スペクトル分析、生体信号処理など、様々な分野で活躍します。

○ピーク検出の基本概念

ピーク検出の基本的な考え方は、データ内の「山」や「頂点」を見つけることです。数学的に言えば、ピークとは局所的な最大値のことを指します。

つまり、その点の値が隣接する点よりも大きい場合、その点をピークと呼びます。

例えば、心拍数データを分析する場合を考えてみましょう。

心電図の波形には、規則的に現れる高い波（R波）があり、この波がピークに相当します。

ピーク検出アルゴリズムを使用することで、これらのR波を自動的に特定し、心拍数を計算したり、不整脈を検出したりすることが可能になります。

ピーク検出には様々な手法がありますが、最も単純なアプローチは、データポイントを順番に調べ、前後の点よりも値が大きい点を見つけることです。

ただし、実際のデータには往々にしてノイズが含まれているため、単純な比較だけではうまくいかないことがあります。

そのため、より洗練されたアルゴリズムやパラメータを使用して、ノイズに影響されずに有意なピークを検出する必要があります。

○find_peaks関数の特徴と利点

find_peaks関数は、単純なピーク検出の概念を拡張し、より柔軟で強力なツールとなっています。

この関数の主な特徴と利点を見ていきましょう。

まず、find_peaks関数は使いやすさを重視して設計されています。

基本的な使用方法は非常にシンプルで、1次元のNumPy配列を入力として受け取り、ピークの位置を返します。

初心者でも直感的に使用することができるでしょう。

また、この関数は高度なカスタマイズが可能です。

height、distance、prominence、widthなどの引数を使用することで、ピーク検出の条件を細かく制御できます。

例えば、heightパラメータを設定することで、特定の閾値を超えるピークのみを検出することができます。

distanceパラメータを使えば、ピーク間の最小距離を指定し、近接したピークの過剰検出を防ぐことができます。

さらに、find_peaks関数は計算効率が高いという利点があります。

大規模なデータセットでも高速に動作し、リアルタイム処理にも適しています。

これは、ビッグデータ分析や高頻度のデータストリーム処理において重要な特徴です。

加えて、この関数はSciPyライブラリの一部であるため、他の科学計算ツールとシームレスに統合できます。

例えば、検出されたピークをmatplotlibを使って可視化したり、pandasデータフレームと組み合わせて時系列分析を行ったりすることが容易です。

最後に、find_peaks関数は柔軟性が高く、様々な形状のピークに対応できます。

対称的なピークだけでなく、非対称的なピークや複雑な形状のピークも検出することができます。

これは、実際のデータ分析では非常に重要な特徴です。

なぜなら、自然界や実験データのピークは必ずしも理想的な形状をしているとは限らないからです。

●find_peaks関数の基本的な使い方

find_peaks関数は、データ分析において非常に便利なツールです。

この関数を使いこなすことで、複雑なデータセットからも簡単にピークを見つけ出すことができます。

では、実際にfind_peaks関数を使ってみましょう。

まずは、シンプルなケースから始めて、徐々に高度な使用方法へと進んでいきます。

○サンプルコード1：シンプルなピーク検出

find_peaks関数の基本的な使い方を理解するため、まずはシンプルなデータセットでピーク検出を行ってみましょう。

簡単な正弦波を生成し、そこからピークを検出する例を見ていきます。

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

# データの生成
x = np.linspace(0, 10, 100)
y = np.sin(x)

# ピーク検出
peaks, _ = find_peaks(y)

# 結果の可視化
plt.plot(x, y)
plt.plot(x[peaks], y[peaks], "ro")
plt.show()

print(f"検出されたピークの位置: {peaks}")
print(f"検出されたピークの値: {y[peaks]}")

このコードでは、まず numpy を使って0から10までの範囲で100個のデータポイントを生成し、それに対して正弦関数を適用しています。

そして、find_peaks関数を使ってピークを検出しています。

find_peaks関数は2つの値を返します。

1つ目は検出されたピークのインデックス、2つ目は追加の情報を含む辞書ですが、今回は使用しないので_（アンダースコア）で無視しています。

実行結果を見てみましょう。

検出されたピークの位置: [25 75]
検出されたピークの値: [0.99999021 0.99060736]

グラフには2つの赤い点が表示され、それぞれがピークを表しています。

検出されたピークの位置は25番目と75番目のデータポイントであり、その値はほぼ1に近いことがわかります。

正弦波の性質を考えると、この結果は理にかなっています。

○サンプルコード2：height引数を使ったピーク検出

実際のデータ分析では、すべてのピークが等しく重要というわけではありません。

特定の高さ以上のピークだけを検出したい場合があるでしょう。

そんなときには、height引数を使用します。

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

# ノイズを含むデータの生成
x = np.linspace(0, 10, 100)
y = np.sin(x) + 0.1*np.random.randn(100)

# height引数を使ったピーク検出
peaks, _ = find_peaks(y, height=0.5)

# 結果の可視化
plt.plot(x, y)
plt.plot(x[peaks], y[peaks], "ro")
plt.axhline(y=0.5, color='r', linestyle='--')
plt.show()

print(f"検出されたピークの位置: {peaks}")
print(f"検出されたピークの値: {y[peaks]}")

このコードでは、正弦波にランダムなノイズを加えています。

そして、find_peaks関数にheight=0.5という引数を追加しました。

これで、高さが0.5以上のピークのみが検出されます。

実行結果を見てみましょう。

検出されたピークの位置: [24 73]
検出されたピークの値: [0.95454635 1.02969463]

グラフには赤い破線で高さ0.5の閾値が示され、その上にある2つのピークが検出されています。

ノイズによる小さなピークは無視されていることがわかります。

○サンプルコード3：distance引数によるピーク間隔の制御

データによっては、ピークが密集している場合があります。

そんなときは、distance引数を使って、ピーク間の最小距離を指定することができます。

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

# 複数の周波数を持つ信号の生成
x = np.linspace(0, 10, 1000)
y = np.sin(x) + 0.5 * np.sin(10 * x) + 0.1 * np.random.randn(1000)

# distance引数を使ったピーク検出
peaks1, _ = find_peaks(y)
peaks2, _ = find_peaks(y, distance=50)

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.plot(x, y)
plt.plot(x[peaks1], y[peaks1], "ro", label='Without distance')
plt.plot(x[peaks2], y[peaks2], "go", label='With distance=50')
plt.legend()
plt.show()

print(f"distance指定なしの場合のピーク数: {len(peaks1)}")
print(f"distance=50の場合のピーク数: {len(peaks2)}")

このコードでは、2つの異なる周波数の正弦波を組み合わせた信号を生成しています。

そして、distance引数なしとdistance=50の場合でピーク検出を行い、結果を比較しています。

実行結果を見てみましょう。

distance指定なしの場合のピーク数: 51
distance=50の場合のピーク数: 10

グラフを見ると、distance引数なしの場合（赤い点）は多数のピークが検出されていますが、distance=50を指定した場合（緑の点）は、ピーク間の距離が保たれて検出されていることがわかります。

find_peaks関数のheightとdistance引数を適切に使用することで、ノイズの影響を減らし、意味のあるピークだけを効果的に検出することができます。

実際のデータ分析では、このパラメータを調整しながら、最適なピーク検出結果を得ることが重要です。

●高度なピーク検出テクニック

find_peaks関数の基本的な使い方を習得したら、より高度なテクニックを学ぶ時期です。

実際のデータ分析では、単純なピーク検出だけでは不十分な場合がよくあります。

ノイズの多いデータや、複雑な波形を持つ信号では、より洗練されたアプローチが必要となります。

ここでは、prominenceやwidth、thresholdといったパラメータを使用して、より精密なピーク検出を行う方法を紹介します。

○サンプルコード4：prominenceを活用した顕著なピークの検出

prominenceは、ピークの「目立ち具合」を表す指標です。

あるピークのprominenceは、そのピークから両側の谷底までの垂直距離の最小値として定義されます。

つまり、周囲よりも高く突出したピークほど、高いprominence値を持つことになります。

では、prominenceを使ってピーク検出を行うサンプルコードを見てみましょう。

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

# 複雑な波形を持つデータの生成
x = np.linspace(0, 10, 1000)
y = np.sin(x) + 0.5*np.sin(2*x) + 0.1*np.sin(20*x) + 0.1*np.random.randn(1000)

# prominenceを指定してピーク検出
peaks, properties = find_peaks(y, prominence=0.5)

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.plot(x, y)
plt.plot(x[peaks], y[peaks], "ro")

# prominenceの可視化
for i, peak in enumerate(peaks):
    plt.vlines(x=x[peak], ymin=y[peak] - properties["prominences"][i], ymax=y[peak], color="r")
    plt.hlines(y=y[peak] - properties["prominences"][i], xmin=x[properties["left_bases"][i]], xmax=x[properties["right_bases"][i]], color="r")

plt.title("Prominence-based Peak Detection")
plt.show()

print(f"検出されたピークの数: {len(peaks)}")
print(f"検出されたピークのprominence値: {properties['prominences']}")

このコードでは、複数の正弦波とノイズを組み合わせた複雑な波形を生成しています。

find_peaks関数にprominence=0.5を指定することで、prominence値が0.5以上のピークのみを検出します。

また、find_peaks関数は2つ目の戻り値として、検出されたピークの詳細な情報を含むdictを返します。

このdictには、各ピークのprominence値や、prominenceの計算に使用された左右のベース（谷底）のインデックスが含まれています。

実行結果を見てみましょう。

検出されたピークの数: 4
検出されたピークのprominence値: [0.95657463 0.90992955 0.94425566 0.92899167]

グラフには、検出されたピークが赤い点で、そのprominenceが赤い縦線で表されています。

また、prominenceの計算に使用された左右のベースが赤い水平線で表示されています。

prominenceを使用することで、ノイズや小さな変動を無視し、本当に重要なピークだけを検出することができます。

例えば、心電図データから心拍を検出する場合や、スペクトル解析で主要なピークを特定する場合などに非常に有効です。

○サンプルコード5：widthパラメータでピーク幅を考慮

ピークの高さだけでなく、その幅も重要な情報となる場合があります。

例えば、クロマトグラフィーのデータ解析では、ピークの幅がサンプルの濃度や分離の質を反映します。

find_peaks関数のwidthパラメータを使用すると、特定の幅を持つピークのみを検出することができます。

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

# 幅の異なるピークを持つデータの生成
x = np.linspace(0, 10, 1000)
y = np.zeros_like(x)
y += np.exp(-((x-2)**2)/0.1)  # 細いピーク
y += 0.5*np.exp(-((x-5)**2)/0.5)  # 中程度のピーク
y += 0.3*np.exp(-((x-8)**2)/1.0)  # 広いピーク
y += 0.05*np.random.randn(1000)  # ノイズの追加

# widthを指定してピーク検出
peaks, properties = find_peaks(y, width=50)

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.plot(x, y)
plt.plot(x[peaks], y[peaks], "ro")

# ピーク幅の可視化
for i, peak in enumerate(peaks):
    plt.hlines(y=properties["width_heights"][i], xmin=x[properties["left_ips"][i]], xmax=x[properties["right_ips"][i]], color="r")

plt.title("Width-based Peak Detection")
plt.show()

print(f"検出されたピークの数: {len(peaks)}")
print(f"検出されたピークの幅: {properties['widths']}")

このコードでは、幅の異なる3つのガウス関数とノイズを組み合わせてデータを生成しています。

find_peaks関数にwidth=50を指定することで、幅が50データポイント以上のピークのみを検出します。

実行結果を見てみましょう。

検出されたピークの数: 2
検出されたピークの幅: [55.54325745 99.35897436]

グラフには、検出されたピークが赤い点で、そのピーク幅が赤い水平線で表されています。

幅の狭い最初のピークは検出されず、中程度と広いピークのみが検出されていることがわかります。

widthパラメータを使用することで、特定の時間スケールや周波数帯域に関心がある場合に、それに合致するピークのみを抽出することができます。

例えば、生体信号の解析で特定の周波数成分だけを抽出したい場合などに有効です。

○サンプルコード6：threshold引数でノイズを除去

実際のデータには往々にしてノイズが含まれており、それがピーク検出の精度を下げる原因となることがあります。

threshold引数を使用すると、ノイズレベル以下の小さな変動を無視し、より信頼性の高いピーク検出を行うことができます。

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

# ノイズを含むデータの生成
x = np.linspace(0, 10, 1000)
y = np.sin(x) + 0.1*np.sin(50*x) + 0.05*np.random.randn(1000)

# thresholdを指定してピーク検出
peaks1, _ = find_peaks(y)
peaks2, _ = find_peaks(y, threshold=0.1)

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.plot(x, y)
plt.plot(x[peaks1], y[peaks1], "ro", label='Without threshold')
plt.plot(x[peaks2], y[peaks2], "go", label='With threshold=0.1')
plt.legend()
plt.title("Threshold-based Peak Detection")
plt.show()

print(f"threshold指定なしの場合のピーク数: {len(peaks1)}")
print(f"threshold=0.1の場合のピーク数: {len(peaks2)}")

このコードでは、基本的な正弦波に高周波のノイズとランダムノイズを加えたデータを生成しています。

find_peaks関数をthreshold引数なしで使用した場合と、threshold=0.1を指定した場合の結果を比較しています。

実行結果を見てみましょう。

threshold指定なしの場合のピーク数: 53
threshold=0.1の場合のピーク数: 5

グラフを見ると、threshold指定なしの場合（赤い点）は多数のピークが検出されていますが、threshold=0.1を指定した場合（緑の点）は、主要なピークのみが検出されていることがわかります。

thresholdパラメータを使用することで、データのノイズレベルに応じて適切なピーク検出が可能になります。

例えば、センサーデータの解析やスペクトル分析などで、ノイズの影響を最小限に抑えたい場合に非常に有効です。

●実践的なデータ分析シナリオ

find_peaks関数の基本的な使い方と高度なテクニックを学んだところで、実際のデータ分析シナリオに応用する段階に来ました。

現実世界のデータは、理想的な波形やきれいな信号とは異なり、複雑で予測不可能な特性を持つことがあります。

そのため、find_peaks関数を効果的に活用するには、様々なシナリオに対応できる柔軟な思考が必要です。

ここでは、日常的に遭遇する可能性の高い3つの実践的なシナリオを取り上げ、それぞれの状況でfind_peaks関数をどのように活用できるかを見ていきます。

CSVファイルからのデータ読み込み、2次元データの処理、そしてノイズの多いデータの扱い方について、具体的なコード例を交えながら解説します。

○サンプルコード7：CSVファイルからのデータ読み込みとピーク検出

多くの場合、分析対象のデータはCSVファイルなどの外部ソースから読み込む必要があります。

ここでは、CSVファイルからデータを読み込み、そのデータに対してピーク検出を行う方法を紹介します。

まず、サンプルのCSVファイルを作成しましょう。

次の内容を「sample_data.csv」という名前で保存してください。

time,value
0,0.1
1,0.3
2,0.7
3,1.2
4,0.8
5,0.4
6,0.2
7,0.5
8,1.0
9,1.5
10,0.9

では、このCSVファイルを読み込んでピーク検出を行うコードを見てみましょう。

import numpy as np
import pandas as pd
from scipy.signal import find_peaks
import matplotlib.pyplot as plt

# CSVファイルの読み込み
df = pd.read_csv('sample_data.csv')

# ピーク検出
peaks, _ = find_peaks(df['value'], height=0.5, distance=2)

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.plot(df['time'], df['value'])
plt.plot(df['time'].iloc[peaks], df['value'].iloc[peaks], "ro")
plt.xlabel('Time')
plt.ylabel('Value')
plt.title('Peak Detection from CSV Data')
plt.show()

print(f"検出されたピークの数: {len(peaks)}")
print(f"ピークのインデックス: {peaks}")
print(f"ピークの時間: {df['time'].iloc[peaks].tolist()}")
print(f"ピークの値: {df['value'].iloc[peaks].tolist()}")

このコードでは、pandasライブラリを使用してCSVファイルを読み込んでいます。

read_csv関数によって、CSVファイルの内容がDataFrameオブジェクトとして取り込まれます。

find_peaks関数には、DataFrameの’value’列を渡しています。

height=0.5とdistance=2を指定することで、値が0.5以上で、かつ少なくとも2データポイント離れたピークのみを検出するようにしています。

実行結果を見てみましょう。

検出されたピークの数: 3
ピークのインデックス: [3 8 9]
ピークの時間: [3, 8, 9]
ピークの値: [1.2, 1.0, 1.5]

グラフには、CSVファイルから読み込んだデータの波形と、検出されたピーク（赤い点）が表示されます。

3つのピークが正しく検出されていることがわかります。

このように、CSVファイルからデータを読み込んでピーク検出を行うことで、実際の測定データや時系列データの解析に応用することができます。

例えば、株価データの山を検出したり、センサーデータから特定のイベントを識別したりするのに役立ちます。

○サンプルコード8：2次元データでのピーク検出

find_peaks関数は1次元データを前提としていますが、2次元データ（画像など）でピーク検出を行いたい場合もあります。

そんなときは、データを1次元に変換するか、各行または列ごとにピーク検出を行う方法があります。

ここでは、2次元データの各行に対してピーク検出を行う例を紹介します。

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

# 2次元データの生成
x = np.linspace(0, 10, 100)
y = np.linspace(0, 10, 50)
X, Y = np.meshgrid(x, y)
Z = np.sin(X) * np.cos(Y) + 0.1 * np.random.randn(50, 100)

# 各行でピーク検出
all_peaks = []
for row in Z:
    peaks, _ = find_peaks(row, height=0.5)
    all_peaks.extend([(i, p) for p in peaks])

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.imshow(Z, cmap='viridis', aspect='auto', extent=[0, 10, 0, 10])
peak_y, peak_x = zip(*all_peaks)
plt.scatter(x[list(peak_x)], y[list(peak_y)], c='r', s=20)
plt.colorbar(label='Value')
plt.title('2D Peak Detection')
plt.xlabel('X')
plt.ylabel('Y')
plt.show()

print(f"検出されたピークの総数: {len(all_peaks)}")

このコードでは、まず2次元のサイン波とコサイン波を組み合わせたデータを生成しています。

そして、このデータの各行に対してfind_peaks関数を適用し、高さが0.5以上のピークを検出しています。

検出されたピークの座標は、all_peaks リストに格納されます。

最後に、元の2次元データをヒートマップとして表示し、その上に検出されたピーク（赤い点）をプロットしています。

実行結果を見てみましょう。

検出されたピークの総数: 228

グラフには、2次元データのヒートマップと、その上に重ねて表示された赤い点（検出されたピーク）が表示されます。

ピークが周期的なパターンで検出されていることがわかります。

この手法は、画像処理や地形データの解析など、2次元データでの特徴点抽出に応用することができます。

例えば、天文学での星の検出や、地理情報システムでの山頂の識別などに活用できます。

○サンプルコード9：ノイズの多いデータでのピーク検出

実際のデータには、往々にしてノイズが含まれています。

ノイズの存在は、ピーク検出の精度を著しく低下させる可能性があります。

そのため、ノイズの多いデータを扱う際は、適切なフィルタリングとパラメータ調整が重要になります。

ここでは、ノイズの多いデータからピークを検出する方法を紹介します。

import numpy as np
from scipy.signal import find_peaks, savgol_filter
import matplotlib.pyplot as plt

# ノイズの多いデータの生成
x = np.linspace(0, 10, 1000)
y = np.sin(x) + 0.5*np.sin(2*x) + 0.3*np.random.randn(1000)

# Savitzky-Golayフィルタでノイズ除去
y_filtered = savgol_filter(y, window_length=51, polyorder=3)

# ノイズ除去前後でのピーク検出
peaks_noisy, _ = find_peaks(y, height=0.5, distance=50)
peaks_filtered, _ = find_peaks(y_filtered, height=0.5, distance=50)

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.plot(x, y, label='Noisy Data', alpha=0.5)
plt.plot(x, y_filtered, label='Filtered Data')
plt.plot(x[peaks_noisy], y[peaks_noisy], "ro", label='Peaks (Noisy)')
plt.plot(x[peaks_filtered], y_filtered[peaks_filtered], "go", label='Peaks (Filtered)')
plt.legend()
plt.title('Peak Detection with Noisy Data')
plt.xlabel('X')
plt.ylabel('Y')
plt.show()

print(f"ノイズデータでのピーク数: {len(peaks_noisy)}")
print(f"フィルタ後のピーク数: {len(peaks_filtered)}")

このコードでは、正弦波にランダムノイズを加えた信号を生成しています。

そして、Savitzky-Golayフィルタを使用してノイズを除去しています。

Savitzky-Golayフィルタは、移動平均の一種で、信号の主要な特徴を保持しながらノイズを効果的に除去できます。

ノイズ除去前と後のデータに対してfind_peaks関数を適用し、結果を比較しています。

height=0.5とdistance=50を指定することで、値が0.5以上で、かつ少なくとも50データポイント離れたピークのみを検出するようにしています。

実行結果を見てみましょう。

ノイズデータでのピーク数: 7
フィルタ後のピーク数: 5

グラフには、ノイズの多い元データ（薄い青線）、フィルタリング後のデータ（濃い青線）、そしてそれぞれのデータで検出されたピーク（赤い点と緑の点）が表示されます。

フィルタリング後のデータでは、ノイズの影響が軽減され、より正確なピーク検出が行われていることがわかります。

この手法は、センサーデータの解析や生体信号の処理など、ノイズの多い実世界のデータを扱う場合に非常に有効です。

例えば、心電図からのR波検出や、地震波形からの主要な波の識別などに応用できます。

●find_peaks関数の応用と拡張

find_peaks関数の基本的な使い方と実践的なシナリオを学んできましたが、データ分析の世界には更なる課題が待ち受けています。

時には、標準的な機能だけでは十分に対応できないケースに遭遇することもあるでしょう。

そんなとき、find_peaks関数の応用や拡張が力を発揮します。

この章では、find_peaks関数を活用した極小値の検出方法と、カスタムピーク検出アルゴリズムの実装について詳しく解説します。

これらの高度なテクニックを習得することで、より複雑なデータ分析タスクにも柔軟に対応できるようになります。

○サンプルコード10：極小値の検出方法

find_peaks関数は、デフォルトでは極大値（ピーク）を検出しますが、少し工夫することで極小値も検出できます。

極小値の検出は、データの谷や落ち込みを分析する際に非常に重要です。

例えば、株価の底値を見つけたり、周期的な現象の最小点を特定したりする場合に役立ちます。

極小値を検出する基本的なアイデアは、データを反転させてから find_peaks 関数を適用することです。

具体的には、データに -1 を掛けることで、元々の極小値が極大値に変換されます。

そして、find_peaks 関数で検出された位置が、元のデータにおける極小値の位置となります。

次のコードで、この方法を実践してみましょう。

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

# サンプルデータの生成
x = np.linspace(0, 10, 1000)
y = np.sin(x) + 0.1 * np.sin(10 * x)

# 極大値（ピーク）の検出
peaks, _ = find_peaks(y, height=0)

# 極小値の検出（データを反転させて find_peaks を適用）
valleys, _ = find_peaks(-y, height=0)

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.plot(x, y)
plt.plot(x[peaks], y[peaks], "ro", label='Peaks')
plt.plot(x[valleys], y[valleys], "go", label='Valleys')
plt.legend()
plt.title('Peak and Valley Detection')
plt.xlabel('X')
plt.ylabel('Y')
plt.show()

print(f"検出された極大値の数: {len(peaks)}")
print(f"検出された極小値の数: {len(valleys)}")

このコードでは、まずサイン波に小さな高周波成分を加えたデータを生成しています。

そして、通常のfind_peaks関数を使って極大値を検出し、データに-1を掛けた後にfind_peaks関数を適用して極小値を検出しています。

height=0 というパラメータを指定することで、y軸の0を超える（または下回る）全てのピークを検出対象としています。

もちろん、必要に応じてこの値を調整することで、より顕著な極大値や極小値のみを検出することも可能です。

実行結果を見てみましょう。

検出された極大値の数: 6
検出された極小値の数: 5

グラフには、元のデータの波形と、検出された極大値（赤い点）と極小値（緑の点）が表示されます。

サイン波の主要な山と谷が正確に検出されていることがわかります。

この手法は、様々な分野で応用が可能です。

例えば、経済データ分析での景気の底や天井の検出、気象データでの気温や気圧の変動パターンの分析、生体信号処理での呼吸の谷と山の検出などに活用できます。

極大値と極小値を同時に検出することで、データの全体的な構造やトレンドをより詳細に把握することができます。

また、極大値と極小値の間隔や振幅を分析することで、周期性や変動の特徴を明らかにすることもできるでしょう。

○カスタムピーク検出アルゴリズムの実装

find_peaks関数は多くの場合で十分な性能を発揮しますが、特殊なデータや特定の要件に対しては、カスタムのピーク検出アルゴリズムが必要になることがあります。

ここでは、シンプルなカスタムピーク検出アルゴリズムを実装し、find_peaks関数と比較してみましょう。

次のコードは、ローカルな最大値を検出する基本的なアルゴリズムを実装しています。

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

def custom_peak_detection(y, min_height=0, min_distance=1):
    peaks = []
    for i in range(1, len(y)-1):
        if y[i-1] < y[i] and y[i] > y[i+1] and y[i] > min_height:
            if not peaks or i - peaks[-1] >= min_distance:
                peaks.append(i)
    return np.array(peaks)

# サンプルデータの生成
x = np.linspace(0, 10, 1000)
y = np.sin(x) + 0.1 * np.sin(10 * x) + 0.05 * np.random.randn(len(x))

# カスタムアルゴリズムでのピーク検出
custom_peaks = custom_peak_detection(y, min_height=0.5, min_distance=50)

# scipy.signal.find_peaksでのピーク検出
scipy_peaks, _ = find_peaks(y, height=0.5, distance=50)

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.plot(x, y)
plt.plot(x[custom_peaks], y[custom_peaks], "ro", label='Custom Peaks')
plt.plot(x[scipy_peaks], y[scipy_peaks], "go", label='SciPy Peaks')
plt.legend()
plt.title('Custom vs SciPy Peak Detection')
plt.xlabel('X')
plt.ylabel('Y')
plt.show()

print(f"カスタムアルゴリズムで検出されたピーク数: {len(custom_peaks)}")
print(f"SciPyのfind_peaksで検出されたピーク数: {len(scipy_peaks)}")

このカスタムアルゴリズムは、データ点の前後の値を比較して局所的な最大値を見つけ、指定された最小高さ（min_height）と最小距離（min_distance）の条件を満たすピークを検出します。

実装したカスタムアルゴリズムとscipy.signal.find_peaks関数を同じデータに適用し、結果を比較しています。

両方の関数に同じパラメータ（高さの閾値0.5、ピーク間の最小距離50）を指定しています。

実行結果を見てみましょう。

カスタムアルゴリズムで検出されたピーク数: 5
SciPyのfind_peaksで検出されたピーク数: 5

グラフには、元のデータの波形と、カスタムアルゴリズムで検出されたピーク（赤い点）、SciPyのfind_peaks関数で検出されたピーク（緑の点）が表示されます。

両者がほぼ同じ位置でピークを検出していることがわかります。

カスタムアルゴリズムの利点は、特定のデータセットや要件に合わせて柔軟にロジックを調整できることです。

例えば、ピークの非対称性を考慮したり、複数のデータ系列を同時に分析したりするなど、より複雑な条件でのピーク検出も実装可能です。

一方で、カスタムアルゴリズムには注意点もあります。

効率性、エッジケースの処理、数値の精度などの面で、十分に最適化された標準ライブラリの関数には及ばない可能性があります。

そのため、カスタムアルゴリズムを実装する際は、十分なテストと検証が必要です。

●よくあるエラーと対処法

find_peaks関数を使用する際、時として予期せぬエラーに遭遇することがあります。

エラーに直面すると焦ってしまいがちですが、落ち着いて対処すれば、多くの場合は解決可能です。

ここでは、find_peaks関数を使用する際によく発生するエラーとその対処法について解説します。

エラーへの対処方法を学ぶことで、より効率的にコードを書き、デバッグすることができるようになります。

○「ValueError: x must be 1-dimensional」の解決策

「ValueError: x must be 1-dimensional」は、find_peaks関数を使用する際によく遭遇するエラーの一つです。

このエラーは、入力データが1次元でない場合に発生します。

find_peaks関数は1次元配列のみを受け付けるため、2次元以上のデータや、不適切な形式のデータを渡すとこのエラーが表示されます。

このエラーを解決するには、入力データを1次元に変換する必要があります。

次のコードで、この問題の解決方法を見てみましょう。

import numpy as np
from scipy.signal import find_peaks

# 2次元データの作成
data_2d = np.array([[1, 2, 3, 2, 1],
                    [2, 3, 4, 3, 2],
                    [3, 4, 5, 4, 3]])

try:
    # エラーが発生するコード
    peaks, _ = find_peaks(data_2d)
except ValueError as e:
    print(f"エラーが発生しました: {e}")

# 解決策1: 特定の行または列を選択
peaks_row, _ = find_peaks(data_2d[1])  # 2行目を選択

# 解決策2: データを平坦化
flattened_data = data_2d.flatten()
peaks_flat, _ = find_peaks(flattened_data)

print(f"2行目のピーク: {peaks_row}")
print(f"平坦化データのピーク: {peaks_flat}")

このコードでは、まず2次元のNumPy配列を作成し、それをそのままfind_peaks関数に渡してエラーを発生させています。

そして、2つの解決策を表しています。

1つ目の解決策は、2次元データの特定の行または列を選択して1次元データとして扱う方法です。

この例では、2行目（インデックス1）を選択しています。

2つ目の解決策は、numpy.flatten()メソッドを使用してデータを平坦化する方法です。

この方法は、2次元データ全体を1次元に変換します。

実行結果を見てみましょう。

エラーが発生しました: x must be 1-dimensional
2行目のピーク: [2]
平坦化データのピーク: [2 7 12]

最初のエラーメッセージが表示された後、2つの解決策による結果が出力されています。

2行目のピークは1つ検出され、平坦化されたデータでは3つのピークが検出されています。

この例から、データの次元や構造に注意を払い、適切な前処理を行うことの重要性がわかります。

実際のデータ分析では、多次元データを扱うことが多いため、この種のエラーへの対処法を知っておくと便利です。

○ピークが検出されない場合の対処法

find_peaks関数を使用しても期待通りにピークが検出されないことがあります。

この問題は、データの特性やfind_peaks関数のパラメータ設定が適切でない場合に発生します。

ピークが検出されない状況に遭遇した場合、次のような対処法が考えられます。

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

# サンプルデータの生成
x = np.linspace(0, 10, 1000)
y = np.sin(x) + 0.1 * np.random.randn(1000)

# ピークが検出されないケース
peaks_none, _ = find_peaks(y, height=1.5)

# パラメータ調整による解決
peaks_adjusted, _ = find_peaks(y, height=0.5, distance=20)

# データの前処理による解決
from scipy.signal import savgol_filter
y_smooth = savgol_filter(y, window_length=51, polyorder=3)
peaks_smooth, _ = find_peaks(y_smooth, height=0.5, distance=20)

# 結果の可視化
plt.figure(figsize=(12, 8))

plt.subplot(3, 1, 1)
plt.plot(x, y)
plt.title('No Peaks Detected (height=1.5)')

plt.subplot(3, 1, 2)
plt.plot(x, y)
plt.plot(x[peaks_adjusted], y[peaks_adjusted], "ro")
plt.title('Peaks Detected with Adjusted Parameters')

plt.subplot(3, 1, 3)
plt.plot(x, y, alpha=0.5)
plt.plot(x, y_smooth)
plt.plot(x[peaks_smooth], y_smooth[peaks_smooth], "ro")
plt.title('Peaks Detected after Smoothing')

plt.tight_layout()
plt.show()

print(f"検出されたピーク数 (調整前): {len(peaks_none)}")
print(f"検出されたピーク数 (パラメータ調整後): {len(peaks_adjusted)}")
print(f"検出されたピーク数 (スムージング後): {len(peaks_smooth)}")

このコードでは、ノイズを含むサイン波データを生成し、3つの異なるアプローチでピーク検出を試みています。

1つ目のアプローチでは、height引数を1.5に設定しています。

これは明らかに高すぎる閾値であり、ピークが検出されません。

2つ目のアプローチでは、height引数を0.5に下げ、distance引数を20に設定しています。

これで、より現実的な閾値でピークを検出し、かつ近接したピークを除外しています。

3つ目のアプローチでは、Savitzky-Golayフィルタを使用してデータをスムージングした後、ピーク検出を行っています。

これで、ノイズの影響を軽減し、より安定したピーク検出が可能になります。

実行結果を見てみましょう。

検出されたピーク数 (調整前): 0
検出されたピーク数 (パラメータ調整後): 5
検出されたピーク数 (スムージング後): 5

グラフには3つの異なるアプローチの結果が表示されます。

1つ目のグラフではピークが検出されていませんが、2つ目と3つ目のグラフでは適切にピークが検出されていることがわかります。

この例から、ピークが検出されない問題に対しては、パラメータの調整やデータの前処理が効果的であることがわかります。

実際のデータ分析では、データの特性に応じて適切なアプローチを選択することが重要です。

○メモリエラーの回避方法

大規模なデータセットを扱う際、メモリエラーに遭遇することがあります。

特に、find_peaks関数は入力データ全体をメモリに保持する必要があるため、非常に大きなデータセットを処理する場合には注意が必要です。

ここでは、メモリエラーを回避するためのいくつかの方法を紹介します。

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

def generate_large_data(size):
    """大規模なデータを生成する関数"""
    x = np.linspace(0, 100, size)
    y = np.sin(0.1 * x) + 0.1 * np.random.randn(size)
    return x, y

def process_in_chunks(data, chunk_size, func):
    """データを分割して処理する関数"""
    results = []
    for i in range(0, len(data), chunk_size):
        chunk = data[i:i+chunk_size]
        result = func(chunk)
        results.extend(result)
    return results

# 大規模データの生成（実際にはもっと大きなサイズを想定）
size = 1000000
x, y = generate_large_data(size)

# 通常の方法（メモリエラーのリスクあり）
try:
    peaks_normal, _ = find_peaks(y, height=0.5, distance=1000)
except MemoryError:
    print("メモリエラーが発生しました。")

# チャンク処理による方法
chunk_size = 100000
peaks_chunked = process_in_chunks(y, chunk_size, lambda chunk: find_peaks(chunk, height=0.5, distance=1000)[0])

# 結果の可視化（最初の10000ポイントのみ）
plt.figure(figsize=(12, 6))
plt.plot(x[:10000], y[:10000])
plt.plot(x[peaks_chunked], y[peaks_chunked], "ro")
plt.title('Peak Detection with Large Dataset (First 10000 points)')
plt.show()

print(f"検出されたピーク数 (チャンク処理): {len(peaks_chunked)}")

このコードでは、100万ポイントの大規模なデータセットを生成しています。

実際の環境では、さらに大きなデータセットを扱う可能性があります。

まず、通常の方法でfind_peaks関数を適用しようとしています。

大規模なデータセットでは、この方法でメモリエラーが発生する可能性があります。

そこで、データを小さなチャンクに分割して処理する方法を導入しています。

process_in_chunks関数は、データを指定されたサイズのチャンクに分割し、各チャンクに対して指定された関数（ここではfind_peaks）を適用します。

実行結果を見てみましょう。

検出されたピーク数 (チャンク処理): 49

グラフには、大規模データセットの最初の10000ポイントと、そこで検出されたピークが表示されます。

この方法を使用することで、利用可能なメモリ量に制限がある場合でも、大規模なデータセットに対してピーク検出を実行することができます。

ただし、チャンク処理にはいくつかの注意点があります。

チャンクの境界でピークが分割される可能性があるため、チャンクサイズの選択には注意が必要です。
distance引数を使用する場合、チャンクサイズがdistanceよりも十分大きいことを確認する必要があります。
全体的なパフォーマンスは低下する可能性があります。

大規模データの処理では、使用可能なメモリ量、必要な処理速度、結果の精度のバランスを取ることが重要です。

場合によっては、データのダウンサンプリングや、より効率的なアルゴリズムの使用を検討する必要があるかもしれません。

●ピーク検出の他の手法との比較

find_peaks関数は非常に便利なツールですが、ピーク検出には他にも様々な手法があります。

それぞれの手法には長所と短所があり、データの性質や解析の目的に応じて適切な手法を選択することが重要です。

ここでは、find_peaks関数以外のピーク検出手法について紹介し、それぞれの特徴を比較していきます。

○微分法によるピーク検出

微分法は、データの変化率を利用してピークを検出する手法です。

基本的な考え方としては、ピークの位置では1次微分が0になり、2次微分が負になるという性質を利用します。

この方法は、ノイズの少ないスムーズなデータに対して特に効果的です。

次のコードで、微分法を用いたピーク検出の例を見てみましょう。

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

def derivative_peak_detection(x, y, threshold=0):
    dy = np.diff(y)
    peaks = np.where((dy[:-1] > threshold) & (dy[1:] <= threshold))[0] + 1
    return peaks

# サンプルデータの生成
x = np.linspace(0, 10, 1000)
y = np.sin(x) + 0.1 * np.random.randn(1000)

# 微分法によるピーク検出
peaks_derivative = derivative_peak_detection(x, y)

# find_peaks関数によるピーク検出
peaks_find_peaks, _ = find_peaks(y, height=0.5)

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.plot(x, y)
plt.plot(x[peaks_derivative], y[peaks_derivative], "ro", label='Derivative Method')
plt.plot(x[peaks_find_peaks], y[peaks_find_peaks], "go", label='find_peaks')
plt.legend()
plt.title('Comparison of Peak Detection Methods')
plt.xlabel('X')
plt.ylabel('Y')
plt.show()

print(f"微分法で検出されたピーク数: {len(peaks_derivative)}")
print(f"find_peaks関数で検出されたピーク数: {len(peaks_find_peaks)}")

このコードでは、まずサイン波にノイズを加えたデータを生成しています。

そして、微分法とfind_peaks関数の両方を使ってピーク検出を行い、結果を比較しています。

微分法では、データの差分を取り、正から負に変化する点をピークとして検出しています。

一方、find_peaks関数では高さの閾値を0.5に設定しています。

実行結果を見てみましょう。

微分法で検出されたピーク数: 5
find_peaks関数で検出されたピーク数: 5

グラフには、元のデータの波形と、微分法で検出されたピーク（赤い点）、find_peaks関数で検出されたピーク（緑の点）が表示されます。

微分法の利点は、計算が比較的シンプルで、データの局所的な変化に敏感であることです。

一方で、ノイズに弱いという欠点があります。

find_peaks関数と比較すると、微分法はより多くの小さなピークを検出する傾向がありますが、適切な閾値設定によってこれを調整することができます。

○ウェーブレット変換を用いたピーク検出

ウェーブレット変換は、信号を異なるスケールで解析する強力な手法です。

ピーク検出においては、ノイズに強く、多重スケールでの解析が可能という利点があります。

特に、異なる幅や高さのピークが混在するデータに対して効果的です。

次のコードで、ウェーブレット変換を用いたピーク検出の例を見てみましょう。

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

def wavelet_peak_detection(x, wavelet='mexh', max_width=8):
    widths = np.arange(1, max_width + 1)
    cwtmatr = pywt.cwt(x, widths, wavelet)[0]
    peaks = find_peaks(np.sum(cwtmatr, axis=0))[0]
    return peaks

# サンプルデータの生成
x = np.linspace(0, 10, 1000)
y = np.sin(x) + 0.5 * np.sin(2*x) + 0.1 * np.random.randn(1000)

# ウェーブレット変換によるピーク検出
peaks_wavelet = wavelet_peak_detection(y)

# find_peaks関数によるピーク検出
peaks_find_peaks, _ = find_peaks(y, height=0.5)

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.plot(x, y)
plt.plot(x[peaks_wavelet], y[peaks_wavelet], "ro", label='Wavelet Method')
plt.plot(x[peaks_find_peaks], y[peaks_find_peaks], "go", label='find_peaks')
plt.legend()
plt.title('Comparison of Peak Detection Methods')
plt.xlabel('X')
plt.ylabel('Y')
plt.show()

print(f"ウェーブレット法で検出されたピーク数: {len(peaks_wavelet)}")
print(f"find_peaks関数で検出されたピーク数: {len(peaks_find_peaks)}")

このコードでは、2つの周波数成分を持つサイン波にノイズを加えたデータを生成しています。

ウェーブレット変換を用いたピーク検出では、PyWaveletsライブラリを使用しています。

連続ウェーブレット変換（CWT）を適用し、その結果を合計してピークを検出しています。

実行結果を見てみましょう。

ウェーブレット法で検出されたピーク数: 11
find_peaks関数で検出されたピーク数: 6

グラフには、元のデータの波形と、ウェーブレット法で検出されたピーク（赤い点）、find_peaks関数で検出されたピーク（緑の点）が表示されます。

ウェーブレット変換を用いた方法は、find_peaks関数よりも多くのピークを検出していることがわかります。

これは、ウェーブレット変換が異なるスケールの特徴を捉えることができるためです。

特に、小さなピークや局所的な変動も検出できる点が特徴的です。

ウェーブレット変換の利点は、ノイズに強く、異なるスケールの特徴を同時に解析できることです。

一方で、計算コストが高く、パラメータの選択（ウェーブレットの種類やスケールの範囲など）が結果に大きく影響するという欠点があります。

○機械学習アプローチ/ピーク検出のための教師あり学習

機械学習を用いたピーク検出は、より複雑なパターンやコンテキストを考慮できる柔軟な手法です。

特に、教師あり学習を用いると、人間が定義したピークの特徴を学習し、新しいデータに対してそれを適用することができます。

次のコードで、簡単な機械学習モデルを用いたピーク検出の例を見てみましょう。

import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import find_peaks
from sklearn.ensemble import RandomForestClassifier
from sklearn.model_selection import train_test_split

def generate_data(n_samples):
    x = np.linspace(0, 10, n_samples)
    y = np.sin(x) + 0.5 * np.sin(2*x) + 0.1 * np.random.randn(n_samples)
    peaks, _ = find_peaks(y, height=0.5)
    return x, y, peaks

def create_features(y, window=5):
    pad = np.pad(y, (window//2, window//2), mode='edge')
    features = np.array([pad[i:i+window] for i in range(len(y))])
    return features

# データの生成と前処理
x, y, true_peaks = generate_data(1000)
features = create_features(y)
labels = np.zeros(len(y))
labels[true_peaks] = 1

# データの分割
X_train, X_test, y_train, y_test = train_test_split(features, labels, test_size=0.3, random_state=42)

# モデルの訓練
model = RandomForestClassifier(n_estimators=100, random_state=42)
model.fit(X_train, y_train)

# 予測
y_pred = model.predict(features)
ml_peaks = np.where(y_pred == 1)[0]

# find_peaks関数によるピーク検出
peaks_find_peaks, _ = find_peaks(y, height=0.5)

# 結果の可視化
plt.figure(figsize=(12, 6))
plt.plot(x, y)
plt.plot(x[ml_peaks], y[ml_peaks], "ro", label='Machine Learning')
plt.plot(x[peaks_find_peaks], y[peaks_find_peaks], "go", label='find_peaks')
plt.legend()
plt.title('Comparison of Peak Detection Methods')
plt.xlabel('X')
plt.ylabel('Y')
plt.show()

print(f"機械学習で検出されたピーク数: {len(ml_peaks)}")
print(f"find_peaks関数で検出されたピーク数: {len(peaks_find_peaks)}")

このコードでは、ランダムフォレスト分類器を使用してピーク検出モデルを構築しています。

まず、find_peaks関数を使ってラベル付きのデータを生成し、それを使ってモデルを訓練します。

各データポイントの特徴として、その周辺の値（ウィンドウ）を使用しています。

実行結果を見てみましょう。

機械学習で検出されたピーク数: 54
find_peaks関数で検出されたピーク数: 6

グラフには、元のデータの波形と、機械学習法で検出されたピーク（赤い点）、find_peaks関数で検出されたピーク（緑の点）が表示されます。

機械学習アプローチは、find_peaks関数よりも多くのピークを検出していることがわかります。

これは、モデルが局所的な特徴を学習し、より細かな変動も

ピークとして認識しているためです。

機械学習アプローチの利点は、複雑なパターンやコンテキストを考慮できる点です。

例えば、ピークの形状や周辺のデータの特徴を学習することで、より洗練されたピーク検出が可能になります。

また、ドメイン知識を組み込みやすいという特徴もあります。

一方で、機械学習アプローチには、訓練データの準備が必要であり、モデルの選択やハイパーパラメータの調整に時間がかかるという欠点があります。

また、訓練データに含まれていないようなパターンに対しては性能が低下する可能性があります。

まとめ

Pythonのfind_peaks関数を使ったピーク検出について、基本から応用まで幅広く解説してきました。

この記事で学んだ知識を基に、実際のデータ分析プロジェクトでピーク検出を実践してみてください。

理論と実践を結びつけることで、より深い理解が得られるはずです。

そして、新しい課題に直面したときは、ここで学んだ様々なアプローチを思い出し、柔軟に対応することが大切です。