統計数理（2016）
第 64 巻 第 1 号 23–38
c 2016 統計数理研究所


特集「生態学における統計モデリング」［研究ノート］

多変量状態空間モデルを用いたリター分解実験のためのパラメータ推定
仁科 一哉†（受付 2015 年 6 月 30 日；改訂 2016 年 1 月 18 日；採択 3 月 1 日）

要

旨

自然生態系において落枝落葉（リター）の分解のプロセスは，光合成による一次生産と同様に，生態系全体の短期的・長期的な炭素収支を決める重要な生態系プロセスである．生態学における研究では，リターバッグ法と呼ばれる野外培養法によって，リターの分解速度の観測を行ってきた．本研究では LIDET と呼ばれる北アメリカを中心にして行われた複数サイトのリター分
解実験の結果を用いて，リター分解速度定数およびその環境応答を推定するモデルを提案する．
リター分解の時系列データに対して多変量状態空間モデルを構築し，MCMC によるパラメータ推定を行った．用いたデータは 26 の異なるサイトから得られたものであり，各サイト内に繰り返し数が 4 つある時系列データである．採取時間間隔は数カ月毎から各年とサイトにより様々であった．本研究では，落葉広葉樹のサトウカエデと，常緑針葉樹のアメリカネズコの観測結果に対して，同じモデルを導入してパラメータ推定を行った．気温応答の樹種間差異は殆ど見られなかったが，降雨応答に関してはアメリカネズコの方が少降雨時においても，高い分解活性を示す傾向にあった．また観測データでは十分に観測できない，リター分解の季節変動を推定することが可能になった．生態学研究としてはサンプルサイズが大きい部類に入るが，利用したデータはバラつきが大きくかつノイジーであるため，いくつかのサイトでは状態空間モデルの導入によってはじめて速度論的な推定および定量的な環境応答が可能となった．またサン
プルサイズが小さいサイトにおいても，ランダム効果を用いてサイト間の情報を縮約させることによって，妥当な推定を行うことができた．
キーワード：状態空間モデル，リター分解，生態系炭素循環，欠損値，階層ベイズモ
デル，縮約推定．
1.

はじめに陸域生態系の炭素循環において，植物の光合成による炭素同化は大気–陸面の CO2 収支を決める上で極めて重要な要素である．一方で，近年の生態系炭素循環モデルによる複数モデル相互比較実験の将来予測では，生態系の炭素収支について，光合成による炭素固定（第一次生産）よりも，光合成生産物のその後の動態，つまり枯死や落葉後の分解過程，により大きな予測の不確実性を含むことが明らかになってきている（e.g., Friend et al., 2014; Nishina et al., 2014）．
陸域生態系の有機物分解は，その炭素賦存量の大きさから大気 CO2 収支予測の精度向上において非常に重要であると同時に，生態系内の窒素やその他の栄養塩の循環を制御して生態系の一†

国立環境研究所 地域環境研究センター：〒 305–8506 茨城県つくば市小野川 16–2

24

統計数理

第 64 巻

第 1 号 2016

次生産にも影響を与えるため（e.g., Handa et al., 2014），定量的把握を行うことが非常に重要なプロセスである．
生態学（特に生態系生態学と呼ばれる分野）では，古くからリターバッグ法と呼ばれる方法で，様々な植物や環境において落枝落葉などの植物遺体（リター）の分解速度の推定を試みてきた（Bocock and Gilbert, 1957; Wider and Lang, 1982; Berg and McClaugherty., 2004）．この方法は，未分解リターをナイロン製などのメッシュで作成された袋に入れた後，それを対象フィー
ルドの土壌表層やリター層に埋め戻して放置し，定期的にその袋を採取してバッグ内の残存リ
ター重量や炭素含量などを測定し，そこで得られる観測値の時系列変化を観察する方法である．
なお，リターバッグは半年から数年にわたって野外設置し，数カ月から各年で継続的に設置したサンプルを回収して分析する方法が一般的である．
リターバッグ法の観測結果や，そこから得られるパラメータは，陸域生態系炭素循環モデル
や IPCC の将来気候予測にも使われる Earth system model（ESM）と呼ばれるモデルなどのリター分解のプロセスに利用されている（e.g., Potter et al., 1993; Bonan et al., 2013; Wieder et al.,
．しかしながら，観測データの解析に使われるリター分解を記述する統計モデルは，強い2014）
定常状態を仮定して，簡単な定式化で済ますことが多い．リター分解の統計モデルは，リター
分解を一次反応式として取扱い，下記のようなモデルを想定して指数関数でフィッティングする事が殆どである．
(1.1)

dMass
= −kMassdt

ここで，Mass はリターバッグ内の残存重量であり，k は分解速度定数（時間−1 ）と定義される．
この分解速度定数の推定値をそのまま，あるいは参考値として，ESM や陸域生態系炭素循環モ
デルの中で利用される．また，この分解定数を様々な樹種や環境条件中において求めることに
よって，分解機能の生態学的な特性を把握することが，既存の生態学研究では一般的な研究アプローチとなっている（e.g., Berg and McClaugherty, 2004; 大園, 2008）．
リターバッグ法で得られる時系列データは，理想的には単調減少傾向を示し，指数関数でフィッティングしても，見かけ上妥当にあてはまっている研究例が多い．それでもなお，このモデルではいくつかの弱点がある．リター分解は主として微生物による働きによるものである
が，生化学的な酵素活性の温度依存性など，設置期間内の環境要因の変動が考慮されていない．
特に短い野外設置期間である場合は，結果の解釈に注意が必要である．例えば，半年間の野外リター分解試験を日本の森林で行った場合を考える．リターバッグを 7 月に設置して定期的に回収した場合，後半は分解活性が落ちることが考えられるが，これが見かけ上の漸近的減少を示す可能性がある．当然，リターの野外設置時期を半年ずらした場合には，設置期間の後半が高温になるため異なる分解速度定数 k が得られる事になる．この場合において，分解速度定数k の定量的意義付けは非常に難しい．単純に設置期間の平均あるいは積算気温で k の変動を説明することは出来ないだろう．また乾燥条件も微生物活性に影響するため，時期に関わらず分解速度が一定というのは仮定が強すぎる．同じように複数年にまたがるサンプル設置を行った
場合の重量変化に対しての，指数関数のあてはめには注意が必要であろう．一方で，多くの生態系炭素循環モデルにあるリター分解サブモデルでは，下記のように，温度や時間変化を伴う定式化が一般的である（e.g., Potter et al., 1993）．
(1.2)

dMass
= −kf (温度)f (水分)Mass
dt

ここで f (温度)，f (水分) はそれぞれ温度と水分条件の分解応答関数である．先に示した指数関数のフィッティングでは問題となるようなケースでも，この定式化では季節変動の表現が可多変量状態空間モデルを用いたリター分解実験のためのパラメータ推定

25

能になり，より詳細なデータの解釈や定量的な情報を与えることができると考えられる．しかしながら，リターバッグで得られる観測値は，多くの場合サンプリングの時間間隔が不等であったり，時系列に対するサンプルサイズが十分でない．また野外観測のバラつきを推定するために繰り返しがあるため複数時系列データになる．このようなリター分解試験データに対して，生態学では十分な統計モデルが提案されていないのが現状である．
そこで本研究では，状態空間モデルを用いて，生態系モデルのリター分解で用いられているような環境要因の時間・空間変動を考慮したモデルを構築し，観測値を活かしたリター分解モ
デルのパラメータ推定を行う方法を提案する．統計モデルとしては欠損値ありの時系列データ
として取り扱った．また複数の時系列データが利用できるために，多変量状態空間モデルを構築して，リター分解の定量的，速度論的評価を試みる．本研究で用いたリター残存重量の時系
列データは，アメリカで広域的に行われた Long-term intersite decomposition experiment team
（LIDET: Harmon, 2013）のデータベースを活用した．本研究では特に落葉広葉樹と常緑針葉樹
について着目して，それぞれに同一モデルを適用して種間のパラメータ比較を行う．針葉樹と広葉樹のリターは物理的・理化学的性質の差異から大きく分解速度が異なることが知られてい，分解の環境応答特性の違いはあまり言及るが（e.g., Prescott et al., 2000; Cornwell et al., 2008）
されてこなかった．より詳細なリター分解特性の把握に対して，状態空間モデルの枠組みが有用であることを本研究で示す．
2.

材料と方法

2.1 リターバッグ実験のデータ
アメリカで行われた LIDET（Harmon, 2013）のデータを使用した．LIDET は米国内の研究機
関，大学による共同プロジェクトであり，1992 年から 2002 年の間にかけて，共通のプロトコルを用いて北米各地でリターバッグ試験を行ったものである（表 1）．このデータベースには，27の異なるサイトにおける，約 30 種の樹木や草本の葉と根のリターバッグ試験結果が掲載されている．具体的には，それらのリターバッグの重量残存量の時系列データ，およびその理化学的な項目が掲載されている．なお反復のために各サイトには 4 つの異なるリターバッグ設置場所を設けてあり，それぞれ独立のブロックのサンプルとして扱われている．
本研究では，北米に一般的なカエデ科の落葉広葉樹である サトウカエデ Acer Saccharinum（以下 ACSA）とヒノキ科の常緑針葉樹であるアメリカネズコ Thuja plicata（以下 THPL）を扱う．
この 2 つの樹種については，スタンダード試料として各地の実験サイトに同一リターサンプルを配布したものを設置培養している（よって本来はこの 2 つの樹種の自然生息域でないサイトにおいても分解実験が行われている）．リター設置期間は，試験サイトによって異なり，最短で 3から 4 年で，最長で 10 年となっている．データ採取間隔は，温度の高い熱帯では 3 ヶ月から半年とし，それ以外では 1 年おきにサンプリングを行っている．それぞれの樹種についての葉リ
ター残存重量の時系列データを図 1 と図 2 に示す．
本データベースには，リター分解に関するデータの他に各サイトの 30 年平均の月毎の平均気温と降雨量が報告されている．本来は設置期間中の観測データの利用が望ましいが，本研究では 30 年平均の各月データを，リターバッグの設置年数に従って繰り返し使用した．本研究で使用した気温，降雨データを図 3 に示す．
LIDET データベースの中で，異常値として報告されている値や，初期値を大きく超えるような値が報告されている場合には，すべて欠損値として扱った．また，気象データとして降雨量が報告されていない一つのサイトについては今回は取り扱わなかった．よって本研究で扱う全サイト数は 26 となった．

26

統計数理

第 64 巻

第 1 号 2016

表 1．本研究で利用した LIDET 実験サイトの基本情報．

2.2 統計モデル
リター分解は一次反応式として定義した．温度と水分条件が微生物活性に大きく影響を与え，分解速度を制御することが知られている．よって，本研究では気温と降雨量の分解速度定数への制限関数として，これらの環境要因の効果を導入した．
残存リター重量は月次時系列データであり，観測のない月のデータはすべて欠損値として扱った．
まず，尤度を決める観測モデルと状態モデルを以下のように定義した．モデルの記述は Gelman（2004）らの記法に従った．
(2.1)

Mo[0,i,n] ∼ N(M t[0,i,n] , σobs0 )(2.2)

Mo[t,i,n] ∼ N(M t[t,i,n] , σobs )(2.3)

M t[t,i,n] ∼ N(M t[t−1,i,n] exp[−keksite[i] f t[t−1,i] f p[t−1,t−2,i] ], σps )ここで t, i, n はそれぞれ，設置月からの経過した月数，i はサイトの識別 ID，n はサイト内の繰り返し数を示す．Mo は観測されたリター残存重量を示し，Mo[0,i,n] は初期値を示している．一方で，M t は真の重量を示す．N(μ, σ) は平均 μ，標準偏差 σ の正規分布を表す．σobs は

多変量状態空間モデルを用いたリター分解実験のためのパラメータ推定27

図 1．各サイトにおける ACSA（サトウカエデ）リターの残存重量の経時変化．

観測誤差を表す標準偏差を示し，σps はプロセス誤差を示す標準偏差である．k は分解速度定数（month−1 ）を示し，ksite[i] は分解速度定数のサイトランダム効果である．f t[t−1,i] と f p[t−1,t−2,i]はそれぞれ温度と降雨量の制限関数であり，本研究では以下のように定義した．

(T[t,i] − 30) log(Q10 )f t[t,i] = exp(2.4)
10
温度制限関数は，改変した van’t Hoﬀ 式を利用した．ここで T は各サイトの月平均気温であは温度依存性を決めるパラメータで，10o C 上昇した時に，分解速度が何倍り，Q10（Q10 > 1）になるかを示している．生態系モデルにおける分解速度の温度依存性に関する定式化についてEa; Ea は活
は，この他には化学反応の温度依存性を扱うアレニウスの定式化（exp − R(T+273.15)
性化エネルギー，R は気体定数）を用いることもあるが（Davidson and Janssens, 2006; Davidson，van’t Hoﬀ 式などの指数関数による定式化は経験的に当てはまりがよく，多くのet al., 2006）
生態系モデルで使用される関数系の為に本研究ではこちらを採用した．
次に降雨量制限関数を定義する．
 1

(2.5)
f p[t,i] = {1 − exp −c(1 − d)P[t,i] } 1−d本関数は降雨量に対して 0 から 1 の値を取る曲線を示すが（式は von Bertalanﬀy の式を利用

28

統計数理

第 64 巻

第 1 号 2016

図 2．各サイトにおける THPL（アメリカネズコ）リターの残存重量の経時変化．

した），ここでパラメータ c, d はそれぞれ降雨応答の曲線のシェイプを決める．本関数は [0,1)の値を通る単調増加の非線形式であり，降雨が 0 mm の時に値が 0 となり分解を全く抑制するように働く．分解の降雨量応答の定式化に関しては，同じような特徴を示すミカエリスメンテン，
式の様な有理関数が使われる事が多いが（e.g., Raich and Potter, 1995; Hararuk and Luo, 2014）
本モデルで採用した関数の方が立ち上がりの形状の自由度が高いために，このような定式化を行った．
分解速度定数 k に対してはサイト間のバラつきを考慮するために，下記のようにランダム効果を導入した．
(2.6)

ksite[i] ∼ N(0, σid )

分解速度定数 k のサイト間のばらつきは，例えば，土壌の窒素などの栄養塩などに由来する
温度や降雨量条件によらない分解速度のサイト間差を表わすと考えられる．
（U（上限, 下限））を各パラメータの事前分布は表 2 にまとめた．個々の事前分布には一様分布利用した．各パラメータの上限，下限は経験的に十分に大きな値を取っているが，幾つかのパラメータに関してはエキスパートジャッジによって決めた．
以上より事後分布は下記の様に定義される．

多変量状態空間モデルを用いたリター分解実験のためのパラメータ推定

29

図 3．LIDET 試験地における月平均気温と月降雨量．折れ線が月平均気温を示し（左軸），棒グラフは月降雨量を示す
（右軸）．

(2.7)

p(θ, σobs0 , σobs , σps , M t|Mo, T, P) ∝
×
×
×


t≥1

i

i

n

t≥1

i



N(M t[0,i,n] |Mo[0,i,n] , σobs0 )

N(Mo[t,i,n] |M t[t,i,n] , σobs )n

N(M t[t,i,n] |M t[t−1,i,n] , σps , θ)nN(ksite[i] |0, σid )i

×U(σobs0 )U(σobs )U(σps )U(σid )U(k)U(Q10 )U(c)U(d)θ はパラメータベクトルを示す（ここでは θ = {k, Q10 , c, d}）．
事後分布のサンプリングには，WinBUGS ver. 1.4.3（Spiegelhalter et al., 1996; Lunn et al.,を使ってマルコフ連鎖モンテカルロ法により行った．サンプリング条件は，独立したチェー
2000）
ン数が 3 つ，試行回数を 50000 とし，サンプリングした前半の半分を棄却した．その内，自己相関を防ぐために 10 ステップ毎にサンプリングを行って事後分布とした．

統計数理

30

3.

第 64 巻

第 1 号 2016

結果と考察

3.1 状態空間モデルの推定結果とその種間差
ACSA，THPL ともに全てのパラメータにおいて，チェーン間の分散の比である R̂（Gelman
が 1.05 以下であり，収束したと判断した．
and Rubin, 1992）
いずれの樹種においても推移に伴うプロセス誤差 σps よりも，観測誤差 σobs の方が大きくなっている（表 2）．しかしながら，両樹種共にプロセス誤差も無視できないほどの大きさを示している．σps から得られる単純な計算では 95%の信用区間で，10 年後（120 ヶ月後）には ACSA で±6.0 g であり，THPL で ±4.5 g と推定される．各サイトの繰り返し間における差は大体この概算値と比較できるが，実際には分解が早いサイトで数年で 0 g 付近になった場合にはこのモデルから推定された概算は一致しない．
次に分解速度定数を比較すると，ACSA の k は THPL に比べて事後平均値が 2 倍程度であ．針葉樹である ACSA のリターは，リグニンなどの難分
り，高い分解速度定数を示した（表 2）
解性の物質が多く，窒素含有量が低い（Harmon, 2013）．一方で，落葉広葉樹である ACSA は，
リグニン含量が低く，窒素含有量が高いために，微生物にとって利用しやすい性質を持っている．こうした理化学的な性質の違いが分解係数の差に反映されていると考えられる．ランダム，σid で推定結果で示されている効果から求められた各サイト毎の分解速度定数を示すと（図 4）
ように，ACSA ではサイト間の差異が大きく，分解速度定数は最大で 6 倍の開きがあるのに対して，THPL ではサイト間差は 3 倍程度であった．この分解速度定数のサイト間のバラつきは
気温や降雨量によって説明できなかった変動である．
．
分解速度への気温に対する依存性を表す Q10 には樹種間差は認められなかった（表 2，図 5 上）
Q10 は将来の気温上昇による生態系に蓄積された炭素の分解促進の程度を定量的に評価する指標として注目されているが，野外で行われるリターバッグ分解試験では直接的にこのパラメータを推定する統計モデルが無かったこともあり，月平均気温（あるいは地温）に対する Q10 として比較可能な文献値は見当たらなかった．しかし Fierer et al.（2005）の 24 樹種のリターの室内の間には，推定された結果は 95%C.I. まで含めて内包
培養実験の結果にある範囲（2.0 から 3.4）
されていた．よって推定結果は妥当であると考えられる．
に関しては，わずかながら種間差が認められた（表 2，図 5 下）．
降雨量応答のパラメータ（c, d）
ACSA では重み付き降雨量が 300 mm を超えたところで，ほぼ飽和するのに対し，THPL では
200 mm 付近である．特により乾燥した範囲（100 mm 以下）では，樹種間の差は顕著で，THPL
は低降雨時にも分解が抑制されにくい．この差は葉リターの形態的な特徴が関係している可能性がある．THPL リターは鱗片状の針葉で，薄い広葉の ACSA よりもリターそのものの水分保表 2．各樹種におけるパラメータの事後推定結果．

多変量状態空間モデルを用いたリター分解実験のためのパラメータ推定31

図 4．各樹種における分解速度定数 “k” の事後分布．全体の分解速度定数（上）およびランダム（下）の事後推定結果．上図におけ効果を用いたサイト毎の分解速度定数 “k + ksite[i] ”
る実線は “k” の，破線は “k + ksite” の事後分布を示す．下図の数値はサイト ID を示しており，各密度分布の事後平均値直上に示した．

図 5．気温制限関数（上図）と降雨制限関数（下図）の各樹種における推定結果．実線は事後平均より求め，点線は 500 回の事後シミュレーションからのサンプリングによる 95%信用区間を示した．ラグは推定に使用したデータを示す．重み付き降雨に対しては，事後平均を用いて計算した値を示している．

持能が高いことが要因の一つと考えられる．
サイト毎の分解速度定数の事後平均推定値（k + ksite[i] ）を比較すると，強い相関が示された（図 6）．この関係は，種差を超えた，つまりリターそのものの理化学的性質ではない，温度や

32

統計数理

第 64 巻

第 1 号 2016

図 6．ACSA と THPL における各サイトの分解速度定数 k の事後平均値の比較．プロットの数値は Site ID を示す．点線は 1：1 のラインを示す．

降雨量以外のリター分解を制御する共通の環境要因が存在することを示していると考えられる．
データは掲載しないが，どちらの樹種の分解速度定数も，温度や降雨量および蒸発散量の年間値との間には，弱い相関関係ですら認められなかった．この結果は，温度や降雨の分解への影
響を上手く推定していると考えることができる．
これ以上の考察は本論文の趣旨から外れるのであまり掘り下げないが，一般的に，リター分解に関わる支配的な制御因子の強さは，気候 > リターの理化学性 > 土壌の生物相と考えられている（Aerts, 1997）．例えば，同じ LIDET データの解析で指数関数を当てはめた分解速度定数の環境因子との比較を行った研究（e.g., Cusack et al., 2009）では，分解速度定数との温度や降雨量との単純な相関解析に終始していた．ある総説の中では，分解速度定数を求めることは有益．しかし本研究では，温度と降雨量な情報をもたらさないとまで言及している（Prescott, 2005）の主要な気候制御要因を適切に統計モデルに導入したことにより，より一歩進んだ生態学的な意義付けを観測データから抽出できるようになったと考えられる．
一方で，温度，降雨量の 2 つの依存性を正しく評価するには，これらの共変量が十分なレンジを持っていることが必要となる．そのためには本研究で提案したように，複数サイトのデータ．同時に，複数サイトデータの利用とランダム効果のを上手く活用する必要がある（図 3，図 5）
利用は，縮約推定量を扱うことになり，特に時系列サンプル数の少ないサイトに関しては，より確からしい推定を可能にしていると考えられる．
最後に，各サイトの分解特性を明らかにするために，得られた事後分布を使って各サイトの各．
月のターンオーバー速度（分解速度定数に温度，降雨量制限関数を乗じたもの）を求めた（図 7）本推定結果から，各サイトの分解の季節性の定量的把握が可能になった．多くのサイトでは，（ID 1），BSF（ID 5）や夏季にあたる 6–8 月にターンオーバー速度が最大となる．一方で，ANDMTV（ID 20）などにおいては少雨によって（図 3），この時期にターンオーバー速度が抑制されては雨季と乾季のある熱帯性気候に属するが，季節変動がいた．また熱帯地方にある BCI（ID 3）
主に降雨によって引き起こされていることが，本モデルから定量的に推定できている．このよ
うに，状態空間モデルを使うことによって不完全な時系列データから，定量的な分解特性を明らかにすることができる．
本研究では，モデルの構造の不確実性を扱うような，複数モデルの提案をしていない．共変多変量状態空間モデルを用いたリター分解実験のためのパラメータ推定

33

図 7．事後分布から計算される各サイトの月別のターンオーバー速度（month−1 ）．灰色実線
は ACSA，黒色点線は THPL の中央値を示し，灰色の塗りつぶしと斜線塗りつぶしは，それぞれ ACSA，THPL の 95% C.I. を示す．

量の選択のレベルでは，利用できるデータが存在すれば，地温や土壌水分などを利用したほうが，よりリター分解に対して直接的な要因を扱うモデルになる．また降雨量応答に関しては，本（実際には現場の土壌にモデルの定式化では観測月の降雨が 0 となると分解も停止してしまうが水分が保持されている），Hashimoto et al.（2015）では降雨量 P を観測月とその前月との重み付けとして利用する事によって，このような極端な応答を回避する定式化を提案している．他にリター分解モデルの改良として考えられるのは，下記の式の様に難分解性と易分解性など複数のコンパートメントに分ける，二項の指数関数を利用したモデルの導入である（e.g., Wider and．
Lang, 1982; Cusack et al., 2009; Currie et al., 2010）
(3.1)

dMass
= −ka Ffast − kb Fslow
dt

ここでは Mass = Ffast + Fslow と定義している．植物遺体はセルロースのような容易に分解される画分（Ffast ）と，リグニンと呼ばれる難分解性の物質を主とした画分（Fslow ）に化学的に分ける事ができる．これらの割合は樹種によっても大きく異なることが知られている．本研究で構

統計数理

34

第 64 巻

第 1 号 2016

築した状態空間モデルにおいても，2 つの画分を潜在変数として考慮したモデル化を行うことは可能である．しかし，各画分の実測値が無い中で複雑なモデルを当てはめることは，節約の原則に反する．いずれにしても，複数モデルセットの比較を行うにあたっては，時系列データ
かつ階層ベイズモデルの枠組みで利用できる情報量規準やクロスバリデーションの利用が必要となる．
3.2 リター分解における状態空間モデルの利用とその有用性
生態学では動植物の個体動態研究などに早くから状態空間モデルの適用が試みられてきたが（Ecosystem ecology）と呼ばれ（e.g., Schnute, 1994; Clark and Bjørnstad, 2004），生態系生態学
る物質循環を扱う分野では，本来は定量という観点からも相性が良いにもかかわらず，状態空間モデルの適用例は近年になっても限られている．少なくともリター分解実験データに対して
状態空間モデルを適用した論文は，著者の調べた限りでは存在しなかった．単調減少である事が殆どであるリター分解実験データには，状態空間モデルの様な枠組みは不要と思われてしま
うのが利用されない一因かもしれない．しかし，動態に影響しない観測誤差と，動態に影響するプロセス誤差を，適切に扱わなければ時系列変動データを定量的に扱うことは難しい（Clark,．例えば，リターバッグからのリターの意図せぬ脱落，逆にリターバッグへの土壌や現場
2007）
リターの混入や藻類の繁茂，またはリターバッグを置いた場所が土壌ではなく採取期間中たまたま成長して露出した樹木根の上に位置して分解は進まなかったなど，野外研究では様々な予，明ら期せぬトラブルがあり得る．本研究で使用した観測データでも確認できるが（図 1，図 2）
かに単調減少ではないデータが散見される．多くの研究では，アンサンブル平均した値に対してモデルの適用を行う事が多く，このようなバラつきは自然と隠されてしまう．より定量的にデータを扱うには，すべてのサンプルの多変量の時系列データを扱って，2 つの誤差を上手く推定することが重要であると考える．また，状態空間モデルがあまり活用されてこなかった他の理由として，時間を離散化することに抵抗があるのかもしれない．一方で，生態系モデルな
どは差分方程式で表されているため，状態空間モデルで得られた分解速度定数や環境応答へのパラメータは親和性が高い．加えて，活用できる観測された環境データの時間解像度に応じて，
上手くモデル化することが出来るため，離散化にはデメリットよりもメリットの方が大きいと感じている．
多くの生態系モデルにおいては，落葉広葉樹や常緑針葉樹といった機能群で生態系を分類し，個別の樹種に対するパラメータ設定は，特に広域を扱うモデルており（e.g., Potter et al., 1993）
では現実的な問題から設定されない．しかし階層ベイズモデルの枠組みを用いれば，パラメー
タにランダム効果を導入して階層化をすることにより，グローバルパラメータとして，より包括的な推定値を観測データから直接的に得ることができる．また LIDET には広葉樹種，針葉
樹種ともに複数あるが，これらを階層化した統計モデルの構築も可能である．もちろん本研究の両樹種についても，より明示的に統合したモデル化が可能であるし，多種間の分解特性に関わるリターの理化学的性質（e.g., Kurokawa and Nakashizuka, 2008）を考慮したモデルの構築も考えられる．こうした推定値は，陸域生態系の炭素循環の将来予測に使われるプロセスベース
モデルのパラメータとして利用でき，より定量的に観測データの情報を活かすことが出来る方法であると考えられる．
4.

おわりに

近年の計算機の発展やアルゴリズムの改良から，より巨大なモデルのデータ同化が可能になっ，光合成から植物成長，分解の系を含んだプロセスベースの生態系モており（e.g., 樋口, 2011）

多変量状態空間モデルを用いたリター分解実験のためのパラメータ推定35

デルにおいても，大気 CO2 濃度や比較的観測が簡単で大規模な衛星データなどを使った逆推定．しかし，更新によるパラメータの更新は一般的なものになりつつある（e.g., Saito et al., 2014）
できるパラメータ数には限りがあり，生態系のすべてのプロセスに対して一括して行うことは現実的に困難である．一方で，生態学的な研究では野外観測を主とした個別生態系要素のデー
タの蓄積と公開が一般的になりつつあり，利用できる観測情報が豊富になってきている．しかしながら本研究で使用したデータのように，生態学で取得される野外観測データはサンプルサ
イズが小さいということに加え，無視できない誤差がそこかしこに紛れ込んでいるデータが多い．このようなデータから定量的で速度論的な情報を推定するには，本研究で提示した様に利用できるデータに合わせたテーラーメイドなデータ統合的な統計モデルが必要であろう．陸域
生態系の物質循環研究分野におけるモデル研究と野外観測研究の連携の必要性が叫ばれて久し
いが，観測の誤差やプロセスの誤差というものを考慮することが殆ど無い．本研究で示したよ
うに，不確実性を考慮してパラメータ推定を行うのに，ベイズ状態空間モデルによるパラメー
タ推定は非常に有用なツールである．
謝辞

本特集号の企画者である北海道大学の久保拓弥先生と統計数理研究所の持橋大地先生には，本特集号への執筆の機会を頂いたことに対し，心より感謝を申し上げます．査読者には本稿を改善する有用なコメントを数多く頂いたことに深謝致します．また，統計数理研究所の島谷健
一郎先生には，セミナーをはじめとする様々な統計を勉強する機会を提供頂いたことに感謝申2014 “アジア地域におけるチャン
し上げます．なお，本研究の一部は，国環研内競争的資金（A）
バー観測ネットワークの活用による森林土壌 CO2 フラックスの定量的評価” の資金によって執り行われた．

参考

文

献Aerts, R. (1997). Climate, leaf litter chemistry and leaf litter decomposition in terrestrial ecosystems:
A triangular relationship, Oikos, 79, 439–449.
Berg, B. and McClaugherty, C. (2004).『森林生態系の落葉分解と腐植形成』（大園享司 訳）, シュプリンガー・フェアラーク東京, 東京.
Bocock, K. L. and Gilbert, O. J. W. (1957). The disappearance of leaf litter under diﬀerent woodland
conditions, Plant and Soil, 9, 179–185.
Bonan, G. B., Hartman, M. D., Parton, W. J. and Wieder, W. R. (2013). Evaluating litter decomposition in earth system models with long-term litterbag experiments: An example using theCommunity Land Model version 4 (CLM4), Global Change Biology, 19, 957–974.
Clark, J. S. (2007). Models for Ecological Data: An Introduction, Princeton University Press, Oxford.
Clark, J. S. and Bjørnstad, O. N. (2004). Population time series: Process variability, observation errors,
missing values, lags, and hidden states, Ecology, 85, 3140–3150.
Cornwell, W. K., Cornelissen, J. H., Amatangelo, K., Dorrepaal, E., Eviner, V. T., Godoy, O. andQuested, H. M. (2008). Plant species traits are the predominant control on litter decomposition rates within biomes worldwide, Ecology Letters, 11, 1065–1071.
Currie, W. S., Harmon, M. E., Burke, I. C., Hart, S. C., Parton, W. J. and Silver, W. L. (2010). Crossbiome transplants of plant litter show decomposition models extend to a broader climatic range
but lose predictability at the decadal time scale, Global Change Biology, 16, 1744–1761.
Cusack, D. F., Chou, W. W., Yang, W. H., Harmon, M. E., Seilver, W. L. and The LIDET Team

36

統計数理

第 64 巻

第 1 号 2016(2009). Controls on long-term root and leaf litter decomposition in neotropical forests, Global
Change Biology, 15, 1339–1355.
Davidson, E. A. and Janssens, I. A. (2006). Temperature sensitivity of soil carbon decomposition andfeedbacks to climate change, Nature, 440, 1339–1355.
Davidson, E. A., Janssens, I. A. and Luo, Y. (2006). On the variability of respiration in terrestrial
ecosystems: Moving beyond Q10, Global Change Biology, 12, 154–164.
Fierer, N., Craine, J. M., McLauchlan, K. and Schimel, J. P. (2005). Litter quality and the temperaturesensitivity of decomposition, Ecology, 86, 320–326.
Friend, A. D., Lucht, W., Rademacher, T. T., Keribin, R., Betts, R., Cadule, P., et al. (2014). Carbon
residence time dominates uncertainty in terrestrial vegetation responses to future climate andatmospheric CO2 , Proceedings of the National Academy of Sciences, 111, 3280–3285.
Gelman, A. and Rubin, D. B. (1992). Inference from iterative simulation using multiple sequences,
Statistical Science, 7, 457–511.
Gelman, A., Carlin, J. B., Stern, H. S. and Rubin, D. B. (2004). Bayesian Data Analysis, 2nd ed.,
Chapman & Hall/CRC, New York.
Handa, I. T., Aerts, R., Berendse, F., Berg, M. P., Bruder, A., Butenschoen, O., Chauvet, E., Gessner,
M. O., Jabiol, J., Makkonen, M., Malmqvist, B., Peeters, E. T. H. M., Scheu, S., Schmid, B.,
van Ruijven, J., Vos, V. C. A. and Hättenschwiler, S. (2014). Consequences of biodiversity loss
for litter decomposition across biomes, Nature, 509, 218–221.
Hararuk, O. and Luo, Y. (2014). Improvement of global litter turnover rate predictions using a Bayesian
MCMC approach, Ecosphere, 5, article 163.
Harmon, M. (2013). LTER Intersite Fine Litter Decomposition Experiment (LIDET), 1990 to 2002,
Long-Term Ecological Research, Forest Science Data Bank, Corvallis, Oregon. [Database].
http://andrewsforest.oregonstate.edu/data/abstract.cfm?dbcode=TD023 (7 May 2015).
Hashimoto, S., Carvalhais, N., Ito, A., Migliavacca, M., Nishina, K. and Reichstein, M. (2015). Global
spatiotemporal distribution of soil respiration modeled using a global database, Biogeosciences,12, 4121–4132.
樋口知之 (2011).『データ同化入門—次世代のシミュレーション技術—』, 朝倉書店, 東京.
Kurokawa, H. and Nakashizuka, T. (2008). Leaf herbivory and decomposability in a Malaysian tropicalrain forest, Ecology, 89, 2645–2656.
Lunn, D. J., Thomas, A., Best, N. and Spiegelhalter, D. (2000). WinBUGS — A Bayesian modelling
framework: Concepts, structure, and extensibility, Statistics and Computing, 10, 325–337.
Nishina, K., Ito, A., Beerling, D. J., Cadule, P., Ciais, P., Clark, D. B., et al. (2014). Quantifying uncertainties in soil carbon responses to changes in global mean temperature and precipitation,
Earth System Dynamics, 5, 197–209.
大園享司 (2008). カナダにおけるリター分解の地域間比較: CIDET プロジェクトの成果と課題, 日本生態学会誌, 58, 87–101.
Potter, C. S., Randerson, J. T., Field, C. B., Matson, P. A., Vitousek, P. M., Mooney, H. A. andKlooster, S. A. (1993). Terrestrial ecosystem production: A process model based on global
satellite and surface data, Global Biogeochemical Cycles, 7, 811–841.
Prescott, C. E. (2005). Do rates of litter decomposition tell us anything we really need to know?, Forest
Ecology and Management, 220, 66–74.
Prescott, C. E., Zabek, L. M., Staley, C. L. and Kabzems, R. (2000). Decomposition of broadleaf andneedle litter in forests of British Columbia: Inﬂuences of litter type, forest type, and litter
mixtures, Canadian Journal of Forest Research, 30, 1742–1750.
Raich, J. W. and Potter, C. S. (1995). Global patterns of carbon dioxide emissions from soils, Global
Biogeochemical Cycles, 9, 23–36.

多変量状態空間モデルを用いたリター分解実験のためのパラメータ推定

37

Saito, M., Ito, A. and Maksyutov, S. (2014). Optimization of a prognostic biosphere model for terrestrial
biomass and atmospheric CO2 variability, Geoscientific Model Development, 7, 1829–1840.
Schnute, J. (1994). A general framework for developing sequential ﬁsheries models, Canadian Journalof Fisheries and Aquatic Sciences, 51, 1676–1688.
Spiegelhalter, D., Thomas, A., Best, N. and Gilks, W. (1996). BUGS: Bayesian Inference Using Gibbs
Sampling, Version 0.50, MRC Biostatistics Unit, Cambridge.
Wider, R. K. and Lang, G. E. (1982). A critique of the analytical methods used in examining decomposition data obtained from litter bags, Ecology, 63, 1636–1642.
Wieder, W. R., Boehnert, J. and Bonan, G. B. (2014). Evaluating soil biogeochemistry parameterizations in Earth system models with observations, Global Biogeochemical Cycles, 28, 211–222.

38

Proceedings of the Institute of Statistical Mathematics Vol. 64, No. 1, 23–38 (2016)

Statistical Modeling for Longitudinal Data in
Multi-sites Litter-bag Experiment Using Multivariate State-space Model
Kazuya Nishina
Regional Environmental Research Center, National Institute for Environmental Studies

The decomposition of litter in natural ecosystems is an important process as is primary production that occurs through photosynthesis, which together determine the shortand long-term C budgets in terrestrial ecosystems. In the study of ecosystems, the litterbag method is commonly used to evaluate the rate of litter decomposition for various
sites and species while using diﬀerent types of experimental manipulation and litter. In
this study, we used Long-term Inter-site Decomposition Experiment (LIDET) data froma litter-bag experiment database in North America to estimate decomposition constantsof two diﬀerent species, a deciduous, broadleaf species: sugar maple, and an evergreen:conifer, red cedar. We applied a multivariate state-space model to evaluate litter decomposition constants and the responses of environmental factors for leaf litter of these two
species. The LIDET database has time-series data of the remaining mass of litter withfour replicates from each of 26 diﬀerent sites. Bayesian estimation of the state space model
revealed the diﬀerences of litter decomposition constants and litter precipitation responses
between the two species, but no obvious diﬀerence in the temperature sensitivity parameter Q10 was observed. In our model, statistical shrinkage using multi-site time series data
enabled the estimation of plausible decomposition constants even in sites with data having
small sample sizes.

Key words: State-space model, litter decomposition, litter-bag method, carbon cycling.