条件
2次元上に \(n\) 個の点$$
(X_1, Y_1),(X_2,Y_2) , \ldots ,(X_{n-1},Y_{n-1}) , (X_n,Y_n)
$$(ただし、\(X_1<X_2<\cdots<X_n\)) が与えられたとき、それらをなるべく 「滑らかに」 つなげるためにはどうしたらよいでしょうか。3次スプライン曲線(関数)は、このような問題の解答のひとつです。
3次スプライン曲線とは、
- 点と点に挟まれる区間は3次関数
- 区間\([-\infty,X_1]\)と\([X_n,\infty]\) は一次関数
- 各点で曲線の値・傾き・曲率が連続
という性質をもつ曲線です。具体的な導出方法を簡単に説明します。
まず、\(m\) 個目の点と、\(m+1\) 個目の点に挟まれる区間\( (X_m ,X_{m+1})\) の関数の形状を$$ F_m(x)=a_mx^3+b_mx^2+c_mx+d_m$$ と定義します。ただし、1つめの点の左側は \( F_0(x)=c_0x+d_0\)、\(n\) 個めの点より右側は \( F_n(x)=c_nx+d_n\) とします。スプライン曲線の要件の最後 (各点で曲線の値・傾き・曲率が連続) を満たすためには、\(n\) 個の点の全てで 0, 1, 2 次微分が等しくないといけません。
まず、\(2 <= m <= n-1\) の点に関しては、以下の4つの条件
- \( F_{m-1}(X_m)= a_{m-1}X_m^3+b_{m-1}X_m^2+c_{m-1}X_m+d_{m-1} =Y_m\)
- \( F_m(X_m)= a_mX_m^3+b_mX_m^2 + c_m X_m + d_m =Y_m\)
- \( F’_{m-1}(X_m) = 3a_{m-1}X_m^2+2b_{m-1}X_m+c_{m-1}= F’_m(X_m) = 3a_mX_m^2+2b_m X_m + c_m\)
- \( F”_{m-1}(X_m) = 6a_{m-1}X_m+2b_{m-1}= F”_m(X_m) = 6a_mX_m+2b_m\)
を満たす必要があります。このような条件は合計 \(4(n-2)\) 個あります。
また、\(m=1\) のときは、
- \(F_0(X_1)= c_0X_0+d_0 =Y_1\)
- \(F_1(X_1)= a_1X_1^3+b_1X_1^2+c_1X_1+d_1 =Y_1\)
- \(F’_{0} (X_1)= c_0 = F’_1(X_1) = 3 a_1 X_1^2 + 2 b_1 X_1 + c_1 \)
- \(F”_{0}(X_1)= 0 = F”_1(X_1)= 6 a_1X_1 + 2 b_1 \)
という4つの条件になります。\(m=n\) のときは、
- \(F_{n-1}(X_n)= a_{n-1}X_n^3+b_{n-1} X_n^2+c_{n-1}X_n+d_{n-1} =Y_n\)
- \(F_n(X_n)= c_n X_n+d_n =Y_n\)
- \(F’_{n-1}(X_n)= 3a_{n-1}X_n^2+2b_{n-1}X^n+c_{n-1} = F’_n(X_n)= c_n\)
- \(F”_{n-1}(X_n)= 6a_{n-1}X_n+2b_{n-1}= F”_n(X_n)=0\)
という4つの条件になります。ここまでをまとめると、結局以下に示す\(4n\) 個の未知パラメータ $$
(c_0, d_0), (a_1,b_1,c_1,d_1), (a_2,b_2,c_2,d_2),\cdots, (a_{n-1}, b_{n-1},c_{n-1},d_{n-1}), (c_n,d_n)
$$を、\(4n\) 個の条件式から求めるという問題に帰結できます。
行列表現
未知パラメータの数と線形の条件式の数が同一ですから、条件式を行列で表現したあと逆行列を解くという方針が最もシンプルでしょう。
まず、\(2 <= m <= n-1\) のときの条件は、
$$ \alpha_m= \begin{pmatrix} X_m^3 & X_m^2 & X_m & 1 & 0 & 0 & 0 & 0 \\ 0 &0 & 0 & 0 & X_m^3 & X_m^2 & X_m &1 \\ 3X_m^2 &2X_m & 1 & 0 &- 3X_m^2 & -2X_m & -1 &0 \\ 6X_m & 2 & 0 & 0 &-6X_m &-2 &0 & 0 \end{pmatrix}$$として、$$
\alpha_m
\begin{pmatrix}
a_{m-1}\\b_{m-1}\\c_{m-1}\\d_{m-1}\\a_m\\b_m\\c_m\\d_m
\end{pmatrix}=
\begin{pmatrix}
Y_m\\Y_m\\0\\0
\end{pmatrix}
$$と書き直すことが出来ます。同様に、\(m=1\) の時の条件は
$$ \alpha_1= \begin{pmatrix}
X_1 & 1 & 0 & 0 & 0 & 0 \\
0 & 0 & X_1^3 & X_1^2 & X_1 &1 \\
1 & 0 & -3X_1^2 & -2X_1 & -1 &0 \\
0 & 0 &-6X_1 &-2 &0 & 0
\end{pmatrix},
\,\,\,\, \alpha_1
\begin{pmatrix}
c_0\\d_0\\a_1\\b_1\\c_1\\d_1
\end{pmatrix}=
\begin{pmatrix}
Y_1\\Y_1\\0\\0
\end{pmatrix}
$$ \(m=n\) の時は$$
\alpha_n= \begin{pmatrix}
X_n^3 & X_n^2 & X_n &1 & 0 & 0\\
0 &0 & 0 & 0 & X_n &1 \\
3X_n^2 & 2X_n & 1 &0 & -1 &0 \\
6X_n & 2 & 0 & 0 & 0 & 0
\end{pmatrix},
\,\,\,\, \alpha_n
\begin{pmatrix}
a_{n-1} \\ b_{n-1} \\ c_{n-1} \\ d_{n-1} \\ c_n \\ d_n
\end{pmatrix}=
\begin{pmatrix}
Y_n\\Y_n\\0\\0
\end{pmatrix}
$$のようにまとめることが出来ます。
さらに、\(\alpha\)を左右に分割した部分行列を以下のように定義します。$$
\alpha^\textrm{L}_m= \begin{pmatrix} X_m^3 & X_m^2 & X_m & 1 \\
0 &0 & 0 & 0 \\
3X_m^2 &2X_m & 1 & 0 \\
6X_m & 2 & 0 & 0 \end{pmatrix},\,\,
\alpha^\textrm{R}_m= \begin{pmatrix}
0 & 0 & 0 & 0 \\
X_m^3 & X_m^2 & X_m &1 \\
-3X_m^2 & -2X_m & -1 &0 \\
-6X_m &-2 &0 & 0 \end{pmatrix}\\ \, \\
\alpha^\textrm{L}_1= \begin{pmatrix} X_1 & 1 \\0 & 0 \\ 1 & 0 \\0 & 0 \end{pmatrix},
\,\,\,\, \alpha^\textrm{R}_1= \begin{pmatrix}
0 & 0 & 0 & 0 \\X_1^3 & X_1^2 & X_1 &1 \\ -3X_1^2 & -2X_1 & -1 &0 \\-6X_1 &-2 &0 & 0
\end{pmatrix}\\ \, \\
\alpha^\textrm{L}_n= \begin{pmatrix}
X_n^3 & X_n^2 & X_n &1\\
0 &0 & 0 & 0 \\
3X_n^2 & 2X_n & 1 &0 \\
6X_n & 2 & 0 & 0
\end{pmatrix}, \,\,\,\,
\alpha^\textrm{R}_n= \begin{pmatrix} 0 & 0\\ X_n &1 \\ -1 &0 \\ 0 & 0 \end{pmatrix}
$$ これで準備は整いました。\(\alpha^\textrm{L}\) や \(\alpha^\textrm{R}\) を適切に並べて一つの大きな行列を作ると、以下のように全ての条件をまとめることができます。$$
\begin{pmatrix}
\alpha^\textrm{L}_1 & \alpha^\textrm{R}_1 & 0 & & & & \\
0 & \alpha^\textrm{L}_2 & \alpha^\textrm{R}_2 & 0 & & \huge{0} & \\
& 0 & \alpha^\textrm{L}_3 & \alpha^\textrm{R}_3 & 0 & & \\
& & \ddots & \ddots & \ddots & \ddots & \\
& \huge{0} & & 0 & \alpha^\textrm{L}_{n-1} & \alpha^\textrm{R}_{n-1} & 0 \\
& & & & 0 & \alpha^\textrm{L}_n & \alpha^\textrm{R}_n \\
\end{pmatrix}
\begin{pmatrix}c_0 \\ d_0 \\ a_1 \\ b_1 \\ c_1 \\ d_1 \\ a_2 \\ b_2 \\ c_2 \\ d_2 \\ \vdots \\ \ a_{n-1} \\ b_{n-1} \\ c_{n-1} \\ d_{n-1} \\ c_n \\ d_n \end{pmatrix}
=
\begin{pmatrix}Y_1 \\ Y_1 \\ 0 \\ 0 \\ Y_2 \\ Y_2 \\ 0 \\ 0 \\ \vdots \\ Y_{n-1} \\ Y_{n-1} \\ 0 \\ 0 \\ Y_n \\ Y_n \\ 0 \\ 0 \end{pmatrix}$$
巨大ですがシンプルな行列演算で表現できました。左辺の \(\alpha\) が作る行列は必ず正則行列ですから、これの逆行列を解いて、右辺に左からかけたら、求めたいパラメータ \(a, b, c, d\) が求まります。あとは簡単なので省略します。