diff_order_noise.py の技術ドキュメント
プログラムの動作
diff_order_noise.py は、ノイズを含むデータに対して異なる次数の数値微分近似（2点、3点、5点、7点）を適用し、その結果を比較評価するためのPythonスクリプトです。
このプログラムの主な機能は以下の通りです。
関数の定義: 被微分関数として \(f(x) = \sin(x)\) を、その解析微分として \(f'(x) = \cos(x)\) を定義します。
ノイズ付きデータ生成: 指定された範囲、ステップサイズで \(f(x)\) の値に一様乱数に基づくノイズを加えたデータセットを生成します。
数値微分: 以下の異なる次数の中央差分近似（2点前方差分を含む）を用いて、ノイズ付きデータポイントの微分値を計算します。
2点前方差分（
diff2forward）3点中心差分（
diff3）5点中心差分（
diff5）7点中心差分（
diff7）
結果の出力: 計算された元の関数値、ノイズ付き関数値、解析微分値、および各数値微分値を標準出力とCSVファイルに書き出します。
グラフ表示: ノイズの有無による関数値の比較、および解析微分値と各数値微分値の比較を2つのサブプロットとしてMatplotlibで視覚化します。
このプログラムは、実際の観測データにノイズが含まれる状況において、どの程度の次数の数値微分法がノイズの影響を受けにくいか、あるいは高次の近似が必ずしも高精度に繋がらないケース（特にノイズが多い場合）を理解するのに役立ちます。
原理
このプログラムは、数値微分における差分近似の原理に基づいています。関数 \(f(x)\) の導関数 \(f'(x)\) を、離散的なデータ点 \(x_i\) とその関数値 \(f(x_i)\) を用いて近似します。データ点の間隔を \(h = x_{i+1} - x_i\) とします。
2点前方差分 (Forward Difference): 最も単純な近似方法の一つで、以下の式で計算されます。 $\( \left(\frac{df}{dx}\right)_i \approx \frac{f(x_{i+1}) - f(x_i)}{h} \)$
3点中心差分 (Central Difference): \(x_i\) の前後2点を用いることで、前方差分よりも高精度な近似が得られます。 $\( \left(\frac{df}{dx}\right)_i \approx \frac{f(x_{i+1}) - f(x_{i-1})}{2h} \)$
5点中心差分: より多くの点を用いることで、さらに高次のテイラー展開に基づく近似が可能になります。 $\( \left(\frac{df}{dx}\right)_i \approx \frac{-f(x_{i+2}) + 8f(x_{i+1}) - 8f(x_{i-1}) + f(x_{i-2})}{12h} \)$
7点中心差分: さらに高次の近似で、より広範囲のデータ点を使用します。 $\( \left(\frac{df}{dx}\right)_i \approx \frac{f(x_{i+3}) - 9f(x_{i+2}) + 45f(x_{i+1}) - 45f(x_{i-1}) + 9f(x_{i-2}) - f(x_{i-3})}{60h} \)$
ノイズのモデル:
データ生成関数 make_data では、真の関数値 \(f(x)\) にノイズを加えます。ノイズは以下の式でモデル化されます。
$\( f_{\text{noisy}}(x) = f(x) + \text{noise} \times (R - 0.5) \times 2.0 \)\(
ここで \)R\( は \)[0, 1)\( の範囲の一様乱数であり、\)(\text{noise} \times (R - 0.5) \times 2.0)\( は \)[- \text{noise}, + \text{noise})$ の範囲で一様分布する乱数となります。
高次の差分近似は、理論的にはテイラー展開のより多くの項を考慮するため、データの精度が高い場合にはより正確な微分値を与えます。しかし、ノイズが含まれるデータの場合、広範囲のデータ点を使用する高次の近似は、より多くのノイズを拾い上げてしまう可能性があり、結果として微分値の変動が大きくなることがあります。このプログラムは、この現象を視覚的に示します。
必要な非標準ライブラリとインストール方法
このプログラムの実行には、以下の非標準ライブラリが必要です。
NumPy: 数値計算（特に乱数生成）に使用されます。
Matplotlib: 結果のグラフ描画に使用されます。
これらのライブラリは pip コマンドでインストールできます。
pip install numpy matplotlib
必要な入力ファイル
diff_order_noise.py は、外部の入力ファイルを必要としません。すべてのデータはプログラム内で生成されます。
生成される出力ファイル
プログラムは以下の出力ファイルを生成します。
diff_order_noise.csv: 計算された各種データがカンマ区切り形式で保存されます。内容: 各行が1つのデータポイントに対応し、以下の列が含まれます。
x: x座標f(x): 真の関数値 \(f(x)\)f+noise: ノイズ付きの関数値 \(f_{\text{noisy}}(x)\)df/dx: 真の導関数 \(f'(x)\) (解析微分値)Df/Dx(2-point): 2点前方差分による数値微分値Df/Dx(3-point): 3点中心差分による数値微分値Df/Dx(5-point): 5点中心差分による数値微分値Df/Dx(7-point): 7点中心差分による数値微分値
グラフウィンドウ: プログラム実行中にMatplotlibによって生成されるグラフが表示されます。このグラフは自動的に保存されませんが、関数値のプロットと微分値のプロットの2つのサブプロットが含まれます。
コマンドラインでの使用例 (Usage)
このプログラムは、コマンドライン引数を必要としません。Pythonインタープリタを使用して直接実行します。
python diff_order_noise.py
コマンドラインでの具体的な使用例
diff_order_noise.py を実行すると、標準出力に計算結果のサマリーが表示され、同時に diff_order_noise.csv ファイルが生成されます。また、計算結果を視覚化したグラフウィンドウが表示されます。グラフはEnterキーを押すまで表示され続けます。
python diff_order_noise.py
実行結果の例（標準出力の一部）:
Numerical differentiation using differnet approximations
Write to [diff_order_noise.csv]
x : f(x) f+noise df/dx Df/Dx(2-point) Df/Dx(3-point) Df/Dx(5-point) Df/Dx(7-point)
0.000: -0.0000 -0.1011 1.0000 0.8690 0.9839 0.9997 1.0000
0.200: 0.1987 -0.0762 0.9801 1.2505 1.0264 0.9818 0.9801
0.400: 0.3894 0.3957 0.9211 0.7107 0.9287 0.9201 0.9211
0.600: 0.5646 0.6974 0.8253 0.8306 0.8291 0.8249 0.8253
0.800: 0.7174 0.8000 0.6967 0.6635 0.6970 0.6967 0.6967
1.000: 0.8415 0.7303 0.5403 0.3995 0.5402 0.5403 0.5403
1.200: 0.9320 0.8441 0.3624 0.2878 0.3624 0.3624 0.3624
1.400: 0.9854 0.8358 0.1700 0.0076 0.1697 0.1700 0.1700
1.600: 0.9996 1.1276 -0.0292 -0.0526 -0.0289 -0.0292 -0.0292
1.800: 0.9738 1.0284 -0.2272 -0.2285 -0.2271 -0.2272 -0.2272
2.000: 0.9093 0.7818 -0.4161 -0.4162 -0.4161 -0.4161 -0.4161
2.200: 0.8085 0.8993 -0.5885 -0.6402 -0.5882 -0.5885 -0.5885
2.400: 0.6755 0.8268 -0.7374 -0.7812 -0.7371 -0.7374 -0.7374
2.600: 0.5155 0.4143 -0.8583 -0.8712 -0.8581 -0.8583 -0.8583
2.800: 0.3349 0.4285 -0.9455 -0.9126 -0.9455 -0.9455 -0.9455
3.000: 0.1411 0.0072 -0.9900 -0.9660 -0.9901 -0.9900 -0.9900
3.200: -0.0583 -0.1243 -0.9983 -0.9705 -0.9985 -0.9983 -0.9983
3.400: -0.2555 -0.4214 -0.9680 -0.9221 -0.9683 -0.9680 -0.9680
3.600: -0.4425 -0.4132 -0.9026 -0.8530 -0.9029 -0.9026 -0.9026
3.800: -0.6119 -0.6120 -0.8066 -0.7584 -0.8068 -0.8066 -0.8066
4.000: -0.7568 -0.7297 -0.6840 -0.6436 -0.6841 -0.6840 -0.6840
4.200: -0.8716 -0.9959 -0.5403 -0.5055 -0.5402 -0.5403 -0.5403
4.400: -0.9516 -1.0249 -0.3809 -0.3400 -0.3807 -0.3809 -0.3809
4.600: -0.9936 -0.9231 -0.2115 -0.1706 -0.2113 -0.2115 -0.2115
4.800: -0.9962 -1.0335 -0.0385 -0.0150 -0.0384 -0.0385 -0.0385
5.000: -0.9589 -1.0560 0.1341 0.1654 0.1342 0.1341 0.1341
5.200: -0.8835 -0.9023 0.3017 0.3340 0.3018 0.3017 0.3017
5.400: -0.7728 -0.7773 0.4578 0.4996 0.4579 0.4578 0.4578
5.600: -0.6300 -0.4851 0.5966 0.6148 0.5967 0.5966 0.5966
5.800: -0.4599 -0.3486 0.7126 0.7225 0.7127 0.7126 0.7126
6.000: -0.2690 -0.1804 0.8011 0.8055 0.8012 0.8011 0.8011
Press ENTER to exit>>
生成される diff_order_noise.csv の内容（一部）:
x,f(x),f+noise,df/dx,Df/Dx(2-point),Df/Dx(3-point),Df/Dx(5-point),Df/Dx(7-point)
0.0,0.0,-0.1010776510167683,1.0,0.8690412555621815,0.9838706346764516,0.9996614451006566,1.0000000000000002
0.2,0.19866933079506122,-0.07621454479577935,0.9800665778412416,1.2505298459527063,1.0264024467776495,0.981759604179782,0.9800665778412416
0.4,0.3894183423086505,0.39570189569733475,0.9210609940028849,0.7107297956485885,0.9287399515599026,0.920067339794025,0.9210609940028849
0.6,0.5646424733950355,0.6973549298282121,0.8253356149096783,0.8306346911677322,0.8291414442220455,0.8249007204965155,0.8253356149096783
0.8,0.7173560908995228,0.8000028458784122,0.6967067093471654,0.663539860086256,0.6970420790250002,0.6967067093471654,0.6967067093471654
1.0,0.8414709848078965,0.7303350383187212,0.5403023058681398,0.3994793836262973,0.5402422799308003,0.5403023058681398,0.5403023058681398
...
グラフは、ノイズ付きの元の関数が真の関数に対してどれだけ変動しているか、また各数値微分法が真の微分値（解析微分）に対してどれだけ近いか（または離れているか）を視覚的に示します。一般に、ノイズが多いデータに対しては、高次の差分近似が必ずしも優れた結果をもたらさないことが観察されるでしょう。