Keyword: 最小二乗, ３次, スプライン, フィット
概要
本サンプルは最小二乗３次スプライン曲線フィットを行うC言語によるサンプルプログラムです。 本サンプルは以下に示されるデータについて最小二乗３次スプライン曲線フィットを行い、スプラインの値を求めて出力します。
※本サンプルはnAG Cライブラリに含まれる関数 nag_1d_spline_fit_knots() のExampleコードです。本サンプル及び関数の詳細情報は nag_1d_spline_fit_knots のマニュアルページをご参照ください。
ご相談やお問い合わせはこちらまで
入力データ
（本関数の詳細はnag_1d_spline_fit_knots のマニュアルページを参照）| このデータをダウンロード |
nag_1d_spline_fit_knots (e02bac) Example Program Data
14
5 2
1.50
2.60
4.00
8.00
0.20 0.00 0.20
0.47 2.00 0.20
0.74 4.00 0.30
1.09 6.00 0.70
1.60 8.00 0.90
1.90 8.62 1.00
2.60 9.10 1.00
3.10 8.90 1.00
4.00 8.15 0.80
5.15 7.00 0.50
6.17 6.00 0.70
8.00 4.54 1.00
10.00 3.39 1.00
12.00 2.56 1.00
- １行目はタイトル行で読み飛ばされます。
- ２行目にデータ点の数(m)を指定しています。
- ３行目にスプラインの区間の数(ncap)、重みづけをするかどうかのフラグ(wght)を指定しています。"1"は重みづけをしないことを意味し、"1"以外は重みづけをすることを意味します。
- ４〜７行目にノット(lamda)を指定しています。
- ８〜２１行目に独立変数xの値、従属変数yの値、重み(weights)を指定しています。
出力結果
（本関数の詳細はnag_1d_spline_fit_knots のマニュアルページを参照）| この出力例をダウンロード |
nag_1d_spline_fit_knots (e02bac) Example Program Results
Number of distinct knots = 6
Distinct knots located at
0.2000 1.5000 2.6000 4.0000 8.0000 12.0000
J B-spline coeff c
1 -0.0465
2 3.6150
3 8.5724
4 9.4261
5 7.2716
6 4.1207
7 3.0822
8 2.5597
Residual sum of squares = 1.78e-03
Cubic spline approximation and residuals
r Abscissa Weight Ordinate Spline Residual
1 0.2000 0.2000 0.0000 -0.0465 -4.7e-02
0.3350 1.0622
2 0.4700 0.2000 2.0000 2.1057 1.1e-01
0.6050 3.0817
3 0.7400 0.3000 4.0000 3.9880 -1.2e-02
0.9150 5.0558
4 1.0900 0.7000 6.0000 5.9983 -1.7e-03
1.3450 7.1376
5 1.6000 0.9000 8.0000 7.9872 -1.3e-02
1.7500 8.3544
6 1.9000 1.0000 8.6200 8.6348 1.5e-02
2.2500 9.0076
7 2.6000 1.0000 9.1000 9.0896 -1.0e-02
2.8500 9.0353
8 3.1000 1.0000 8.9000 8.9125 1.2e-02
3.5500 8.5660
9 4.0000 0.8000 8.1500 8.1321 -1.8e-02
4.5750 7.5592
10 5.1500 0.5000 7.0000 6.9925 -7.5e-03
5.6600 6.5010
11 6.1700 0.7000 6.0000 6.0255 2.6e-02
7.0850 5.2292
12 8.0000 1.0000 4.5400 4.5315 -8.5e-03
9.0000 3.9045
13 10.0000 1.0000 3.3900 3.3928 2.8e-03
11.0000 2.9574
14 12.0000 1.0000 2.5600 2.5597 -3.5e-04
- ３行目にノットの数が出力されています。
- ５〜７行目にノットの位置が出力されています。
- １０〜１９行目にBスプライン係数が出力されています。
- ２１行目に残差平方和が出力されています。
- ２３〜５２行目にインデックス、独立変数xの値、重み、従属変数yの値、３次スプラインフィットの値、、残差が出力されています。
ソースコード
（本関数の詳細はnag_1d_spline_fit_knots のマニュアルページを参照）
※本サンプルソースコードはnAG数値計算ライブラリ（Windows, Linux, MAC等に対応）の関数を呼び出します。
サンプルのコンパイル及び実行方法
| このソースコードをダウンロード |
/* nag_1d_spline_fit_knots (e02bac) Example Program.
*
* CLL6I261D/CLL6I261DL Version.
*
* Copyright 2017 Numerical Algorithms Group.
*
* Mark 26.1, 2017.
*
*/
#include <nag.h>
#include <stdio.h>
#include <nag_stdlib.h>
#include <nage02.h>
int main(void)
{
Integer exit_status = 0, j, m, ncap, ncap7, r, wght;
NagError fail;
Nag_Spline spline;
double fit, ss, *weights = 0, *x = 0, xarg, *y = 0;
INIT_FAIL(fail);
/* Initialize spline */
spline.lamda = 0;
spline.c = 0;
printf("nag_1d_spline_fit_knots (e02bac) Example Program Results\n");
scanf("%*[^\n]"); /* Skip heading in data file */
while (scanf("%ld", &m) != EOF)
{
if (m >= 4) {
if (!(weights = nAG_ALLOC(m, double)) ||
!(x = nAG_ALLOC(m, double)) || !(y = nAG_ALLOC(m, double))
)
{
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
}
else {
printf("Invalid m.\n");
exit_status = 1;
goto END;
}
scanf("%ld%ld", &ncap, &wght);
if (ncap > 0) {
ncap7 = ncap + 7;
spline.n = ncap7;
if (!(spline.lamda = nAG_ALLOC(ncap7, double)))
{
printf("Allocation failure\n");
exit_status = -1;
goto END;
}
}
else {
printf("Invalid ncap.\n");
exit_status = 1;
goto END;
}
for (j = 4; j < ncap + 3; ++j)
scanf("%lf", &(spline.lamda[j]));
for (r = 0; r < m; ++r) {
if (wght == 1) {
scanf("%lf%lf", &x[r], &y[r]);
weights[r] = 1.0;
}
else
scanf("%lf%lf%lf", &x[r], &y[r], &weights[r]);
}
/* nag_1d_spline_fit_knots (e02bac).
* Least squares curve cubic spline fit (including
* interpolation), one variable
*/
nag_1d_spline_fit_knots(m, x, y, weights, &ss, &spline, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_1d_spline_fit_knots (e02bac).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
printf("\nNumber of distinct knots = %ld\n\n", ncap + 1);
printf("Distinct knots located at \n\n");
for (j = 3; j < ncap + 4; j++)
printf("%8.4f%s", spline.lamda[j],
(j - 3) % 6 == 5 || j == ncap + 3 ? "\n" : " ");
printf("\n\n J B-spline coeff c\n\n");
for (j = 0; j < ncap + 3; ++j)
printf(" %ld %13.4f\n", j + 1, spline.c[j]);
printf("\nResidual sum of squares = ");
printf("%11.2e\n\n", ss);
printf("Cubic spline approximation and residuals\n");
printf(" r Abscissa Weight Ordinate"
" Spline Residual\n\n");
for (r = 0; r < m; ++r) {
/* nag_1d_spline_evaluate (e02bbc).
* Evaluation of fitted cubic spline, function only
*/
nag_1d_spline_evaluate(x[r], &fit, &spline, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_1d_spline_evaluate (e02bbc).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
printf("%3ld %11.4f %11.4f %11.4f %11.4f"
" %10.1e\n", r + 1, x[r], weights[r], y[r], fit, fit - y[r]);
if (r < m - 1) {
xarg = (x[r] + x[r + 1]) * 0.5;
/* nag_1d_spline_evaluate (e02bbc), see above. */
nag_1d_spline_evaluate(xarg, &fit, &spline, &fail);
if (fail.code != NE_NOERROR) {
printf("Error from nag_1d_spline_evaluate (e02bbc).\n%s\n",
fail.message);
exit_status = 1;
goto END;
}
printf(" %14.4f %33.4f\n", xarg, fit);
}
}
/* Free memory used by spline */
nAG_FREE(spline.lamda);
nAG_FREE(spline.c);
END:
nAG_FREE(weights);
nAG_FREE(x);
nAG_FREE(y);
}
return exit_status;
}