470 likes | 917 Views
データ 同化 ワークショップ 2014/01/08. 4次元変分法を 用いた 初期値 とパラメータの同時最適化. 伊藤耕 介 海洋 研究開発 機構. ■謝辞:石川洋一先生、淡路敏之先生、 杉浦 望実様 、 川畑拓矢様、本田有機様、加藤輝之様. 数値予報の誤差要因. 初期値の誤差 境界値 の誤差 モデルの不完全性 ( 必要な要素の欠如, パラメタリゼーションスキーム , 固定パラメータの値など ). 今日のメッセージ. 重要な固定パラメータに不確定性があるのなら、制御すべき。. 発表 内容. 4 次元変分法における初期値の最適化
E N D
データ同化ワークショップ • 2014/01/08 4次元変分法を用いた初期値とパラメータの同時最適化 伊藤耕介 海洋研究開発機構 ■謝辞:石川洋一先生、淡路敏之先生、杉浦望実様、 川畑拓矢様、本田有機様、加藤輝之様
数値予報の誤差要因 • 初期値の誤差 • 境界値の誤差 • モデルの不完全性(必要な要素の欠如,パラメタリゼーションスキーム,固定パラメータの値など)
今日のメッセージ 重要な固定パラメータに不確定性があるのなら、制御すべき。
発表内容 • 4次元変分法における初期値の最適化 • 初期値と固定パラメータの同時最適化 • 理論的な話 • 台風状況下の海面交換係数の最適化 • Frequently Asked Question • (余裕があれば)Hybrid EnKF-4DVAR ※簡単のため、今日の話は擾乱の時間発展に線形を仮 定し、観測演算子は線形であるとします。また、 内積が要素の単純積の和である標準内積を考えます。
4次元変分法における初期値の最適化(まずは固定パラメータに誤差がないと仮定する)4次元変分法における初期値の最適化(まずは固定パラメータに誤差がないと仮定する)
典型的な4次元変分法 観測値がない所で同化結果が 第一推定値から離れないように • 評価関数Jを定義 • 解いている問題 • ある固定区間(同化ウィンドウ)で、xがモデルの時間発展に従うという条件のもと、Jを最小化する初期値x0を求める + 観測値と同化結果が 近づくように
評価関数Jを最小にするx0を探索 • 評価関数J の分布は事前に分からない • 方法1: x0の値を変えて予報を繰り返す. ⇒単純明快だが変数が多いと実行不可能 • 方法2: ∂J/∂x0を求める ⇒効率的に最適なx0にたどり着く J x0 ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ xb xupdate xoptimal
∂J/∂x0の計算方法 • xの摂動の時間発展をdδx/dt=Mδxと書く • 勾配に関し以下の式が得られる(導出は略) • 同化ウィンドウをt=0からt=T(>0)まで, ∂J/∂x=0(t=T)として、逆向きに時間積分 • これで∂J/∂x0が計算できる! 観測値とモデル結果 のミスフィット等 アジョイント行列(転置行列に等しい)
x 非線形モデル δxo t 接線形モデル アジョイントモデル
第一推定値からちゃんと真値に近づく 真値・ 同化結果 第一推定値 〇 真の初期値 〇 第一推定値の初期値 〇 同化結果の初期値 同化ウィンドウ:t=0-1, ⊿tobs=0.1 σ=10.0, r=28.0, b=8/3 x0t=10.0, y0t=14.0, z0t=24.0, x0f=11.0, y0f=13.0, z0f=25.0
固定パラメータにも誤差 • ここまでは、固定パラメータに誤差がないとしていた。 • Lorenz63モデルで、固定パラメータrにも誤差があるのに、初期値だけを最適化するとどうなる?(他の固定パラメータは真値とする)
真値 同化結果 第一推定値 〇 真の初期値 〇 第一推定値の初期値 〇 同化結果の初期値
うまくいかなかった • 初期値と固定パラメータの両方に誤差があるのに、初期値だけを制御してもうまくいかない。 • 固定パラメータにも誤差があることをちゃんと考慮しよう。 • J を定義し直して、∂J/∂rを計算。 +
接線形モデル アジョイントモデル ∂J/∂x,∂J/∂y, ∂J/∂zは変更なし&∂J/∂rは既知の∂J/∂yから計算可
今度はうまくいった! 〇 真の初期値 〇 第一推定値の初期値 〇 同化結果の初期値
位相空間での概念図 Jの極小値 ★ 固定パラメータの値 ★ ★ ★ 初期値
現実問題における固定パラメータ • 固定パラメータと言ってもいろいろ • 重力加速度のように、精密に値が観測されており不確定性が小さいもの。 • パラメタリゼーションを行うにあたり、素過程のモデル化に現れる係数 • 後者の場合、観測と照らし合わせることになるが、そのような観測に基づいた値は往々にして、ユニバーサルに適用可能でなかったり、観測の困難さから不確定性が高くなっている。 • もちろん、良いパラメタリゼーションと十分な観測が望ましいが、それが完璧でないときには不確定性があることも考慮に入れておいたほうがよい。
NationalHurricaneCenterの予報結果 120時間予報 120時間予報 96時間予報 96時間予報 72時間予報 72時間予報 最大風速予報の誤差(kt) 進路予報の誤差 48時間予報 48時間予報 24時間予報 24時間予報 1990年 2008年 1990年 2008年 RSMC Tokyo– Typhoon Centerの予報結果 72時間予報 進路予報の誤差 最大風速予報の誤差(m/s) 48時間予報 24時間予報 (沢田雅洋博士提供) 1982年 2010年
WISHE: 台風の最大風速の理論式(Emanuel, 1986) k ~ ■ Ck: 熱交換係数 (~CE: 潜熱交換係数),CD: 摩擦係数 壁雲 眼 CD摩擦係数 CE潜熱交換係数 波浪や波しぶきに依存する係数
2003年以降に提案された摩擦係数と水蒸気交換係数の風速依存性2003年以降に提案された摩擦係数と水蒸気交換係数の風速依存性 ・依然として、強風速域の海面交換係数は 観測が難しく、不確定性が大きい
軸対称台風モデルを用いた最大風速の強風速域におけるCD, CEに対する依存性 100 ⊿CE 0.0 20 0.0 ⊿CD (伊藤, 2011)
水蒸気交換係数に対する勾配 • モデル最下層の比湿の時間発展 • 水蒸気交換係数に対する勾配 水蒸気交換係数 モデル最下層の比湿
物理空間における模式図4次元変分法による海面交換係数の最適化(簡単のため,定常な移流過程だけを考えています)物理空間における模式図4次元変分法による海面交換係数の最適化(簡単のため,定常な移流過程だけを考えています)
軸対称台風モデルを用いた双子実験 ■ある数値実験の結果を「真値」とみなし,「誤った初期値・交換係数」から開始した数値計算に対して観測値を同化し、真値に近い値を復元できるか確認する. ■本研究で行った双子実験の構成要素 DA DA
4 r(km) z(km) NOAA/HRDによる台風の集中観測の一例 (Montgomery et al., 2006) 0 60 r(km) 0
データ同化 ~ 「真の」CDとCE CDとCEの最適化の結果 NoAsm:「誤った」CD, CEによる実験 CDとCEの復元値 擬似的に生成したドロップゾンデ観測 誤ったCD 誤ったCD 真のCE 真のCE データ同化後のCD 真のCD 真のCD データ同化後のCE 誤ったCE 誤ったCE (Ito et al. 2010)
最大接線風速の再現性 • 真値からの平均二乗誤差 “NoAsm”: 7.3 m/s “Asm_NoCoef” 6.3 m/s “Asm_Coef” 2.1 m/s • “NoAsm”や”Asm_NoCoef” • では,台風が系統的に弱い。 • “Asm_Coef” では真値に近い。潜熱交換係数CEが大きくなり,摩擦係数CDが小さい場合に,台風が強くなるというWISHEメカニズムと整合的. 観測値 (Ito et al. 2010)
風速・温位・水蒸気場の誤差(同化期間平均)風速・温位・水蒸気場の誤差(同化期間平均) NoAsm AsmNoCoef AsmCoef 陰影⊿v, ベクトル(⊿u,⊿w) 陰影⊿v, ベクトル(⊿u,⊿w) 陰影⊿v, ベクトル(⊿u,⊿w) 接線風速 の弱まり 接線風速 の弱まり コンター⊿θ, 陰影⊿q コンター⊿θ, 陰影⊿q コンター⊿θ, 陰影⊿q ⊿θ<-2K ⊿θ<-4K ⊿θ<-4K 壁雲域での 水蒸気不足 壁雲域での 水蒸気不足 ⊿θ<-2K ⊿θ<-2K
非静力学メソ4次元変分法(JNoVA)への適用 • 2010年台風第14号(Chaba) • 計算機資源の制約から最盛期の1日間のみ
海面交換係数の補正量(第一推定値=1) • 摩擦係数は強風速域で小さく、水蒸気交換係数で大きくするような修正。 水蒸気交換係数 摩擦係数 Ito et al. 2013
海面交換係数の風速依存性 • NoCoef:初期値のみを最適化 • Coef:初期値及び海面交換係数を最適化
Jobs: 第一推定値と観測値のミスフィット (主要項で全体の90%以上を占める) 評価関数 バルク係数を制御変数に加えることで, ミスフィットは5-20%程度低下した. 1 2 3 4 5 6 7 サイクル数 ミスフィット項の変化率(%) Coef NoCoef Jobs Jobs NoCoef Jobs
TC Chabaの構造の変化:軸対称モデルと整合的 Wind field [Coef] Wind field [Coef-NoCoef] 海面付近での接線風速に大きな差が出た Psrf, (Usrf, Vsrf) ⊿Psrf, (⊿Usrf, ⊿Vsrf) [Coef] [Coef - NoCoef] 気圧が低下し、 中心位置はベストトラックに近く 200 km 200 km
強度と中心位置のずれの改善 • 解析場における台風強度の再現性が向上したほか、解析時刻の中心位置のずれも改善。 中心位置のずれ(km) 最大風速(m/s) 初期値のみ最適化 ベストトラック 初期値と海面交換 係数の同時最適化 初期値のみ最適化 初期値と海面交換 係数の同時最適化
台風の進路予報精度向上 • 初期値と海面交換係数の同時最適化を施した結果、進路予報の精度も向上した。
FAQ 1:チューニングと同じ? • 非常に効率的なチューニングと言える。→5個の固定パラメータを10種類ずつ試すと10^5回の実験が必要。本手法は既存の4DVARシステムとほぼ同じ計算量 J x ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ xa xb xupdate
FAQ 2: 固定パラメータの不確定性をどのようにして与えればよいか? • 物理的にあり得ない固定パラメータにならないように、評価関数の定義で制約をかけておくべきである。 • 台風状況下の海面交換係数については、2003年以降に提案された値を集め、そのばらつきを不確定性を代表する値として用いた。
FAQ 3: 固定パラメータの値にその他の誤差を押し付ける可能性は? • 現実的なデータ同化において、その可能性は十分にある。 • 提案されたモデルと手に入る観測値の範囲内で最適な固定パラメータを探そうとしている。 • ただし、固定パラメータに誤差があることを全く想定せずに初期値に押し付けてしまうよりは、一歩前進している。
FAQ 4: 初期値と固定パラメータどちらの制御が大事か? 問題によるが、一般的に☆相対的に長い時間スケールを扱う☆固定パラメータの不確定性が大☆固定パラメータの場に対する影響が大きいときには固定パラメータの制御が重要になるだろう。
FAQ 5: アジョイントモデルは難しそう • 示した通り、固定パラメータに対する勾配は、アジョイントモデルMTの出力結果だけを使って計算できる。 • グリーン関数法やEnKFを通した固定パラメータ推定法など、アジョイントモデルを使わないパラメータ推定もある。(Menemenlis, 2005;Ruiz et al. 2013)
まとめ (Ito et al. 2010, 2013) • 気象学におけるデータ同化では、初期値を制御するのが一般的であった。 • しかし、不確定性が高く、場に大きな影響を与える固定パラメータについては、最適化を図るのが自然。 • その例として、台風状況下の海面交換係数を制御変数に加えることにより、同化システムの性能が改善することを示した。
4DVARでは、第一推定値の誤差共分散Bを事前に仮定しておく必要がある。4DVARでは、第一推定値の誤差共分散Bを事前に仮定しておく必要がある。 • 通常は、気候値的なBが使われるが、アンサンブル計算の出力を使ってBを作ることで、よりもっともらしくなる。 第一推定値と同化結果のずれ 観測値と同化結果のずれ It is used for B. Model Analysis Obs. EnKFCycle 4D-Var
2011年台風12号における1点疑似観測実験 ・同化ウィンドウ:t=0-1h, Obs: U(t=0,z=20m) ・ベクトル:水平風の修正量 (⊿U, ⊿V) (t=0h) NMC-based B (Typical) Hybrid B 基本場 東風 U, Vの修正量 基本場 西風 渦を南偏させる水平風の修正 渦の構造とあまり関係のない修正 ■EnKFの結果は気象研究所国井勝氏に提供していただいた。