lsq_func_tkFit.py 技術ドキュメント
プログラムの動作
lsq_func_tkFit.py は、Pythonで記述された汎用的な非線形最小二乗フィッティングプログラムです。このプログラムは、与えられたデータセットに対し、ユーザーが定義した関数モデルをフィッティングし、モデルの最適パラメータを決定することを目的としています。
主な機能は以下の通りです。
データフィッティング: 指定されたExcelファイルからデータを読み込み、コマンドライン引数で定義された数式にフィッティングします。
最適化アルゴリズムの選択: BFGS法、Nelder-Mead法、Powell法、共役勾配法 (CG) など、複数の最適化アルゴリズムから選択して使用できます。
初期パラメータとフィッティング範囲の指定: ユーザーはフィッティングの初期パラメータと、各独立変数のフィッティング範囲を細かく設定できます。
結果の出力: フィッティングされたパラメータ、適合度、収束過程、およびフィッティング結果をExcelファイルやパラメータファイルとして出力します。
グラフ表示: フィッティングの進捗と最終結果をリアルタイムでグラフ表示し、視覚的に確認できます。
このプログラムは、実験データやシミュレーションデータに対して、特定の物理モデルや数学モデルを当てはめ、そのモデルが持つ未知のパラメータを客観的に決定するという課題を解決します。例えば、スペクトルデータへのガウス関数フィッティングや、時間減衰データへの指数関数フィッティングなどに利用できます。
原理
lsq_func_tkFit.py は、非線形最小二乗法に基づいています。この方法の目的は、観測データとモデル関数との間の残差の二乗和を最小化するようなモデルパラメータ \(\mathbf{p}\) を見つけることです。
目的関数: プログラムは、与えられたモデル関数 \(f(\mathbf{x}; \mathbf{p})\) と観測データ \(( \mathbf{x}_i, y_i )\) に対し、以下の残差二乗和 (Sum of Squared Errors, SSE) を最小化します。
\[S^2 = \sum_{i=1}^N (y_i - f(\mathbf{x}_i; \mathbf{p}))^2\]ここで、\(N\) はデータ点の総数、\(y_i\) は \(i\) 番目の測定された従属変数の値、\(f(\mathbf{x}_i; \mathbf{p})\) は \(i\) 番目の独立変数 \(\mathbf{x}_i\) とパラメータ \(\mathbf{p}\) を用いたモデル関数の予測値です。
最適化アルゴリズム: 目的関数 \(S^2\) を最小化するために、
scipy.optimize.minimize関数が使用されます。このプログラムは、以下の最適化メソッドをサポートしており、コマンドライン引数で選択可能です。BFGS (Broyden–Fletcher–Goldfarb–Shanno) 法: 準ニュートン法の一種で、ヘッセ行列を近似的に計算し、その逆行列を用いて探索方向を決定します。勾配情報に基づいて効率的に最適解に収束することが特徴です。デフォルトのメソッドとして設定されています。
Nelder-Mead (Downhill Simplex) 法: 導関数を必要としない直接探索法です。シンプレックス（単体）と呼ばれる幾何学的な図形を動かしながら最適解を探します。
Powell 法: 導関数を必要としない最適化アルゴリズムです。共役方向探索に基づき、各軸に沿って最小化を繰り返します。
CG (Conjugate Gradient) 法 (Polak-Ribiere method): 勾配情報を使用し、共役な方向へ探索を進めることで、効率的な最適化を実現します。
モデル関数の評価: モデル関数は、
cparams.funcとして文字列で与えられ、Pythonのeval()関数を使用して実行時に評価されます。この関数は、独立変数リストxとパラメータリストpを入力として受け取ります。例えば、ガウス関数はp[0]*exp(-((x[0]-p[1])/p[2])**2)のように記述されます。データ処理: 入力データは
pandasライブラリを使用してExcelファイルから読み込まれ、numpy配列として処理されます。指定されたフィッティング範囲外のデータは、フィッティング計算から除外されます。
必要な非標準ライブラリとインストール方法
このプログラムを実行するには、以下の非標準ライブラリが必要です。
numpy: 数値計算のための基本的なライブラリです。
pandas: データ構造とデータ解析ツールを提供するライブラリです。特にExcelファイルの読み書きに使用されます。
scipy: 科学技術計算のためのライブラリで、このプログラムでは特に最適化 (
scipy.optimize.minimize) に使用されます。matplotlib: グラフ描画のためのライブラリです。
tklib: このプログラムの内部で利用されているカスタムユーティリティモジュール群です (
tkutils,tkapplication,tkparams,tksci.tksci,tksci.tkFit,tksci.tkFit_mなど)。これらは通常、標準的なpipコマンドではインストールできません。本プログラムと同時に提供されるか、別途入手・設定が必要なライブラリです。
tklib を除く必要なライブラリは、以下のコマンドでインストールできます。
pip install numpy pandas scipy matplotlib
tklib については、本プログラムの実行環境に事前に配置されている必要があります。
必要な入力ファイル
プログラムは以下の形式の入力ファイルを期待します。
ファイル形式: Excelファイル (
.xlsx)ファイル名: コマンドライン引数 (
infile) で指定します。デフォルトはpeak.xlsxです。データ構造: データは複数列で構成され、少なくとも1列の独立変数（x軸）と1列の従属変数（y軸）を含んでいる必要があります。列のヘッダーは、コマンドライン引数
xlabel,ylabel,zlabel,olabelで指定されたラベルと一致する必要があります。例:
peak.xlsxファイルの想定される内容（CSV形式で記述しますが、Excelファイルとして保存してください）X_col,Y_col -2.0,0.1 -1.5,0.3 -1.0,0.8 -0.5,1.2 0.0,1.3 0.5,1.2 1.0,0.8 1.5,0.3 2.0,0.1
この例では、
X_colが独立変数、Y_colが従属変数として使用されます。
生成される出力ファイル
プログラムの実行が完了すると、以下のファイルが生成されます。ファイル名は、入力ファイル名に基づいて自動的に生成されます（例: peak.xlsx の場合、peak-fitting.xlsx となります）。
{infile_body}-fitting.xlsx:内容: フィッティング結果の詳細をまとめたExcelファイル。元の入力データ、初期パラメータで計算された関数値 (
initial)、最終パラメータで計算された関数値 (final)、および各データ点の残差などが含まれます。用途: フィッティングの品質を詳細に分析し、元のデータとフィッティング曲線がどれだけ一致しているかを数値的に確認するために使用します。
{infile_body}-convergence.xlsx:内容: 最適化プロセスの収束過程を記録したExcelファイル。イテレーション数 (
iter) ごとの目的関数値（残差二乗和MSE）が含まれます。用途: 最適化がスムーズに収束したか、あるいは振動したり停滞したりしたかを追跡し、最適化アルゴリズムや初期パラメータの選択が適切であったかを評価するために使用します。
{infile_body}.prm:内容: フィッティングに使用されたパラメータと最終結果を記述したテキストファイル。初期パラメータ値、最終的に決定されたパラメータ値、最終的なMAE (Mean Absolute Error、ここでは最終的な残差二乗和
res.funのことを指す)、および最適化のイテレーション数などが保存されます。用途: フィッティング条件と結果の主要な数値を簡潔に記録し、他の解析ツールでの再利用や記録管理に役立てます。
グラフウィンドウ:
内容: プログラム実行中にmatplotlibによって生成されるグラフウィンドウ。初期データと初期パラメータによるモデル曲線、最適化過程の収束状況、そして最終的なデータとフィッティング曲線が視覚的に表示されます。
用途: フィッティングがデータに適切にフィットしているかを直感的に判断し、必要に応じてパラメータや関数モデルの調整を検討するために使用します。プログラムの実行終了時に一時停止し、ユーザーがウィンドウを閉じるまで表示されます。
コマンドラインでの使用例 (Usage)
lsq_func_tkFit.py は、コマンドライン引数を通じて実行されます。基本的なコマンドラインの書式は以下の通りです。
python lsq_func_tkFit.py [mode] [infile] [method] [func] [p0s] [fit_range] [maxiter] [tol] [h] [olabel] [xlabel] [ylabel] [zlabel]
引数説明
mode: プログラムの動作モード。fit(デフォルト): データをフィッティングします。sim: シミュレーションモード。初期関数のみをプロットし、フィッティングは行いません。
infile: 入力データファイル名（例:peak.xlsx）。method: 最適化アルゴリズム。bfgs(デフォルト): BFGS法。nelder-mead: Nelder-Mead法。powell: Powell法。cg: 共役勾配法 (Polak-Ribiere)。
func: フィッティング関数を表す文字列。xは入力データのリスト（例:x[0]は1番目の独立変数）、pはパラメータのリスト（例:p[0]は1番目のパラメータ）として使用します。例:
"p[0]*exp(-((x[0]-p[1])/p[2])**2)"(ガウス関数)
p0s: 初期パラメータのカンマ区切り文字列。例:
"1.3,0.6,0.1"
fit_range: フィッティング範囲のカンマ区切り文字列。各独立変数の範囲をxmin:xmaxの形式で指定します。複数変数の場合はカンマで区切ります。無制限の場合は-1e100:1e100を使用します。例:
"-1e100:1e100"(単一の独立変数で全範囲)例:
"-10:10,-5:5"(2つの独立変数でそれぞれ範囲を指定)
maxiter: 最大イテレーション数。デフォルトは100。tol: 収束判定の許容誤差。目的関数の変化がこの値以下になると収束と見なされます。デフォルトは1.0e-5。h: 数値微分の微小変位。tkFit_mクラス内部で勾配計算に使用される可能性があります。デフォルトは0.01。olabel: 従属変数（目的変数、y軸データ）のデータ列ラベル。xlabel: 最初の独立変数（x軸データ）のデータ列ラベル。ylabel: 2番目の独立変数（オプション）のデータ列ラベル。使用しない場合は---を指定。zlabel: 3番目の独立変数（オプション）のデータ列ラベル。使用しない場合は---を指定。
コマンドラインでの具体的な使用例
ここでは、前述の「必要な入力ファイル」セクションで示した peak.xlsx ファイルが存在することを前提としたガウス関数フィッティングの例を示します。
peak.xlsx の内容（例）
X_col,Y_col
-2.0,0.1
-1.5,0.3
-1.0,0.8
-0.5,1.2
0.0,1.3
0.5,1.2
1.0,0.8
1.5,0.3
2.0,0.1
実行コマンド
python lsq_func_tkFit.py fit peak.xlsx bfgs "p[0]*exp(-((x[0]-p[1])/p[2])**2)" "1.3,0.0,0.5" "-1e100:1e100" 1000 1e-6 0.01 Y_col X_col --- ---
コマンド引数の説明
fit: フィッティングモードで実行します。peak.xlsx: 入力ファイル名です。bfgs: 最適化メソッドとしてBFGS法を選択します。"p[0]*exp(-((x[0]-p[1])/p[2])**2)": フィッティング関数としてガウス関数A * exp(-((x - mu) / sigma)^2)を指定しています。ここで、\(A=p[0]\)、\(mu=p[1]\)、\(sigma=p[2]\)、\(x=x[0]\) です。"1.3,0.0,0.5": パラメータ \(p[0], p[1], p[2]\) の初期値として、\(A=1.3\)、\(mu=0.0\)、\(sigma=0.5\) を与えます。"-1e100:1e100": 独立変数x[0]のフィッティング範囲を全体（-無限大から+無限大）に設定します。1000: 最大イテレーション数を1000回に設定します。1e-6: 収束判定の許容誤差を \(1.0 \times 10^{-6}\) に設定します。0.01: 数値微分の微小変位を \(0.01\) に設定します。Y_col: 従属変数のデータ列ラベルとしてY_colを指定します。X_col: 最初の独立変数のデータ列ラベルとしてX_colを指定します。---: 2番目の独立変数 (ylabel) は使用しないため---を指定します。---: 3番目の独立変数 (zlabel) は使用しないため---を指定します。
実行結果
プログラムを実行すると、標準出力にフィッティングの進捗状況（イテレーションごとの目的関数値など）が表示されます。
Initialize parameters
Update parameters by command-line arguments
Fitting program to given function
mode : fit
method : bfgs
input file : peak.xlsx
fit_range : [[-1e+100, 1e+100]]
output file : peak-fitting.xlsx
convergence file: peak-convergence.xlsx
parameter file : peak.prm
func : p[0]*exp(-((x[0]-p[1])/p[2])**2)
Target : Y_col
x : X_col
y : ---
z : ---
tol : 1e-06
maxiter : 1000
Read [peak.xlsx]
For ranges: [[-1.e+100 -1.e+100]
[ 1.e+100 1.e+100]]
x ranges:
x[0]: -2.0 - 2.0
fit in: -1e+100 - 1e+100
Calculate initial function:
Initial functions
x(0) Y_col initial | diff(initial)
-2.000000 0.100000 0.000000 | -0.100000
-1.500000 0.300000 0.000000 | -0.300000
-1.000000 0.800000 0.000000 | -0.800000
-0.500000 1.200000 0.000000 | -1.200000
0.000000 1.300000 1.300000 | 0.000000
0.500000 1.200000 0.000000 | -1.200000
1.000000 0.800000 0.000000 | -0.800000
1.500000 0.300000 0.000000 | -0.300000
2.000000 0.100000 0.000000 | -0.100000
plot
Nonlinear least-squares fitting by [bfgs]
Variables for fitting
p(0) : 1.3
p(1) : 0.0
p(2) : 0.5
tol=1e-06
nmaxiter=1000
Iteration 0: SSE = 5.80800000e+00, Parameters: [1.3, 0.0, 0.5]
Iteration 1: SSE = 1.32839938e-01, Parameters: [1.29999999, 0.0001, 0.5]
...
Iteration 20: SSE = 3.25055276e-03, Parameters: [1.29988581, -0.00001859, 0.99958376]
Iteration 21: SSE = 3.25055276e-03, Parameters: [1.29988581, -0.00001859, 0.99958376]
Converged at iteration: 21
Final parameters
p(0) : 1.3
p(1) : -1.859e-05
p(2) : 0.9996
f=0.00325055
Optimized at S2=0.00325055:
Variables for fitting
p(0) : 1.3
p(1) : -1.859e-05
p(2) : 0.99958376
Scores between y(input) and y(fit)
R2 score : 0.99912
MAE : 0.0151
RMSE : 0.0201
Save results to [peak-fitting.xlsx]
Save convergence process to [peak-convergence.xlsx]
Save parameters to [peak.prm]
Plot optimized
また、同時に以下のファイルが生成されます。
peak-fitting.xlsx:X_col,Y_col(入力データ),initial(初期パラメータによる関数値),final(最終パラメータによる関数値) などが含まれます。peak-convergence.xlsx:iter(イテレーション数) とMSE(目的関数値) が記録されます。peak.prm: フィッティングの初期・最終パラメータ、MAE、イテレーション数が記録されます。フィッティング結果を示すグラフウィンドウが表示されます。このウィンドウは、ユーザーが閉じるまで表示され続けます。