全球 気SPEEDYモデルの実習 (SPEEDY

全球⼤気SPEEDYモデルの実習
(SPEEDY-LETKF)
⼩槻 峻司
理化学研究所 計算科学研究機構
With many thanks to
RIKEN Data Assimilation Research Team
(特に, 近藤さんと岡崎さん!!!)
iTHES データ同化スクール, Sep. 12, 2016 @ 理研IIB (神⼾)
資料1
• スケジュール
• 配布資料の確認
• 注意事項
• 実習班と課題
実習の⽬標
• OSSEを理解する
• データ同化の基本を理解する
• 重要なパラメータの役割を理解する
• SPEEDY-LETKFを⾃分で回せるようになる
OSSE (観測システムシミュレーション実験)
別に,双⼦実験 (idealized twin experiment)とも呼ばれる
③ データ同化実験
② 観測
の⽣成
④ 検証
① 真値の作成
淡路(2009)
OSSEのメリット
• 真値を知っている
• 観測誤差が既知である
• モデルが完全である
• 架空の状況を設定可能
データ同化のコンセプト
データ同化サイクル
誤差共分散⾏列
E[ ] : 統計的期待値
n
平均
E x    xi
分散
2

VAR  x   E  x  E  x  


(mean)
(variance)
共分散
(covariance)
相関
(correlation)
i 1
COV  x, y   E  x  E  x   y  E  y   
CORR  x, y   COV  x, y  / VAR  x   VAR  y 
誤差共分散⾏列
 x  x  x truth
P  E  x  x  


T
対⾓成分: variance
推定変数の確からしさ
⾮対⾓成分: covariance
変数間の関係
Kalman Filter (wo/ model error)
x  Mx
f
t
2. 誤差の発展
Pt f  MPta 1MT
K t  Pt H  HPt H  R 
f
T
x: 状態変数
M: モデル
1. 時間発展
a
t1
f
T
1
P: 誤差共分散
3. カルマンゲインKの取得
H: 観測演算⼦
R: 観測誤差共分散
x ta  x tf  K t  y ot  Hx tf 
4. 状態の解析
Pta  [I  K t H]Pt f
5. 誤差の解析
三好(2005) : アンサンブル・カルマンフィルタ―データ同化とアンサンブル予報の接点
解析の式
1
o
f



x  x  P H  HP H  R   y  Hx 
a
f
x: 推定値
y: 観測値
f
T
f
T
P: 推定誤差共分散⾏列
R: 観測誤差共分散⾏列
H: 観測演算⼦
f: first guess
a: analysis
o: observation
簡単な例
もし⼀変数の問題だったら︖ (Hは変数を直接観測と想定)
p
o
f
 yt  xt 
x  xt 
pr
a
t
f
xtf
yto
a
t
x
p
r
結局は,⾏列で重み付け平均を計算している
 ⼤事なことは,p と r を確からしく設定すること
実習課題
SPEEDY-LETKFによる
データ同化感度実験
SPEEDY-LETKF
SPEEDY: 簡易な全球⼤気モデル
LETKF: EnKFのひとつであるデータ同化
• 感度実験
– 観測分布
– 観測誤差
– インフレーション
– アンサンブルサイズ
– 局所化スケール
観測分布を変えた感度実験
実験する観測分布
観測誤差
OSSEでは, ⾃分で作成しているためRは既知.
実際の問題では, Rは未知である.
R  R
truth
γ: 感度実験のパラメータ
このパラメータを変えて実験する
1
x  x  P H  HP H  R   y  Hx 
a
f
x: 推定値
y: 観測値
f
T
f
T
P: 推定誤差共分散⾏列
R: 観測誤差共分散⾏列
H: 観測演算⼦
o
f
f: first guess
a: analysis
o: observation
インフレーション
EnKFでは,背景誤差共分散⾏列が過⼩評価される
(モデルの⾮線形性, アンサンブルメンバーの不⾜により)
P  P
f
inf
f
Δ: inflation parameter
このパラメータを変えて実験する
1
x  x  P H  HP H  R   y  Hx 
a
f
x: 推定値
y: 観測値
f
T
f
T
P: 推定誤差共分散⾏列
R: 観測誤差共分散⾏列
H: 観測演算⼦
o
f
f: first guess
a: analysis
o: observation
アンサンブルサイズ
背景誤差共分散⾏列をアンサンブルメンバーで近似
スプレッドの⼤きさで確からしさを測る





, xm  x 
X   x1  x , x 2  x , ... 

m: アンサンブルメンバー
1
T
P
X X
m 1
1
x  x  P H  HP H  R   y  Hx 
a
f
x: 推定値
y: 観測値
f
T
f
T
P: 推定誤差共分散⾏列
R: 観測誤差共分散⾏列
H: 観測演算⼦
o
f
f: first guess
a: analysis
o: observation
局所化 (localization)
遠くの観測の情報が⼊りすぎないようにする
(サンプリングエラーの緩和)
1
T
P
X X
m 1
Gaspari and Cohn
(1999)
S(d)
P  S d   P
d
1
x  x  P H  HP H  R   y  Hx 
a
f
x: 推定値
y: 観測値
f
T
f
T
P: 推定誤差共分散⾏列
R: 観測誤差共分散⾏列
H: 観測演算⼦
o
f
f: first guess
a: analysis
o: observation
実習の流れ
• FX-10 (or ⾃前サーバー)にログイン
• SPEEDY-LETKFを⼀度実⾏ (OSSE)
• 感度実験
• 感度実験結果の解析
結果の解析
• 3ツールを準備
– 時系列のRMS Error とアンサンブルスプレッド
– RMS Errorの空間分布
– Adaptive inflation の時系列
ある実験の解析例
まとめ
• OSSEを理解する
• データ同化の基本を理解する
• 重要なパラメータの役割を理解する
• SPEEDY-LETKFを⾃分で回せるようになる
Adaptive inflation methods
iTHES Data Assimilation School
(Sep 13, 2016)
Shunji KOTSUKI
Data Assimilation Research Team, RIKEN-AICS
Today’s goals
• Review: covariance inflation methods
• Understand innovation statistics
• Understand adaptive inflation methods
Kalman Filter
x  Mx
f
t
a
t1
1. 時間発展
η
2. 誤差の発展
Pt f  MPta 1MT  Q
K t  Pt H  HPt H  R 
f
T
x: 状態変数
M: モデル
η: ランダム誤差
f
T
1
P: 誤差共分散
Q: モデル誤差共分散
3. カルマンゲインKの取得
H: 観測演算⼦
R: 観測誤差共分散
x ta  x tf  K t  y ot  Hx tf 
4. 状態の解析
Pta  [I  K t H]Pt f
5. 誤差の解析
Multiplicative inflation
Pb '    Pb
Multiplicative inflation
FCST
b
Ptmp
inflation
P
b
inf
DA
P
a
Posterior
Prior
Prior’
FCST
Prior’
Covariance inflations in Kalman filters
x  Mx
f
t
a
t 1
1. 時間発展
η
2. 誤差の発展
Pt f  MPta1MT  Q
Multiplicative inflation
K t  Pt H  HPt H  R 
f
T
x: 状態変数
M: モデル
η: ランダム誤差
f
T
1
P: 誤差共分散
Q: モデル誤差共分散
3. カルマンゲインKの取得
H: 観測演算⼦
R: 観測誤差共分散
x ta  x tf  K t  y ot  Hx tf 
4. 状態の解析
Pta  [I  K t H]Pt f
5. 誤差の解析
Covariance inflation factor
Tunable parameters of KFs/EnKFs
How can we optimize
these parameters ?
Localization length scale
Miyoshi (2005), Fig. 2.7
Innovation statistics
d o b  y o  Hx b
b: background
d o a  y o  Hx a
a: analysis
d
a b
 Hx  Hx
a
b
o: observation
Desroziers’ statistics (Desroziers et al. 2005)
d
o b
d
oa
(d
o b T
(d
)
o b T
)
 HP H  R
d a b (d o b ) T  HPb H T
R
d a b (d o a ) T  HP a H T
b
T
Derivation
 : 統計的期待値
Definition
ε o  y o  y t , ε b  x b  x t, ε a  x a  x t
R  ε o (ε o ) T , Pb  ε b (ε b ) T , P a  ε a (ε a ) T
Derivation
d o b  y o  Hx b  y o  y t  Hx t  Hx b   o  H b
d o b (d o b ) T   o ( o ) T  H  b ( b ) T H T
  o ( H b ) T  H  b ( o ) T
 R  HPb H T
Adaptive multiplicative inflation
Miyoshi (2011)
d o b (d o b ) T  HPb H T  R
tr  d o b (d o b ) T  R   tr  HPb H T 

tr  d o b (d o b ) T  R 
b
tr  HPest
H T 
Δ: Inflation Factor
Estimated
この式で,Rが不完全(真でない)場合はどうなるでしょうか︖
Kalman filter
Desroziers’ statistics
x tf  Mx ta1  η
d o b (d o b ) T  HBH T  R
Pt f  MPta1MT  Q
Multiplicative inflation
K t  Pt H  HPt H  R 
f
T
f
T
x ta  x tf  K t  y ot  Hx tf 
Pta  [I  K t H]Pt f
1
d
a b
(d
o b T
)
 HBH
T
Implementation: global weather forecasts
Radar
Aircraft
Satellite
Weather balloon
Ship
Buoy
Surface station
Observation data (6-h period)
World’s effort!
(no border in the atmosphere)
(Courtesy of JMA)
A case for Relaxation w/ NICAM-LETKF
Adaptively
estimated α
adaptive
Change in Error
manual
Summary
• Review: covariance inflation methods
• Understand innovation statistics
– Desroziers et al. (2005)
• Understand adaptive inflation methods
– Multiplicative: Miyoshi (2011), Li et al. (2009)
– Relaxation to prior spread: Ying and Zhang (2015)
– Relaxation to prior perturbation: Kotsuki et al.
(prep)