3次エルミートスプライン

数値解析において、3次エルミートスプラインまたは3次エルミート補間は、各部分がエルミート形式で指定される3次多項式、つまり対応する区間の両端の値と1次微分によって指定されるスプラインである 。

3次エルミートスプラインは、通常、与えられる引数値 $x_1,x_2,\dots ,x_n$ で指定される数値データを補間して連続関数を得るために使用される。 データは、各 $x_k$ についての所望の関数値と微分から構成される必要がある（値のみが与えられている場合、微分はそれから推定する必要がある）。 エルミート補間多項式は各区間 $(x_k,x_{k+1})$ に個別に適用される。 結果のスプラインは連続、かつ、一次微分連続となる。

3次多項式スプラインは他の方法でも指定でき、ベジェ形式が最も一般的である。 しかし、これら2つの方法で得られるスプラインセットは同一であり、ベジェ形式とエルミート形式の間でデータの変換も容易であるため、これらの名称はしばしば同義語のように使用される。

3次多項式スプラインは、コンピュータグラフィックスやテンプレート:日本語版にない記事リンクにおいて、平面または3次元空間上の指定された点を通過する曲線や運動軌跡を得るために広く用いられている。 これらのアプリケーションでは、平面または空間の各座標は、個別のパラメータ t を持つ3次スプライン関数によって個別に補間される。 3次多項式スプラインは、テンプレート:日本語版にない記事リンクなどの構造解析アプリケーションでも広く用いられている。

3次多項式スプラインは、死亡率分析 や死亡率予測 にも応用されている。

3次スプラインは、いくつかの方法で2つ以上のパラメータを持つ関数に拡張できる。 バイキュービックスプライン（バイキュービック補間）は、デジタル画像のピクセル値や地形の高度データなど、規則的な矩形グリッド上のデータを補間するためによく使用される。 バイキュービック曲面パッチ は、コンピュータグラフィックスにおいて不可欠なツールである。

3次スプラインは、特にコンピュータグラフィックスではしばしばcsplineと呼ばれる。エルミートスプラインは、シャルル・エルミートにちなんで名付けられた。

単一区間の補間

単位区間 [0, 1]

4つのエルミート基底関数。各部分区間における補間関数は、これら4つの関数の線形結合である。

単位間隔 $[0,1]$ について、 $t=0$ における始点 ${𝒑}_0$、$t=1$ における終点 ${𝒑}_1$、および、 $t=0$ における開始接ベクトル ${𝒎}_0$、$t=1$ における終了接ベクトル ${𝒎}_1$ が与えられると、 多項式は次のように定義される。

$𝒑(t)=(2t^3-3t^2+1){𝒑}_0+(t^3-2t^2+t){𝒎}_0+(-2t^3+3t^2){𝒑}_1+(t^3-t^2){𝒎}_1$

ただし t ∈ [0, 1] である。

任意区間の補間

任意の区間 $(x_k,x_{k+1})$ 内の $x$ での補間は、アフィン（1次）変数変換により区間を $[0,1]$ にマッピングすることで行われる。 式は次のようになる。

$𝒑(x)=h_{00}(t){𝒑}_k+h_{10}(t)(x_{k+1}-x_k){𝒎}_k+h_{01}(t){𝒑}_{k+1}+h_{11}(t)(x_{k+1}-x_k){𝒎}_{k+1}$

ただし $t=(x-x_k)/(x_{k+1}-x_k)$ と $h$ は、以下に定義される基底関数を参照。 なお接ベクトルの長さは、単位区間の方程式に対し $x_{k+1}-x_k$ 倍にスケーリングされる。

一意性

上記の式は、指定された接ベクトルを持つ2点間の、一意の3次多項式の経路を示す。

証明

$P,Q$ は、与えられた境界条件を満たす2つの3次多項式とする。$R=Q-P$ とすると、

$R(0)=Q(0)-P(0)=0$

$R(1)=Q(1)-P(1)=0$

である。 $Q$ と $P$ はいずれも3次多項式であるので、$R$ は最大で3次の多項式である。 そのため、$R$ は以下の式でなければならない。

$R(x)=ax(x-1)(x-r)$

導関数を計算すると、

$R^{\prime }(x)=ax(x-1)+ax(x-r)+a(x-1)(x-r)$

さらに、以下が明らかである。

$R^{\prime }(0)=Q^{\prime }(0)-P^{\prime }(0)=0,$

$R^{\prime }(1)=Q^{\prime }(1)-P^{\prime }(1)=0,$

式 (テンプレート:EquationNote) と (テンプレート:EquationNote) を合わせると、$a=0$ が導かれ、したがって $R=0$、よって $P=Q$ となる。

表現

単位区間上の補間多項式は次のように書ける（任意区間については、上記再スケール版を参照）。 $$𝒑(t)=h_{00}(t){𝒑}_0+h_{10}(t){𝒎}_0+h_{01}(t){𝒑}_1+h_{11}(t){𝒎}_1$$ ただし $h_{00}$, $h_{10}$, $h_{01}$, $h_{11}$ はエルミート基底関数である。 これらは様々な方法で記述することができ、それぞれ異なる性質を示す。

	展開	因数分解	バーンスタイン
$h_{00}(t)$	$2t^3-3t^2+1$	$(1+2t)(1-t{)}^2$	$B_0(t)+B_1(t)$
$h_{10}(t)$	$t^3-2t^2+t$	$t(1-t{)}^2$	$\frac{1}{3}{B}_1(t)$
$h_{01}(t)$	$-2t^3+3t^2$	$t^2(3-2t)$	$B_3(t)+B_2(t)$
$h_{11}(t)$	$t^3-t^2$	$t^2(t-1)$	$-\frac{1}{3}{B}_2(t)$

「展開」列は上記の定義で使用された表現を示している。 「因数分解」列は、$h_{10}$ と $h_{11}$ が端点でゼロであることを直接示している。 さらに、$h_{01}$ と $h_{11}$ は0で重複度2の零点を持ち、$h_{00}$ と $h_{10}$ は1で重複度2の零点を持つため、それら端点で傾きが0になると言える。 「バーンスタイン」列は、エルミート基底関数を3次のバーンスタイン多項式に分解したものを示しており、

$B_k(t)=(\genfrac{}{}{0ex}{}{3}{k})\cdot t^k\cdot (1-t{)}^{3-k}$

である。

この関係を利用すると、3次エルミート補間を 3次ベジェ曲線の形式で、次の4つの値 ${𝒑}_0,{𝒑}_0+\frac{1}{3}{𝒎}_0,{𝒑}_1-\frac{1}{3}{𝒎}_1,{𝒑}_1$ によって表現でき、ド・カステリョのアルゴリズムを用いてエルミート補間を行うことができる。 これは、3次ベジェ曲線 において、補間曲線の両端点における接ベクトルが、中間の2つの制御点により決定されることを示している。

この多項式は次の標準形式で記述することもできる： $$𝒑(t)=(2{𝒑}_0+{𝒎}_0-2{𝒑}_1+{𝒎}_1)t^3+(-3{𝒑}_0+3{𝒑}_1-2{𝒎}_0-{𝒎}_1)t^2+{𝒎}_{0}t+{𝒑}_0$$ ただし、制御点と接ベクトルが係数である。 この形式では係数が一定となり、一度計算して再利用できるため、t の様々な値について多項式を効率的に評価できる。

行列形式では $$𝒑(t)=\left[\begin{array}{cccc}{t}^3& t^2& t& 1\end{array}\right]\left[\begin{array}{cccc}2& 1& -2& 1\\ -3& -2& 3& -1\\ 0& 1& 0& 0\\ 1& 0& 0& 0\end{array}\right]\left[\begin{array}{c}{𝒑}_0\\ {𝒎}_0\\ {𝒑}_1\\ {𝒎}_1\end{array}\right]$$ となる。

ベジェ曲線への変換

3次ベジェ曲線の行列形式は $$𝒑(t)=\left[\begin{array}{cccc}{t}^3& t^2& t& 1\end{array}\right]\left[\begin{array}{cccc}-1& 3& -3& 1\\ 3& -6& 3& 0\\ -3& 3& 0& 0\\ 1& 0& 0& 0\end{array}\right]\left[\begin{array}{c}{{𝒑}_{𝒃𝒛}}_0\\ {{𝒑}_{𝒃𝒛}}_1\\ {{𝒑}_{𝒃𝒛}}_2\\ {{𝒑}_{𝒃𝒛}}_3\end{array}\right]$$ であることから、3次エルミートスプラインと等しくなる3次ベジェ曲線の制御点は $$\left[\begin{array}{c}{{𝒑}_{𝒃𝒛}}_0\\ {{𝒑}_{𝒃𝒛}}_1\\ {{𝒑}_{𝒃𝒛}}_2\\ {{𝒑}_{𝒃𝒛}}_3\end{array}\right]={\left[\begin{array}{cccc}-1& 3& -3& 1\\ 3& -6& 3& 0\\ -3& 3& 0& 0\\ 1& 0& 0& 0\end{array}\right]}^{-1}\left[\begin{array}{cccc}2& 1& -2& 1\\ -3& -2& 3& -1\\ 0& 1& 0& 0\\ 1& 0& 0& 0\end{array}\right]\left[\begin{array}{c}{𝒑}_0\\ {𝒎}_0\\ {𝒑}_1\\ {𝒎}_1\end{array}\right]=\left[\begin{array}{c}{𝒑}_0\\ {𝒑}_0+\frac{{𝒎}_0}{3}\\ {𝒑}_1-\frac{{𝒎}_1}{3}\\ {𝒑}_1\end{array}\right]$$ である。

データセットの補間

あるデータセット $(x_k,{𝒑}_k)$（$k=1,\dots ,n$）は、各区間に上記手順を適用することで補間できる。ただし、接ベクトルを適切な方法で選択する、つまり端点を共有する区間どうしで接線を等しくする。 補間された曲線は区分3次エルミートスプラインで構成され、$(x_1,x_n)$ において大域的に連続して微分可能である。

接ベクトルの選択は一意ではなく、いくつかのオプションが利用可能である。

有限差分

有限差分の例

最も単純な選択肢は3点差分であり、一定間隔のデータ点を必要としない。 内部の点 $k=2,\dots ,n-1$ に対しては、

${𝒎}_k=\frac{1}{2}(\frac{{𝒑}_{k+1}-{𝒑}_k}{x_{k+1}-x_k}+\frac{{𝒑}_k-{𝒑}_{k-1}}{x_k-x_{k-1}})$

とし、データ集合の端点では片側差分とする。

カーディナルスプライン

2Dにおけるカーディナルスプラインの例。 線は曲線を表し、四角形は制御点 ${𝒑}_k$ を表す。 曲線は最初の点と最後の点には達していないことに注意。 ただし、これらの点は曲線の形状に影響を与える。テンションパラメータは0.1を用いた。

${𝒎}_k=(1-c)\frac{{𝒑}_{k+1}-{𝒑}_{k-1}}{x_{k+1}-x_{k-1}}$

を用いて接ベクトルを計算すると、カーディナルスプライン（英語: cardinal spline）が得られる 。なお、カーディナルスプラインは正準スプライン（英語: canonical spline）と呼ばれることもある 。 パラメータ テンプレート:Mvar はテンションパラメータであり、テンプレート:Math の範囲でなければならない。 ある意味では、これは接ベクトルの「長さ」と解釈できる。 テンプレート:Math を選択するとすべての接ベクトルがゼロになり、テンプレート:Math を選択すると一様パラメータ化（英語: uniform parameterization）Catmull-Romスプラインになる。

Catmull–Romスプライン

等間隔の横軸における黒点のCatmull–Rom 3次補間の幾何学的解釈

接ベクトルとして

${𝒎}_k=\frac{{𝒑}_{k+1}-{𝒑}_{k-1}}{2}$

を選択することによりカーディナルスプラインの特殊な場合のひとつであるCatmull-Romスプラインが得られる。 「一様な」パラメータ間隔を前提とする。

Kochanek–Bartelsスプライン

テンプレート:Main Kochanek–Bartels スプラインは、データ ポイント ${𝒑}_{k-1}$, ${𝒑}_k$, ${𝒑}_{k+1}$ が与えられた場合に接ベクトルを選択する方法をさらに一般化したものであり、テンション、バイアス、連続性の3つのパラメータを指定できる。

単調3次補間

テンプレート:Main 上記のいずれかのタイプの3次エルミートスプラインを単調データセットの補間に使用する場合、補間された関数は必ずしも単調にはならないが、接ベクトルを調整することで単調性を維持できる。

端点の微分が一致する単位区間の補間

点 ${𝒑}_{n-1},{𝒑}_n,{𝒑}_{n+1},{𝒑}_{n+2}$ のひとつの座標成分を、 整数座標 x = n − 1, n, n + 1, n + 2 で関数 f(x) が取る値とすると、

$p_n=f(n)\mathrm{\forall }n\in \mathbb{Z}$

と表すことができる。

さらに、端点における接ベクトルは、隣接点の中心差分 $$m_n=\frac{f(n+1)-f(n-1)}{2}=\frac{p_{n+1}-p_{n-1}}{2}\mathrm{\forall }n\in \mathbb{Z}.$$ で定義されると仮定する。

実数値 x における補間値 f(x) を評価するには、まず x を整数部 n と小数部 u に分離する。

$x=n+u,$

$n=\lfloor x\rfloor =\mathrm{floor}(x),$

$u=x-n=x-\lfloor x\rfloor ,$

$0\le u<1,$

ここで $\lfloor x\rfloor$ は床関数を表し、x 以下の最大の整数を返すものとする。

Catmull-Romスプラインは $$\begin{array}{cc}f(x)=f(n+u)& ={\text{CINT}}_u(p_{n-1},p_n,p_{n+1},p_{n+2})\\ & =\left[\begin{array}{cccc}1& u& u^2& u^3\end{array}\right]\left[\begin{array}{cccc}0& 1& 0& 0\\ -\frac{1}{2}& 0& \frac{1}{2}& 0\\ 1& -\frac{5}{2}& 2& -\frac{1}{2}\\ -\frac{1}{2}& \frac{3}{2}& -\frac{3}{2}& \frac{1}{2}\end{array}\right]\left[\begin{array}{c}{p}_{n-1}\\ p_n\\ p_{n+1}\\ p_{n+2}\end{array}\right]\\ & =\frac{1}{2}{\left[\begin{array}{c}-u^3+2u^2-u\\ 3u^3-5u^2+2\\ -3u^3+4u^2+u\\ u^3-u^2\end{array}\right]}^{\mathrm{T}}\left[\begin{array}{c}{p}_{n-1}\\ p_n\\ p_{n+1}\\ p_{n+2}\end{array}\right]\\ & =\frac{1}{2}{\left[\begin{array}{c}u((2-u)u-1)\\ u^2(3u-5)+2\\ u((4-3u)u+1)\\ u^2(u-1)\end{array}\right]}^{\mathrm{T}}\left[\begin{array}{c}{p}_{n-1}\\ p_n\\ p_{n+1}\\ p_{n+2}\end{array}\right]\\ & =\frac{1}{2}\left(\right(u^2(2-u)-u)p_{n-1}+(u^2(3u-5)+2)p_n+(u^2(4-3u)+u)p_{n+1}+u^2(u-1)p_{n+2})\\ & =\frac{1}{2}((-u^3+2u^2-u)p_{n-1}+(3u^3-5u^2+2)p_n+(-3u^3+4u^2+u)p_{n+1}+(u^3-u^2)p_{n+2})\\ & =\frac{1}{2}((-p_{n-1}+3p_n-3p_{n+1}+p_{n+2})u^3+(2p_{n-1}-5p_n+4p_{n+1}-p_{n+2})u^2+(-p_{n-1}+p_{n+1})u+2p_n)\\ & =\frac{1}{2}\left(\right((-p_{n-1}+3p_n-3p_{n+1}+p_{n+2})u+(2p_{n-1}-5p_n+4p_{n+1}-p_{n+2}))u+(-p_{n-1}+p_{n+1}))u+p_n,\end{array}$$ ただし、$\mathrm{T}$ は行列の転置を表す。 最後の式はホーナー法を用いた計算手法を示す。

この記述はトリキュービック補間（1回の最適化につき、同じ u と異なる p を使用して CINT_u を16回計算する必要がある）にも関連する。

関連項目

• バイキュービック補間、2次元への一般化
• テンプレート:日本語版にない記事リンク、3次元への一般化
• テンプレート:日本語版にない記事リンク
• テンプレート:日本語版にない記事リンク
• テンプレート:日本語版にない記事リンク
• テンプレート:日本語版にない記事リンク

外部リンク

• Spline Curves, Prof. Donald H. House Clemson University
• Multi-dimensional Hermite Interpolation and Approximation, Prof. Chandrajit Bajaj, Purdue University
• Introduction to Catmull–Rom Splines, MVPs.org
• Interpolating Cardinal and Catmull–Rom splines
• Interpolation methods: linear, cosine, cubic and hermite (with C sources)
• Common Spline Equations