井上 亮 rinoue@plan.civil.tohoku.ac.jp 2011年1月28日 第4章　空間解析 8.　空間補間（応用） 井上 亮 rinoue@plan.civil.tohoku.ac.jp.

Presentation on theme: "井上 亮 rinoue@plan.civil.tohoku.ac.jp 2011年1月28日 第4章　空間解析 8.　空間補間（応用） 井上 亮 rinoue@plan.civil.tohoku.ac.jp."— Presentation transcript:

1 井上 亮 rinoue@plan.civil.tohoku.ac.jp
2011年1月28日 第4章　空間解析 8.　空間補間（応用） 井上 亮

2 ここで学ぶこと 空間補間（空間内挿） について学びます． 最近隣の観測点の値を利用した補間 近隣の観測点の値の平均を利用した補間
近隣点の観測値の単純平均 三角網分割(Triangulation)を用いた加重平均 逆距離による加重平均（IDW） 放射基底関数(Radial-basis function)を利用した補間 ・スプライン補間 バリオグラムを利用した補間（クリギング） 　　について学びます． ここでは，空間補間について学びます．空間補間は空間内挿と呼ばれることもあります． 単純な補間法である，最近隣の観測点の値を利用した補間や 近隣の観測点の値の平均を利用した補間， 放射基底関数を利用した補間， バリオグラムを利用した補間， について学びます ．

3 空間補間とは 限られた地点（観測点）において観測された 空間的事象の値から，観測されていない 任意地点（補間点）における値を推定すること．
（例）アメダスの気温から，自宅の気温を知りたい． 観測点C 19℃ 観測点A 18℃ ?℃ 空間補間とは 空間的事象の限られた地点における観測値を用いて，観測されていない任意の地点における値を推定することです． なお，この講義では，値が観測されている地点のことを「観測点」と呼び，値が観測されておらず空間補間を行いたい地点のことを補間点と呼ぶことにします． 空間補間の例としては，「気象庁がアメダスで観測した気温データを用いて，自宅や職場の気温を推測する」「10mメッシュのDEMの情報から，1mメッシュの標高情報を知りたい」などが挙げられます． 観測点D 22℃ 観測点B 20℃

4 空間補間を行う上での前提 前提 空間補間を行う対象の空間情報には， 空間的自己相関が存在する → 「距離の近い点は，類似性が高い」 より近い点における情報が 空間補間を行う時により参考になる この空間補間を行う上での前提は，4-7で学習した「空間的自己相関」が空間補間をしようとする対象である空間情報に存在するということです． すなわち，「距離の近い点の情報は類似度が高い」ということをいっており， 空間補間を行う際には，より近い点において観測された情報が，より参考になるということを意味します． そこで，まず最も単純な空間補間法として，最近隣点の情報をそのまま用いる方法を紹介します．

5 最近隣点データを利用した空間補間 za zb zd zc 最近隣点法：最も近い観測点の値を用いて 補間点における補間値とする方法．
4章3の「ボロノイ分割」を活用して領域分割しておけば，簡単に任意地点の補間値を知ることができる． 観測点A: za 観測点C: zc 観測点B: zb 観測点D: zd za zb 最も単純な空間補間は．最近隣の観測点における値を補間点の値とする方法です．観測点を母点としたボロノイ分割(4章3)を用いて対象領域を分割しておけば，対象領域内の任意の地点に対して，データを補間することが可能です． zd zc

6 最近隣点データを利用した空間補間 za zb zd zc
<性質> ボロノイ領域の境界で補間値は不連続． 複数の観測点の情報を平均して 空間補間したほうがよいだろう． この辺りでは (za+ zc)/2 ぐらいが 補間値として 自然では? za zb zd zc しかし，最近隣法では，ボロノイ領域の境界で，補間値は急激に変化してしまいます．境界付近では，2点の観測点の平均値ぐらいの値になっていると考えるのが自然でしょう． そこで，1点だけの観測点の情報を使うのではなく，複数の観測点の情報を平均して用いることを考えましょう．

7 近隣点の単純平均 n番目までの近隣観測点の平均値で空間補間 例えば，2番目に近い点までの情報を用いるとすると… B D C A
A-B C-A C-D D-C D-B A-C B-A B-C B-D C-B それでは，複数の近隣観測点の観測値の平均を計算して，空間補間することを考えましょう． 例えば，最近隣とその次に近い観測点の情報を用いて空間補間してみよう．先ほどと同じように2番目に近い点までを求めるようにボロノイ領域を求めます．

8 近隣点の単純平均 <性質> 観測点の近傍の補間点でも観測値と異なる値に なってしまう． →単純平均による内挿には限界 B
<性質> 観測点の近傍の補間点でも観測値と異なる値に なってしまう．　　→単純平均による内挿には限界 B (za+ zb)/2 (za+ zc)/2 (zc+ zd)/2 (zb+ zc)/2 (zb+ zd)/2 A 求まったボロノイ領域に対して，観測値の平均を計算して与えると図のようになります．最近隣点だけを使う場合よりも，多少スムーズに値が変化しているようにも見えますが，やはり補間値は空間的に不連続になります． また，観測点の近傍の補間点の値を推定する場合を考えると，例えば観測値zaの観測点Aの近傍では補間値は(za+ zc)/2となってしまい，一致しません．観測点の位置における補間値が観測値と一致しないのは 補間方法としてあまり適切だとはいえないでしょう．このように近隣点の単純平均による内挿には限界があります． そこで，単純平均ではなく，観測点の近さに応じて重み付けをした平均値を用いて空間補間する方法を考えましょう． D C

9 補間値 z*=za+α (zb-za)+β(zc-za)
三角網分割を用いた加重平均 補間点の”近隣”の3観測点を用いて線形補間． 観測点を母点として作成された 不整形三角網(TIN: Triangular Irregular Network)で領域を分割， 三角形毎に異なる線形の式で補間する． A za C zc B zb D zd 1 α β 補間値 z*=za+α (zb-za)+β(zc-za) 近隣の観測点重み付け平均を行う方法の一つは三角網分割をベースにした補間法です． 補間点の近隣の３つの観測点における観測値を用いて線形補間をします．観測点を母点として作成した不整形三角網(TIN)を用いて補間対象の領域を分割し，三角網毎に線形の式を用いて補間します． 例えば，観測点ABCから構成される三角形の内側にある補間点の値を補間する場合，観測点ABCの情報のみを用いて線形の式で補間します．この式で補間する場合，観測点の位置における補間値は必ず観測値と一致し， 単純平均による補間のように，観測値近傍の補間値が観測値から大きくずれることはありません． また，三角網毎に異なる観測点の情報を用いて補間しますが，その領域の境界上(例えば図の線分BC上)ではどちらの三角網を用いても同じ補間値が得られるため，補間値が空間的に不連続に変化することはないという点も， 単純平均による補間よりも優れていると言えます 観測値の値をz軸にとって三角網分割による空間補間のイメージを示すのが右図です．三角網の頂点が観測点，そのz座標の値が観測値を表しています．観測点以外の補間値は，観測点から構成される三角形の面で表されます． ところで，この補間法で用いる三角網はどのように作成するのでしょうか．

10 三角網分割を用いた加重平均 空間補間を行うのに，どちらの三角網がふさわしいか?
三角網は最小角最大の基準で作られた Delaunay三角形分割で作成 A za C zc B zb D zd A za C zc B zb D zd 赤点の位置の補間値を求める上で，どちらの三角網がふさわしいでしょうか？右図は遠く離れた観測点Dの情報を利用しており，空間相関の性質「より近くの点との類似度が高い」を考えると適切とは言えません． すなわち，扁平な三角形を避ける，すなわち，正三角形に近い三角形を利用することが重要です． Delaunay三角形分割は，正三角形に近い三角形から構成される三角網の作成方法です．三角網の最も小さい内角を最大にするという基準で作成する三角網分割であり，扁平な三角形を避けた 三角網分割が可能です．通常，三角網分割を用いた空間補間を行う場合，Delaunay三角形分割を用いて三角網を作成します．

11 三角網分割を用いた加重平均 補間点の座標(x, y)→補間値 z*を求める
補間点を含む三角形 　　頂点座標(x1, y1), (x2, y2), (x3, y3) 　　頂点での値　z1, z2, z3 頂点と補間点が(x, y, z)空間において同一平面上にあることから 平面線形方程式 を解けばz*が求まる <性質> 観測点では観測値に一致 三角形の境界部分でも連続 三角形毎に異なる線形式で補間するため， 三角形の境界ではスムーズにならない（１階微分は不連続） それでは，より一般的な形で補間値を示しましょう． いま，座標(ｘ，ｙ)の地点における値を補間することを考え，z*を得たい． 補間点を含む三角形の頂点座標が，それぞれ(x1, y1), (x2, y2), (x3, y3)，頂点での値がz1, z2, z3だとします． 頂点と補間点が同一平面上にあることから平面線形方程式を解けば，補間点での補間値z*を求めることができます． さて，三角網分割を用いた補間法では，観測点では補間値は観測値と一致し，かつ，三角形の境界部分でも補間値は不連続ではないことを述べました． しかし，三角形毎に異なる線形式を用いて補間を行うため，図から明らかなように，補間値は境界部分ではスムーズに変化せず，一回微分は不連続となってしまいます． この課題を解決し，スムーズな補間を行う方法として，逆距離を用いた加重平均法（IDW）を次に紹介します．

12 zada-2 +zbdb-2 +zcdc-2 +zddd-2
逆距離を用いた加重平均 逆距離加重平均法 (IDW: Inverse Distance Weighted Average) 補間点から観測点までの距離に関する重みを付けて， 加重平均する方法． (例) 距離2の逆数で重み付け　 zada-2 +zbdb-2 +zcdc-2 +zddd-2 da-2 +db-2 +dc-2 +dd-2 補間値z*= A C B D 観測点までの距離の逆数(逆距離)を 重み付けした加重平均を用いることにより， 近くの観測値を大きく評価した補間が可能． db da 逆距離加重平均法とは，補間点から観測点までの距離の逆数，あるいはそのべき乗などを重みとして加重平均する方法です． 例えば，観測点から補間点までの距離の2乗の逆数で重み付けする場合，補間値z*はこの式のように表されます． 距離の逆数で重み付けすることにより，近くの観測点における観測値を大きく重み付けして平均することになります． dd dc

13 位置ベクトルsiの観測点iにおける観測値ziを用いて (ただしi=1…N)
逆距離を用いた加重平均 位置ベクトルsiの観測点iにおける観測値ziを用いて (ただしi=1…N) 位置ベクトルsの点の補間値 z*を求める ただし， ：補間点と観測点kの距離 ：重みの関数　 通常　　　　　　　が用いられる <性質> 観測点では観測値に一致 観測点以外では１階微分も連続 観測点では1階微分は不連続 適切な重みを定める基準が存在しない (距離減衰パラメータαの設定で補間値が大きく変わる) より一般的な形で記述すると，補間点と観測点の距離から求められる重みによって加重平均するという式で表されます． なお，dは補間点と観測点の距離を求める関数，wは距離から重みを算出する関数で，通常，重みは距離の逆数の定数乗や指数関数で定義されます． 補間値の性質として，観測点では観測値に一致します．また，これまでの方法と異なり，補間する際に領域分割を行わず距離を変数として補間を行うため，連続的な補間値が得られます． ただし，観測点では１階微分は連続にはならず，スムーズな補間値は得られません． また，この方法を用いて補間をする際の問題として，適切な重みを定める基準が存在せず，重みwの関数形を定めたとしても距離減衰パラメータαの設定によって補間結果が大きく変わることがあります． 次に，観測点における１階微分が不連続となる問題を解消する方法である，Radial-basis fuctionを利用したスプライン補間について説明します．

14 放射基底関数を利用した補間 放射基底関数(Radial-basis function)とは 中心 c からのユークリッド距離のみに依存する関数 放射基底関数の足し合わせ で位置ベクトルsiの観測点iで得られた観測値ziを用いて表現し， 観測点以外における補間値を算出する． 放射基底関数の1つである，スプライン関数を使った補間法を紹介 放射基底関数とは，ある中心cからのユークリッド距離のみに依存する関数のことです． この放射基底関数の足し合わせで，位置ベクトルsiの観測点iで得られた観測値ziを表すことにより，観測点以外における補間値を算出します． ここでは，放射基底関数の1つである，スプライン関数を使った補間法を紹介します．

15 スプライン関数による補間 「スプライン関数」とは 連続条件を満たすように多項式を接続した区分多項式 m次のスプライン曲線は1, …, m-1次微分が連続であり 観測点上のスプライン関数の値が観測点の値と 等しいという特性 de Boorの方法 矩形領域・格子上に配置された 観測点におけるデータに対する双3次スプライン内挿 スプライン関数とは，連続条件を満たすように多項式を接続した区分多項式です． それぞれの観測点間で異なる多項式を用いますが，その境界となる観測点で連続的に繋がるようにします． m次のスプライン曲線では，1からm－1次微分が連続であり，かつ観測点状のスプライン関数の値が，観測点の値と等しくなります． ここでは，領域が矩形で，観測点が格子状に配置されている場合に，3次のスプライン関数，すなわち観測点において2階微分まで連続である関数を用いた 双3次スプライン内挿(2次元の位置における属性データに対する3次スプライン内挿)に関して，de Boorの方法を紹介します．

16 スプライン関数による補間 yJ 矩形領域R: x0≦x≦xI; y0≦y≦yJ y0 x0 xI
格子点(xi, yj)の観測値 zij ＝f(xi, yj) (i = 0, 1,…, I; j = 0, 1,…, J) y0 x0 xI 格子点(xi, yj)の観測値 zij ＝f(xi, yj) 領域の境界点における法線方向の1次微分係数 　pij = fx(xi, yj) (i = 0, I; j = 0, 1,…, J) 　qij = fy(xi, yj) (i = 0, 1,…, I; j = 0, J) 領域Rijの4頂点での2次微分係数 　rij = fxy(xi, yj) (i = 0, I; j = 0, J) が与えられている時， これらを満たす双3次スプライン関数S(x, y)は ただ一つだけ存在することが証明されている． x0≦x≦xI; y0≦y≦yJの矩形領域Rを考え，領域内の格子点(xi, yj)の観測値を zij ＝f(xi, yj) とする． また，領域の境界上にある点における法線方向の1次微分係数pij，qijおよび領域の4頂点での2次微分係数rijが与えられている時，これらを満たす双3次スプライン関数はただ一つだけ存在することが証明されている．

17 スプライン関数による補間 まずpij, qij, rijを2I + J - 5個の線形系によって決定 j = 0, 1,…, J に関して
i = 0, 1,…, I に関して i = 0, I に関して 双3次スプライン関数の求め方は，まず，pij, qij, rijをこれらの式を用いて算出．

18 スプライン関数による補間 yj 領域Rij yj-1 xi-1 xi 双3次スプライン多項式 係数γij,mnは で求められる
補間を行う際には，補間点が含まれている格子領域を探し，その領域の多項式を用いて補間値を算出する． 補間点(x, y)の補間値z*の算出 補間点が含まれる領域Rijを探してi, jを決定→z* = Sij(x, y)

19 バリオグラムを利用（クリギング） 空間的自己相関を距離の関数としてモデル化した バリオグラム（あるいは，コバリオグラム（共分散関数））を 利用した空間補間法．(モデル化の前提は4章7を参照 ) 基本的には観測点からの距離に応じて 観測点の値を加重平均する方法． 他の補間法に比べて， 補間値の統計的な性質が明快，かつ， 統計的に優れた性質を持つ． 普遍クリギング(Universal kriging)を用いると 観測値以外の属性情報を利用した空間補間が可能． 次に，4-7で学習した空間的自己相関を距離の関数としてモデル化したバリオグラムを利用した空間補間法について説明します．この方法はクリギングと呼ばれます． 基本的には観測点からの距離に応じて観測点の値を加重平均をする方法です．ただし，逆距離加算法とは異なり，得られた補間値の統計的な性質が明快で，また，統計的に優れた性質を持っています． さらに，普遍クリギングと呼ばれる手法を用いると，観測値以外の属性情報を利用した空間補間が可能です．例えば，気温の内挿を行う際に，観測点における気温情報だけではなく，観測点や補間点の標高値など補間したい情報と関連がある他の情報を活用することが可能です． 以降では，普遍クリギングに関して説明をします．

20 普遍クリギングによる空間補間 補間対象の変数を被説明変数，観測点の属性を説明変数とした線形回帰モデルをたてる．
補間対象の変数を被説明変数，観測点の属性を説明変数とした線形回帰モデルをたてる．　 ただし，線形回帰モデルだけでは， 補間対象の変数を十分に表現ができない． →　線形回帰モデルの攪乱項に空間的自己相関あり． （例）補間対象の変数：気温 　　　観測点の属性：標高 　　　考慮されていない要因：観測時の風向き，植生，土地被覆 そこで，線形回帰モデルの攪乱項の空間的自己相関に対して二次定常性を仮定し，共分散関数を用いて構造化する． 普遍クリギングでは，まず補間対象の変数の観測点における値を被説明変数，その変数と関連性のある観測点の属性を説明変数とした線形回帰モデルをたてます． ただし，線形回帰モデルだけでは補間対象の変数を十分に表現することができないため，線形回帰モデルの攪乱項には空間的自己相関が残ってしまうとします． 例えば，気温を補間する際に，気温を標高で説明するモデルを立てますが，気温観測時の風向きや植生・土地被覆などモデルに反映されていない要因が気温に影響を与えているため， 攪乱項には空間的自己相関が残ってしまうという場合を考えて下さい． そこで，線形回帰モデルの攪乱項の空間的自己相関に対して二次定常性を仮定し，共分散関数を用いて構造化することとします．

21 普遍クリギングによる空間補間 観測点のモデル
y: 被説明変数ベクトル(内挿対象)，X: 説明変数行列， β: パラメータ，u: 攪乱項ベクトル，V: 攪乱項の分散共分散行列， C: 共分散関数，dij: 観測点ij間の距離 一般化最小二乗法(GLS: Generalized Least Squares)より，パラメータの推定値は 観測点のデータからこのようなモデルを作ったとしましょう．攪乱項の分散共分散行列Vの各要素は，バリオグラムを用いて推定した共分散関数を利用して与えられています． このように分散共分散行列の形が分かっている基では，この線形回帰モデルのパラメータ推定は一般化最小二乗法により求められます．

22 普遍クリギングによる空間補間 補間点のモデル このとき，補間点iの被説明変数の線形推定量 を求める．
yi: 補間点の値（確率変数），xi: 補間点の説明変数ベクトル， ui:補間点の攪乱項，σi2: 攪乱項の分散， c:補間点と観測点の攪乱項の共分散ベクトル このとき，補間点iの被説明変数の線形推定量 を求める． 次に，補間点におけるモデルを考えます．観測点のモデルと同様に，補間点における被説明変数の値は，補間点の属性データベクトルとパラメータベクトルの積，足す，補間点の攪乱項という形でモデルが立てられます． 更に，補間点と観測点の攪乱項の共分散も観測点間と同じバリオグラムに従っていると考えるのが自然であるので，攪乱項の共分散構造に関して観測点と同様のモデル化を行います． このとき，補間点における被説明変数，すなわち補間値の線形推定量を求めてみましょう．

23 普遍クリギングによる空間補間 望ましい補間値の推定量は，不偏性（期待値が真の値と等しい）を持ち，かつ，分散が最小
を解くと，補間値が得られる． <性質> Vやcで表される攪乱項の空間的自己相関の構造が真であれば，最良線形不偏予測量 望ましい補間値の推定量は，期待値が真の値と等しくなる不偏性を持ち，かつ分散が最小であるものであるので，その条件から，補間値の式を得ることができます． この補間値の性質としては，バリオグラムを基にVやcで与えられる攪乱項の空間的自己相関の構造が真であれば，最良線形不偏予測量，すなわち期待値が真の値になる不偏性を持った線形の推定量のなかで最も分散が小さいという統計的に優れた性質を持っている．他の補間方法による補間値は統計的な意味が明らかではないのに対して，クリギングによる補間値がこのような統計的性質を持っていることは大きなアドバンテージといえます．

24 参考文献 市田浩三 et al.：スプライン関数とその応用，教育出版，pp.62-74, 高阪宏行：地理情報技術ハンドブック, 朝倉書店, pp.27-58, 張 長平：地理情報システムを用いた空間データ分析, 古今書院, pp , 間瀬 茂・武田 純：空間データモデリング, 共立出版, pp , Cressie, N. A. C.：Statistics for Spatial Data, John Wiley & Sons, pp , Weibel, R. et al.：GIS原典 地理情報システムの原理と応用 [Ⅰ], 古今書院，pp , 1998.