lsq_peakfit_test.py
このドキュメントは、Pythonプログラム lsq_peakfit_test.py のソースコードを解析し、Sphinx (MyST) で問題なくビルド可能なMarkdown形式で記述したものです。
概要
lsq_peakfit_test.py は、擬Voigt関数と線形ベースラインを組み合わせたモデルを使用して、データに対するピークフィットを実行するPythonスクリプトです。データファイルの読み込み、モデルのシミュレーション、線形最小二乗フィット、非線形最小二乗フィット、およびモデル選択（ eta パラメータの固定値比較）の各モードを提供します。このプログラムは、主に numpy, pandas, matplotlib, scipy と、カスタムライブラリである tklsq を利用しています。
非標準ライブラリ
本プログラムは、以下の非標準ライブラリに依存しています。
tklsq: データ入出力 (tkdataio)、パラメータ管理 (tkparamio)、最小二乗フィットのヘルパー関数 (tkminfit,tknlsq)、フィット診断 (tkfitdiag) など、このプログラムの中核となる機能を提供するカスタムパッケージです。
入出力
入力
データファイル (
--infile):XY形式のデータを含むファイル。
tklsq.tkdataio.read_xy()関数によって読み込まれるため、CSVファイルやExcelファイルなど、対応する形式が利用可能です。読み込み時に
--xcol,--ycol,--sheet,--xmin,--xmaxオプションでデータの範囲や列を指定できます。
パラメータファイル (
--paramfile):フィットパラメータの初期値や制約などを定義するCSVファイル。
指定がない場合、
<infileのステム>_params.csvという名前で自動的に生成されます。プログラム内部では、
DEFAULT_PARAMSで定義された構造とデフォルト値が使用されます。
出力
画像ファイル:
--saveオプションが1の場合、フィット結果のプロットや最適化過程のアニメーションが画像ファイルとして出力されます。ファイル名は、
--outfileオプションで指定されたプレフィックスとモードに応じたサフィックス (例:<outfile>_read.png,<outfile>_sim.png,<outfile>_lfit.png,<outfile>_fit.png,<outfile>_fit_animation.gif) が結合されたものになります。
Excelファイル:
--saveオプションが1の場合、フィットの詳細結果や診断情報がExcelファイルとして出力されます。ファイル名は、
--outfileオプションで指定されたプレフィックスとモードに応じたサフィックス (例:<outfile>_fit.xlsx,<outfile>_model_select.xlsx) が結合されたものになります。
パラメータファイル (
--paramfile):run_lfit()またはrun_fit()モードで--saveオプションが1の場合、フィット後のパラメータがこのCSVファイルに更新されて保存されます。
標準出力:
プログラムの実行状況、フィットの進捗（コールバックによる）、最終的なフィット結果、パラメータ値とその標準誤差、共分散行列の診断情報などが表示されます。
コマンドライン引数
プログラムは argparse モジュールを使用してコマンドライン引数を解析します。
--mode <MODE>必須。実行するモードを指定します。
選択肢:
read,sim,lfit,fit,model_select
--infile <FILE>必須。入力データファイルへのパスを指定します。
--paramfile <FILE>パラメータCSVファイルのパスを指定します。
デフォルト:
<infileのステム>_params.csv
--outfile <PREFIX>出力ファイル名のプレフィックスを指定します。
デフォルト:
<infileのステム>
--save <0|1>(
int型) 結果をファイルに保存するかどうかを指定します。0: 保存しない,1: 保存する。デフォルト:
1
--show <0|1>(
int型) 結果のプロットを表示するかどうかを指定します。0: 表示しない,1: 表示する。デフォルト:
1
--xcol <COLUMN>入力データファイルからX軸のデータを読み込む列名またはインデックスを指定します。
デフォルト:
0
--ycol <COLUMN>入力データファイルからY軸のデータを読み込む列名またはインデックスを指定します。
デフォルト:
1
--sheet <SHEET>Excelファイルの場合、データを読み込むシート名またはインデックスを指定します。
デフォルト:
0
--xmin <VALUE>(
float型) データを読み込むX軸の下限値を指定します。デフォルト:
-1.0e100
--xmax <VALUE>(
float型) データを読み込むX軸の上限値を指定します。デフォルト:
1.0e100
--nplot_interval <INTERVAL>(
int型)fitモードで、最適化過程のアニメーションフレームを保存する間隔（コールバック呼び出し回数）を指定します。デフォルト:
10
--method <METHOD>scipy.optimize.minimize()で使用する最適化手法を指定します。デフォルト:
Nelder-Mead
--maxiter <ITER>(
int型) 最適化の最大イテレーション数を指定します。デフォルト:
1000
--xatol <TOL>(
float型) 最適化におけるxの変化に関する許容誤差を指定します。デフォルト:
1.0e-9
--fatol <TOL>(
float型) 最適化における目的関数の変化に関する許容誤差を指定します。デフォルト:
1.0e-9
--criterion <CRITERION>model_selectモードで、モデル選択の基準を指定します。選択肢:
AIC,BICデフォルト:
BIC
プログラムの機能
主要な定数
DEFAULT_PARAMS擬Voigtピークモデル (
amp,x0,fwhm,eta) および線形ベースラインモデル (b0,b1) の初期パラメータを定義するリストです。各パラメータは辞書形式で、以下のキーを持ちます。
キー
説明
varnameパラメータ名。
optid非線形最適化において、このパラメータを最適化対象とするか (
1) しないか (0) を指定します。optid_lin線形最適化において、このパラメータを最適化対象とするか (
1) しないか (0) を指定します。p0パラメータの初期値。
pscaleパラメータの最適化スケール (
linearまたはlog)。dp初期シンプレックスの構築や数値微分に使用されるステップサイズ。
pminパラメータの下限値。
pmaxパラメータの上限値。
kpenalty境界条件からの逸脱に対するペナルティの係数。
主要な関数
auto_stem(infile)入力ファイルパス (
infile) から拡張子を除いたファイル名を返します。
pseudo_voigt_unit(x, x0, fwhm, eta)擬Voigt関数の正規化された形状を計算します。これはガウス関数とローレンツ関数の加重和として定義されます。
fwhmは半値幅、etaはローレンツ成分の割合 (0.0: ガウス, 1.0: ローレンツ) を示します。内部で
fwhmは最小値1.0e-300に、etaはfloat型にそれぞれ変換されます。
model(x, p)擬Voigtピークと線形ベースラインの合計を計算するモデル関数です。
ピーク部分は
p["amp"] * pseudo_voigt_unit(x, p["x0"], p["fwhm"], p["eta"])で計算されます。ベースラインは
p.get("b0", 0.0) + p.get("b1", 0.0) * xで計算されます。
read_params(args)--paramfileで指定されたパラメータファイルを読み込み、DEFAULT_PARAMSをデフォルト値として使用します。ファイルが存在しない場合は新規作成されます。
load_data(args)--infileで指定されたデータファイルを読み込み、X軸、Y軸のデータ、およびラベルを返します。--xcol,--ycol,--sheet,--xmin,--xmaxオプションが適用されます。
free_names_all_optid1(params)渡されたパラメータ辞書
paramsから、optid=1が設定されているすべてのパラメータ名をリストとして抽出して返します。
make_progress_figure(),update_progress_plot(),save_final_plot(),save_animation()matplotlibを用いたプロットおよびアニメーション生成に関連する関数群です。update_progress_plot()は、フィットの進行状況をリアルタイムで更新するためのプロット処理を行います。save_final_plot()は、最終的なフィット結果のプロットを生成・保存・表示します。save_animation()は、run_fit()モードでの最適化過程をアニメーションGIFとして保存します。
compute_fit_band(x, p_fit, free_names, cov_free)フィットされたモデル (
p_fit) と、最適化された自由パラメータの共分散行列 (cov_free) を用いて、フィット曲線の ±1σ 信頼帯を計算します。デルタ法が使用されます。
モードごとの詳細
run_read(args) モード
機能:
--infileで指定されたデータファイルを読み込み、散布図として表示します。出力:
プロット表示 (
--show 1の場合)。<outfile>_read.png(--save 1の場合)。
run_sim(args) モード
機能:
--paramfileからパラメータを読み込み、そのパラメータでmodel()関数を使用してシミュレーション曲線を生成します。データとシミュレーション曲線を比較するプロットを生成します。出力:
プロット表示 (
--show 1の場合)。<outfile>_sim.png(--save 1の場合)。
run_lfit(args) モード
機能: 擬Voigt形状 (
x0,fwhm,eta) を固定し、DEFAULT_PARAMSでoptid_lin=1が設定されているパラメータ (amp,b0,b1) のみを対象とした線形最小二乗フィットを実行します。tklsq.tkminfit.solve_linear_block()を使用します。出力:
フィット結果のパラメータと標準誤差が標準出力に表示されます。
プロット表示 (
--show 1の場合)。<outfile>_lfit.png(--save 1の場合)。フィット後のパラメータが
--paramfileに更新されて保存されます (--save 1の場合)。
run_fit(args) モード
機能:
DEFAULT_PARAMSでoptid=1が設定されているすべてのパラメータを対象とした非線形最小二乗フィットを実行します。scipy.optimize.minimize()を使用します。最適化の進行状況は
callback()関数を通じてリアルタイムでプロット表示され (--show 1の場合)、アニメーションフレームが収集されます。フィット完了後、残差二乗和 (RSS)、目的関数値、自由度、分散、最適化されたパラメータとその標準誤差、共分散行列の診断情報 (相関行列、JTJ行列の条件数と固有値) が標準出力に表示されます。
出力:
プロット表示 (
--show 1の場合)。<outfile>_fit.png(--save 1の場合): 最終フィット結果（データ、フィット曲線、±1σバンド、残差）のプロット。<outfile>_fit_animation.gif(--save 1の場合): 最適化過程のアニメーション。<outfile>_fit.xlsx(--save 1の場合): フィットされたデータ、パラメータ、相関行列、JTJ固有値を含むExcelファイル。フィット後のパラメータが
--paramfileに更新されて保存されます (--save 1の場合)。
run_model_select(args) モード
機能: 擬Voigt関数の
etaパラメータを[0.0, 0.25, 0.5, 0.75, 1.0]の候補値にそれぞれ固定し、他のoptid=1のパラメータに対して非線形フィットを実行します。各フィットの結果に基づいて、残差二乗和 (RSS)、AIC (Akaike Information Criterion)、BIC (Bayesian Information Criterion) を計算し、比較表を標準出力に表示します。表は
--criterionオプションで指定された基準（AICまたはBIC）に基づいてソートされます。出力:
モデル選択結果の表が標準出力に表示されます。
<outfile>_model_select.xlsx(--save 1の場合): モデル選択の比較表を含むExcelファイル。
プログラムの実行方法
以下の例は、 lsq_peakfit_test.py を実行する方法を示しています。
# readモードでデータを表示し、画像を保存
python lsq_peakfit_test.py --mode read --infile my_data.csv --save 1 --show 1
# simモードでパラメータファイルからモデルをシミュレーションし、データと比較
python lsq_peakfit_test.py --mode sim --infile my_data.csv --paramfile my_params.csv --save 1
# lfitモードで線形最小二乗フィットを実行し、結果を保存・表示
python lsq_peakfit_test.py --mode lfit --infile my_data.csv --paramfile my_params.csv --save 1 --show 1
# fitモードで非線形最小二乗フィットを実行し、結果を保存・表示し、アニメーションも生成
python lsq_peakfit_test.py --mode fit --infile my_data.csv --paramfile my_params.csv --save 1 --show 1 --nplot_interval 5
# model_selectモードで異なるeta値でのフィットを比較し、結果を保存
python lsq_peakfit_test.py --mode model_select --infile my_data.csv --paramfile my_params.csv --save 1 --criterion AIC