読み込み中...

Pythonを使った簡単なピークフィッティングの手法・活用15選

ピークフィッティング 徹底解説 Python
この記事は約34分で読めます。

【サイト内のコードはご自由に個人利用・商用利用いただけます】

この記事では、プログラムの基礎知識を前提に話を進めています。

説明のためのコードや、サンプルコードもありますので、もちろん初心者でも理解できるように表現してあります。

本記事のサンプルコードを活用して機能追加、目的を達成できるように作ってありますので、是非ご活用ください。

※この記事は、一般的にプロフェッショナルの指標とされる『実務経験10,000時間以上』を満たす現役のプログラマチームによって監修されています。

※Japanシーモアは、常に解説内容のわかりやすさや記事の品質に注力しております。不具合、分かりにくい説明や不適切な表現、動かないコードなど気になることがございましたら、記事の品質向上の為にお問い合わせフォームにてご共有いただけますと幸いです。
(送信された情報は、プライバシーポリシーのもと、厳正に取扱い、処分させていただきます。)

●Pythonピークフィッティングとは？

データ解析において、ピークフィッティングは重要な手法です。

実験結果から意味ある情報を抽出する際に欠かせません。

Pythonを使うと、この複雑な作業が驚くほど簡単になります。

○ピークフィッティングの基本概念と重要性

ピークフィッティングとは、データ中の山なりの形状（ピーク）を数学的な関数で近似する技術です。

実験データや測定結果には、しばしばこのようなピークが現れます。

ピークの位置、高さ、幅を正確に把握することで、物質の特性や反応の挙動を深く理解できます。

例えば、分光分析では、特定の波長における光の吸収や放出がピークとして現れます。

このピークを解析することで、サンプル中の物質を同定したり、濃度を定量したりできるのです。

○Pythonを選ぶ5つの理由

Pythonがピークフィッティングに適している理由は多岐にわたります。

  1. 豊富なライブラリ -> NumPy、SciPy、Matplotlib など、データ処理に特化したツールが揃っています。
  2. 高い可読性 -> 初心者でも理解しやすい明快な構文を持っています。
  3. 無料で利用可能 -> オープンソースソフトウェアなので、コストを抑えられます。
  4. 大規模なコミュニティ -> 問題解決のためのリソースが豊富です。
  5. 柔軟性 -> 様々な分野のデータ解析に応用できます。

○必要なライブラリとセットアップ方法

Pythonでピークフィッティングを始めるには、いくつかのライブラリを準備する必要があります。

次のコマンドをコマンドプロンプトやターミナルで実行しましょう。

pip install numpy scipy matplotlib

このコマンドで、基本的なライブラリがインストールされます。

環境構築が完了したら、いよいよコーディングの開始です。

●ガウス関数によるピークフィッティング

ガウス関数は、ピークフィッティングでよく使われる関数の一つです。

その対称的な釣り鐘型の形状は、多くの自然現象や実験データのピークを表現するのに適しています。

○サンプルコード1：シンプルなガウスフィッティング

まずは、単一のガウス関数によるフィッティングを行ってみましょう。

次のコードは、ノイズを含むデータにガウス関数をフィッティングする例です。

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

def gaussian(x, amp, cen, wid):
    return amp * np.exp(-(x - cen)**2 / (2 * wid**2))

# データの生成
x = np.linspace(0, 10, 100)
y = gaussian(x, 1, 5, 2) + np.random.normal(0, 0.1, x.shape)

# フィッティング
popt, _ = curve_fit(gaussian, x, y)

# 結果のプロット
plt.scatter(x, y, label='データ')
plt.plot(x, gaussian(x, *popt), 'r-', label='フィッティング')
plt.legend()
plt.show()

print(f"振幅: {popt[0]:.2f}, 中心: {popt[1]:.2f}, 幅: {popt[2]:.2f}")

このコードを実行すると、ノイズを含むデータ点と、それにフィッティングされたガウス関数のグラフが表示されます。

また、フィッティングで得られたパラメータ（振幅、中心、幅）が出力されます。

○サンプルコード2：複数ピークのガウスフィッティング

実際のデータでは、複数のピークが重なっていることがよくあります。

そのような場合、複数のガウス関数の和でフィッティングを行います。

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

def double_gaussian(x, amp1, cen1, wid1, amp2, cen2, wid2):
    return (amp1 * np.exp(-(x - cen1)**2 / (2 * wid1**2)) +
            amp2 * np.exp(-(x - cen2)**2 / (2 * wid2**2)))

# データの生成
x = np.linspace(0, 10, 200)
y = double_gaussian(x, 1, 4, 0.5, 0.7, 6, 0.3) + np.random.normal(0, 0.05, x.shape)

# フィッティング
popt, _ = curve_fit(double_gaussian, x, y)

# 結果のプロット
plt.scatter(x, y, label='データ')
plt.plot(x, double_gaussian(x, *popt), 'r-', label='フィッティング')
plt.legend()
plt.show()

print(f"ピーク1 - 振幅: {popt[0]:.2f}, 中心: {popt[1]:.2f}, 幅: {popt[2]:.2f}")
print(f"ピーク2 - 振幅: {popt[3]:.2f}, 中心: {popt[4]:.2f}, 幅: {popt[5]:.2f}")

このコードでは、2つのガウス関数の和でデータをフィッティングしています。

結果として、2つのピークの位置と形状を同時に推定できます。

○サンプルコード3：ベースライン補正付きガウスフィッティング

実験データには、ピークの他にバックグラウンドノイズやベースラインが含まれることがあります。

正確なフィッティングのためには、そのようなベースラインも考慮する必要があります。

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

def gaussian_with_baseline(x, amp, cen, wid, a, b):
    return amp * np.exp(-(x - cen)**2 / (2 * wid**2)) + a * x + b

# データの生成
x = np.linspace(0, 10, 100)
y = gaussian_with_baseline(x, 1, 5, 1, 0.1, 0.5) + np.random.normal(0, 0.05, x.shape)

# フィッティング
popt, _ = curve_fit(gaussian_with_baseline, x, y)

# 結果のプロット
plt.scatter(x, y, label='データ')
plt.plot(x, gaussian_with_baseline(x, *popt), 'r-', label='フィッティング')
plt.plot(x, popt[3] * x + popt[4], 'g--', label='ベースライン')
plt.legend()
plt.show()

print(f"振幅: {popt[0]:.2f}, 中心: {popt[1]:.2f}, 幅: {popt[2]:.2f}")
print(f"ベースライン傾き: {popt[3]:.2f}, 切片: {popt[4]:.2f}")

このコードでは、ガウス関数に線形のベースラインを加えてフィッティングを行っています。

結果として、ピークの形状だけでなく、バックグラウンドの傾向も同時に推定できます。

●高度なフィッティング手法

ピークフィッティングの世界は奥深く、ガウス関数だけでは捉えきれない複雑なデータに遭遇することがあります。

そんな時こそ、より高度なフィッティング手法が威力を発揮します。

ローレンツ関数、Voigt関数、非線形最小二乗法など、様々な手法を駆使することで、より精密な解析が可能になります。

○サンプルコード4：ローレンツ関数によるフィッティング

ローレンツ関数は、ガウス関数よりも裾野が広がった形状を持ち、特定の物理現象や分光データの解析に適しています。

Pythonを使ってローレンツ関数によるフィッティングを行う方法を見てみましょう。

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

def lorentzian(x, amp, cen, wid):
    return amp * wid**2 / ((x - cen)**2 + wid**2)

# データの生成
x = np.linspace(0, 10, 100)
y = lorentzian(x, 1, 5, 0.5) + np.random.normal(0, 0.05, x.shape)

# フィッティング
popt, _ = curve_fit(lorentzian, x, y)

# 結果のプロット
plt.scatter(x, y, label='データ')
plt.plot(x, lorentzian(x, *popt), 'r-', label='フィッティング')
plt.legend()
plt.show()

print(f"振幅: {popt[0]:.2f}, 中心: {popt[1]:.2f}, 幅: {popt[2]:.2f}")

実行結果を見ると、ローレンツ関数特有の裾野の広がりがよくわかります。

振幅、中心、幅のパラメータが出力され、データの特性を数値で把握できます。

○サンプルコード5：Voigt関数を用いたピーク分離

Voigt関数は、ガウス関数とローレンツ関数の畳み込みで表される関数です。

複雑なスペクトルデータの解析に威力を発揮します。

Pythonでのモナコ公国をめぐるVoigt関数フィッティングの実装を見てみましょう。

import numpy as np
import matplotlib.pyplot as plt
from scipy.special import voigt_profile
from scipy.optimize import curve_fit

def voigt(x, amp, cen, sigma, gamma):
    return amp * voigt_profile(x - cen, sigma, gamma)

# データの生成
x = np.linspace(-10, 10, 200)
y = voigt(x, 1, 0, 0.5, 0.5) + np.random.normal(0, 0.05, x.shape)

# フィッティング
popt, _ = curve_fit(voigt, x, y, p0=[1, 0, 0.5, 0.5])

# 結果のプロット
plt.scatter(x, y, label='データ')
plt.plot(x, voigt(x, *popt), 'r-', label='フィッティング')
plt.legend()
plt.show()

print(f"振幅: {popt[0]:.2f}, 中心: {popt[1]:.2f}, シグマ: {popt[2]:.2f}, ガンマ: {popt[3]:.2f}")

Voigt関数を使うと、ガウス成分とローレンツ成分を同時に考慮できます。

実行結果では、複雑な形状のピークがうまくフィッティングされているのがわかります。

○サンプルコード6：非線形最小二乗法の実装

非線形最小二乗法は、より一般的なフィッティング手法です。

複雑な関数形でも柔軟に対応できる優れものです。

Pythonで非線形最小二乗法を実装してみましょう。

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import least_squares

def model(x, params):
    a, b, c = params
    return a * np.exp(-b * x) + c

def residuals(params, x, y):
    return model(x, params) - y

# データの生成
x = np.linspace(0, 10, 100)
true_params = [2.5, 0.5, 0.5]
y_true = model(x, true_params)
y = y_true + np.random.normal(0, 0.2, x.shape)

# フィッティング
initial_guess = [1.0, 1.0, 1.0]
res = least_squares(residuals, initial_guess, args=(x, y))

# 結果のプロット
plt.scatter(x, y, label='データ')
plt.plot(x, model(x, res.x), 'r-', label='フィッティング')
plt.legend()
plt.show()

print(f"最適化されたパラメータ: a={res.x[0]:.2f}, b={res.x[1]:.2f}, c={res.x[2]:.2f}")

非線形最小二乗法を使うと、複雑な関数形でもフィッティングが可能です。

実行結果を見ると、指数関数的な挙動をうまく捉えていることがわかります。

●scipyライブラリを活用したプロ級テクニック

scipyライブラリは、科学技術計算のための強力なツールボックスです。

ピークフィッティングにおいても、scipyの機能を活用することで、より高度な解析が可能になります。

○サンプルコード7：scipy.optimizeでのカスタムフィッティング

scipy.optimizeモジュールを使えば、より複雑なカスタム関数でのフィッティングが可能です。

例えば、非対称なピーク形状をフィッティングしてみましょう。

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

def asymmetric_peak(x, amp, cen, wid, asym):
    return amp * np.exp(-(x - cen)**2 / (2 * wid**2 * (1 + asym * (x - cen))))

# データの生成
x = np.linspace(-5, 5, 100)
y = asymmetric_peak(x, 1, 0, 1, 0.5) + np.random.normal(0, 0.05, x.shape)

# フィッティング
popt, _ = curve_fit(asymmetric_peak, x, y)

# 結果のプロット
plt.scatter(x, y, label='データ')
plt.plot(x, asymmetric_peak(x, *popt), 'r-', label='フィッティング')
plt.legend()
plt.show()

print(f"振幅: {popt[0]:.2f}, 中心: {popt[1]:.2f}, 幅: {popt[2]:.2f}, 非対称性: {popt[3]:.2f}")

実行結果を見ると、左右非対称なピーク形状がうまくフィッティングできています。

非対称性パラメータを導入することで、より柔軟なフィッティングが可能になりました。

○サンプルコード8：scipyのピーク検出機能の活用

scipyには、ピーク検出のための便利な関数が用意されています。

複数のピークがある場合、まずピークを検出してから個別にフィッティングする方法が有効です。

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

def gaussian(x, amp, cen, wid):
    return amp * np.exp(-(x - cen)**2 / (2 * wid**2))

# データの生成
x = np.linspace(0, 10, 1000)
y = gaussian(x, 1, 2, 0.5) + gaussian(x, 0.7, 5, 0.3) + gaussian(x, 0.5, 8, 0.2) + np.random.normal(0, 0.05, x.shape)

# ピーク検出
peaks, _ = find_peaks(y, height=0.1, distance=100)

# 各ピークをフィッティング
fit_results = []
for peak in peaks:
    peak_range = slice(max(0, peak-50), min(len(x), peak+50))
    popt, _ = curve_fit(gaussian, x[peak_range], y[peak_range], p0=[y[peak], x[peak], 0.5])
    fit_results.append(popt)

# 結果のプロット
plt.plot(x, y, label='データ')
plt.plot(x[peaks], y[peaks], "x", label='検出されたピーク')
for popt in fit_results:
    plt.plot(x, gaussian(x, *popt), '--', label=f'ピーク (中心: {popt[1]:.2f})')
plt.legend()
plt.show()

for i, popt in enumerate(fit_results):
    print(f"ピーク {i+1}: 振幅: {popt[0]:.2f}, 中心: {popt[1]:.2f}, 幅: {popt[2]:.2f}")

実行結果では、複数のピークが自動的に検出され、個別にフィッティングされているのがわかります。

各ピークの特性を簡単に把握できるため、複雑なスペクトルデータの解析に非常に便利です。

○サンプルコード9：2次元データのフィッティング

実際の実験データでは、2次元のピークフィッティングが必要になることがあります。

例えば、画像データや地形データなどがそうです。

scipyを使って2次元ガウス関数でのフィッティングを行ってみましょう。

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

def gaussian_2d(xy, amp, x0, y0, sigma_x, sigma_y):
    x, y = xy
    return amp * np.exp(-((x-x0)**2/(2*sigma_x**2) + (y-y0)**2/(2*sigma_y**2)))

# データの生成
x = np.linspace(-5, 5, 100)
y = np.linspace(-5, 5, 100)
x, y = np.meshgrid(x, y)

data = gaussian_2d((x, y), 1, 0, 0, 2, 1) + np.random.normal(0, 0.1, x.shape)

# フィッティング
params, _ = curve_fit(gaussian_2d, (x.ravel(), y.ravel()), data.ravel(), p0=[1, 0, 0, 1, 1])

# 結果のプロット
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 5))
ax1.imshow(data, cmap='viridis')
ax1.set_title('元データ')
ax2.imshow(gaussian_2d((x, y), *params).reshape(x.shape), cmap='viridis')
ax2.set_title('フィッティング結果')
plt.show()

print(f"振幅: {params[0]:.2f}, x0: {params[1]:.2f}, y0: {params[2]:.2f}, sigma_x: {params[3]:.2f}, sigma_y: {params[4]:.2f}")

実行結果では、2次元のガウス分布がうまくフィッティングされているのがわかります。

左側が元のノイズを含むデータ、右側がフィッティング結果です。

2次元フィッティングは、画像処理や空間データ解析など、様々な分野で活用できる重要な技術です。

●実データの処理と応用

研究室で得られた実験データ、山積みになった測定結果。

Pythonを駆使して、眠っているデータから宝物を掘り出す時が来ました。

実データの処理と応用は、ピークフィッティングの真価が発揮される場面です。

エクセルファイルやFITSファイル、ノイズだらけのデータも、Pythonなら軽々と解析できます。

さあ、データ解析の冒険に出発しましょう。

○サンプルコード10：エクセルデータの読み込みとフィッティング

研究室でよく使われるエクセルファイル。

Pythonを使えば、エクセルデータを簡単に読み込んでフィッティングできます。

pandasライブラリを使って、エクセルファイルからデータを読み込み、ガウス関数でフィッティングする方法を見てみましょう。

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

# エクセルファイルの読み込み
df = pd.read_excel('data.xlsx')
x = df['X'].values
y = df['Y'].values

# ガウス関数の定義
def gaussian(x, amp, cen, wid):
    return amp * np.exp(-(x - cen)**2 / (2 * wid**2))

# フィッティング
popt, _ = curve_fit(gaussian, x, y)

# 結果のプロット
plt.scatter(x, y, label='エクセルデータ')
plt.plot(x, gaussian(x, *popt), 'r-', label='フィッティング')
plt.legend()
plt.show()

print(f"振幅: {popt[0]:.2f}, 中心: {popt[1]:.2f}, 幅: {popt[2]:.2f}")

pandasを使うと、エクセルファイルがすんなり読み込めます。

データを読み込んだら、あとはお馴染みのガウスフィッティング。

実行結果を見ると、エクセルデータがきれいにフィッティングされているのがわかります。

振幅、中心、幅のパラメータも即座に算出されます。

研究データの解析が、グッと身近になりましたね。

○サンプルコード11：FITSファイルの解析とピーク検出

天文学や分光学でよく使われるFITSファイル。

Pythonを使えば、FITSファイルの解析も朝飯前です。

astropyライブラリを使って、FITSファイルからデータを読み込み、ピーク検出を行う方法を見てみましょう。

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

# FITSファイルの読み込み
hdul = fits.open('spectrum.fits')
data = hdul[1].data
wavelength = data['WAVELENGTH']
flux = data['FLUX']
hdul.close()

# ピーク検出
peaks, _ = find_peaks(flux, height=0.1*np.max(flux), distance=50)

# 結果のプロット
plt.plot(wavelength, flux, label='スペクトル')
plt.plot(wavelength[peaks], flux[peaks], "x", label='検出されたピーク')
plt.xlabel('波長')
plt.ylabel('フラックス')
plt.legend()
plt.show()

print(f"検出されたピークの数: {len(peaks)}")
for i, peak in enumerate(peaks):
    print(f"ピーク {i+1}: 波長 = {wavelength[peak]:.2f}, フラックス = {flux[peak]:.2f}")

astropyを使えば、FITSファイルがスルスルと読み込めます。

ピーク検出には、scipyのfind_peaks関数を使いました。

実行結果を見ると、スペクトル中のピークが自動的に検出されているのがわかります。

各ピークの波長とフラックスも即座に算出されます。天文学的な発見も、もう目の前です。

○サンプルコード12：ノイズの多いデータのロバストフィッティング

実験データには、往々にしてノイズがつきものです。

ノイズに惑わされず、真のシグナルを捉えるには、ロバストなフィッティング手法が必要です。

scipyのoptimizeモジュールを使って、ノイズに強いロバストフィッティングを行う方法を見てみましょう。

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import least_squares

def model(x, params):
    return params[0] * np.exp(-((x - params[1])**2) / (2 * params[2]**2))

def residuals(params, x, y):
    return (y - model(x, params)) / (np.abs(y) + 0.1)

# ノイズの多いデータの生成
x = np.linspace(0, 10, 100)
y_true = model(x, [1, 5, 1])
noise = np.random.normal(0, 0.2, x.shape) + np.random.normal(0, 1, x.shape) * (np.random.rand(100) < 0.1)
y = y_true + noise

# ロバストフィッティング
res = least_squares(residuals, [1, 5, 1], args=(x, y), loss='soft_l1')

# 結果のプロット
plt.scatter(x, y, label='ノイズデータ')
plt.plot(x, y_true, 'g-', label='真の関数')
plt.plot(x, model(x, res.x), 'r-', label='ロバストフィッティング')
plt.legend()
plt.show()

print(f"振幅: {res.x[0]:.2f}, 中心: {res.x[1]:.2f}, 幅: {res.x[2]:.2f}")

least_squares関数を使うと、ノイズに強いロバストフィッティングが可能です。

lossパラメータを’soft_l1’に設定することで、外れ値の影響を軽減しています。

実行結果を見ると、ノイズだらけのデータの中から、真の関数がきれいに抽出されているのがわかります。

振幅、中心、幅のパラメータも、ノイズにほとんど影響されずに算出されています。ノイズに負けない、タフなデータ解析の極意です。

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

Pythonでピークフィッティングを行っていると、時として思わぬエラーに遭遇します。

でも、心配ありません。エラーは、プログラミングの道を歩む者にとって、成長のチャンス。

よくあるエラーとその対処法を学べば、データ解析の腕前がグッと上がります。

○「ValueError: math domain error」の解決策

「ValueError: math domain error」。

数学の領域外、つまり計算不能な状態に陥った時に出るエラーです。

よくあるケースは、負の数の平方根や、対数の引数が0以下の場合です。

対処法を見てみましょう。

import numpy as np
from scipy.optimize import curve_fit

def model(x, a, b):
    return a * np.sqrt(x) + b

x = np.array([-1, 0, 1, 2, 3, 4])
y = np.array([1, 1, 2, 3, 4, 5])

try:
    popt, _ = curve_fit(model, x, y)
except ValueError as e:
    print(f"エラーが発生しました: {e}")
    x = np.maximum(x, 0)  # 負の値を0に置き換え
    popt, _ = curve_fit(model, x, y)

print(f"フィッティング結果: a = {popt[0]:.2f}, b = {popt[1]:.2f}")

エラーの原因は、sqrt関数に負の値が入力されたこと。対処法は、負の値を0に置き換えるなど、入力データの前処理を行うことです。

実行結果を見ると、エラーをキャッチして適切に対処し、フィッティングが成功しているのがわかります。

エラーも、上手く付き合えば怖くありません。

○フィッティングが収束しない場合の対処法

フィッティングが収束しない。データ解析者にとって、頭を抱えるような問題です。

でも、諦めるのはまだ早い。

初期値の調整や、制約条件の追加で、収束しないフィッティングを救う方法があります。

import numpy as np
from scipy.optimize import curve_fit

def model(x, a, b, c):
    return a * np.exp(-b * x) + c

x = np.linspace(0, 10, 100)
y = model(x, 1, 0.1, 0.5) + np.random.normal(0, 0.1, x.shape)

# 不適切な初期値でのフィッティング
try:
    popt, _ = curve_fit(model, x, y, p0=[100, 100, 100], maxfev=1000)
except RuntimeError as e:
    print(f"エラーが発生しました: {e}")

# 適切な初期値と制約条件を設定したフィッティング
popt, _ = curve_fit(model, x, y, p0=[1, 0.1, 0.5], bounds=([0, 0, 0], [10, 1, 1]))

print(f"フィッティング結果: a = {popt[0]:.2f}, b = {popt[1]:.2f}, c = {popt[2]:.2f}")

初期値が不適切だと、フィッティングが収束しないことがあります。

また、パラメータに物理的な制約がある場合、bounds引数で範囲を指定すると良いでしょう。

実行結果を見ると、適切な初期値と制約条件を設定することで、フィッティングが成功しているのがわかります。

粘り強さが、データ解析の成功を左右します。

○外れ値がある場合のデータクリーニング技術

外れ値。

データ解析の大敵です。

でも、適切なデータクリーニング技術を使えば、外れ値にも怯まない頑健な解析ができます。

外れ値を検出し、除去する方法を見てみましょう。

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

# 外れ値を含むデータの生成
np.random.seed(0)
x = np.linspace(0, 10, 100)
y = 2 * x + 1 + np.random.normal(0, 1, 100)
y[80] = 30  # 外れ値の追加

# 外れ値の検出と除去
z_scores = np.abs(stats.zscore(y))
y_cleaned = y[z_scores < 3]
x_cleaned = x[z_scores < 3]

# 線形回帰
slope, intercept, _, _, _ = stats.linregress(x_cleaned, y_cleaned)

# 結果のプロット
plt.scatter(x, y, label='元のデータ')
plt.scatter(x_cleaned, y_cleaned, label='クリーニング後のデータ')
plt.plot(x, slope * x + intercept, 'r-', label='回帰直線')
plt.legend()
plt.show()

print(f"回帰直線: y = {slope:.2f}x + {intercept:.2f}")

z-scoreを使って外れ値を検出し、除去しています。

実行結果を見ると、外れ値が適切に除去され、より正確な回帰直線が得られているのがわかります。

データクリーニングは、信頼性の高い解析結果を得るための重要なステップです。

外れ値に惑わされない、冷静な目を持つことが大切です。

●ピークフィッティングの応用例

ピークフィッティングは、様々な分野で活用される強力な解析手法です。

スペクトル解析から画像処理、さらには機械学習との融合まで、応用範囲は広大です。

実践的な例を通じて、ピークフィッティングの可能性を探ってみましょう。

○サンプルコード13：スペクトル解析での活用

分光分析は、ピークフィッティングが真価を発揮する分野です。

複雑なスペクトルから、個々の成分を分離し、定量分析を行うことができます。

Pythonを使って、多成分スペクトルの解析を行う方法を見てみましょう。

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

def multi_peak(x, *params):
    y = np.zeros_like(x)
    for i in range(0, len(params), 3):
        amp, cen, wid = params[i:i+3]
        y += amp * np.exp(-(x - cen)**2 / (2 * wid**2))
    return y

# スペクトルデータの生成
x = np.linspace(0, 10, 1000)
true_params = [1, 2, 0.5, 0.7, 5, 0.3, 0.5, 8, 0.2]
y = multi_peak(x, *true_params) + np.random.normal(0, 0.05, x.shape)

# フィッティング
popt, _ = curve_fit(multi_peak, x, y, p0=true_params)

# 結果のプロット
plt.plot(x, y, 'b.', label='観測データ')
plt.plot(x, multi_peak(x, *popt), 'r-', label='フィッティング')
for i in range(0, len(popt), 3):
    plt.plot(x, multi_peak(x, *popt[i:i+3]), '--', label=f'ピーク {i//3 + 1}')
plt.legend()
plt.show()

for i in range(0, len(popt), 3):
    print(f"ピーク {i//3 + 1}: 振幅 = {popt[i]:.2f}, 中心 = {popt[i+1]:.2f}, 幅 = {popt[i+2]:.2f}")

多成分スペクトルの解析は、まるで宝探しのようです。

複数のガウス関数を組み合わせて、複雑なスペクトルを再現しています。

実行結果を見ると、重なり合ったピークがきれいに分離されているのがわかります。

各ピークの振幅、中心、幅が数値として出力され、定量分析の基礎となります。

○サンプルコード14：画像処理におけるピーク検出

画像処理の分野でも、ピーク検出は重要な役割を果たします。

例えば、天体画像から星を検出したり、顕微鏡画像から細胞核を識別したりするのに使われます。

2次元画像でのピーク検出を、Pythonで実装してみましょう。

import numpy as np
import matplotlib.pyplot as plt
from scipy import ndimage as ndi
from skimage.feature import peak_local_max

# 2D画像データの生成
x, y = np.meshgrid(np.linspace(-4, 4, 100), np.linspace(-4, 4, 100))
d = np.sqrt(x*x + y*y)
sigma, mu = 0.5, 0.0
g = np.exp(-( (d-mu)**2 / ( 2.0 * sigma**2 ) ) )
g += np.exp(-( ((x-2)-mu)**2 + ((y-2)-mu)**2 ) / ( 2.0 * sigma**2 ) )
g += np.random.random(g.shape) * 0.1

# ピーク検出
coordinates = peak_local_max(g, min_distance=20)

# 結果のプロット
fig, ax = plt.subplots()
ax.imshow(g, cmap=plt.cm.gray)
ax.plot(coordinates[:, 1], coordinates[:, 0], 'r.')
plt.show()

print(f"検出されたピークの数: {len(coordinates)}")
for i, coord in enumerate(coordinates):
    print(f"ピーク {i+1}: y = {coord[0]}, x = {coord[1]}")

2次元画像でのピーク検出は、まるで星座を見つけるようです。

skimageライブラリのpeak_local_max関数を使うと、局所的な最大値を簡単に検出できます。

実行結果を見ると、画像中の明るい部分（ピーク）が自動的に検出されているのがわかります。

天体画像解析や細胞画像解析など、画像処理の可能性が広がります。

○サンプルコード15：機械学習との組み合わせ

ピークフィッティングと機械学習を組み合わせると、より高度なデータ解析が可能になります。

例えば、スペクトルデータの自動分類を行ってみましょう。

Pythonを使って、ピークフィッティングと機械学習を融合させる方法を見てみます。

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
from sklearn.cluster import KMeans

def gaussian(x, amp, cen, wid):
    return amp * np.exp(-(x - cen)**2 / (2 * wid**2))

# スペクトルデータの生成
x = np.linspace(0, 10, 100)
spectra = []
labels = []
for _ in range(50):
    if np.random.rand() < 0.5:
        y = gaussian(x, 1, 3, 0.5) + np.random.normal(0, 0.1, x.shape)
        labels.append(0)
    else:
        y = gaussian(x, 0.7, 7, 0.3) + np.random.normal(0, 0.1, x.shape)
        labels.append(1)
    spectra.append(y)

# ピークフィッティングと特徴抽出
features = []
for spectrum in spectra:
    popt, _ = curve_fit(gaussian, x, spectrum, p0=[1, 5, 1])
    features.append(popt)

# K-means クラスタリング
kmeans = KMeans(n_clusters=2, random_state=42)
predicted_labels = kmeans.fit_predict(features)

# 結果のプロット
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 5))
for spectrum, label in zip(spectra, predicted_labels):
    ax1.plot(x, spectrum, color='red' if label == 0 else 'blue', alpha=0.1)
ax1.set_title('クラスタリングされたスペクトル')

ax2.scatter([f[1] for f in features], [f[0] for f in features], c=predicted_labels)
ax2.set_xlabel('ピーク位置')
ax2.set_ylabel('ピーク強度')
ax2.set_title('特徴空間でのクラスタリング')

plt.show()

accuracy = np.mean(np.array(labels) == predicted_labels)
print(f"分類精度: {accuracy:.2f}")

機械学習とピークフィッティングの融合は、まるで錬金術のようです。

まず、ガウスフィッティングで各スペクトルの特徴（振幅、中心、幅）を抽出し、それをK-meansクラスタリングの入力として使用しています。

実行結果を見ると、2種類のスペクトルが自動的に分類されているのがわかります。

左側のグラフは分類されたスペクトル、右側は特徴空間でのクラスタリング結果です。

分類精度も高く、スペクトルの自動分類に成功しています。

ピークフィッティングと機械学習の組み合わせは、データ解析の新たな地平を切り開きます。

まとめ

ピークフィッティングは、データ解析の強力な武器です。

Pythonを使えば、複雑なピークフィッティングも簡単に実装できます。

本記事で学んだテクニックを活用して、皆さんの研究やプロジェクトのレベルを引き上げてみましょう。