JP2017187371A - 気象予測装置及び風力発電所 - Google Patents

Abstract

【課題】数時間先の風速を高い精度で予測する。
【解決手段】本発明の気象予測装置１は、対象地点における過去の気象観測値、及び、対象地点以外の１又は複数の参照地点における過去の気象観測値が格納される記憶部１５と、対象地点及び参照地点の過去の気象観測値を入力変数とし、対象地点における将来の気象予測値を出力変数とする予測モデルに対して、対象地点及び参照地点の過去の気象観測値を使用して学習することにより、予測モデルのパラメータを最適化するモデル最適化部２２と、パラメータが最適化された予測モデルを使用し、対象地点における将来の気象予測値を算出する予測値算出部２３と、を備えることを特徴とする。
【選択図】図５

Description

本発明は、気象予測装置及び風力発電所に関する。

風力発電機は、枯渇の可能性がある化石エネルギーではなく、自然界に存在する再生可能な風力エネルギーを電力エネルギーに変換する。風力発電機は、地球温暖化の主因となる二酸化炭素を殆ど発生させないので、地球規模の環境問題を解決する手段として世界中で導入されつつある。日本においても、電力会社が有する既存の電力系統に対して、電力会社以外の事業者が有する風力発電機を連系する（接続する）試みが既に始まっている。

事業者と電力会社との間の取り決めの一形態として、電力会社は、事業者に対して“出力一定化制御”を要求する場合がある。例えば、ある電力会社は、風力発電機の事業者に対して以下を要求している。
・風力発電機とともに蓄電池を併設し、蓄電池出力を制御することによって、風力発電に起因する出力変動をほぼゼロとすること。
・風力発電機と蓄電池との合成出力の目標値からの偏差を、合成出力の最大値の２％以内とすること。

一般に、既存の電力系統の規模が小さい地域（人口過疎地など）においては、電力の需給をバランスさせるために既存の電力系統の供給力を機動的に調整すること（いわゆる“下げ代”の確保）が難しい。そのような地域において、前記の要求はより切実になる。このような要求を満足させるためには、蓄電池を増設すればよい。しかしながら、蓄電池の増設は大きな費用及びスペースを必要とする。そこで、風力発電機の発電量を正確に予測することが肝要となってくる。

特許文献１の発電量予測装置は、過去のデータにおける異なる時刻間の統計的相関に基づいて、又は、異なる発電機間の統計的相関に基づいて発電機の発電量を予測する。このとき、発電量予測装置は、分散共分散を作成する。
特許文献２の風車発電システムは、過去の風向に基づいて将来の風向を予測し、その風向に風車を向ける。このとき、風車発電システムは、埋め込み定理に基づき状態空間内にアトラクタ（風向の軌跡）を再構成し、再構成されたアトラクタに基づいて、決定論的非線形手法によって風向を予測する。

特開２０１１−２０００４０号公報 特開２００５−９８１８１号公報

特許文献１の発電量予測装置は、過去のデータがガウシアン分布に従うとの仮定のもとで乱数を発生させて風速を予測するが、予測の精度が低い。特許文献２の風車発電システムは、過去データの平均的な風向を予測値として出力するに過ぎず、やはり精度が低いうえ、そもそも風速を予測することを前提としていない。そして、特許文献１の発電量予測装置及び特許文献２の風車発電システムの両者ともに、数分程度の将来の風速・風向を予測する場合であれば精度が高くても、１０分を超える将来の風速・風速を予測する場合は、精度が大幅に低下する。

前記した出力一定化制御において電力が供給過剰となった場合、電力会社は、風力発電機を連係する事業者に対し、出力を抑制する旨の遠方指令を送信する実験が行なわれている。このようなレベルでの出力制御を行なうためには、数時間先の風速を正確に予測しなければならない。また、風力発電の供給過剰分に見合う火力発電等の下げ代の準備をするためにも、数時間先の風速を正確に予測しなければならない。予測の精度が向上すれば、蓄電池の放充電サイクル数も減少し、蓄電池容量も低減できる。さらに、風力発電の不足分を補うために火力等の化石エネルギーを待機させておく量が減少し、二酸化炭素の排出量は減少する。
そこで、本発明は、数時間先の風速を高い精度で予測することを目的とする。

本発明の気象予測装置は、対象地点における過去の気象観測値、及び、対象地点以外の１又は複数の参照地点における過去の気象観測値が格納される記憶部と、対象地点及び参照地点の過去の気象観測値を入力変数とし、対象地点における将来の気象予測値を出力変数とする予測モデルに対して、対象地点及び参照地点の過去の気象観測値を使用して学習することにより、予測モデルのパラメータを最適化するモデル最適化部と、パラメータが最適化された予測モデルを使用し、対象地点における将来の気象予測値を算出する予測値算出部と、を備えることを特徴とする。
その他の手段については、発明を実施するための形態のなかで説明する。

本発明によれば、数時間先の風速を高い精度で予測することができる。

気象データの構造を説明する図である。 パラメータの最適化を説明する図である。 地点ごとの重みを説明する図である。 時点ごとの重みを説明する図である。 気象予測装置の構成を説明する図である。 処理手順のフローチャートである。 風速の時系列の波形である。 風速及び風向の時系列の波形である。 各地点の位置関係を示す図である。 振舞が異なる風速の時系列の波形である。 ケースごとに予測誤差等を比較した表である。 (ａ)は、時系列変動が大きい場合の風速の波形である。(ｂ)は、時系列変動が小さい場合の風速の波形である。 (ａ)は、風速の時系列変動が大きい場合の最適化されたパラメータである。(ｂ)は、風速の時系列変動が小さい場合の最適化されたパラメータである。 （ａ）及び（ｂ）は、他地点カバー率のグラフである。 ケースごとに予測値と観測値とを比較した表である。

以降、本発明を実施するための形態（“本実施形態”という）を、図等を参照しながら詳細に説明する。
３つの地点を想定し、それぞれを地点Ａ、地点Ｂ及び地点Ｃとする。地点Ａには、風力発電機が配置されている。風力発電機は、自然の風力エネルギーを電力エネルギーに変換する設備であり、通常屋外に配置された風車と風車の回転力を利用する発電機とを有する。地点Ｂ及び地点Ｃのそれぞれには、パイロットプラントが配置されている。パイロットプラントとは、実用プラント（実際の風力発電所）を建設するに先立ち、その設計に必要なデータを収集するための試験的プラント（例えば鉄塔）である。なお、“対象地点”には地点Ａが相当し、“参照地点”には地点Ａ以外の地点が相当する。

（予測モデル）
本実施形態の予測モデル（式１）は、地点Ａ、地点Ｂ及び地点Ｃにおいて過去に観測された気象データ（観測値）を入力変数として、地点Ａにおける将来の気象データ（予測値）を出力変数とする数理モデルである。
Ｘ_ｔ＋１＝
ｗ_ｘ（ｓ_ｔａ_１Ｘ_ｔ＋ｓ_ｔ−１ａ_２Ｘ_ｔ−１＋ｓ_ｔ−２ａ_３Ｘ_ｔ−２＋・・・）
＋ｗ_ｙ（ｓ_ｔｂ_１Ｙ_ｔ＋ｓ_ｔ−１ｂ_２Ｙ_ｔ−１＋ｓ_ｔ−２ｂ_３Ｙ_ｔ−２＋・・・）
＋ｗ_ｚ（ｓ_ｔｃ_１Ｚ_ｔ＋ｓ_ｔ−１ｃ_２Ｚ_ｔ−１＋ｓ_ｔ−２ｃ_３Ｚ_ｔ−２＋・・・）
（式１）

式１において、ｔ＋１、ｔ、ｔ−１、ｔ−２、・・・は時点である。時点の数は特に限定されていない。時点又は時点の数を“期”と呼ぶ場合がある。
Ｘ_ｉ、Ｙ_ｉ及びＺ_ｉ（ｉ＜ｔ＋１）は、それぞれ、地点Ａ、地点Ｂ及び地点Ｃにおける気象データの観測値である。Ｘ_ｉ（ｉ＝ｔ＋１）は、地点Ａにおける気象データの予測値である。気象データとは、センサによって観測できる自然界のあらゆる物理量であり、典型的には、風速（ｍ／秒）である。
ａ_１、ａ_２、ａ_３、・・・、ｂ_１、ｂ_２、ｂ_３、・・・、ｃ_１、ｃ_２、ｃ_３、・・・は、パラメータである。詳細は後記するが、本実施形態の気象予測装置は、過去の観測値を使用して機械学習することによって、パラメータの値を最適化する。

ｗ_ｘ、ｗ_ｙ及びｗ_ｚは、地点ごとの重みである。０≦ｗ_ｘ≦１、０≦ｗ_ｙ≦１、０≦ｗ_ｚ≦１、及び、ｗ_ｘ＋ｗ_ｙ＋ｗ_ｚ＝１が成立するものとする。
ｓ_ｔ、ｓ_ｔ−１、ｓ_ｔ−２、・・・は、過去の時点ごとの重みである。０≦ｓ_ｔ≦１、０≦ｓ_ｔ−１≦１、０≦ｓ_ｔ−２≦１、・・・、及び、ｓ_ｔ＋ｓ_ｔ−１＋ｓ_ｔ−２＋・・・＝１が成立するものとする。
つまり、過去のすべての観測値に対して、３つの値（パラメータ、場所ごとの重み及び時点ごとの重み）が乗算される。

（ベクトル自己回帰モデル）
仮に式１の右辺が“ｗ_ｘ（ｓ_ｔａ_１Ｘ_ｔ＋ｓ_ｔ−１ａ_２Ｘ_ｔ−１＋ｓ_ｔ−２ａ_３Ｘ_ｔ−２＋・・・）”だけである場合、この予測モデルを自己回帰（ＡＲ：Auto-regressive）タイプと呼ぶ。ＡＲタイプの予測モデルに対しては、単一の入力変数が時系列で入力される。つまり、ＡＲタイプの予測モデルは、例えば、ある地点の過去の風速を使用して当該地点の将来の風速を予測するためのものである。これに対して、本実施形態の予測モデルは、多変数を入力変数とするベクトル自己回帰（ＶＡＲ：Vector Auto-regressive）タイプである。ＶＡＲタイプの予測モデルは、風速を含む多くの種類の気象データを入力でき、また、これらを地点ごとにも入力できる。

（パラメータ）
そのうち、パラメータａ_１、ａ_２、ａ_３、・・・、ｂ_１、ｂ_２、ｂ_３、・・・、ｃ_１、ｃ_２、ｃ_３、・・・は、機械学習の結果、様々な値に最適化される。例えば、地点Ｂにおける時点ｔの気象データ“Ｙ_ｔ”が、地点Ａにおける時点ｔ＋１の気象データ“Ｘ_ｔ＋１”に大きな影響を与える場合、パラメータ“ｂ_１”の値の絶対値が大きくなり、他のパラメータの値は“０”に近くなる。例えば、過去の時点の数が４つであり、地点の数が３つである予測モデルにおいては、パラメータは４×３＝１２個存在し、その１２個のそれぞれが機械学習の結果異なる値に最適化される。

（地点ごとの重み及び時点ごとの重み）
パラメータが機械学習によって最適化される対象であるのに対して、地点ごとの重み及び時点ごとの重みは機械学習とは直接関係がなく、ユーザ又は気象予測装置が設定する値である。時点の数がどのように変化しても、地点ごとの重み“ｗ”の値は、地点ごと固定される。地点の数がどのように変化しても、時点ごとの重み“ｓ”の値は、時点ごとに固定される。

（パラメータと重み）
重みもまたパラメータと同様、入力変数に対して乗算される“係数”であることには違いない。したがって、重みもまた機械学習の対象となり得る。しかしながら、本実施形態においては、重みを機械学習の対象とすることは想定していない。むしろ“地点ごと又は時点ごとに値が固定される”という縦横の拘束を受ける係数を“重み”とし、当該拘束を受けない係数を“パラメータ”として区別している。

前記した式１の予測モデルは、単純な線形のモデルであるが、予測モデルは、非線形を含むあらゆる形態を取ることが可能であり、ホワイトノイズのような定数項を含んでもかまわない。入力変数の地点の数が“３”でなくとも、“２”以上であればよい。要するに、予測モデルは、以下の条件をすべて満たしておればどのようなものでもよい。
・風力発電機の地点（対象地点）及びパイロットプラントの地点（参照地点）の過去の気象観測値を入力変数とし、対象地点における将来の気象予測値を出力変数とすること。
・入力変数のそれぞれに対するパラメータが機械学習の対象となり、入力変数のそれぞれに対して地点ごとかつ時点ごとに重みが定義できること。

（気象データ）
図１は、気象データの構造を説明する図である。地点Ａ、地点Ｂ及び地点Ｃのそれぞれにおいて、風速、風向、気温、相対湿度、気圧、日照時間及び降水量が観測される。このうち、風速は、２次元風速計によって計測される。つまり風速は、南北方向の風速成分と、東西方向の風速成分とが合成されるベクトルの長さ（ｍ／秒）である。風速の実際の時系列の波形は、１地点の場合図７に示したようになり、４地点の場合、図８の符号５２、５３、５４及び５５のようになる。

風向は、２次元風向計によって計測される。風向の実際の時系列の波形は、図８の符号５１）に示したようになる。通常、北に０度（３６０度）を割り当て、東に９０度を割り当て、南に１８０度を割り当て、西に２７０度を割り当てる。すると、風向（風上の方角）は、０度〜３６０度の範囲で定義される。但しこのままでは、風速及び風向を同じ時系列のグラフとして表現する場合、縦軸の“ｍ／秒”及び“度”が、しばしば重なる。そこで、風向（０度〜３６０度）を“３６”で除算し、さらに除算結果に例えば“２０”を加算する。すると、風向は、２０〜３０の範囲に正規化され重なりは殆ど解消される。図８では、風向の右目盛０度〜３６０度が、風速の左目盛２０〜３０（ｍ／秒）に対応している。図１の気温以下の気象データについては、詳しい説明を省略する。

図１の気象データの横軸は、時点である。気象予測装置のユーザは、時点の刻み幅を、１０分、１時間、１日等任意に設定することができる。いま、地点Ａ、地点Ｂ及び地点Ｃのそれぞれにおいて、時点７までの各種気象データが既に観測されているとする。そして、時点８における風速を予測する必要があるとする。このとき、気象予測装置は、予測モデルに対して、最大７×７×３＝１４７個の気象データの値（網掛けした“・・”）を入力し、気象データの値（例えば時点８の地点Ａの風速）を出力させる。もちろん、気象予測装置は、例えば、各地点の風速（７×３＝２１個）のみを入力してもよい。

（パラメータの機械学習）
図２は、パラメータの機械学習を説明する図である。説明の単純化のため、気象データとして代表的に“風速”のみを用いる予測モデルを想定する。図２では、地点Ａ、地点Ｂ及び地点Ｃのそれぞれにおいて、少なくとも時点１から時点１４までの風速が既に観測されているとする。気象予測装置は、次の処理を実行する。

（１）学習機会１において、すべての地点ごとの重み“ｗ”及びすべての時点ごとの重み“ｓ”に“１”を代入する。
（２）風速の観測値Ｘ_１、Ｘ_２、・・・、Ｘ_７、Ｙ_１、Ｙ_２、・・・、Ｙ_７、Ｚ_１、Ｚ_２、・・・、Ｚ_７を予測モデルに入力する。
（３）パラメータａ_１、ａ_２、・・・、ａ_７、ｂ_１、ｂ_２、・・・、ｂ_７、ｃ_１、ｃ_２、・・・、ｃ_７の値をランダムに発生させ、パラメータ値の組合せとする。
（４）パラメータ値の組合せを予測モデルに代入し、Ｘ_８を出力させる。
（５）予測値として出力されたＸ_８と観測値としてのＸ_８との差異を算出する。
（６）差異を最小にするパラメータ値の組合せを特定し、特定した組合せを暫定候補１とする。

（７）学習機会２において、前記（１）〜（６）を繰り返す。
ただし、繰り返し処理中の（２）においては、風速の観測値Ｘ_２、Ｘ_３、・・・、Ｘ_８、Ｙ_２、Ｙ_３、・・・、Ｙ_８、Ｚ_２、Ｚ_３、・・・、Ｚ_８を予測モデルに入力する。（４）においては、Ｘ_９を出力させる。（５）においては、予測値として出力されたＸ_９と観測値としてのＸ_９との差異を算出する。（６）においては、特定した組合せを暫定候補２とする。学習機会３〜７においても同様の処理を繰り返す。
（８）暫定候補１〜７のうち、差異が最も小さいものを最終候補とする。

（地点ごとの重み）
図３は、地点ごとの重みを説明する図である。いま、気象予測装置が、時点１〜時点７における風速その他の観測値を予測モデルに入力し、時点８における風速の予測値Ｘ_８を出力させる例を考える。気象予測装置は、地点Ａ以外の各地点の時点１〜時点７における観測値と、地点Ａの時点１〜時点７における観測値との類似度を算出する。

類似度を算出する最も簡単な例は、例えば、直近の時点７の風向同士の差分を類似度（値が小さいほど類似している）とする方法である。地点Ｂの風向と地点Ａの風向との差分は、２０−２０＝０である（北風で一致）。地点Ｃの風向と地点Ａの風向との差分は、２５−２０＝５である（南風と北風で相反）。このとき、気象予測装置は、地点Ｂの地点ごとの重み“ｗ_ｙ”を“０.８”とする。そして、地点Ａ及び地点Ｃの地点ごとの重み“ｗ_ｘ”及び“ｗ_ｚ”を両者とも“０.１”とする。地点Ｂは、地点Ａの風上になっている可能性が高い。この例では、ｗ_ｘ＝ｗ_ｚとしたが、必ずそうでなくてもよい。ｗ_ｘ＝ｗ_ｚでなくても、ｗ_ｙ＞ｗ_ｚであればよい。前記以外のより複雑な例として、気象予測装置は、風向以外の１又は複数の気象データ同士を比較してもよいし、時点１〜時点７における気象データの値の時系列変化（波形）同士を比較してもよい。

（時点ごとの重み）
図４は、時点ごとの重みを説明する図である。今度もまた、気象予測装置が、時点１〜時点７における風速その他の観測値を予測モデルに入力し、時点８における風速の予測値Ｘ_８を出力させる例を考える。気象予測装置は、地点Ａの時点７における観測値と、地点Ａの時点１〜時点６における各観測値との類似度を算出する。前記の風向同士の差分の例をそのまま利用すれば、時点４及び時点７において、差分（類似度）は“０”となり、その他の時点において、差分は“５”となる。このとき、気象予測装置は、時点４及び時点７の時点ごとの重み“ｓ_４”及び“ｓ_７”を両者とも“０.２５”とする。そして、他の時点の時点ごとの重みをすべて“０.１”とする。

この例では、ｓ_４＝ｓ_７としたが、必ずそうでなくてもよい。ｓ_４＝ｓ_７でなくても、ｓ_４＞ｓ_１＝ｓ_２＝ｓ_３＝ｓ_５＝ｓ_６であればよい。前記以外のより複雑な例として、気象予測装置は、時点ごとに１又は複数の気象データがなすベクトル同士を比較してもよい。さらに、風上の地点Ｑが特定できている場合は、気象予測装置は、地点Ｑの時点７における観測値と、地点Ｑの時点１〜時点６における観測値との類似度を算出してもよい。

（気象予測装置の構成）
図５に沿って、気象予測装置１の構成を説明する。気象予測装置１は、一般的なコンピュータであり、中央制御装置１１、キーボード等の入力装置１２、ディスプレイ等の出力装置１３、主記憶装置１４及び補助記憶装置１５を有する。これらの装置は、バスで相互に接続されている。補助記憶装置１５は、気象データ３１、パラメータ３２、地点ごとの重み３３、時点ごとの重み３４及び予測モデル３５を格納している。

主記憶装置１４における、観測値受信部２１、モデル最適化部２２及び予測値算出部２３は、プログラムである。以降、“○○部は”と主体を記した場合は、中央制御装置１１が、補助記憶装置１５から各プログラムを読み出し、主記憶装置１４にロードしたうえで、各プログラムの機能（詳細後記）を実現するものとする。
風力発電所ＷＦは、風力発電機ＷＧ及び気象予測装置１を有する。風力発電所ＷＦは、気象予測装置１を使用し気象予測装置１から将来の気象データの予測値を受け取る。風力発電所ＷＦは、予測された気象データに基づいて、風力発電機ＷＧを制御する。風力発電機ＷＧには蓄電設備（図示せず）が接続されている。

（処理手順）
図６に沿って、処理手順を説明する。ステップＳ１０１が開始される前提として、気象予測装置１の観測値受信部２１は、定期的に観測値を受信しているものとする。具体的には、観測値受信部２１は、風力発電所ＷＦ及びパイロットプラント（図示せず）からネットワーク（図示せず）等を介してリアルタイムで気象データを受信する。そして、受信の都度、観測値受信部２１は、受信した気象データを地点ごと、気象データの種類（風速、風向、・・・）ごとに時系列で記憶するものとする。

ステップＳ１０１において、気象予測装置１のモデル最適化部２２は、パラメータを更新する必要があるか否かを判断する。具体的には、モデル最適化部２２は、前回ステップＳ１０１の処理を実行した時点以降において少なくとも１度観測値受信部２１が観測値を受信したか否かを判断する。少なくとも１度受信した場合（ステップＳ１０１“ＹＥＳ”）、ステップＳ１０３に進み、それ以外の場合（ステップＳ１０１“ＮＯ”）、ステップＳ１０２に進む。

ステップＳ１０２において、モデル最適化部２２は、ユーザから開始指示を受け付けたか否かを判断する。具体的には、モデル最適化部２２は、ユーザが入力装置１２を介して“予測処理を開始する”旨の指示を入力するのを受け付けた場合（ステップＳ１０２“ＹＥＳ”）、ステップＳ１０３に進む。それ以外の場合は、そのまま待機する。

ステップＳ１０３において、モデル最適化部２２は、パラメータを最適化する。具体的には、第１に、モデル最適化部２２は、補助記憶装置１５から気象データ３１を読み出す。
第２に、モデル最適化部２２は、前記した図２の方法で、パラメータの値の組合せ（前記の最終候補）を特定し、特定した組合せをパラメータ３２（図５）として補助記憶装置１５に記憶する。

ステップＳ１０４において、モデル最適化部２２は、地点ごとの重みを決定する。具体的には、モデル最適化部２２は、前記した図３の方法で、地点ごとの重みを算出し、算出した地点ごとの重みを、地点に関連付けたうえで、地点ごとの重み３３（図５）として補助記憶装置１５に記憶する。なお、補助記憶装置１５は、前記の類似度に関連付けて地点ごとの重みの値が記憶されたテーブルを記憶しているものとする。なお、モデル最適化部２２は、出力装置１３に気象データ３１（図１）を表示したうえで、ユーザが地点ごとの重みを、入力装置１２を介して入力するのを受け付けてもよい。

ステップＳ１０５において、モデル最適化部２２は、時点ごとの重みを決定する。具体的には、モデル最適化部２２は、前記した図４の方法で、時点ごとの重みを算出し、算出した時点ごとの重みを、時点に関連付けたうえで、時点ごとの重み３４（図５）として補助記憶装置１５に記憶する。なお、補助記憶装置１５は、前記の類似度に関連付けて時点ごとの重みの値が記憶されたテーブルを記憶しているものとする。なお、モデル最適化部２２は、出力装置１３に気象データ３１（図１）を表示したうえで、ユーザが時点ごとの重みを、入力装置１２を介して入力するのを受け付けてもよい。

ステップＳ１０６において、気象予測装置１の予測値算出部２３は、予測値を算出する。具体的には、第１に、予測値算出部２３は、予測モデルのパラメータに対してステップＳ１０３の“第２”において最適化したパラメータ値を代入する。
第２に、予測値算出部２３は、予測モデルの地点ごとの重みに対してステップＳ１０４において決定した地点ごとの重みの値を代入する。
第３に、予測値算出部２３は、予測モデルの時点ごとの重みに対してステップＳ１０５において決定した時点ごとの重みの値を代入する。

第４に、予測値算出部２３は、予測モデルに気象データを入力する。ここで入力する気象データは、例えば、地点Ａ、地点Ｂ及び地点Ｃの直近の過去７期における風速である。
第５に、予測値算出部２３は、予測モデルに予測値を出力させる。ここで出力させるデータは、例えば、地点Ａにおける将来（次の期）の風速の予測値である。
第６に、予測値算出部２３は、“第５”において出力させた予測値を出力装置１３に表示する。
第７に、予測値算出部２３は、発電量換算式に対して風速の予測値を適用することによって発電量を算出し、出力装置１３に表示する。なお、補助記憶装置１５は、風速と発電量との関係を定義する発電量換算式を記憶しているものとする。
その後処理手順を終了する。

以降、検定の手法と結果を説明する。検定の目的は、主として他地点のデータを使用することの有効性を確認し、予測に向かない振舞（詳細後記）をするデータの影響を確認することである。

（検定１）
検定１では、ＶＡＲタイプの予測モデルを使用した気象予測装置の予測精度を見て行く。地点として、４つの地点Ａ、地点Ｂ、地点Ｃ及び地点Ｄを採用した。このうち、地点Ａには風力発電機が配置されている。地点Ｂ、地点Ｃ及び地点Ｄには、パイロットプラントが配置されている。これらの実際の位置関係は図９の通りである。そして、このうち地点Ａ、地点Ｃ及び地点Ｄにおける風速の時系列の波形は、図１０の通りである。なお、図１０の風速の目盛は、風向も勘案した正負値になっている。

図１０において、枠６１で囲んだ期間に注目すると、地点Ａの波形と地点Ｄの波形とは、よく似ている。それに対して、地点Ａの波形と地点Ｃの波形とは、全く異なる。気象分野において風の動きを一般的に“振舞”と呼ぶことがある。つまり、地点Ａと地点Ｃとでは、枠６１の期間中、振舞が明らかに異なっている。そこで、枠６１の期間中において予測誤差（予測値と観測値との差分）を地点ごとに比較した。

枠６１の期間の長さは５００分である。５００分の期間において、風速は１０分ごとに５０期分観測されている。検定においては、地点Ａにおける３期先（３０分先）の風速を予測することとした。そして、３期先の風速を予測するために、５０期分の風速のうち直前の連続する２０期分を使用する（予測モデルの入力変数とする）こととした。つまり、検定１においては、第１〜第２０期分を使用して予測した場合、第２〜第２１期分を使用して予測した場合、・・・、第３１〜第５０期分を使用して予測した場合というように、予測誤差が計３１回分算出されることになる。予測誤差は、平均誤差（３１個の予測誤差の平均値）及び標準偏差（３１個の予測誤差の標準偏差）の２つの指標で評価した。

比較対象となるケースは、以下の６つである（図１１も参照）。
（ケ―ス１）ＡＲタイプの予測モデルを使用して、地点Ａの風速に基づき地点Ａの風速を予測する。
（ケース２）ＶＡＲタイプの予測モデルを使用して、地点Ａ及び地点Ｃの風速に基づき地点Ａの風速を予測する。
（ケース３）ＶＡＲタイプの予測モデルを使用して、地点Ａ及び地点Ｄの風速に基づき地点Ａの風速を予測する。
（ケース４）ＶＡＲタイプの予測モデルを使用して、地点Ａ、地点Ｃ及び地点Ｄの風速に基づき地点Ａの風速を予測する。
（ケース５）ＶＡＲタイプの予測モデルを使用して、地点Ａ、地点Ｃ及び地点Ｄの風速に基づき地点Ａの風速を予測する。このとき、図２に示したように観測値を更新し、時点ごとの重みを使用した。
（ケース６）ＶＡＲタイプの予測モデルを使用して、地点Ａ、地点Ｃ及び地点Ｄの風速に基づき地点Ａの風速を予測する。このとき、地点ごとの重みを使用した。

図１１は、検定１の結果である。これを見ると以下のことがわかる。
・ケース１とケース２との比較
平均誤差は、ケース１の方が小さい。標準誤差も、ケース１の方が小さい。ＶＡＲタイプの予測モデルの予測誤差がより大きいということになり、一見意外ではあるが、これは、地点Ｃの振舞が、地点Ａの振舞とは明らかに異なることに起因すると思われる。
・ケース１とケース３との比較
平均誤差は、ケース３の方が小さい。標準偏差は、ケース１の方が小さい。この結果だけからは何とも言えない。

・ケース１とケース４との比較
平均誤差は、ケース４の方が小さい。標準誤差も、ケース４の方が小さい。したがって、ＡＲタイプの予測モデルを使用する場合に比して、ＶＡＲタイプの予測モデルを使用する場合のほうが予測誤差は小さいことがわかる。地点Ａとは振舞の異なる地点Ｃの影響が、振舞が類似する地点Ｄによって治癒されている。このように、予測に向かないデータが存在していても、予測に向くデータが存在していれば、パラメータが自動的に最適化されて、予測に向くデータ（地点Ｄ）の影響が大きくなると考えられる。つまり、地点Ａと振舞が類似する地点を特定できることが重要である。

・ケース４とケース５との比較
標準偏差は、ケース５の方が小さい。したがって、観測値を更新したうえでパラメータを最適化すること、及び、時点ごとの重みを使用すると予測誤差が小さくなることがわかる。つまり、観測値を取得する時点を増加させることが重要である。
・ケース４とケース６との比較
標準偏差は、ケース６の方が小さい。したがって、地点ごとの重みを使用すると予測誤差が小さくなることがわかる。なお、“＜１.５７”は、その値が１.５７を下回ることを意味する。

以上を総合すると、ＶＡＲタイプの予測モデル及び複数地点のデータを使用すると予測誤差が小さくなること、複数地点のなかに地点Ａとは振舞の異なる地点が含まれていてもかまわないことがわかる。さらに、観測値を更新し、地点ごとの重み及び時点ごとの重みを使用すると予測誤差が小さくなることがわかる。

（検定２）
検定２では、実際に最適化されたパラメータを見て行く。地点Ａにおいて風速の時系列変動が大きい場合、風速の波形は図１２（ａ）に示されたようになる。地点Ａにおいて風速の時系列変動が小さい場合、風速の波形は図１２（ｂ）に示されたようになる。気象予測装置１が、ＶＡＲタイプの予測モデル並びに地点Ａ、地点Ｂ、地点Ｃ及び地点Ｄの風速のデータを用いて予測モデルのパラメータを最適化させた結果が図１３（ａ）及び図１３（ｂ）である。予測モデルは、地点Ａ、地点Ｂ、地点Ｃ及び地点Ｄの過去の風速の観測値（直前の連続する２０期分）を入力とし、地点Ａの将来の風速の予測値を出力する。

図１３（ａ）は、地点Ａにおいて風速の時系列変動が大きい場合の結果（最適化されたパラメータ値）であり。図１３（ｂ）は、地点Ａにおいて風速の時系列変動が小さい場合の結果である。
・図１３（ａ）のパラメータ７１ａが“０.４２６７”であるのに対して、図１３（ｂ）のパラメータ７１ｂは、“０.８３４８”である。このことから、地点Ａにおける風速の変動が小さい場合に、地点Ａの予測値に対する地点Ａの直近の観測値による影響がより大きくなることがわかる。“０.４２６７”及び“０.８３４８”は、他と比して際立って高い値であり、この例では、仮にＡＲタイプの予測モデルが使用されたとしても予測精度は充分高かったと思われる。

図１３（ａ）において、４期前の地点Ｂのパラメータ７２は“０.１４９９”となっており、この値は比較的大きい。この時点において、地点Ｂから地点Ａに向かって風が吹いていた（地点Ｂが地点Ａの風上であった）ことが考えられる。

（検定３）
検定３では、ある地点の風上の他の地点の風速が当該地点の風速に与える影響、及び、当該地点の風上ではない他の地点の風速が当該地点に与える影響を見て行く。検定２と同様に、ＶＡＲタイプの予測モデルに対して４地点（地点Ａ、地点Ｂ、地点Ｃ及び地点Ｄ）の風速の観測値を２０期分入力する例を考える。

次に記すのは、“他地点カバー率算出アルゴリズム”（単に“アルゴリズム”ともいう）である。
for(tin tstart:tend){
Count=0
if(max(｜a₁｜,…, ｜a_n｜)＜α((max(｜b₁｜,…, ｜b_n｜)
+ max(｜c₁｜,…, ｜c_n｜) + max(｜d₁｜,…, ｜d_n｜))){
Count= Count + 1
}
}

ここで、ａ_１等は地点Ａの各時点のパラメータであり、ｂ_１等は地点Ｂの各時点のパラメータであり、ｃ_１等は地点Ｃの各時点のパラメータであり、ｄ_１等は地点Ｄの各時点のパラメータである。“ｍａｘ（ ）”は、（ ）内の複数の値のうちの最大値を取得する関数である。当該関数は、過去２０期のパラメータのうち、絶対値が最大であるものを、地点ごとに取得する。

“ｔｓｔａｒｔ”は、最初の期であり、“ｔｅｎｄ”は、最後の期である。“ｔｅｎｄ−ｔｓｔａｒｔ”の値は、“２０”に比して充分に大きくなっている。いま仮に１００期分の観測値があるとすると、連続する２０期分の観測値が８１通り存在することになる。したがって、アルゴリズムは、８１回繰り返されることになる。この繰り返し処理の各回において、ｉｆ文の中の不等式が満たされれば“Ｃｏｕｎｔ”の値（初期値は“０”である）が“１”ずつインクリメントされ、ｉｆ文の中の不等式が満たされなければ“Ｃｏｕｎｔ”の値は変化しない。

いま、倍率係数αに“１.００”を代入し、アルゴリズムの繰り返し処理が終了した時点での“Ｃｏｕｎｔ”の値を取得する。次に、倍率係数αを“０.０１”だけ増加させ、“１.０１”とする。そして再度、倍率係数αに“１.０１”を代入し、アルゴリズムの繰り返し処理が終了した時点での“Ｃｏｕｎｔ”の値を取得する。このような処理を、倍率係数αを増加させながら繰り返す。例えば“Ｃｏｕｎｔ”の値が“３５”であるとする。このとき、“３５／８１×１００＝４３.２％”を他地点カバー率と呼ぶ。

他地点カバー率は、αの値（１００個ある）のそれぞれについて算出される。換言すれば、他地点カバー率は、“地点Ａ以外の地点のパラメータの絶対値の最大値の和のα倍が、地点Ａのパラメータの絶対値の最大値を超えるようなケースが、すべてのケースのうちの何％を占めるか”を示している。つまり他地点カバー率が高いほど、地点Ａ自身の単独の観測値よりも、地点Ａ以外の複数の地点の観測値の方が、地点Ａの風速の予測に対して大きく寄与していることを意味している。

図１４（ａ）及び図１４（ｂ）は、他地点カバー率算出アルゴリズムを使用した結果である。図１４（ａ）及び図１４（ｂ）の横軸は倍率係数であり、縦軸は他地点カバー率である。図１４（ａ）は、地点Ａ及び地点Ｂの観測値を使用した結果であり、アルゴリズムのｉｆ文の中の不等式の右辺には、地点Ｃ及び地点Ｄについての項がない。図１４（ｂ）は、地点Ａ、地点Ｂ、地点Ｃ及び地点Ｄの観測値を使用した結果であり、アルゴリズムのｉｆ文の中の不等式は、前記した通りのものである。

図１４（ａ）に注目する。他地点カバー率のグラフ（連続ドット）８１は、地点Ｂが地点Ａの風上にあり、かつ、地点Ａの風速の変動が大きい場合のグラフである。グラフ８２は、地点Ｂが地点Ａの風上になく、かつ、地点Ａの風速の変動が小さい場合のグラフである。グラフ８２は、他地点カバー率がほぼ“０”％付近で推移している。つまり、地点Ａの風速は、地点Ｂの影響を殆ど受けていない。風力発電機がある地点Ａの風速に対し、地点Ａの風上にない地点Ｂは殆ど影響を与えない。グラフ８１は、緩やかな右上りの形状を示し、地点Ａの風速に対し、地点Ａの風上の地点Ｂは相当の影響を与えている。

図１４（ｂ）に注目する。グラフ８３は、地点Ｂが地点Ａの風上にあり、かつ、地点Ａの風速の変動が大きい場合のグラフである。グラフ８４は、地点Ｂが地点Ａの風上になく、かつ、地点Ａの風速の変動が小さい場合のグラフである。一見して明らかなように、図１４（ｂ）の他地点カバー率は、図１４（ａ）の他地点カバー率に比して飛躍的に向上している。特に倍率係数αが“１”である場合であっても、他地点カバー率は相当に高い。このことは、地点Ｂが地点Ａの風上である場合でも風上ではない場合でも、地点Ａの風速の予測に対し、地点Ｃ及び地点Ｄは相当の影響を与えていることを意味する。
図１４（ａ）及び図１４（ｂ）のいずれもが、過去において地点Ａの風上になった地点の重みの値を大きくすることに意味があることを示している。

（検定４）
検定４では、重みを使用することによる効果を再度見て行く。既存の風力発電所の風速を予測する検定を行なった。当該風力発電所は、我が国の日本海側の地点Ｅに立地し、定格出力が２ＭＷである風力発電機を７基有し、合計最大出力は、１４ＭＷである。その結果は、図１５に示す通りである。比較対象となるケースは、以下の４つである（図１５も参照）。

（ケース１１）風速の観測値である。
（ケース１２）ＶＡＲタイプの予測モデルを使用して、地点Ｂ、地点Ｃ、地点Ｄ及び地点Ｅの風速に基づき地点Ｅの３０分先の風速を予測する。但し、場所ごとの重みも時点ごとの重みも使用しない。
（ケース１３）ＶＡＲタイプの予測モデルを使用して、地点Ｂ、地点Ｃ、地点Ｄ及び地点Ｅの風速に基づき地点Ｅの３０分先の風速を予測する。但し、場所ごとの重み及び時点ごとの重みを使用している。
（ケース１４）気象庁の“全球数値予報モデル ＧＳＭ（登録商標）”を使用して地点Ｅの風速を予測する。同モデルは、http://www.jmbsc.or.jp/hp/online/f-online0a.htmlから配信されている。

網掛けの数値は、ケース１２、ケース１３及びケース１４の予測値のうち、ケース１１の観測値との差分が最も小さいものである。さらに予測精度が±２０％である出現頻度割合も比較した。この出現頻割合は、予測値と観測値との差分が観測値の−２０％〜＋２０％の範囲に収まる時点が、全体の何％を占めるかを示している。ケース１２、ケース１３及びケース１４の出現頻度割合のうちで、ケース１３の出現頻度割合が“９３％”と最大であった。以上から、ケース１３（重み使用）が最も高精度で風速を予測できることがわかった。地点ごとの重み及び時点ごとの重みの値をさらに見直すことで、ケース１３の結果は、さらに高精度なものになり得ると思われる。
なお、図にはないが、パラメータをリアルタイムに更新していくと、そうでない場合に比して、予測精度が５％以上向上することが確認されている。

（重み付けの工夫）
以下のことが経験的に分っている。
・我が国の日本海側においては、過去において西風（西から吹く風）が吹いた時点の重みの値を大きくすると予測の精度が向上する。
・我が国の太平洋側においては、過去において東風が吹いた時点の重みの値を大きくすると予測の精度が向上する。
・過去において風速及び風向の両者が類似している時点の重みの値を大きくすると予測の精度が向上する。
・海岸沿いの地点の風力を予測する場合は、過去において海側から風力発電機のある地点に向かって風が吹いた時点の重みの値を大きくし、海側の地点の重みの値を大きくすると予測の精度が向上する。このことは、パイロットプラントが風力発電機の風上にあり、かつ、風力発電機との間に障害物が少ないことが理想であることを意味する。

（最適化の拡張）
前記では、気象予測装置１が最適化する対象はパラメータである例を説明した。しかしながら、気象予測装置１は、地点ごとの重み、時点ごとの重み及びパラメータを含むすべての“係数”をランダムに又は所定の規則に則って発生させ、過去の観測値を使用し機械学習を行うことによって、すべての係数を最適化してもよい。すると、例えば風上の地点の地点ごとの重みの値が大きくなり、振舞が類似する過去の時点の時点ごとの重みの値が大きくなる。

（実施形態の効果）
本実施形態の気象予測装置は、以下の効果を奏する。
（１）ある地点の風速等を、それ以外の地点の風速等を使用して高い精度で予測することができる。
（２）リアルタイムで観測値が更新されて行く場合、風速等をより高い精度で予測することができる。
（３）風力発電機の発電量を高い精度で予測することができる。
（４）高い精度での予測に相応しい地点及び時点を特定することができる。
（５）風速、風向、気圧等の入手しやすいデータに基づいて、高い精度での予測に相応しい地点及び時点を自動的に特定することができる
（６）風力発電機とともに蓄電設備を備える事業者が予測結果を使用することができる。
（７）特に巨大な経営規模で運営される風力発電所が予測結果を使用することができる。

なお、本発明は前記した実施例に限定されるものではなく、様々な変形例が含まれる。例えば、前記した実施例は、本発明を分かり易く説明するために詳細に説明したものであり、必ずしも説明したすべての構成を備えるものに限定されるものではない。また、ある実施例の構成の一部を他の実施例の構成に置き換えることが可能であり、また、ある実施例の構成に他の実施例の構成を加えることも可能である。また、各実施例の構成の一部について、他の構成の追加・削除・置換をすることが可能である。

また、前記の各構成、機能、処理部、処理手段等は、それらの一部又は全部を、例えば集積回路で設計する等によりハードウエアで実現してもよい。また、前記の各構成、機能等は、プロセッサがそれぞれの機能を実現するプログラムを解釈し、実行することによりソフトウエアで実現してもよい。各機能を実現するプログラム、テーブル、ファイル等の情報は、メモリや、ハードディスク、ＳＳＤ（Solid State Drive）等の記録装置、又は、ＩＣカード、ＳＤカード、ＤＶＤ等の記録媒体に置くことができる。
また、制御線や情報線は説明上必要と考えられるものを示しており、製品上必ずしもすべての制御線や情報線を示しているとは限らない。実際には殆どすべての構成が相互に接続されていると考えてもよい。

１ 気象予測装置
１１ 中央制御装置
１２ 入力装置
１３ 出力装置
１４ 主記憶装置（記憶部）
１５ 補助記憶装置（記憶部）
２１ 観測値受信部
２２ モデル最適化部
２３ 予測値算出部
３１ 気象データ
３２ パラメータ
３３ 地点ごとの重み
３４ 時点ごとの重み
３５ 予測モデル
ＷＦ 風力発電所

Claims (7)

1. 対象地点における過去の気象観測値、及び、前記対象地点以外の１又は複数の参照地点における過去の気象観測値が格納される記憶部と、
   前記対象地点及び前記参照地点の前記過去の気象観測値を入力変数とし、前記対象地点における将来の気象予測値を出力変数とする予測モデルに対して、前記対象地点及び前記参照地点の前記過去の気象観測値を使用して学習することにより、前記予測モデルのパラメータを最適化するモデル最適化部と、
   前記パラメータが最適化された予測モデルを使用し、前記対象地点における将来の気象予測値を算出する予測値算出部と、
   を備えることを特徴とする気象予測装置。
2. 前記気象予測装置は、
   前記対象地点及び前記参照地点の気象観測値をリアルタイムで受信する観測値受信部を備え、
   前記モデル最適化部は、
   前記対象地点及び前記参照地点の気象観測値がリアルタイムで受信される都度、前記予測モデルのパラメータを更新すること、
   を特徴とする請求項１に記載の気象予測装置。
3. 前記対象地点は、
   風力発電機が設置される地点であり、
   前記気象観測値及び前記気象予測値は、
   風速であり、
   前記予測値算出部は、
   所定の換算ルールを使用し、前記算出した気象予測値を前記風力発電機の発電量に換算すること、
   を特徴とする請求項２に記載の気象予測装置。
4. 前記予測モデルは、
   前記対象地点及び前記参照地点ごとに又は過去の時点ごとに値が変化する重みを、前記入力変数に対して乗算するものであること、
   を特徴とする請求項３に記載の気象予測装置。
5. 前記モデル最適化部は、
   風速、風向、気圧、降水量、気温、相対湿度及び日照時間のうちの少なくとも１つに応じて、前記重みの値を算出すること、
   を特徴とする請求項４に記載の気象予測装置。
6. 前記風力発電機は、
   蓄電設備を備えていること、
   を特徴とする請求項３ないし請求項５のいずれか１項に記載の気象予測装置。
7. 請求項１ないし請求項６のいずれか１項に記載の気象予想装置を使用する風力発電所。

Images