水資源工学 ( 立川担当分 1 回目 :12 月 21 日資料 ) 流況評価の方法 (1) 講義の目標 :1 年を通じた河川流量の特徴を河川流況という 流況を評価することは 河川による水資源を把握することである 今日の講義では 流況を評価する基本的な手法と流況の量的な特性を表す物理量 ( 水文量 ) を理解する 次に 水資源確保の基本量である 平均して 10 年に 1 回程度の割合で発生する渇水流量 ( 基準渇水流量 ) を推定する手法を理解する 1) 河川流量の時間変化 ( 流量ハイドログラフ ) 河川流量の時間変化を図示したものを流量ハイドログラフという 流量データは水位流量曲線を用いて水位データを変換して得られることが大半である 連続的な流量データから時間平均したデータを時間流量 日平均したデータを日流量 月平均した流量を月流量という 水位 流量曲線の例 2) 様々な河川のハイドログラフ 最上川の年間の日流量変化 吉野川年間の日流量変化 3) ハイドログラフから得られる河川流況の情報 4) 流況曲線一年間の日流量ハイドログラフを流量の多い順に並べ替えて 大きい順に左から図示したものを流況 1
曲線という 流況曲線から得られる流況の特徴量として 豊水流量 ( 一年間を通じて 95 日はこれを下 らない流量 ) 平水流量 ( 一年を通じて 185 日はこれを下らない流量 ) 低水流量 ( 一年通じて 275 日は これを下らない流量 ) 渇水流量 ( 一年を通じて 355 日はこれを下らない流量 ) が得られる 流況曲線の例 5) 基準渇水流量平均的に 10 年に 1 回の頻度で発生する渇水流量を基準渇水流量という 新たに取水を予定する場合は 基準渇水流量 から 河川維持流量と他の水利使用者の双方を満足する水量 ( 正常流量 ) を差し引いた水量の範囲内として決定される 河川維持流量 とは 塩害防止 各種排水の希釈浄化 河道の維持 河口埋塞防止 水生動植物の生存繁殖等 河川に関する公利の確保 公害の除去若しくは軽減のため流水の果す機能を確保するための流量をいう 6) 水文頻度解析手法を用いた基準渇水流量の設定 10 年確率渇水流量 ( 基準渇水流量 ) の求め方は以下の通りである 1 数十年以上の日流量データを収集する 2 毎年の渇水流量を抽出する 3 水文頻度解析手法 * により 渇水流量が適合する確率分布モデルとその母数を推定する 確率分布モデルは 渇水流量についてはワイブル分布を用いることが多い 4 非超過確率が 1/10 となる渇水流量を基準渇水流量とする * 参考図書 : 工学のための確率 統計 : 北村 堀 ( 編著 ) 朝倉書店 2006. 例題で学ぶ水文学 : 椎葉 立川 市川 森北出版 2010. 2
水資源工学 ( 立川担当分 2 回目 :1 月 4 日資料 ) 流況評価の方法 (2) 講義の目標 : 流況の特性として 量的な特性と時間的な特性がある 量的な特性として 流況曲線や基準渇水流量がある 時間的な特性としては 降水量や水位 流量の時系列的な特性がある 今日の講義では流量データの時系列特性と時系列モデルの基礎を理解することを目的とする 1) 水文時系列時間の経過とともに変動する現象の記録を時系列 (time series) といい 時間変化する降雨強度や河川の水位 流量などの水文量を総称して水文時系列という 2) 水文時系列解析の目的と手順 1 水文時系列データの確率 統計的な特性を分析する 2 時系列データの確率 統計的な特性を再現する時系列モデルを構成する 3 構成した時系列モデルを用いて 観測時系列と同じ確率 統計的な特性を有する長期間の時系列を模擬発生させ 水工施設の設計に利用する あるいは 時系列モデルを用い 現在までに得られた観測時系列から将来の水文量の変動を予測する 3) 定常時系列と非定常時系列時系列は その時間的は変動によって定常時系列と非定常時系列に分類される 時間とともに時系列の確率的特性が変化しない時系列を定常時系列という 一方 ある期間で見た時系列の平均値が時間とともに変動したり 平均値の周りのばらつきが時間とともに変化したりする場合など 時間とともに時系列の確率 統計的特性が変化する時系列を非定常時系列という 4) 定常時系列とその特性値離散的な時系列を,,, とし これらの変量の同時確率分布関数を考える 時系列を表す確率分布関数が 時間をずらしても変化しない場合 すなわち,,, を時系列の同時確率分布関数とし を任意の時間間隔を表すインデックスとして,,,,,, が成り立つ場合 その時系列は強定常であるという 特に 二次までのモーメントについて μe E Var Var Cov, E Cov, E が成り立つ場合 弱定常あるいは二次定常という E は期待値 Var は分散 Cov は共分散を求める演算記号である これらをもとに自己共分散関数 自己相関関数 が定義される Cov, E Cov, Var Var 3
5) 時系列モデル (AR モデル ) 定常時系列 を過去の時系列の値と白色ノイズの線形和で表現したモデル を次数 の自己回帰モデルあるいは AR モデル (AutoRegressive model) といい AR() と表す また AR モデルに従う時系列を AR 過程という 1,2,, は自己回帰係数 は平均値 0 で分散 σ の正規分布に従う白色雑音であり はそれより前の時系列 とは無相関とする もっとも簡単な 1 次の自己回帰モデル AR(1) は である AR(1) のモデルパラメータは,, となる の観測時系列が得られており それが AR(1) 過程で表現できるとして モデルパラメータを求める まず μe である また の分散を とすれば 2 となるので 11 という条件のもとに 1 が得られる 次に 時間差が 1 の自己共分散関数 は 時間差が 2 の自己共分散関数 は となる 一般に時間差 k > 0 について の関係が得られるので 時間差 k > 0 の自己相関関数 は / となる 観測時系列データから標本自己相関関数を求めれば の推定値 を定めることができる が 定まれば 観測時系列データから得られる標本分散 用いて 1 となり 1 次の自己回帰モデル AR(1) のパラメータが定まる 4
3 : 1 11 AR III AR(1) AR(1) AR(1) y n 3 { 1 f Y (y n ) = (y n c)ζ 2π exp 1 2 cλζ 3 y n ( ) } 2 ln(yn c) λ ζ (1) x n = ln(y n c) λ ζ (2) x n N(0, 1) y n 1 ρ y x n 1 ρ x ρ y = exp(ζ2 ρ x ) 1 exp(ζ 2 ) 1 (3) 1) 1 ρ y (3) ρ x 2) (2) y n x n σ 2 v = (1 a 2 )σ 2 x v n a = ρ x 3) AR(1) x n = ax n 1 + v n x n 4) (2) x n y n (3) X 1 N(m X1, σ 2 X 1 ), X 2 N(m X2, σ 2 X 2 ) Y 1 = e X1 Y 2 = e X2 C Y1,Y 2 Y 1,2 m Y1,2 σ 2 Y 1,2 m Yi ( ) = exp m Xi + σ2 X i, i = 1, 2 (4) 2 σ 2 Y i = m 2 Y i {exp(σ 2 X i ) 1}, i = 1, 2 (5) 5
X 1 X 2 C X1,X 2 X 1 + X 2 N(m X1 + m X2, σ 2 X 1 + σ 2 X 2 + 2C X1,X 2 ) (6) C Y1,Y 2 = E[Y 1 Y 2 ] E[Y 1 ]E[Y 2 ] = exp{m X1 + m X2 + 1 2 (σ2 X1 + σx2 2 + 2C X1,X 2 )} m Y1 m Y2 = m Y1 m Y2 {exp(c X1,X 2 ) 1} (7) ρ(x 1, X 2 ) = C X1,X 2 /(σ X1 σ X2 )ρ(y 1, Y 2 ) = C Y1,Y 2 /(σ Y1 σ Y2 ) ρ(y 1, Y 2 ) = C Y 1,Y 2 = m Y 1 m Y2 (exp(c X1,X 2 ) 1) σ Y1 σ Y2 σ Y1 σ Y2 exp(c X1,X = 2 ) 1 = exp(σ X 1 σ X2 ρ(x 1, X 2 )) 1 exp(σx 2 1 ) 1 exp(σx 2 2 ) 1 exp(σx 2 1 ) 1 exp(σx 2 2 ) 1 (8) (3) AR(1) n y n x n = (y n µ n )/σ n (9) x n µ n σ n n x n N(0,1) AR(1) x n = a n,n 1 x n 1 + v n (10) a n,n 1 n n 1 x n σ 2 x 1 σ 2 v = (1 a n,n 1 2 )σ 2 x v n N(0, 1 a 2 n,n 1) (10) x n (9) y n a n,n 1 v n Thomas Fiering Thomas-Fiering n y n 3 { 1 f Y (y n ) = exp 1 ( ) } 2 ln(yn c n ) λ n (y n c n )ζ n 2π 2 ζ n y n (11) x n = ln(y n c n ) λ n ζ n (12) (11) (8) y n ρ y,n,n 1 x n ρ x,n,n 1 exp[ζ n 1 ζ n ρ x,n,n 1 ] 1 ρ y,n,n 1 = exp(ζn 1 2 ) 1 exp(ζn) 2 1 (13) 6
水資源工学 ( 立川担当分 4 回目 :1 月 18 日資料 ) 流況予測のためのモデル (2) 講義の目標 : 河川流量をシミュレーション発生させる手法は 前回の講義で示した確率 統計的な時系列モデルによる方法だけでなく 物理的なモデルを用いて降雨を発生させ それを流出モデルを介して河川流量に変換することが考えられる 地球温暖化によって将来の降水量が予測され 流出モデルや河川流モデルを介して河川流量の変化が分析されている 今日の講義では 温暖化時の流量変化予測を例として 河川流況のシミュレーション手法を理解することを目標とする 時系列モデルは 観測された時系列データに内在する確率 統計的な特性を分析して その特性を再現する時系列を模擬発生させることに主眼が置かれる 時系列データが生み出される物理的背景よりも 結果としての観測データを分析して 現象の確率的な変動を再現しようとする 一方で 自然現象を物理的な方程式で表現するモデルも時系列データを発生させるために用いられる たとえば 気候変動にともなう将来の気候推計データは 物理法則にしたがう支配方程式を解いて 将来の気温や降水量などの時系列データを生成し それを流出モデルを介して河川流量に変換して将来の流況が評価されている 1) 全球大気大循環モデル将来の気候を推計する大気大循環モデル (GCMs, General Circulation Models) の空間分解能が著しく高解像度化し その出力データはわが国の河川流域を対象とする水工シミュレーションモデルの入力データとしてそのまま利用できる解像度を有するようになってきている 気象庁気象研究所の全球大気モデルの空間分解能は約 20km であり アメダス観測網に匹敵する空間分解能となっている 全球を対象とする気候モデルとしては 現時点では世界最高の空間分解能を有する気候モデルである 将来気候の推計計算は 温室効果ガス排出量の変化シナリオに従って温室効果ガスの濃度の時間変化が設定され その条件のもとで将来の気候が予測される 気象庁気象研究所の全球 20km 格子大気モデルは A1B シナリオに従って温室効果気体の濃度変化が設定されており 以下の期間のデータが提供されている 1979 年 1 月 ~2003 年 12 月 : 現在気候実験 2015 年 1 月 ~2039 年 12 月 : 近未来気候実験 2075 年 1 月 ~2099 年 12 月 : 21 世紀末気候実験 図 1: 全球大気大循環気候モデルの模式図 ( 左 ) とそれが出力するデータ ( 右 ) 7
2) 流域モデルと流出モデル地球上の任意全域を対象とする 1km 空間分解能のグリッド型の分布型流出モデルを構成し, キネマティックウェーブモデルを用いて大気大循環モデルによる出力データを河川流量に変換する 標高データを用いて グリッドごとに周り 8 方向のうちの最急勾配方向を流水方向と定め それに従って一次元的に流れを追跡する 図 2( 左 ) は流水方向データを標高データから定める模式図であり 右図はこれに従って関東地方の流域モデルを示したものである 図 2: 標高データを用いた流下方向の決定と流域モデル ( 左 ) と関東地方の流域モデル ( 右 ) 3) 大気大循環モデルの出力を流出モデルを介して河川流量に変換した結果最上川 ( 砂越地点 6500km2) での計算流量の時系列データを以下に示す 左から順に現在気候実験 近未来気候実験 将来気候実験の結果であり 一つの図に 25 年間の河川流量を重ねて示している 4) 分析結果の例 図 3:100 年確率年最大流量の変化 ( 左 ) と基準渇水流量 (10 年確率渇水流量 ) の変化 ( 右 ) : 現在気候に対する 21 世紀末気候の比率を示している 参考文献 : 立川康人, 滝野晶平, 藤岡優子, 萬和明, キムスンミン, 椎葉充晴 : 気候変化が日本の河川流量に及ぼす影響の予測, 土木学会論文集, 67(1), pp. 1-15, 2011. 8