群知能最適化手法を用いた分布型流出モデルのパラメーター同定

2012年3月8日 第56回水工学講演会○小槻峻司*・田中賢治**・小尻利治**・浜口俊雄*** 京都大学大学院 工学研究科** 京都大学防災研究所※ JSSwarm

群知能最適化手法を用いた分布型流出モデルのパラメーター同定
研究の背景
1. 何故，分布型モデルなのか
2 何故，群知能最適化なのか
2.
何故 群知能最適化なのか
3. パラメーター同定は必要なのか

開発中のモデル
人間活動を組み込んだ広域水文モデルを開発中
鉛直方向水収支解析
土壌水分量

陸面過程モデル

灌漑スキームSiBUC (Tanaka,(k 2005))

流出 or
蒸発散

SiBUC (Tanaka,(k 2005))

灌漑取・排水流出 デ
流出モデル
Kinematic wave

放流要求

全球解析例

貯水池操作
流入・放流

ROM(Hanasaki, 2006)

河川流量

水平方向水移動解析

日本域解析例

広域分布型モデルの利点
2011年7月新潟‐福島豪雨事例解析
MAOROSHI
16000

0ave
ensemble
prec

14000

馬下

1010000
8000
6000

204000
2000
030
7/28

7/29

7/30

7/31

8/1

Date (JST)CReSiBUC Rainfall [mm/hr]

River Discharge [m3/s]

Ensemble Hydrograph [m3/s]

複数流域を同時に扱う事が出来る
複数流域を同時に扱う事が出来る水文量を分布情報として得ることができる

Rainfall [mm/m/hr]

Discharge[m3m3/s]

12000

陸面・流出モデルの結合方法
流出解析は物理or概念モデルの選択が可能
水収支解析
水熱収支

流出解析

河道流下解析

物理モデル or 概念モデル

Kinematic wave

どういった解析に概念モデルを用いるか
特定の流域を対象とした場合‐ ピーク流量,タイミングを再現できる
ピグキャリブレーションを要するのが強みかつ弱みキ
リブレ シ ンを要するのが強みかつ弱み‐ トライアンドエラー方式で同定されることが多い‐ 客観的指標により同定する事が望まれる

パラメーター同定
 目的関数を設定する事で最適化問題と認識可能

minimize f ( )  (1 , ,  n )s.t . g j ( )  0, j  1,...q

θ：パラメーターベクトル
f : 目的関数
g : 制約条件

 最適化手法を用いたパラメーター同定
集中型モデルに対する適用例‐ GA, PSO, 進化戦略 等
分布型では 計算時間の制約も大きい
分布型では，計算時間の制約も大きい

 本研究では粒子群最適化手法（PSO）を用いる
PSO手法はGA等従来の最適化手法に比して収束が早い法
等従来 最法束 早

PSO (Particle Swarm Optimization)
PSO: 繰り返し計算により，解を収束させる手法
鳥群や魚群の， 匹の発見した良経路に倣う現象を模倣‐鳥群や魚群の，一匹の発見した良経路に倣う現象を模倣θkθjθi
収束のイメージ

PSO (Particle Swarm Optimization)
PSO: 繰り返し計算により，解を収束させる手法
鳥群や魚群の， 匹の発見した良経路に倣う現象を模倣‐鳥群や魚群の，一匹の発見した良経路に倣う現象を模倣解析

pbesti  min f ( i ) 1, ,t

各個体の最適解

t→t+1

gbesti  min pbestii

グループの最適解

移動ベクトルt t 1
* t
* tv i  wv i  c1rand ( i   i )  c2 rand ( G   i )更新
 t 1itt  i  vi

概念図c,w:収束速度の支配パラメーター（提唱値を用いる）
rand: (0:1)の乱数

目的PSO手法は分布型流出モデルのパラメーター同定に適用可能であるか局所解のリスクと計算時間のトレ ド オフ
局所解のリスクと計算時間のトレード・オフ‐ 局所解リスク ：粒子数
‐ 計算時間：粒子数と計算時間
粒子数と計算時間f(θ)θk

最適解局所解θ
計算時間 =1回の解析時間×粒子数×繰り返し計算数

θjθi
PSO手法に於いて粒子数は局所解に陥らない安定性を表す

計算設定 (阿武隈川流域)□分布型流出モデル
Hydro‐BEAM
Hydro
BEAM (小尻ら(小尻ら, 1998)
基底流出は貯留関数法で表現□目的関数 (流量同定)
CVRMSE 

阿久津観測点
集水面積：1858km2

 t  Qobs (t )  Qsim (t ) n2/ Qobs

□落水線情報GDBD(250m, Masutomi)からアップスケール□気象データ
白河観測点
集水面積：159km2

AMEDAS及び地上気象観測データから作成
蒸発散量は陸面過程モデルSiBUC出力値

Hydro-BEAM(分布型流出モデル)

input rain (rainfall+snowmelt)

表面流出：kinematic wave 法

water level

hSurface runoff
Subsurface runoffA Qrt xQ  A mevapo‐transpiration
recoveryflow

precipitation

surfacerunoff
waste water

基底流出：貯留関数法
Subsurface runoff
InfiltrationdS I Odt

O  (k1  k2 )  S

infiltration groundwaterrunoff

基底流出は概念モデル
2～3段のタンクで表現パラメータ同定値による再現性
 2段タンクの同定結果（横流出2項，縦流出1項）
同定期間(CVRMSE=0 32)
同定期間(CVRMSE=0.32)

検証期間(CVRMSE=0 35)
検証期間(CVRMSE=0.35)

700
600
500
4003002001000AKUTSU (RMSE, 2002)

2discharge [m /s]

obs
Pbest

7/1

8/1

9/1

10/1

10/31

1800
1600
1400
1200
1000
800
600
400
2000obs
Pbest

6/1

7/18/1date○ ：青線 ： 解析流量[m3/s]

9/1

10/1

10/31

date

検証期間(CVRMSE=0.41)

観測流量[ 3/s]
観測流量[m/]AKUTSU (RMSE, 2003)
600

obs
Pbest

500

26/1

discharge [m /s]2discharge [m /s]

AKUTSU (RMSE, 2001)

400
300
200
100
06/1

7/1

8/1

9/1

date

PSO手法は分布型流出モデルのパラメーター同定に適用可能

10/110/31

粒子数の設定
 パラメータの感度を変え4実験を行った
流出率 粗度の寄与
流出率・粗度の寄与
パラメータ感度の違いが収束に与える影響
タンク
3 param.

流出率・粗度

5 param.

3 param.

evapo‐
transpiration
recoveryflow

precipitation

surfacerunoffwaste water

5 param.

パラメーター数が収束に与える影響

i filt ti
infiltration
groundwater
runoff

タンクの寄与 各粒子数で4回収束計算を行った‐ 粒子数：10,粒 数
25, 50, 100, 250, 500 の6通り
通
局所解に陥らない粒子数を把握
局所解に陥らな粒子数を把握※収束の定義：繰り返し計算5回の間，最良解が更新されない事

繰り返し計算数
 収束に要した繰り返し計算数を比較
流出率・粗度の同定

30252015
10503param
5param102550100

Number of Particle

250500

Numbber of Repeat

Numbe
ber of Repeat

タンクパラメータの同定

3025201510503param
5param102550100250500

Number of Particle

赤線： 3パラメーター同定
青線： 5パラメ
5パラメーター同定
タ 同定

収束に必要となる計算数には
１．同定パラメータ数の影響は小さい
２ 感度が高い程多くの計算を要する
２．感度が高い程多くの計算を要する

粒子数の設定
 5パラメータ同定の場合に，局所解に陥る粒子があった
タンク同定の例

RMSE MORE THAN 100BIRD

Pdiv
Pcov

0.60.40.20khb

khc

kvb
Parameter

khd

kvc

Normarized Parameterr

Normarized Parameterr

RMSE LESS THAN 50BIRD
10.8

10.8

Pcov
Pdiv

0.60.40.20

粒子数10 25 50
粒子数10，25，50khb

khc

kvb
Parameter

khdkvc

粒子数100 250 500粒子数100，250，500

 各粒子数4組の計算から，局所解に陥った割合を算出
タンク
粒子数102550100
250500

流出率・粗度

3param.

5 param.

3 param.

5 param.

5075250252575502502505000000255025000１．感度が高いモデルでは局所解に陥るリスクが高い
２．粒子数を増やしても局所解に陥るリスクは存在するため複数の計算が理想．

※局所解の定義：目的関数値が全24回の平均値+σ/2より大きい場合

結論
 PSOアルゴリズムは容易であるが，効果は大きい
 新たにPSOを適用する際に，どう粒子数の設定するか‐ かけられる計算時間からの逆算‐ 計算時間 = 1回の解析時間×粒子数×繰り返し計算数
1週間

1分

700個15回程度‐ 3～5個程度の同定であれば粒子数100は確保したい
3 5個程度の同定であれば粒子数100は確保したい
ただし，同定数の増加は必要粒子数を乗数で増やす見込み

 PSOに，局所解から逃れるスキームを導入する事も可能‐ しかし，結局計算数が増える
しかし 結局計算数が増える‐ 実用的には，複数計算による担保が良い

稲成長モデルへの応用 田植日から，気温を入力条件に稲の成長を予報

DVI (t )   DVR
Development Index

全国一律のパラメーターで各県の出穂・収穫日を予報可能か？

Development Rate

SIMRIW 2000

稲成長を気温（T）と日長（L）の関数で記述
同定パラメ タ同定パラメータ11f j (T ) G j 1  expp  Aj (T  Th j )g j ( L)  1  exp  B j ( L  Lc )

300

Heading(DVI = 2)

Harvest(DVI = 3)

320

Harvest

300

280

280

260

260

240

240

220

220

Heading180

160160

harvestheading±5days

320

200

Planting(DVI = 0)

SIMRIW 2001
340

harvestheading±5daysDOY (Simulated)
DOd)

DVR  fi (T )  gi ( L)

340

2000
180

200

220

240

260

280

300

DOY (Statistical)

320

HarvestHeading

200

180

160
340 160

2001
180200

220

240

260

280

300

DOY (Statistical)

320

340

ご清聴ありがとうございました．
ご清聴ありがとうございました