特集「時空間統計解析：新たなる分野横断的展開」［原著論文］

統計数理（2012）
第 60 巻 第 1 号 109–119
c 2012 統計数理研究所


森林被覆率の非線形回帰モデリング1

2

宮田 大毅 ・西井 龍映 ・田中 章司郎

3（受付 2011 年 7 月 12 日；改訂 11 月 11 日；採択 12 月 27 日）

要旨

森林減少は人間の活動が主要因であるが，同時に地形的な制約も受けている．ここでは森林減少の定量的評価のため，森林被覆率を人口密度および土地の起伏量により説明する精密な回帰モデルを考察する．Tanaka and Nishii（2009, IEEE Transactions on Geoscience and Remote
は森林被覆率のロジット変換を目的変数，人口密度と起伏量を説明変数としたパラメ
Sensing ）トリックな非線形回帰モデルを探索した．ここでは平均構造を加法的な自然スプライン関数で表現した非線形回帰モデルを探索する．また誤差として 1 次近傍，2 次近傍での空間依存性も考慮した精密なモデル化を行う．広島県で観測した実データでの検証により，提案モデルが従来モデルよりはるかに優れていることが示された．
キーワード： ロジスティック回帰分析，自然スプライン，森林減少，人口密度，起
伏量．
1.

はじめに

オゾン層濃度，土壌水分，地表面温度，土地被覆状況など，テラバイト・オーダーの情報が日々，人工衛星・自動化測定機器から集められている．しかし，それらのデータ相互の関係について，特に社会経済的要因がどのように地表面に影響しているかについて，数理モデルを通した研究が不足していると Lambin（1997）が指摘している．森林減少もその例外ではない．
Myers（1990）
は森林減少の主要因は発展途上国における換金作物栽培，牛などの放牧，過度の焼き畑，材木の採取や燃料需要の増加，先進国における空気中汚染物質や酸性雨であると述べている．つまり，これらの要因はすべて人為的な活動と関連している（Holdren and Ehrlich, 1974;
Chown et al., 2003 を参照せよ）．森林減少の定量的評価の研究として，時系列モデルにおけるロジスティック曲線（Esser, 1989），線形回帰モデル（Grainger, 1999; Mather and Needle, 2000;，またロジスティック型曲線（Tanaka and Nishii, 2005, 2009）が議論された．まScrieciu, 2007）
た森林の成長モデルとして，Candy（1997）and Nathoo（2010）, Kneib et al.（2009）がある．
本研究は 田中・西井（2003）, Tanaka and Nishii（2005, 2009）の研究に続くものである．観測した各メッシュの森林被覆率を，当該メッシュにおける人口密度，および地形的要因により森林被覆率を説明することを試みる．約 1km2 の観測区画の座標 s = (u, v) において，森林被覆率（森林面積比率）を Fs（0 ≤ Fs ≤ 1，生態学的な疎林を含む），人口密度を Ns ，Rs をそのメッ
シュでの起伏量，すなわち最高標高と最低標高の差とする．図 1，図 2 を参照せよ．各観測区
123

九州大学大学院 数理学府：〒819–0395 福岡市西区元岡 744
九州大学 マス・フォア・インダストリ研究所：〒819–0395 福岡市西区元岡 744
島根大学 総合理工学部：〒113–8654 島根県松江市西川津町 1060

110

統計数理

図 1.

第 60 巻

第 1 号 2012

観測区画の座標と観測変数．

図 2.

起伏量 R．

画に 1, 2, . . . , n と通し番号をつけ，番号 s により平面上の位置を表すことにする．また全区画の集合を D = {1, 2, . . . , n} で表す．
Tanaka and Nishii（2009）は森林被覆率 F をロジスティック変換を含む巾変換後の値を目的変数とし，人口密度と起伏量を説明変数としたパラメトリックな非線形回帰モデルを探索した．
ここでは各説明変数の効果を加法的な自然スプライン関数で表現した非線形回帰モデルを探索する．また周辺からの影響を， 1 次近傍，2 次近傍それぞれからの影響を取り入れた正規モデルを推定する．ただし 1 次近傍 Us1 は s に隣接した区画番号の集合， 2 次近傍 Us2 はななめに位置する区画番号の集合を表す．図 3 を参照せよ．
2.

パラメトリックな非線形回帰関数

区画 s における森林被覆率 Fs は 0，1，またはその中間の値をとる．そこで全観測区画の集合 D を次のように 3 つに分割する．
(2.1)

D0 = {s ∈ D | Fs = 0}, D1 = {s ∈ D | Fs = 1}, D2 = {s ∈ D | 0 < Fs < 1}.

D0 は森林が全くない領域，D1 は森林に完全に覆われている領域を表す．
Tanaka and Nishii（2005）
は広島県のメッシュデータ（4 節の数値解析例でデータ諸元を紹介）

森林被覆率の非線形回帰モデリング

図 3.

表 1.

表 2.

111

観測区画 s とその 1 次近傍 Us1 と 2 次近傍 Us2 ．

パラメトリックな非線形回帰関数（Tanaka and Nishii, 2005）．

AIC によるモデルの比較（太字は上位 3 モデル，Tanaka and Nishii, 2005）．

で 0 < Fs < 1 を満たす被覆率のロジット変換が，平均が非線形関数 β + g(Ns ) + h(Rs ) で与えられ，誤差が空間依存性を持つ正規回帰モデルを適応した．最初に g(Ns ), h(Rs ) それぞれを階段関数で推定することにより人口密度 Ns ，起伏量 Rs の効果の概形を知り，その階段関数に近いものとして， 表 1 の連続関数の組み合わせを評価した．表 2 は非線形回帰関数の組み合わせを AIC で評価したものであり，（spatial）
は 1 次近傍の影響を考慮したモデルである．その結果，平均構造としては β + g2b (N ) + hs2 (R) が最適であり，1 次近傍での空間従属性を持つ
正規誤差を考慮すると AIC を大きく改良できることが示された．最適モデルは同論文で，日本全体の広域データ，Harbin，Wuhan，Tientsin のメッシュデータでも最適，あるいは準最適と

統計数理

112

第 60 巻

第 1 号 2012

して選ばれた．ここでは森林被覆率 Fs（s ∈ D2 ）のさらなる精密なモデル化を試みる．
2.1 自然スプラインによる平均構造の推定
表 1 で議論したパラメトリックな回帰関数は，推定した母数が不安定であり，また細かい変動を取り込めない欠点がある．そこで回帰関数 g(Ns )，h(Rs ) を 自然 3 次スプラインで近似して，回帰関数を柔軟に推定する．Hastie et al.（2009）の 5 章を参照せよ．
一般に自然 3 次スプラインは基底関数 b1 (x), . . . , bK (x) の 1 次結合で関数 f (x) を次のように
近似する．
(2.2)

f (x) =

K


γi bk (x)k=1

ただし K 個の基底関数は x の 3 次以下の区分多項式：(2.3)

b1 (x) = 1,b2 (x) = x,bk+2 (x) = dk (x) − dK−1 (x)(k = 1, 2, . . . , K − 2)で与えられ，dk (x) は次で定義される．
(2.4)

dk (x) =(x − ξk )3+ − (x − ξK )3+
ξK − ξk

ただし (x)+ = x (x > 0 のとき); = 0 (x ≤ 0 のとき) を表し，ξ1 < ξ2 < · · · < ξK は変数 x の定義域を分割した節点である．自然 3 次スプラインは 2 回微分が連続となる区分多項式であり，両端の区間 (−∞, ξ1 ], [ξK , ∞) では 1 次式となる．基底関数を固定すると式（2.2）は線形モデルとなるため，回帰係数 γ1 , . . . , γK の推定は容易である．
また森林被覆率のロジット変換を巾変換に拡張したモデルも考察する．すなわち 0 < x < 1を満たす値 x に対して巾変換型のロジット変換（p-logit 変換）を定義する．
⎧
θ = 0 のとき,
⎨ log{x/(1 − x)}(2.5)ψθ (x) =θ⎩ {x/(1 − x)} − 1 θ = 0 のとき.
θp-logit 変換は，全ての θ に対して x の単調増加関数となっていることに注意せよ．
3.

森林被覆率のモデル化とモデル評価

森林被覆率 Fs（s ∈ D2 ）に対して 1）目的変数の変換，2）平均構造の再検討，3）誤差の空間依
存性，の観点から精密なモデル化を試みる．
3.1 平均構造の推定の精密化
区画 s における説明変数ベクトルを xs = (Ns , Rs ) とおく．さて s での森林被覆率 Fs は近傍の区画における被覆率，人口密度，起伏量すべてに影響されていると考えられる．そのため本稿では，森林被覆率の p-logit 変換は，xs や s の近傍の目的変数 Ft ，説明変数 xt が与えられたとき，次の正規分布で与えられると仮定する．


2



(3.1) ψθ (Fs ){xs , (Ft , xt ); t ∈ Us1 ∪ Us2 } ∼ N β(s) + g(Ns ) + h(Rs ) +ρret , σ 2r=1t∈Usr ∩D2

ここで定数 β(s) は β(s) = β + β0 c0 (s) + β1 c1 (s) (ci (s) = |Di ∩ Us1 | : 1 次近傍における Di の区画
数，i = 0, 1) と分解され，s の近傍内で D0 , D1 に含まれる区画による影響を考慮した定数である．
また g(N ) は人口密度による影響を kN 個の自然 3 次スプライン g1 (N ), g2 (N ), . . . , gkN (N ) の 1

森林被覆率の非線形回帰モデリング

113

次結合で，h(R) は起伏量による影響を kR 個の自然 3 次スプライン h1 (N ), h2 (N ), . . . , hkR (R)の 1 次結合で近似した関数を表す． Us1 ，Us2 はそれぞれ s の 1 次近傍，2 次近傍，et = ψθ (Ft ) −β(t) − g(Nt ) − h(Rt ) (t ∈ D2 ) は誤差，ρr は ψθ (F ) と近傍 Usr (r = 1, 2) での説明変数等との空間依存性を示すパラメータである．
式（3.1）は中心区画 s および近傍区画における残差 es , et が空間的な相関を持つモデルであるため，2 次近傍までの説明変数も中心区画の森林被覆率に影響するモデルとなっている．な
お Tanaka and Nishii（2005）でのモデル化では 1 次近傍における D2 の区画だけを用いた．
3.2 同時分布と母数推定
森林被覆率の p-logit 変換 ψθ (Fs ) (s ∈ D2 ) の同時分布を導くため，D2 における n 個の目的
変数ベクトルを z θ = (ψθ (F1 ), . . . , ψθ (Fn ))T とおく．このとき，森林被覆率の同時分布は z θ を目的変数とした線形回帰モデルで表現できる．ただし β(s) の加法的分解には定数項 β が含まれているため，自然 3 次スプラインの基底関数のうち g1 (N ) = h1 (R) ≡ 1 は不要である．その結果，計画行列 X : n × (kN + kR + 1) と回帰係数ベクトル β は次で与えられる．
⎛⎞
1 c0 (1) c1 (1) g2 (N1 ) . . . gkN (N1 ) h2 (R1 ) . . . hkR (R1 )
⎜. .
⎟
..
..
..
..
..
..
..
⎟.
(3.2)
X =⎜
.
.
.
.
.
.
.
⎝ .. ..
⎠
1 c0 (n) c1 (n) g2 (Nn ) . . . gkN (Nn ) h2 (Rn ) . . . hkR (Rn )
β = β, β0 , β1 , γ2N , . . . , γkNN , γ2R , . . . , γkRR .
(3.3)

ただし ci (s) = |Us1 ∩ Di |(i = 0, 1) を表す．このとき，式（3.1）より 同時分布z θ ∼ Nn (Xβ, σ 2 (I − ρ1 B1 − ρ2 B2 )−1 )(3.4)

が導かれる．ただし B1 の (s, t) 成分は 区画 t が s の 1 次近傍にあれば 1，その他は 0 の n 次対称行列であり，B2 は 2 次近傍から同様に得られる対称行列である．
母数 (θ, ρ1 , ρ2 ) が既知のとき，式（3.4）から回帰係数ベクトルと分散の最尤推定量は(3.5)

β̂(θ, ρ1 , ρ2 ) = {X  (I − C)X}−1 X  (I − C)z θ ,(3.6)

σ̂ 2 (θ, ρ1 , ρ2 ) = {z θ − X β̂(θ, ρ1 , ρ2 )} (I − C){z θ − X β̂(θ, ρ1 , ρ2 )}/n

（C ≡ ρ1 B1 + ρ2 B2 ）と求められる．これを式（3.4）から与えられる同時密度関数に代入した(θ, ρ1 , ρ2 ) に関するプロファイル尤度関数は次となる．
1
 n|I − ρ1 B1 − ρ2 B2 | 2L(θ, ρ1 , ρ2 ) =n exp −
22{2πσ̂ (θ, ρ1 , ρ2 )} 2
これを用いると β, σ 2 の最尤推定量 β̂ = β̂(θ̂, ρ̂1 , ρ̂2 )，σ̂ 2 = σ̂ 2 (θ̂, ρ̂1 , ρ̂2 ) が得られ，AIC は次で与
えられる．
(3.7)

AIC = n log(2πe) − log(|I − ρ̂1 B1 − ρ̂2 B2 |) + n log σ̂ 2



Fs
+2
log{Fs (1 − Fs )} − θ̂ log+ 2(kN + kR + 1 + 4)
1 − Fs
s∈D
2

ここで
（3.7）
式 2 行目の第 1 項は，確率変数 Fs を巾変換したことによる変数変換のヤコビアンに由来する．またペナルティー項のうち kN + kR + 1 は未知回帰係数ベクトル（3.3）の次元数，
4 は未知母数 {θ, ρ1 , ρ2 , σ 2 } の数を表す．
注意：式（3.7）はフルモデルに対する AIC であり，θ = 0（ロジット変換），ρ2 = 0（1 次近傍だ
けの空間相関を考慮）， ρ1 = ρ2 = 0（空間相関のない独立モデル）等の部分モデルを考える場合
は，未知母数の数を減らす必要がある．

統計数理

114

第 60 巻

第 1 号 2012

図 4. 広島県の人口密度と起伏量の頻度分布．

4.

解析例：広島県（a）
人口密度 N （，b）
起伏量 R の頻度分布である．観測メッシュは
図 4 は広島県の観測領域の1km2 であり，被覆率で 3 グループに分けたときの D0 , D1 , D2 の標本数はそれぞれ 159，1713，
6825 である．地形的に見ると，沿岸部の平坦な領域に広島市や福山市等の都市部があり，県東北部や県西部は山間部となっている．
4.1 平均構造のスプラインによる推定
0 < Fs < 1 を満たす
（s ∈ D2 ）被覆率を p-logit 変換し，その平均を 2 つの説明変数の加法モデルで表現する．回帰関数として自然 3 次スプラインを用いるが，その引数としては人口密度や起伏量を対数変換した log(Ns + 1)，log(Rs + 1) を採用する．
まず人口密度に関する回帰関数 g(Ns ) を log(Ns + 1) の多項式を用いて次のように近似する．
log(Ns + 1) の変域を [0, ξ1 ], [ξ1 , ξ2 ], . . . , [ξkN −1 , ξkN ], [ξkN , ∞] と等分割し，対応する自然 3 次ス
プラインの kN 個の基底関数 gi (N ), i = 1, . . . , kN により g(Ns ) を近似した．同様に h(Rs ) に対して，log(Rs + 1) の多項式 hi (N ), i = 1, . . . , kR で近似した．
次に節点の個数 kN , kR を決定する．まず被覆率のロジット変換 z 0 = (ψ0 (F1 ), . . . , ψ0 (Fn )) に対して，等分散で独立な誤差を持つ正規線形回帰モデル：
z 0 = Xβ + eを考え，次の手順で kN , kR を決定した．
（3.2）
を用いた z 0 の回帰モデルを考（1） 3 以上 10 以下の節点数 kN , kR に対して，計画行列
える．
（2） 回帰モデル（1）について，AIC に基づく自然 3 次スプラインの係数に対する変数減少法
で，より小さいモデルが選ばれなかった（全ての係数が選ばれた）(kN , kR ) を記録する．
（3）（2）で選ばれた節点数の組み合わせのなかで，AIC 最小のものを最適な節点数とする．
この探索により kN = 5，kR = 7 が選ばれた．そこで今回考察する全てのモデルで，この節点数を用いた．もちろん回帰係数はモデルごとに再推定している．また ρ1 , ρ2 の最尤推定量は
分散共分散行列 I − ρ1 B1 − ρ2 B2 が正定値行列となる条件下でグリッドサーチにより求めた．
表 3 は 1）
ロジット変換と p-logit 変換，2）独立モデルと空間従属モデルを比較している．表 2

森林被覆率の非線形回帰モデリング
表 3.

115

自然スプラインに基づく加法モデルによる母数推定とモデル比較（太字は上位 3 モデル）．

表 4.

自然 3 次スプラインに基づく加法モデルによる回帰係数の推定値．

から 1 次近傍までの従属性を考えた森林被覆率をロジット変換した最適パラメトリックモデルの AIC は −10559.6 であった．同じ設定で平均構造をスプラインで近似すると −11689.2（表 3の 3 行目）となり改善された．これは平均構造を柔軟な回帰関数に変えたことによる効果である．さらに p-logit 変換したモデルは単なるロジット変換したモデルより AIC が小さいことがわかる．また空間依存性を考慮すると，AIC は大幅に改善されている．
最適モデルでは空間相関が 0 < ρ2 < ρ1 と推定されていることや，β0 < 0（森林被覆率ゼロの近傍区画の影響），β1 > 0（森林で覆われた区画の影響）と推定されていることから，モデルが有効に機能していることがわかる．
表 4 はスプラインにかかる回帰係数をまとめたものである．独立モデルと従属モデルには違いがあるが，その他のモデルでは，係数に大きな差は見られない．表 5 は 2 次近傍での被覆率0 または 1 の区画からの影響を考慮したモデルの評価表であり， 1 次近傍で同様の考察を行った表 3 の AIC をわずかに改良している．
さてロジット変換かつ独立モデルの場合と，p-logit かつ 2 次近傍までを考慮した場合のスプラインによる平均構造を比較する（表 3 における最悪モデルと最良モデル）．
図 5 は被覆率に対する人口密度の影響をみたものである（横軸の節点の設定位置は表 6 の 2行目を参照）．左図は表 4 の 2 行目，右図は表 4 の 10 行目で与えられている回帰係数を用いてグラフ化したものである．人口密度の増加に伴い，森林被覆率を変換した値が減少することが

統計数理

116
表 5.

第 60 巻

第 1 号 2012

空間相関が 2 次近傍までの自然スプラインによる 6 つの推定モデルと AIC（周辺画素の平均への影響を定数にモデル化，太字は上位 3 モデル）．

わかる．両図の違いは，人口密度の高いときの挙動にある．
図 6 は被覆率に対する起伏量の影響をみたものである（横軸の節点の設定位置は表 6 の 3 行
目を参照）．左図は表 4 の 6 行目，右図は表 4 の最下行に対応する．両グラフは似ていて， Rが大きいところでは直線的に増加している．これは Tanaka and Nishii（2005）で選ばれた R の関数形は 表 1 の hs2 であり，対数関数となっていることと合致している．
一方 R が小さい（横軸の値が 4 以下）では図 6 は複雑な変動をしている．表 1 の hs2 では立
ち上がりの θ（巾指数ではない）は 田中・西井（2003）により 25.84（m）と推定されている．なお
図 6 の横軸の目盛り 4 は R = exp(3.1824) − 1 = 23.1045（m）であるため，本研究は，R < 25 の領域における R の森林被覆率への影響をより細かく観察したといえる．ただし，R ≤ 23.1045を満たす標本数は 271 と他区間の標本数と比べて少なく，この区間での関数の変動を過度に重要視する必要はないであろう．
5.

ディスカッション

森林被覆率が 0 < F < 1 の場合，人口密度と起伏量で説明するロジット型回帰問題を次のようにモデルを精密にした．1）F = 0, 1 となる周辺からの区画の影響をモデル化し，2）スプラ
インによる平均構造を導入し，3）2 次近傍までの空間依存性を考慮し，4）森林被覆率の p-logit
変換を導入し， 5）同時分布から母数推定法を確立した．以上を広島県のデータに応用して下記が判明した．
（1） 回帰関数をパラメトリックな非線形関数から自然 3 次スプラインに変え，また D0 や D1の近傍区画の影響をモデル化することにより，AIC が 1000 程度改善した．また自然ス
プラインにより，平均構造の柔軟な変動を検出できた．
（2） 被覆率のロジット変換から p-logit 変換に変えた変数を目的変数とすることにより，AIC
が 50 − 100 程度改善した．
（3） 空間従属性について，独立モデルから 1 次近傍の従属モデルに変えることにより，AIC
が 400 程度改善した．さらに 2 次近傍までの従属性を考えると， AIC が 50 程度改善

森林被覆率の非線形回帰モデリング

図 5.

変換した森林被覆率に対する人口密度 N の回帰関数（左図：ロジット変換 (θ = 0)，独立モデル，右図：p-logit 変換 (θ = −0.03297)， 2 次近傍まで）．

図 6.

変換した森林被覆率に対する起伏量 R の回帰関数（左図：ロジット変換 (θ = 0)，独立
モデル，右図：p-logit 変換 (θ = −0.03297)， 2 次近傍まで）．

表 6.

117

図 5，図 6 における横軸の節点の設定位置．

した．
母数は 0 < ρ2 < ρ1 ，β0 < 0，β1 > 0 と期待通りの関係を満たす値として推定され，モデルが有効に機能していることがわかる．
本稿のモデリング手法は，空間従属性を持つ目的変数，従属変数の関係を精密に探索する新
たな可能性を提示したといえる．平均構造を自然スプラインで近似する場合，もとより区間数が重要である．また生の変数の多項式を用いるのか，あるいは単調変換した変数の多項式を用いるのかも重要なポイントである．区間幅の取り方にも注意が必要である．
次の課題として，近傍の説明変数も平均モデルに取り込んだ平均構造のモデル化が考えられる．また同時に他地域への応用と先行研究との比較が必要である．
謝

辞

2 人の査読者には拙論を丁寧に読んでいただき，本論文を改訂する上で建設的なコメントを数多く頂戴しました．また編集者には誤りの指摘や編集上のご配慮をいただきました．あわせ

118

統計数理

第 60 巻

第 1 号 2012てお礼申し上げます．本研究は，科学研究費基盤研究（A（）19200020）時空間現象データに対す
る統計科学モデルの構築及び解析に関する組織的研究（研究代表者：矢島美寛（東京大学）
）
，及
び基盤研究
（B（
）19300096）超高次元データの分類手法の導出とその理論的性質の解明，および
実データへの応用の研究（研究代表者：西井龍映（九州大学））の援助をうけました．

参 考 文 献Candy, S.（1997）. Estimation in forest yield models using composite link functions with random（1）, 146–160.
eﬀects, Biometrics, 53
Chown, S. L., van Rensburg, B. J., Gaston, K. J., Rodrigues, A. S. L. and van Jaarsveld, A. S.

（2003）. Energy, species richness, and human population size: Conservation implications at a（5）
, 1233–1241.
national scale, Ecological Applications, 13
Esser, G.（1989）. Global land-use changes from 1860 to 1980 and future projections to 2500, EcologicalModelling, 44, 307–316.
. Constraints on modelling the deforestation and degradation of tropical open
Grainger, A.（1999）
woodlands, Global Ecology and Biogeography, 8, 179–190.
. The Elements of Statistical Learning, 2nd ed.,Hastie, T., Tibshirani, R. and Friedman, J.（2009）
Springer-Verlag, New York.
Holdren, J. P. and Ehrlich, P. R.（1974）. Human population and the global environment, American
Scientist, 62, 282–292.
. Variable selection and model choice in geoadditive
Kneib, T., Hothorn, T. and Tutz, G.（2009）regression models, Biometrics, 65（2）, 626–634.
. Modelling and monitoring land-cover change processes in tropical regions,
Lambin, F. E.（1997）
Progress in Physical Geography, 21（3）, 375–393.
. The relationships of population and forest trends, The
Mather, A. S. and Needle, C. L.（2000）
Geographical Journal, 166（1）, 2–13.
Myers, N.（1990）. The world’s forests and human populations: The environmental interconnections,（supplement）, 1–15.
Population and Development Review, 16
Nathoo, F. S.（2010）. Space-time regression modeling of tree growth using the skew-t distribution,
Environmetrics, 21（7–8）, 817–833.
. Can economic causes of tropical deforestation be identiﬁed at a global level?,
Scrieciu, S.（2007）
Ecological Economics, 62, 603–612.
田中章司郎，西井龍映（2003）. 人口増加に伴う森林減少の空間モデル，応用統計学，32（1）, 1–5.
Tanaka, S. and Nishii, R.（2005）. Veriﬁcation of deforestation in East Asia by spatial logit models due
to population and relief energy, Proceedings of 12th SPIE Europe International Symposium
on Remote Sensing 5976, No. 59760W-1, 1–10.
Tanaka, S. and Nishii, R.（2009）. Non-linear regression models to identify functional forms of deforestation in East Asia, IEEE Transactions on Geoscience and Remote Sensing, 47（8）, 2617–2626.

Proceedings of the Institute of Statistical Mathematics Vol. 60, No. 1, 109–119 (2012)

119

Nonlinear Regression Modeling of Forest-coverage Ratios
Daiki Miyata1 , Ryuei Nishii2 and Shojiro Tanaka3
1 Graduate School of Mathematics, Kyushu University
2 Institute of Mathematics for Industry, Kyushu University
3 Interdisciplinary Faculty of Science and Engineering, Shimane University

The worldwide decrease of forest area is chieﬂy due to the negative impact of humanactivities, but this is restricted to topographical circumstances. In this research, detailed
regression models are proposed to elucidate forest area ratio observed with a grid-cellsystem by two covariates: population density and relief energy, which is the diﬀerence
between the maximum and minimum altitudes at each site.
Tanaka and Nishii (2009, IEEE Transaction on Geoscience and Remote Sensing)
explored a regression model of logit-transformed forest area ratios, where the mean function is given by the sum of two non-linear parametric functions of the covariates. Basedon their results, we consider the mean structure represented by the sum of two natural
spline functions of the respective covariates. The spatial dependency from the ﬁrst- andsecond-order neighborhoods is also considered. The proposed models applied to real data
in Hiroshima prefecture demonstrate outstanding superiority to the previous models interms of AIC.

Key words: Logistic regression, natural spline, deforestation, human population density, relief energy.