四点((−9, 5) , (−4, 2) , (−1, −2) , (7, 9) ) に対する三次補間多項式 L (x ) (黒破線) は基底多項式 y0 ℓ 0 (x ) , y1 ℓ 1 (x ) , y2 ℓ 2 (x ) , y3 ℓ 3 (x ) のスケールを変えたものの和である。この補間多項式は与えられた四つの制御点を通り、各スケールされた基底多項式は対応する制御点を通り、それ以外の三つの制御点に対応する x のところでは 0 になっている。
数値解析 におけるラグランジュ補間 (ラグランジュほかん、英 : Lagrange interpolation )は、多項式補間 に用いられる。相異なる点の集合 xj および数値 yj に対し、そのラグランジュ補間多項式は、各 xj において対応する値として yj をとるような次数 最小の多項式である。このように次数最小の多項式は一意に決まるが、決定する方法は複数存在するため、「ラグランジュ補間多項式」という名称をその一意な多項式の「ラグランジュ形」というふうに言及するのは正確でない。
名称はジョゼフ=ルイ・ラグランジュ に因んだものだが、ラグランジュの発表する1795年よりも以前に、この方法を初めて発見したのは1779年のエドワード・ワーリング である。ラグランジュの結果はレオンハルト・オイラー が1783年に発表したより複雑な形の公式の簡単な帰結となるものであった[1]
ラグランジュ補間多項式は数値積分 法の一種ニュートン–コーツ法 でも用いられ、また有限体 上で計算されたラグランジュ補間多項式は暗号理論 におけるシャミアの秘密分散法 (英語版 ) でも用いられる。
ラグランジュ補間は巨大振幅に関するルンゲ現象 の影響を受けやすい。また評価点 xj の変更に関して補間の計算を全くやり直す必要があるから、そのような目的では変更が容易にできるニュートン補間 がしばしば用いられる。
定義
k + 1 個のデータ点集合
(
x
0
,
y
0
)
,
… … -->
,
(
x
j
,
y
j
)
,
… … -->
,
(
x
k
,
y
k
)
{\textstyle (x_{0},y_{0}),\ldots ,(x_{j},y_{j}),\ldots ,(x_{k},y_{k})}
はどの二つの xj も同じでないとする。これらのデータのラグランジュ型の補間多項式 とは、ラグランジュ基底多項式
ℓ ℓ -->
j
(
x
)
:=
∏ ∏ -->
0
≤ ≤ -->
m
≤ ≤ -->
k
m
≠ ≠ -->
j
x
− − -->
x
m
x
j
− − -->
x
m
=
(
x
− − -->
x
0
)
(
x
j
− − -->
x
0
)
⋯ ⋯ -->
(
x
− − -->
x
j
− − -->
1
)
(
x
j
− − -->
x
j
− − -->
1
)
(
x
− − -->
x
j
+
1
)
(
x
j
− − -->
x
j
+
1
)
⋯ ⋯ -->
(
x
− − -->
x
k
)
(
x
j
− − -->
x
k
)
{\displaystyle \ell _{j}(x):=\prod _{0\leq m\leq k \atop m\neq j}{\frac {x-x_{m}}{x_{j}-x_{m}}}={\frac {(x-x_{0})}{(x_{j}-x_{0})}}\cdots {\frac {(x-x_{j-1})}{(x_{j}-x_{j-1})}}{\frac {(x-x_{j+1})}{(x_{j}-x_{j+1})}}\cdots {\frac {(x-x_{k})}{(x_{j}-x_{k})}}}
の線型結合
L
(
x
)
:=
∑ ∑ -->
j
=
0
k
y
j
ℓ ℓ -->
j
(
x
)
{\displaystyle L(x):=\sum _{j=0}^{k}y_{j}\ell _{j}(x)}
を言う。
まず最初の仮定により
x
j
− − -->
x
m
≠ ≠ -->
0
{\textstyle x_{j}-x_{m}\neq 0}
であるから、この式は常に矛盾なく定まる 。
x
i
=
x
j
{\textstyle x_{i}=x_{j}}
かつ
y
i
≠ ≠ -->
y
j
{\displaystyle y_{i}\neq y_{j}}
となるような対を許さない理由は、必要な補間多項式 L (
y
i
=
L
(
x
i
)
{\textstyle y_{i}=L(x_{i})}
) が存在しないからである(各引数 xi に対する値は一つでなければならない)。他方、
y
i
=
y
j
y_{i}=y_{j}
もそれら二点が実際には単一の点となる。
任意の i ≠ j に対して、基底多項式 ℓj (x ) は分子に (x − xi ) を因子として持つから、x = xi のとき積全体も零となる。他方、i = j のときは明らかに ℓi (xi ) = 1 が成り立つ。すなわち、
ℓ ℓ -->
j
(
x
i
)
=
δ δ -->
i
j
{\textstyle \ell _{j}(x_{i})=\delta _{ij}}
である。ここに右辺はクロネッカーのデルタ 。したがって各点 xi において
L
(
x
i
)
=
y
i
+
0
+
0
+
⋯ ⋯ -->
+
0
=
y
i
{\textstyle L(x_{i})=y_{i}+0+0+\dotsb +0=y_{i}}
となり L は所期の函数の補間を与えるものとなる。
注意
ラグランジュ多項式の集合に対する補間の発散の例
ラグランジュ型の補間では、多項式補間の線型性と補間多項式の一意性があるため、証明や理論的な議論では都合がよい。この一意性はヴァンデルモンド行列 の可逆性(同じことだがヴァンデルモンド行列式 が至る所消えないこと)からも導けることである。
しかしながら、構成からわかる通り、節点 xk を変更するごとに、ラグランジュ基底多項式をすべて計算し直さなければならない。そこで実用上(あるいは計算上)の目的では、重心形のラグランジュ補間(後述 )やニュートン補間 を用いるほうが良い。
等間隔に取った節点に関するラグランジュ補間やその他の補間は、真の函数の上下に振動する多項式を生じるものである。この振る舞いは節点の数を多くするにつれてルンゲ現象 と呼ばれる発散に逢着しやすくなっていく。この問題は、補間に用いる点をチェビシェフ節点 (英語版 ) に選ぶことで解消できる[2] 。
線型代数学的説明
補間問題を解くことは、逆行列を計算する線型代数学 の問題につながる。標準的な単項式基底 を用いた補間多項式
L
(
x
)
=
∑ ∑ -->
j
=
0
k
x
j
m
j
{\textstyle L(x)=\sum \limits _{j=0}^{k}x^{j}m_{j}}
では、L (xi ) = yi を L の係数 mj に関してヴァンデルモンド行列
(
(
x
i
)
j
)
{\textstyle ((x_{i})^{j})}
を逆に解かなければならない。対して、ラグランジュ基底を用いて
L
(
x
)
=
∑ ∑ -->
j
=
0
k
l
j
(
x
)
y
j
{\textstyle L(x)=\sum \limits _{j=0}^{k}l_{j}(x)y_{j}}
を作れば、この場合先ほどはヴァンデルモンド行列が現れた部分には単位行列
(
δ δ -->
i
j
)
{\textstyle (\delta _{ij})}
が現れ、逆行列は単位行列自身であるから、ラグランジュ基底は自動的に逆に解かれていることになる。
この構成は中国の剰余定理 と類似対応している(つまり、素数を法とする整数の剰余を調べる代わりに、一次式を法とする多項式の剰余を見ている)。
重心形
ラグランジュ基底多項式を
ℓ ℓ -->
(
x
)
=
(
x
− − -->
x
0
)
(
x
− − -->
x
1
)
⋯ ⋯ -->
(
x
− − -->
x
k
)
ℓ ℓ -->
′
(
x
j
)
=
d
ℓ ℓ -->
(
x
)
d
x
|
x
=
x
j
=
∏ ∏ -->
i
=
0
,
i
≠ ≠ -->
j
k
(
x
j
− − -->
x
i
)
{\displaystyle {\begin{aligned}\ell (x)&=(x-x_{0})(x-x_{1})\cdots (x-x_{k})\\[5pt]\ell '(x_{j})&={\frac {d\ell (x)}{dx}}{\Big |}_{x=x_{j}}=\prod _{i=0,i\neq j}^{k}(x_{j}-x_{i})\end{aligned}}}
を用いて、
ℓ ℓ -->
j
(
x
)
=
ℓ ℓ -->
(
x
)
ℓ ℓ -->
′
(
x
j
)
(
x
− − -->
x
j
)
{\textstyle \ell _{j}(x)={\frac {\ell (x)}{\ell '(x_{j})(x-x_{j})}}}
と書き直す。これは重心重み付け (barycentric weight )[3] を
w
j
=
1
ℓ ℓ -->
′
(
x
j
)
{\textstyle w_{j}={\frac {1}{\ell '(x_{j})}}}
と定めれば簡潔に
ℓ ℓ -->
j
(
x
)
=
ℓ ℓ -->
(
x
)
w
j
x
− − -->
x
j
{\displaystyle \ell _{j}(x)=\ell (x){\frac {w_{j}}{x-x_{j}}}}
と書くことができる、これを重心ラグランジュ補間の「第一形」と呼ぶ。
この形の多項式補間を考えることの利点は、補間多項式
L
(
x
)
=
ℓ ℓ -->
(
x
)
∑ ∑ -->
j
=
0
k
w
j
x
− − -->
x
j
y
j
{\textstyle L(x)=\ell (x)\sum \limits _{j=0}^{k}{\frac {w_{j}}{x-x_{j}}}y_{j}}
を評価するときに、重み wj が事前に分かっていれば、O(n ) で計算できることである(ラグランジュ基底 ℓj (x ) を個別に計算するのは O(n 2 ) 掛かる)。もうひとつ重心形補間の利点として、新しい節点 x k +1 の追加も各 wj を
(
x
j
− − -->
x
k
+
1
)
{\textstyle (x_{j}-x_{k+1})}
で割って、新しい w j +1 を計算するだけで容易にできる点が挙げられる。
さらに第一形を単純化することもできて、まず定数函数 g (x ) ≡ 1 の重心補間
g
(
x
)
=
ℓ ℓ -->
(
x
)
∑ ∑ -->
j
=
0
k
w
j
x
− − -->
x
j
{\textstyle g(x)=\ell (x)\sum \limits _{j=0}^{k}{\frac {w_{j}}{x-x_{j}}}}
を計算してから L を g で割れば、得られる
L
(
x
)
=
∑ ∑ -->
j
=
0
k
w
j
x
− − -->
x
j
y
j
∑ ∑ -->
j
=
0
k
w
j
x
− − -->
x
j
{\displaystyle L(x)={\frac {\textstyle \sum \limits _{j=0}^{k}{\frac {w_{j}}{x-x_{j}}}y_{j}}{\sum \limits _{j=0}^{k}{\frac {w_{j}}{x-x_{j}}}}}}
は与えられた節点における補間性を失わない。この補間多項式を重心補間多項式の「第二形」あるいは「真の形」という。真の重心補間多項式は、L の各節点における評価に際してラグランジュ基底 ℓ を評価しなくてよいという点で有利である。
誤差項
与えられた函数 f を、節点 x 0 , …, xn において次数 n の多項式で補完するとき、誤差
R
(
x
)
=
f
(
x
)
− − -->
L
(
x
)
{\displaystyle R(x)=f(x)-L(x)}
は
R
(
x
)
=
f
[
x
0
,
… … -->
,
x
n
,
x
]
ℓ ℓ -->
(
x
)
=
ℓ ℓ -->
(
x
)
f
(
n
+
1
)
(
ξ ξ -->
)
(
n
+
1
)
!
,
x
0
<
ξ ξ -->
<
x
n
,
{\displaystyle R(x)=f[x_{0},\ldots ,x_{n},x]\ell (x)=\ell (x){\frac {f^{(n+1)}(\xi )}{(n+1)!}},\quad \quad x_{0}<\xi <x_{n},}
と表せる[4] 。ただし、
f
[
x
0
,
… … -->
,
x
n
,
x
]
{\textstyle f[x_{0},\ldots ,x_{n},x]}
は差商 である[注釈 1] . またこの誤差項を複素領域における周回積分
R
(
z
)
=
ℓ ℓ -->
(
z
)
2
π π -->
i
∫ ∫ -->
C
f
(
t
)
(
t
− − -->
z
)
(
t
− − -->
z
0
)
⋯ ⋯ -->
(
t
− − -->
z
n
)
d
t
=
ℓ ℓ -->
(
z
)
2
π π -->
i
∫ ∫ -->
C
f
(
t
)
(
t
− − -->
z
)
ℓ ℓ -->
(
t
)
d
t
{\displaystyle R(z)={\frac {\ell (z)}{2\pi i}}\int _{C}{\frac {f(t)}{(t-z)(t-z_{0})\cdots (t-z_{n})}}dt={\frac {\ell (z)}{2\pi i}}\int _{C}{\frac {f(t)}{(t-z)\ell (t)}}dt}
と書くこともできる。
この誤差項によって、誤差の範囲を
|
R
(
x
)
|
≤ ≤ -->
(
x
n
− − -->
x
0
)
n
+
1
(
n
+
1
)
!
max
x
0
≤ ≤ -->
ξ ξ -->
≤ ≤ -->
x
n
|
f
(
n
+
1
)
(
ξ ξ -->
)
|
{\displaystyle |R(x)|\leq {\frac {(x_{n}-x_{0})^{n+1}}{(n+1)!}}\max _{x_{0}\leq \xi \leq x_{n}}|f^{(n+1)}(\xi )|}
と見積もることができる。
脚注
^ 差分商とも呼ぶ。近い概念として、有限差分 の商 Δy/Δx を差分商 あるいは差商と呼ぶ場合もあるので混同してはならない。
参考文献
^
Meijering, Erik (2002), “A chronology of interpolation: from ancient astronomy to modern signal and image processing”, Proceedings of the IEEE 90 (3): 319–342, doi :10.1109/5.993400 .
^ Quarteroni, Alfio; Saleri, Fausto (2003), Scientific Computing with MATLAB , Texts in computational science and engineering, 2 , Springer, p. 66, ISBN 9783540443636 , https://books.google.com/books?id=fE1W5jsU4zoC&pg=PA66 .
^ Jean-Paul Berrut & Lloyd N. Trefethen (2004). “Barycentric Lagrange Interpolation”. SIAM Review 46 (3): 501-517. doi :10.1137/S0036144502417715 .
^ Abramowitz and Stegun, "Handbook of Mathematical Functions," p.878
関連項目
外部リンク