スライド 1

第1回データ同化ワークショップ
(2011年4月22日、気象研究所)
データ同化と数値天気予報
1.データ同化法のレビュー
2.マルチスケール・データ同化
気象研究所予報研究部
露木 義
データ同化
観測と数値シミュレーションの融合手段
境
界
条
件
観
測
デ
ー
タ
包括的な4次元データセット
数値モデル
データ同化の目的

対象
時間発展する大規模かつ複雑なシステム
(大気、海洋、・・・)

目的
数値モデルの初期値の推定
数値天気予報、数値シミュレーション
包括的な4次元データセットの作成
研究の基盤的データ
数値モデルのパラメータの推定
経験的定数、モデル誤差、発生・吸収源
観測システムの評価・設計 ・・・
データ同化・予測問題
• 状態変数 x(t ) の時間発展( dβ t :White noise )
dx  F(x, t )dt  G(x, t )dβt ,
dβt dβt   Q(t ) (t  t)dtdt
T
• 確率密度関数 p(x, t ) の時間発展:Fokker-Planck 方程式

p n 
1 n 2

( pFi )  
p GQ GT
t i 1 xi
2 i , j 1 xi x j

ij
(伊藤積分の場合)
• 観測データ y による p(x, t )の収縮:Bayesの定理
p(x(tk ) | y 0 ,, y k ) 
p(y k | x(tk )) p(x(tk ) | y 0 ,, y k 1 )
n
p
(
y
|
x
(
t
))
p
(
x
(
t
)
|
y
,

,
y
)
d
x(tk )
k
k
k
0
k

1

この問題を逐次モンテカルロ法で解く
粒子フィルタ
(大自由度系への適用可能性が大きな課題)
大自由度系のデータ同化法
確率密度関数 p(x, t ) の平均値またはモードの
直接推定、及びそれらの周りの共分散の推定
• 平均値

x a   x :  x p(x | x f , y o ) d n x

カルマンフィルタ (Kalman, 1960)
線形最小分散推定


x a  x f  K y o  H (x f ) ,

K : BHT R  HBHT

1
アンサンブルカルマンフィルタ (EnKF) (Evensen, 1994)
• モード
x a  x m : arg max p(x | x f , y o )
x
MAP推定または最尤推定
変分法解析 (Sasaki, 1958)
Minimize: J (x)   log p(x f | x)  log p(y o | x)  log p(x)
Adjoint法による4次元変分法 (4DVar)
(Lewis and Derber, 1985; LeDimet and Talagrand, 1986)
状態変数の時間発展
•初期値:平均値
d  x
1
 F(x, t )  F( x, t )  TrFP     F( x, t )
dt
2
モデルが非線形だと、予測値は平均値でない
アンサンブル予報
•初期値:モード(ランダム誤差がない場合)
dx m
 F(x m , t )   p( p) 1   F x x ,
p(x, t ) : 確率密度関数
dt
m
モデルが非線形でも   F   0を満たせば、予測値もモード
決定論的予報
•正準Hamilton力学系:   F  0
•Lorenz-63モデル、 Lorenz-96モデル、順圧渦度方程式
のスペクトルモデル、・・・
4DVarとEnKFの比較
4DVar
統計的推定法
解析値
時間発展
EnKF
MAP推定(最尤推定)
線形最小分散推定
モード
平均値
決定論的予報
アンサンブル予報
・数値モデルと観測演算
子のアジョイント・コード
・サンプリングエラーや
システム構築の手間 の作成
フィルタの発散などを抑
・背景誤差共分散行列の えるためのチューニング
設計
解析誤差共分散
非線形の影響
別途、要計算
自然に得られる
最適値探索に困難
平均値でなくなる
弱拘束4DVar
モデル誤差の考慮
数値モデルの方程式
dx  F(x, t )dt  εt
1.制御変数: 初期値とモデル誤差
J (x 0 , ε t ) 




1
f T
f
x 0  x 0 B 1 x 0  x 0
2
T
1
1 T T
 H ( M (x 0 , ε t ))  y o R 1 H ( M (x 0 , ε t ))  y o   ε t Q 1ε t dt
2
2 0


2.制御変数: 状態変数
J (x t ) 


1
f
x0  x0
2
 B x
T
1

0  x0

f




1
1 T
T
o T
1
o
 H ( M (x t ))  y R H ( M (x t ))  y   x t  ~
xt  Q 1 x t  ~
xt dt
2
2 0
(~
xt :モデル誤差を無視して計算した状態変数)
現業数値天気予報では、2009年にECMWFで初めて導入された。
Long-window 弱拘束4DVar
弱拘束にすることによって、長い同化ウィンドウを採用
強拘束4DVar
x1
x4
x2
x3
Tremolet (2006)
評価関数
J (x1, , xn )  J1 (x1 )  J n (xn )  J q (x1, x2 )   J q (xn1, xn )
J q (x k , x k 1 ) :
1
( M (x k )  x k 1 ) T Q 1 ( M (x k )  x k 1 ),
2
• J k (xk ), J k (xk ) (k  1,, n)
(Q : モデル誤差共分散)
は並列的に計算できる
• 同化ウィンドウが十分長ければ、背景誤差共分散は必要ない
• 最適値探索のためのPre-conditioningが課題
Ensemble 4DVar (En4DVar)
Adjoint モデルが不要
・ アンサンブル予報(m メンバー)から背景誤差共分散を近似計算

,


1
x1f  x f ,, x mf  x f
m 1
・ 状態変数 x をアンサンブル空間で近似(w :m 次元ベクトル)
B  X f X f
X f :
T
x  x f  X f w
・ 評価関数とその勾配ベクトル(制御変数:w )



T
1 T
1
w w  H ( M (x f  X f w ))  y o R 1 H ( M (x f  X f w ))  y o
2
2
1
J (w )  w  (HMX f ) T R 1 H ( M ( x f  X f w ))  y o
2
J (w ) 



・ 勾配ベクトルをアンサンブル予報を使って近似計算
HMX f 


1
H ( M (x1f ))  H ( M (x f )),, H ( M (x mf ))  H ( M (x f ))
m 1
Liu et al (2008), Zupanski, (2005), …
En4DVarの共分散局所化
誤差共分散のサンプリングエラー を抑える
C:相関行列
Lorenc (2003)
En4DVarでは・・・


B  C  X f X f


T
(  :Schur 積)

f
f
x  x f  C  X1 ,  , C  Xm w
r個
1
xkf  x f ,  , xkf  x f 
C  CCT ,
X k f :
m 1
C : n  r 行列 (r  n)
w : m  r ベクトル
高度なデータ同化法の現状
- 変分法
> 4DVar (strong constraint, weak constraint)
> En4DVar
- カルマンフィルタ
> Perturbed observation EnKF
> Square root EnKF (EnSRF, EAKF, LETKF)
- ハイブリッド
> 3DVar + EnKF
> 4DVar + EnKF
- 粒子フィルタ
Background
Observation
EnKF (PO)
EnKF (SRF)
Lawson and Hansen (2004)
データ同化法の比較
カナダ気象局の全球モデル(400×200×58L):2007年2月
NH-Z500
SH-Z500
NH-Z500
4DVar/Benkf
4DVar
SH-Z500
4DVar/Benkf
3DVar/Benkf
En4DVar
EnKF
• EnKFと4DVarは同程度の精度
• 4DVarや3DVarにEnKFで計算したB を使うと精度向上(ハイブリッド)
• En4DVarは、EnKFや4DVarとほぼ同程度の精度
Buehner et al. (2010)
数値予報センターの予報誤差の比較
3DVar
4DVar
3DVar
3DVar 4DVar
4DVar
全球数値予報の精度の推移
データ同化法の課題
1.非線形
• データ同化における非線形の強さは、観測
データの時空間密度に依存する。
• 非線形は非ガウス分布をもたらす。
MTSAT-1R
・Long window 弱拘束4DVar
・粒子フィルタとのハイブリッド、・・・
2.マルチスケール
• 大気海洋結合データ同化
• 全球雲解像データ同化、・・・
NICAM-3.5km
• 階層的アプローチ、・・・
H Miura et al. Science 2007;318:1763-1765
次世代スパコンで初めてできること
アンサンブル手法に基づく
広域雲解像データ同化システム
を構築
高精度な雲解像
マルチスケール解析値
積雲対流と大規模環境場
を整合的に解析
① 集中豪雨などの1日程度先
までの確率的予測の初期値
② 積雲対流を含むスケール
間相互作用の解明に貢献
後半は途中結果の発表でしたので、誤解
を避けるために、今回は割愛しました。