(57)【特許請求の範囲】
【請求項１】
第１気象変数の複数の予測値を含む複数のアンサンブルと、第２気象変数の観測値を含む複数の観測データとに基づき、前記第１気象変数の予測値と前記第２気象変数の観測値間の同時確率密度分布を生成する第１分布生成部と、
前記同時確率密度分布に基づき、予測対象時刻における前記第１気象変数の複数の予測値が与えられた場合の前記第２気象変数の観測値の確率密度分布である第１確率密度分布を生成し、複数の前記第１確率密度分布を合成することにより、第２確率密度分布を生成する、第２分布生成部と、
を備えた気象予測装置。
【請求項２】
前記第１分布生成部は、前記第１気象変数の複数の予測値と、前記第２気象変数の観測値との複数の組を生成し、前記組ごとに第３確率密度分布を対応付け、前記組に対応づけた第３確率密度分布を足し合わせることにより、前記同時確率密度分布を生成する
請求項１に記載の気象予測装置。
【請求項３】
前記第１分布生成部は、前記第１気象変数の前記複数の予測値の代表値と前記第２気象変数の観測値との組を生成し、前記組ごとに第３確率密度関数を対応付け、前記組に対応づけた第３確率密度関数を足し合わせることにより、前記同時確率密度分布を生成する
請求項１に記載の気象予測装置。
【請求項４】
前記代表値は、前記第１気象変数の前記複数の予測値の平均値である
請求項３に記載の気象予測装置。
【請求項５】
前記代表値は、前記第１気象変数の前記複数の予測値のうち、前記第２気象変数の観測値との差分が最も小さい前記予測値である
請求項３に記載の気象予測装置。
【請求項６】
前記第３確率密度関数は、カーネル密度関数である
請求項３ないし５のいずれか一項に記載の気象予測装置。
【請求項７】
前記第２分布生成部は、前記予測対象時刻における前記第１気象変数の複数の予測値を適用した前記同時確率密度分布を、前記予測対象時刻における前記第１気象変数の複数の予測値を適用した前記同時確率密度分布を前記第２気象変数が取り得る観測値の範囲で積分したもので除算することにより、前記第２確率密度分布を生成する、
請求項１ないし６のいずれか一項に記載の気象予測装置。
【請求項８】
前記第１気象変数は、前記第２気象変数と同じ気象変数である
請求項１ないし７のいずれか一項に記載の気象予測装置。
【請求項９】
前記第１気象変数は、前記第２気象変数と異なる気象変数である
請求項１ないし７のいずれか一項に記載の気象予測装置。
【請求項１０】
前記第１分布生成部は、第３気象変数の複数の予測値を含む複数のアンサンブルをさらに用いて、前記第１気象変数の予測値と前記第３気象変数の予測値と前記第２気象変数の観測値間の同時確率密度分布を生成し、
前記第２分布生成部は、前記同時確率密度分布に基づき、前記予測対象時刻における前記第１気象変数の複数の予測値と前記第３気象変数の複数の予測値との組毎に、前記組が与えられた場合の前記第２気象変数の観測値の確率密度分布である前記第１確率密度分布を生成する
請求項１ないし９のいずれか一項に記載の気象予測装置。
【請求項１１】
第１気象変数の複数の予測値を含む複数のアンサンブルと、第２気象変数の観測値を含む複数の観測データとに基づき、前記第１気象変数の予測値と前記第２気象変数の観測値間の同時確率密度分布を生成するステップと、
前記同時確率密度分布に基づき、予測対象時刻における前記第１気象変数の複数の予測値が与えられた場合の前記第２気象変数の観測値の確率密度分布である第１確率密度分布を生成し、複数の前記第１確率密度分布を合成することにより、第２確率密度分布を生成するステップと
をコンピューターが実行する気象予測方法。
【請求項１２】
第１気象変数の複数の予測値を含む複数のアンサンブルと、第２気象変数の観測値を含む複数の観測データとに基づき、前記第１気象変数の予測値と前記第２気象変数の観測値間の同時確率密度分布を生成するステップと、
前記同時確率密度分布に基づき、予測対象時刻における前記第１気象変数の複数の予測値が与えられた場合の前記第２気象変数の観測値の確率密度分布である第１確率密度分布を生成し、複数の前記第１確率密度分布を合成することにより、第２確率密度分布を生成するステップと
をコンピューターに実行させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明の実施形態は、確率的気象予測装置、確率的気象予測方法およびプログラムに関する。
【背景技術】
【０００２】
気象予測は、未来の気象を予測するものである。気象予測は必ず誤差を伴うため、気象予測に基づき意思決定を行う場面では、多くの場合、それがどの程度信頼できるかを考慮する。
【０００３】
予測の分散性を表現する技術として、統計分布のパラメトリックな表現を利用するものがある。しかし、複雑な気象条件をごく限られたパラメーターの組み合わせで表現するのは難しいと考えられている。今日、予測の不確定性を評価する標準的な手法は、アンサンブル予測と呼ばれる手法である。
【０００４】
しかしながら、アンサンブル予測で考慮できる予測の不確定性は限られている。現実的には数値計算の初期値の誤差のみが考慮されることが多い。このため、不確定性を低く評価してしまう（分散性が不足する）場合が多いという問題がある。
【先行技術文献】
【特許文献】
【０００５】
【特許文献１】特開第２０１４−２１５５５号公報
【非特許文献】
【０００６】
【非特許文献１】T. N. Palmer. "The economic value of ensemble forecasts as a tool for risk assessment:from days to decades". Q. J. R. Meteorol. Soc., Vol. 128, No. 581, pp. 747-774, april 2002.
【発明の概要】
【発明が解決しようとする課題】
【０００７】
本発明の実施形態は、気象予測をその不確定性を適正に評価しつつ行うことを可能にしようとするものである。
【課題を解決するための手段】
【０００８】
本実施形態としての確率的気象予測装置は第１分布生成部と、第２分布生成部とを備える。第１分布生成部は、第１気象変数の複数の予測値を含む複数のアンサンブルと、第２気象変数の観測値を含む複数の観測データとに基づき、前記第１気象変数の予測値と前記第２気象変数の観測値間の同時確率密度分布を生成する。第２分布生成部は、前記同時確率密度分布に基づき、予測対象時刻における前記第１気象変数の複数の予測値のそれぞれが与えられた場合の前記第２気象変数の観測値の確率密度分布である第１確率密度分布を生成し、複数の前記第１確率密度分布を合成することにより、第２確率密度分布を生成する。
【図面の簡単な説明】
【０００９】
【図１】第１実施形態に係る確率的気象予測装置のブロック図。
【図２】観測履歴記憶部のデータ構造例を示す図。
【図３】アンサンブル記憶部のデータ構造例を示す図。
【図４】対象位置およびその周辺位置の配置例と、対象位置および周辺位置の予測値の分布の例を示す図。
【図５】同時確率密度分布の推定方法を説明するための図。
【図６】同時確率密度分布の推定方法を説明するための図。
【図７】第１実施形態の動作を説明するための図。
【図８】第１実施形態の動作を説明するための図。
【図９】同時確率密度分布生成フェーズの動作のフローチャートを示す図。
【図１０】確率的予測フェーズの動作のフローチャートを示す図。
【図１１】第１実施形態に係る確率的気象予測装置を備えたシステムを示す図。
【図１２】第１実施形態に係る確率的気象予測装置のハードウェア構成の一例を示す図。
【図１３】第２実施形態の動作を説明するための図。
【図１４】第３実施形態の動作を説明するための図。
【発明を実施するための形態】
【００１０】
以下、図面を参照しながら、本発明の実施形態について説明する。
【００１１】
＜第１実施形態＞
図１に本発明の実施形態に係る確率的気象予測装置を示す。この確率的気象予測装置は、予測対象となる気象変数について、予測対象時刻における当該気象変数の観測値の確率密度分布を推定するものである。本実施形態では、気象変数として、太陽光発電にとって重要な日射強度または当該日射強度に基づく指標（快晴指数）を例にして説明する。ただし、予測する気象変数は特定の気象変数に限定されず、他の気象変数の場合も同様にして行うことができる。
【００１２】
図１の確率的気象予測装置は、アンサンブル記憶部１０１と、観測履歴記憶部１０２と、本実施形態に係る第１分布生成部である同時確率密度分布生成部１０３と、本実施形態に係る第２分布生成部である条件付確率密度分布生成部１０４を備える。条件付確率密度分布生成部１０４は、分布生成部１０5と、合成部１０６とを備える。
【００１３】
観測履歴記憶部１０２は、センサーにより測定された気象変数の観測値の観測履歴を記憶している。観測履歴は、観測値、観測位置および観測時刻を含む観測データを複数含む。気象変数が日射強度の場合、センサーとして日射計を用いることができる。観測位置ごとにセンサーが配置され、観測値が測定される。各観測位置で測定された観測値を観測位置および観測時刻の情報とともに外部の装置から受信し、観測履歴記憶部１０２に格納する。
【００１４】
図２に観測履歴記憶部１０２のデータ構造の例を示す。観測位置は、緯度と経度によって表されている。一分間隔で日射強度が観測位置（緯度、経度）とともに記憶されている。図２では、１つの位置（３５．５、１３９．５）のみ示されているが、実際には、複数の位置が記憶されている。
【００１５】
観測値は、例えば一定時間毎に各観測位置に設置されたセンサーで測定され、本装置に観測位置および観測時刻とともに送信されて、観測履歴記憶部１０２に格納されてもよい。センサーの観測値を収集する収集装置が存在し、当該収集装置からセンサー毎の観測値を観測位置および観測時刻とともに受信して、観測履歴記憶部１０２に格納してもよい。
【００１６】
アンサンブル記憶部１０１は、アンサンブル予測データを記憶している。アンサンブル予測データは、時刻と、位置と、気象変数のＭ（Ｍは２以上の整数）個の予測値とのデータを複数含む。Ｍ個の予測値の組をアンサンブルと呼ぶ。また、Ｍ個の予測値の各々を、アンサンブルのメンバーと呼ぶ。アンサンブル予測データは、後述する数値気象予測（数値気象シミュレーション）により得られる。アンサンブル予測データに含まれる時刻および位置、すなわち、数値気象シミュレーションで算出する時刻および位置を、観測データに含まれる観測位置および観測時刻と区別して、予測時刻および予測位置と呼ぶこともある。
【００１７】
図４にアンサンブル記憶部１０１のデータ構造の例を示す。時刻と、位置と、気象変数のＭ個の予測値（Ｍ個のメンバー）を含むアンサンブルと、が記録されている。当該気象変数は、日射強度であり、本実施形態に係る第１気象変数に対応する位置は、緯度と経度によって表されている。前述した観測履歴記憶部１０２に記憶される観測履歴の気象変数は、本実施形態に係る第２気象変数に対応する。
【００１８】
アンサンブルにおける各メンバーの意味は、数値気象シミュレーションの方法によって異なってよい。一例として、少しずつ条件を変えて並行してシミュレーションを行う方法の場合、各メンバーの値は、シミュレーションにおいて条件毎に計算される、同じ予測位置の日射強度を意味する。すなわち、この場合、メンバーのすべての予測位置は一致している。別の例として、数値気象シミュレーションでは一つしか条件を設定せず、メンバーの中のひとつはこのシミュレーションでの予測位置の日射強度を表し、他のメンバーは、シミュレーションの格子点のうち、予測位置の周辺の位置の日射強度を表す（詳細は後述するが、周辺位置の予測値が予測位置の予測値に一致または近い場合があることを考慮したものである）。本実施形態では、これらのいずれの場合も可能である。いずれにしても、アンサンブルの個々のメンバーは、予測時刻および予測位置での気象変数の予測値を表している。
【００１９】
図３では、１つの位置（３５．０、１４０．０）のみ示されているが、実際には、複数の位置が記録されている。また、図３では、気象変数として日射強度しか示されていないが、その他の気象変数が記録されていてもよい。例えば、大気温度、圧力、風速、風向などの気象変数が記録されていて良い。また、これらの大気温度、圧力、風速、風向などの気象変数の場合、位置を表す値として、緯度および経度に加えて、高度を含めてもよい。すなわち、これらの大気変数の場合、緯度、経度および高度を含む、三次元的なデータが記録される。
【００２０】
通常、気象予測は、主要な国の機関（日本の気象庁、ヨーロッパのＥＣＭＷＦなど）が定期的に発行する地球規模の予測がベースになっている。これらの機関の出す予測の時間スケジュールで定期的にアンサンブルを取得して、アンサンブル記憶部１０１に格納してもよい。
【００２１】
なお、アンサンブル記憶部１０１に記録される予測時刻と予測位置は、観測履歴記憶部１０２に記録されている観測時刻と観測位置と必ずしも一致している必要はない。一致していない場合、アンサンブル記憶部１０１のデータを、観測履歴記憶部１０２の観測時刻と観測位置に対して、適切に補間して用いればよい。以下では、アンサンブル記憶部１０１のデータを適切に補間したものを用いる、もしくは両者が一致している場合を想定して説明を行う。
【００２２】
ここで、数値気象予測（数値気象シミュレーション）、およびアンサンブルの詳細について説明する。数値気象予測は、大気の運動や時間発展の仕組み（ダイナミクス）をコンピューターで数値計算するものである。通常、このダイナミクスは決定論的なもので確率性は組み込まれていない。しかしながら、系内に不確定な要素がある場合は、確率性を表現することが重要になる。数値気象予測の場合、初期値の不確定性、およびダイナミクスの不確定性といった不確定性が存在する。
【００２３】
数値気象予測の場合、大気の状態（各位置の気圧、風速・風向、気温、水蒸気量など）のある初期値から出発して、その後の時間発展を計算する。しかしながら、この初期値を求めるための観測データは、観測位置および観測対象の種類ともに限られている。また観測誤差も発生し得る。このため、各位置の大気の状態の初期値を一意に決定できない。観測値をベースにして、統計的に、各位置の大気の状態の初期値を推定（各位置での気象変数の値を推定）することになる。このため、各位置における初期値を確率的に表現する必要がある。
【００２４】
また気象のダイナミクスは、流体力学および熱力学などの物理現象をベースに組み立てられており、今日、現象の理解は進んでいる。しかし、まだ理解されていない部分も多くある。とくに、数値計算の必要上、空間を離散的な格子点に切ることで、連続的な状態を表現することになるため、ある種の近似を余儀なくされる場合が多い。このためダイナミクスにも不確定性があり、ダイナミクスを確率的に表現する必要がある。
【００２５】
確率性を表現するために用いられるものとしてアンサンブルと呼ばれるものがある。ひとつの気象数値計算の結果は、ひとつの初期値、ひとつの決められたダイナミクスに基づく時間発展（これをパスと呼ぶことにする）を表現するものでしかない。このため、確率性を表現するために複数のパスを作成することが行われる。このような複数のパスの集合をアンサンブルと呼ぶ。
【００２６】
初期値の不確定性をあらわすために、乱数などを用いて初期値に摂動を加えて複数の初期値を生成し、それぞれに対するパスを造ることを行う。このようにして作られたアンサンブルを、初期値アンサンブルと呼ぶ。また、ダイナミクスの不確定性をあらわすためには、数値気象予測のモデルの中で、ダイナミクスの性質を現すパラメーターに対して、乱数などを用いて摂動を加えて、複数の数値気象モデルを作り、それぞれに対応するパスを作る。これをモデルアンサンブルと呼ぶ。
【００２７】
上記のように数値気象予測において、複数のパスを作成する方法は、代表的なアンサンブル作成法である。アンサンブルは、広く解釈すれば、予測対象（気象変数）の各時刻における複数の予測値の集まり、ということもできる。このように考えると、以下のようにアンサンブルを生成することも可能である。
【００２８】
数値気象予測の誤差として、位置ズレ誤差がよく発生することが知られている。この位置ズレ誤差とは、気象パターンの予測はおおむね正しいものの、細かく見ると位置ずれが生じていることである。数値気象予測に発生するこの誤差の性質から考えて、予測対象位置の周辺のいずれかの位置の予測が、対象位置の気象を正確に表している可能性がある。このことから考えて、対象位置の周辺の予測値を集めることでアンサンブルを作ることができる。図４（Ａ）に、対象位置およびその周辺位置の例を示し、図４（Ｂ）に対象位置および周辺位置の予測値の分布の例を示す。対象位置およびその周辺が広域に同様の気象（すなわち予測値が同じまたは近い）であれば、分布はその気象を中心とした鋭いピークを持つようになり、その対象位置の予測値が示す気象が実現する可能性が高いといえる。逆に対象位置およびその周辺の気象が空間的に大きく変動していれば、分布の幅も広くなり、対象位置に対して与えられた予測の不確定性が高いといえる。このように、対象位置の周辺の予測値から作られるアンサンブルは理にかなっている。具体的には、予測対象称点を中心として、あらかじめ決められた数Ｎがあって、Ｎ×Ｎ個の数値気象シミュレーションの格子点の予測値を取り、これらをまとめてアンサンブルとすることができる。
【００２９】
本実施形態では、図３を用いて説明したように、これらの対象位置および周辺位置（複数の格子点）の予測値をまとめたアンサンブルでもよいし、前述した初期値アンサンブルおよびモデルアンサンブルでもよい。要するに、予測時刻および予測時刻に対して、気象変数の複数の予測値をあつめたものであればよい。
【００３０】
同時確率密度分布生成部（第１分布生成部）１０３は、アンサンブル記憶部１０１に記憶されたアンサンブル予測データと、観測履歴記憶部１０２に記憶された観測履歴に基づき、気象変数の予測値と、気象変数の観測値間の同時確率密度分布を生成する。本実施形態では、予測値の気象変数（第１気象変数）と、観測値の気象変数（第２気象変数）は同じであるとするが、後述する他の実施形態のように、これらの気象変数が異なってもよい。以下、同時確率密度分布の生成について詳細に説明する。
【００３１】
同時確率密度分布とは、二つの事象Ａ、Ｂが同時に発生する確率分布のことである。ここでは、同時確率密度分布は、予測値が示す日射強度（予測日射強度）と、観測値が示す日射強度（観測日射強度）の同時確率分布を指す。以下では、日射強度の変わりに、次式で与えられる快晴指数を考える。
【数１】
0006552928.tif 000002
【００３２】
快晴指数は、０以上１以下（すなわち[0,1]）の範囲の値をとる。観測日射強度と、予測日射強度を区別するために、観測日射強度をR、予測日射強度をrと表す。また、これら二つの強度の同時確率密度分布をp(r,R)と表す。
【００３３】
以下、カーネル密度推定を用いた同時確率密度分布の推定手法を説明する。図５はその方法を説明するための図である。ある観測点で過去のある時刻tに観測値として、快晴指数（観測快晴指数）
【数２】
0006552928.tif 000003
が観測され、その観測点でその時刻tにおける数値気象予測により得られる快晴指数の予測値（予測快晴指数）を
【数３】
0006552928.tif 000004
と記載する。過去に行われた予測の予測値、過去に行われた観測の観測値を今後の予測のために利用する時、上に“〜”が付いた文字で表す。観測快晴指数
【数４】
0006552928.tif 000005
と予測快晴指数
【数５】
0006552928.tif 000006
を、図５（Ａ）に示すように、予測快晴指数（予測値）を横軸に、観測快晴指数（観測値）を縦軸にとった平面上の点
【数６】
0006552928.tif 000007
で表すことができる。
【００３４】
カーネル密度推定では、この平面上の点に、図５（Ｂ）に示すように、カーネル密度関数を対応付ける。カーネル密度関数は、式（１）のような多重正規分布を用いるのが一般的である。カーネル密度関数は、本実施形態に係る第３確率密度関数に対応する。第３確率密度関数として、カーネル密度関数以外の関数または確率密度分布を用いてもかまわない。
【数７】
0006552928.tif 000008
ここで
【数８】
0006552928.tif 000009
である。
また
【数９】
0006552928.tif 000010
である。
【数１０】
0006552928.tif 000011
はそれぞれ予測快晴指数軸、観測快晴指数軸方向のカーネル密度関数の広がりの大きさを表し、それぞれ標準偏差に対応する。|Σ|はΣの判別式を表す。
【数１１】
0006552928.tif 000012
である。
【００３５】
ある予測対象位置について、ある時刻における複数のメンバーの予測値（予測快晴指数）のそれぞれを、当該時刻の観測値（観測快晴指数）と組み合わせて （予測快晴指数、観測快晴指数）の組を複数生成する。これを複数の時刻について行う。メンバー毎に、複数の時刻分の組が得られる。これらを
【数１２】
0006552928.tif 000013
と表すことができる。Ｎを時刻の個数を表すとすると、Ｎ＝ｔ１−ｔ０＋１である。Ｍはアンサンブルのメンバー数を表す。予測値（予測快晴指数）と観測値（観測快晴指数）とを軸とする座標系に、各組を配置した例を図６（Ａ）に示す。実際にはより多数の組が配置されるが、ここでは表記の簡潔のため概略的に示している。図６（Ａ）の各位置に、カーネル密度関数を対応づけた様子を図６（Ｂ）に示す。メンバー毎の複数の時刻分の組のそれぞれにカーネル密度関数を対応付け、互いに足し合わせることにより、以下の式（２）のように、同時確率密度関数を求めることができる。
【数１３】
0006552928.tif 000014
【００３６】
上述した同時確率密度分布の生成処理のイメージを図７の左に示す。図７の左では、ある時刻における複数の予測快晴指数（メンバー１〜Ｍ）を、
【数１４】
0006552928.tif 000015
とした場合において、予測快晴指数のそれぞれを、観測快晴指数（観測値）と組み合わせて、複数の組を生成する例が示される（なお図では“(t)”の表記は省略している）。
【００３７】
条件付確率密度分布生成部（第２分布生成部）１０４は、同時確率密度分布生成部１０３が生成した同時確率密度分布と、将来のある時刻（予測対象時刻）における複数の予測快晴指数に基づき、当該複数の予測快晴指数が与えられた場合の予測対象時刻における観測快晴指数の確率密度分布（条件付確率密度分布）を生成する。条件付確率密度分布生成部１０４は、分布生成部１０５と、合成部１０６とを備えている。
【００３８】
分布生成部１０５は、上記同時確率密度分布に基づき、予測対象時刻におけるアンサンブルのメンバー毎の予測快晴指数が与えられた場合の、観測快晴指数の条件付確率密度分布を求める。当該メンバー毎の条件付確率密度分布は、本実施形態に係る第１確率密度分布に対応する。予測対象時刻の日射量は、アンサンブル記憶部に格納されており、したがって、予測対象時刻の各メンバーの快晴指数は、予測対象時刻の各メンバーの日射強度から算出すればよい。
【００３９】
ここで、上記予測対象時刻における快晴指数の観測値の条件付確率密度分布は、条件付確率密度分布を与える以下の一般的な式（３）を用いて求めることができる。式（３）において、Rは快晴指数であり、Rは[0,1]の範囲の値をとるため、式（３）の分母におけるRの積分範囲は[0,1]である。
【数１５】
0006552928.tif 000016
【００４０】
式（３）に従って、予測対象時刻におけるメンバー毎の予測快晴指数を与えることで、メンバー毎に条件付確率密度分布を算出する。この処理のイメージを、図７の右、および図８に示す。すなわち、予測対象時刻tにおける予測値
【数１６】
0006552928.tif 000017
が与えられた場合、式（３）に当該予測値を代入することで、観測快晴指数の条件付確率密度分布
【数１７】
0006552928.tif 000018
が計算できる（なお図７および図８では“(t)”の表記は省略している）。なお、図８の左下の観測値（観測快晴指数）と予測値（予測快晴指数）間の分布は、図６（Ｂ）に対応している。
【００４１】
合成部１０６は、予測対象時刻tにおけるアンサンブルの全メンバーの条件付確率密度分布を合成することにより、アンサンブルの全メンバー１〜Ｍの予測値（予測快晴指数）が与えられた場合の観測快晴指数の条件付確率密度分布を算出する。当該全メンバーの予測値が与えられた場合の条件付確率密度分布は、本実施形態に係る第２確率密度分布に対応する。具体的に、アンサンブルの全メンバー１〜Ｍの予測快晴指数が与えられた場合の観測快晴指数の条件付確率密度分布を、アンサンブルの各メンバーの実現確率が等しく、それぞれ1/M であるとして、以下の式（４）により計算する。
【数１８】
0006552928.tif 000019
【００４２】
ここでは、各メンバーの重みを同じにしたが、メンバーに応じて重みを変えることも可能である。図８、および図７の右下に、各メンバーの観測快晴指数の条件付確率密度分布を合成する様子が示されている（なお図では“(t)”の表記は省略している）。これにより、複数の予測快晴指数が与えられた場合の予測対象時刻における観測快晴指数の条件付確率密度分布が得られる。
【００４３】
このようにして求めた予測対象時刻における観測快晴指数の条件付確率密度分布は、種々の場合に応用可能である。例えば、確率的な気象予測が必要な例として、災害に対する対策がある。災害が予想される場合、重要な資源を保護する、あるいは避難させるなどの対策が考えられる。しかし、この対策には一定の費用（コスト）がかかるのが普通である。したがって、災害が起こった場合の被害の大きさ、対策のコスト、および災害の発生する確率を考慮して、対策を実施するか否かの判断を下す必要がある。太陽光発電の場合だと、日射量の不足が予想される場合には、あらかじめ火力発電を立ち上げておくなどの対策が考えられる。しかし、日射量不足の予測が大きく外れた場合、火力発電立ち上げコストは無駄になる。このため、対策を採ることが経済的に見合うことなのか、確率的な試算により判断する必要がある。このような場合に日射量またはこれに基づく指標（観測快晴指数）の確率的な予測が必要になる。このように確率的な予測が与えられた場合に、それに基づいた対応をとることができる。
【００４４】
以下、図９および図１０のフローチャートを用いて、本実施形態の動作を説明する。本動作は大きく同時確率密度分布生成フェーズと、確率的予測フェーズとを備える。気象予測は翌日の計画を立てるために用いられることが多いため、典型例として、前日のしかるべき時間（例えば前日の１２：００）に、一日分の予測を作るために、同時確率密度分布生成フェーズと確率的予測フェーズを連続して行う場合を示す。ただし、それぞれのフェーズは独立に実行することができるため、必ずしも連続して実行される必要はない。
【００４５】
図９に、同時確率密度分布生成フェーズの動作のフローチャートを示す。
【００４６】
ステップＳ１０１では、同時確率密度分布生成部１０３が、観測履歴記憶部１０２から同時確率密度分布の生成に必要となる観測データを取得する。本例では、ある日における気象変数の観測値の確率密度分布の予測のために、その前Ｌ日分の期間のデータを用いるものとする。Ｌの値が小さいと同時確率密度の精度が低くなり、一方、Ｌの値が大きすぎると、季節変化などの影響で不適切なデータを用いることになる可能性がある。したがって、Ｌの値は適切に調整することが望ましい。典型的には一月分（Ｌ＝３０）である。
【００４７】
ステップＳ１０２では、同時確率密度分布生成部１０３が、アンサンブル記憶部１０１から同時確率密度作成に必要となるアンサンブル予測データを取得する。取得するアンサンブル予測データの期間の範囲は、一例として、観測データと同じにする。
【００４８】
ステップＳ１０３では、同時確率密度分布生成部１０３が、同時刻の予測値rと観測値Rの組ρ=（r,R）を生成する。より詳細には、同時刻の予測値が、アンサンブルのメンバーの数M存在し、これらを{r1,r2,…,rM}とする。全てのメンバー１〜Ｍの予測値と、観測値Rとの組を以下のように生成する。これらの組を上記Ｌ日分の各時刻について生成する。
【数１９】
0006552928.tif 000020
【００４９】
ステップＳ１０４では、同時確率密度分布生成部１０３が、ステップＳ１０３で取得した予測値と観測値の組を用いて、カーネル密度推定により、同時確率密度分布p(r,R)を作成する（式（２）参照）。
【００５０】
図１０に、確率的予測フェーズの動作のフローチャートを示す。
【００５１】
ステップＳ２０１では、条件付確率密度分布生成部１０４の分布生成部１０５が、アンサンブル記憶部１０１から、将来の時刻（予測対象時刻）tのアンサンブルの予測値
【数２０】
0006552928.tif 000021
を取り出す。
【００５２】
ステップＳ２０２では、分布生成部１０５が、同時確率密度分布生成部１０３が生成した同時確率密度分布に基づき、前述した式（３）に従って、メンバー毎に条件付確率密度分布を計算する。すなわち、予測対象時刻におけるアンサンブルの全メンバーの予測値
【数２１】
0006552928.tif 000022
のそれぞれに対して、当該同時確率密度分布から条件付確率密度分布を得る。すなわち、
【数２２】
0006552928.tif 000023
(m=1,2,…,M)を得る。
【００５３】
ステップＳ２０３では、条件付確率密度分布生成部１０４の合成部１０６が、前述した式（４）に基づいて、各メンバーの条件付確率分布を合成する。これにより予測対象時刻tにおける観測値Rの条件付確率密度分布を得る。
【００５４】
図１１は、本実施形態に係る確率的気象予測装置を備えたシステムの全体構成を示す。このシステムは、図１の確率的気象予測装置と、アンサンブル生成装置２０１と、観測装置２０２と、応用装置２０３と、出力装置２０４とを備える。
【００５５】
アンサンブル生成装置２０１は、図３に示したようなアンサンブル予測データを生成し、確率的気象予測装置に送信する。確率的気象予測装置は、アンサンブル生成装置２０１から受信したアンサンブル予測データを、アンサンブル記憶部１０１に格納する。アンサンブル予測データは、気象庁等の公的機関、またはその他の民間の機関が作成したデータでもよく、このときアンサンブル生成装置２０１は当該機関またはその委託先の機関の装置でもよい。アンサンブル予測データは一定時間間隔または任意のタイミングでアンサンブル生成装置２０１から受信してもよい。アンサンブル記憶部１０１に記憶されるアンサンブル予測データは、現在時刻に対して、過去の時刻のデータのみならず、将来の時刻のデータも含まれる。
【００５６】
観測装置２０２は、複数の観測値に配置されたセンサーを用いて、気象変数の値を測定し、観測値、観測位置および観測時刻を含む観測データを、確率的気象予測装置に送信する。確率的気象予測装置は、観測装置２０２から受信した観測データを、観測履歴記憶部１０２に格納する。観測データは一定時間間隔または任意のタイミングで観測装置２０２から受信してもよい。
【００５７】
応用装置２０３は、合成部１０６により生成された条件付確率密度分布を受信し、受信した分布に基づき、種々の判断処理を行う。例えば、当該分布に基づき、日射量の不足の可能性を計算し、火力発電立ち上げコストの観点から、あらかじめ火力発電を立ち上げておくか否かの判断を行ってもよい。応用装置２０３で行う処理を何でもよく、特定のものに限定されない。
【００５８】
出力装置２０４は、確率的気象予測装置の合成部１０６により生成された条件付確率密度分布を受信し、受信した分布を画面に表示する。ユーザは、表示装置に表示された条件付確率密度分布を参照することで、各種の判断を行うことができる。応用装置２０３の判断の結果を出力装置２０４に送信して、その結果を画面に表示してもよい。ユーザは、その結果を参照して、各種の判断を行ってもよい。
【００５９】
図１２は、本実施形態に係る確率的気象予測装置のハードウェア構成の一例を示す。図１２の確率的気象予測装置は、プロセッサ４０１、メモリ４０２、ネットワークインタフェース４０３、外部インタフェース４０４、およびストレージ４０６を備え、これらがバス４０５を介して接続されている。
【００６０】
ネットワークインタフェース４０３は、有線または無線のネットワークに接続され、当該ネットワークを介して、外部の装置と通信する。ここではネットワークインタフェースを１つのみ示しているが、複数のネットワークインタフェースが搭載されていてもよい。ネットワークインタフェース４０３は、アンサンブル生成装置２０１、観測装置２０２、応用装置２０３および出力装置２０４のうちの少なくとも１つと、ネットワークを介して通信してもよい。
【００６１】
メモリ４０２は、プロセッサ４０１が実行する命令、およびプロセッサ４０１が利用する各種データ等を一時的に記憶する。メモリ４０２は、ＳＲＡＭ、ＤＲＡＭ等の揮発性メモリでも、ＮＡＮＤ、ＭＲＡＭ等の不揮発性メモリでもよい。ストレージ４０６は、ハードディスク、光ディスク、ＳＳＤ等のデータを永続的に記憶する記憶装置である。メモリ４０２およびストレージ４０６の少なくとも一方により、図１のアンサンブル記憶部１０１および観測履歴記憶部１０２が実現される。
【００６２】
プロセッサ４０１は、メモリ４０２からプログラムを読み出して、実行することで、図１の同時確率密度分布生成部１０３、条件付確率密度分布生成部１０４、分布生成部１０５、および合成部１０６の機能が実現される。
【００６３】
外部インタフェース４０４は、外部装置に接続するためのインタフェースである。ここでは外部インタフェースを１つのみ示しているが、複数の外部インタフェースが搭載されていてもよい。外部インタフェース４０４に、アンサンブル生成装置２０１、観測装置２０２，応用装置２０３および出力装置２０４のうちの少なくとも１つが接続されてもよい。
【００６４】
以上、本実施形態によれば、複数の時刻のそれぞれにおける複数の予測値と観測値とに基づき、予測値と観測値間の同時確率密度分布を生成する。この分布と、予測対象時刻の複数の予測値とに基づき、予測対象時刻における観測値の条件付確率密度分布を生成する。この条件付確率密度分布を利用することで、アンサンブル予測に発生しやすい分散性の不足を補って、予測の不確定性を正しく評価できる。すなわち、気象予測の不確定性を適正に反映した予測が可能となる。
【００６５】
＜第２実施形態＞
第２実施形態は、同時確率密度分布生成と、条件付確率密度分布の推定および合成の動作が第１実施形態と異なる。図１３は、第２実施形態に係る同時確率密度分布生成と、条件付確率密度分布の推定および合成の動作説明図である。
【００６６】
同時確率密度分布生成部１０３は、アンサンブルにおける予測値の組{r1,r2,…,rM}に対して、これらの予測値の代表値として、以下の式（５）により、平均（平均予測値と呼ぶ）を計算する。
【数２３】
0006552928.tif 000024
【００６７】
同時確率密度分布生成部１０３は、複数の時刻のそれぞれについて、平均予測値と観測値との組
【数２４】
0006552928.tif 000025
を生成する。そして、第１実施形態と同様に、カーネル密度推定を用いて、同時確率密度分布p(r,R)を求める。
【００６８】
条件付確率密度分布生成部１０４の分布生成部１０５は、同時確率密度分布生成部１０３により生成された同時確率密度分布p (r,R)に基づき、式（３）に従って、予測対象時刻tのアンサンブルの全てのメンバー１〜Ｍの予測値
【数２５】
0006552928.tif 000026
に対して、条件付確率密度分布
【数２６】
0006552928.tif 000027
(m=1,2,…,M)を計算する（なお図１３では“(t)”の表記は省略している）。
【００６９】
合成部１０６は、分布生成部１０５により生成された予測対象時刻tにおけるメンバー毎の条件付確率密度分布を合成する（（式（４）参照））。これにより、予測対象時刻における観測値Rの条件付確率密度分布を得る。
【００７０】
＜第３実施形態＞
第３実施形態は、同時確率密度分布生成と、条件付確率密度分布の推定および合成の動作が第１および第２実施形態と異なる。図１４は、第３実施形態に係る同時確率密度分布生成と、条件付確率密度分布の推定および合成の動作説明図である。
【００７１】
同時確率密度分布生成部１０３は、アンサンブルにおける予測値の組{r1,r2,…,rM}に対して、これらの予測値の代表値として、観測値Rに最も近い予測値を取り出す。これを最良メンバー(Best Member)と呼び、それを
【数２７】
0006552928.tif 000028
と表す。
【００７２】
同時確率密度分布生成部１０３は、複数の時刻のそれぞれについて、最良メンバーと観測値との組
【数２８】
0006552928.tif 000029
を生成する。そして、第１実施形態と同様に、カーネル密度推定を用いて、同時確率密度分布p(r,R)を求める。
【００７３】
条件付確率密度分布生成部１０４の分布生成部１０５は、同時確率密度分布生成部１０３により生成された同時確率密度分布p (r,R)に基づき、式（３）に従って、予測対象時刻tにおける全てのメンバー１〜Ｍの予測値
【数２９】
0006552928.tif 000030
に対して、条件付確率密度分布
【数３０】
0006552928.tif 000031
(m=1,2,…,M)を計算する（なお図１１では“(t)”の表記は省略している）。
【００７４】
合成部１０６は、分布生成部１０５により生成された予測対象時刻tにおけるメンバー毎の条件付確率密度分布を合成する（式（４）参照）。これにより、予測対象時刻における観測値Rの条件付確率密度分布を得る。
【００７５】
第３実施形態では、複数の予測値のうち観測値Rに最も近い予測値を、複数の予測値の代表値として選択し、第２実施形態では複数の予測値の平均値を、当該代表値として選択したが、代表値はこれらの値に限定されない。例えば最頻値、中央値、など、別の値でもよい。
【００７６】
＜第４実施形態＞
第１〜第３実施形態では、同じ気象変数（日射量または快晴指数）に関する予測値と観測値とに基づき同時確率密度分布を生成し、そこからメンバー毎に条件付確率密度分布を推定していた。本実施形態では、予測値と観測値間で異なる気象変数を用いる例（例１）と、２つ以上の気象変数を用いる例（例２）を示す。
【００７７】
（例１）同時確率密度分布を生成するのに、必ずしも観測値と予測値とで同じ気象変数を利用する必要は無い。例えば、日射強度を予測したいとする。この場合、観測値は日射強度であるが、予測値は必ずしも、日射強度である必要はない。数値気象予測の場合、多様な気象変数が利用可能であり、これらから適切なものを利用できる。例えば、日射強度を観測値、気温を予測値として用い、気温と日射強度と組にして同時確率密度分布を生成することも可能である。
【００７８】
（例２）同時確率密度分布を構成する気象変数は２つ以上でも可能である。一例として、快晴指数rと気温Tを気象変数（それぞれ第１気象変数、第３気象変数に対応）として用い、快晴指数を予測したいとする。快晴指数の観測値をRとする。同時確率密度分布生成部１０３は、時刻tにおけるアンサンブルの気象変数の予測値の組
【数３１】
0006552928.tif 000032
を取り出す。また、この時の快晴指数の観測値を
【数３２】
0006552928.tif 000033
とすると、次のように予測値および観測値の組を作る。
【数３３】
0006552928.tif 000034
【００７９】
ここでρ=(r,T,R)に対して、以下のようなカーネル関数を定義する。
【数３４】
0006552928.tif 000035
ここで
【数３５】
0006552928.tif 000036
である。同時確率密度分布は式（２）と同様に、
【数３６】
0006552928.tif 000037
として計算できる。Ｎは、過去の時刻t0から、時刻t1まで範囲の時刻数（データ数）を表す。
【００８０】
将来の時刻（予測対象時刻）tの日射強度を予測したい場合、その予測対象時刻tにおけるアンサンブルの予測値を
【数３７】
0006552928.tif 000038
とすると
【数３８】
0006552928.tif 000039
により、予測対象時刻tにおける観測値Rの条件付確率密度分布を計算できる。
【００８１】
なお、本発明は上記実施形態そのままに限定されるものではなく、実施段階ではその要旨を逸脱しない範囲で構成要素を変形して具体化できる。また、上記実施形態に開示されている複数の構成要素の適宜な組み合わせにより、種々の発明を形成できる。例えば、実施形態に示される全構成要素から幾つかの構成要素を削除してもよい。さらに、異なる実施形態にわたる構成要素を適宜組み合わせてもよい。
【符号の説明】
【００８２】
１０１：アンサンブル記憶部
１０２：観測履歴記憶部
１０３：同時確率密度分布生成部
１０４：条件付確率密度分布生成部
１０５：分布生成部
１０６：合成部
２０１：アンサンブル生成装置
２０２：観測装置
２０３：応用装置
２０４：出力装置
【図１】
0006552928.tif 000040
【図２】
0006552928.tif 000041
【図３】
0006552928.tif 000042
【図４】
0006552928.tif 000043
【図５】
0006552928.tif 000044
【図６】
0006552928.tif 000045
【図７】
0006552928.tif 000046
【図８】
0006552928.tif 000047
【図９】
0006552928.tif 000048
【図１０】
0006552928.tif 000049
【図１１】
0006552928.tif 000050
【図１２】
0006552928.tif 000051
【図１３】
0006552928.tif 000052
【図１４】
0006552928.tif 000053