(57)【要約】【課題】調査ほ場数をできる限り抑制しながら、空撮画像と農作物の特定の生育ステージにおける時系列気象データを用いて、収量を高精度に推定する。【解決手段】過去に蓄積された調査地域の画像やほ場の属性情報と、調査対象であるほ場が撮影された空撮画像から、調査すべき実測ほ場を選定するための判断基準となるパラメータ群を決定し、前記パラメータ群ができる限り分散を持つように、実測ほ場を選定し、また調査の負担をできる限り減らすために、実測ほ場の候補ができる限り位置的に集中するように選定する。また気象データの時系列パターンを各生育ステージごとに解析することにより、生育状況に相関のあるパラメータ群を算出し、画像特徴量、ほ場の属性情報と前記パラメータ群を説明変数とした収量推計を実施する。【選択図】図１
2015000049.tif 000002
【特許請求の範囲】
【請求項１】
ほ場を含む空撮画像を受信できる受信部と、
当該ほ場を含む複数のほ場の形状情報と、位置情報とをそれぞれ記憶する記憶部と、
受信した前記ほ場を含む空撮画像と、当該ほ場の形状情報と位置情報とから、当該ほ場の画像特徴量を算出する画像解析部と、
算出された前記画像特徴量から前記ほ場の予測収量を算出する収量推計部とを有し、
前記収量推計部は、
予め記憶された気象データの時系列パターンを基に、前記ほ場で栽培される作物の特定の生育ステージごとに、前記作物の生育状況と相関を有する第一パラメータ群を出力する時系列パターン解析部を備え、
前記画像特徴量と前記第一パラメータ群を用いて前記ほ場の予測収量を算出することを特徴とする収量予測システム。
【請求項２】
請求項１に記載の収量予測システムにおいて、
前記パラメータ群は、
高周波成分を除去した気象データの時系列パターンから抽出したピークの値、位置、数、前記時系列パターンの積算値のうち少なくとも一つを用いて算出されること、
を特徴とする収量予測システム。
【請求項３】
請求項１に記載の収量予測システムにおいて、
前記記憶部に記憶された当年度以前のほ場ごとの属性情報から、実測ほ場を選定するためのパラメータの優先度を算出するパラメータ優先度算出部と、
前記優先度と、前記画像特徴量と、前記ほ場の位置情報から、実測ほ場を選定する実測ほ場選択部を有することを特徴とする収量予測システム。
【請求項４】
請求項１に記載の収量予測システムにおいて、
前記空撮画像の取得時における前記作物の生育ステージを推定する生育ステージ分類部をさらに有し、
前記生育ステージの推定は、当年度の時系列空撮画像と標高データの少なくとも一つと、当年度の時系列気象データから算出され、
前記生育ステージと、前記画像特徴量と、前記ほ場の位置情報から、実測ほ場を選定する実測ほ場選択部を有することを特徴とする収量予測システム。
【請求項５】
請求項３または４に記載の収量予測システムにおいて、
選定された前記実測ほ場の実測収量を入力できる実測データ入力部を有することを特徴とする収量予測システム。
【請求項６】
請求項１に記載の収量予測システムにおいて、
前記時系列パターン解析部は、
前記気象データの時系列パターンを取得し、前記ほ場で栽培される作物の特定の生育ステージごとに、また前記ほ場の撮影された空撮画像の撮影日前後ごとに、前記作物の生育状況と相関を有する第二パラメータ群を出力することを特徴とする収量予測システム。
【請求項７】
請求項３に記載の収量予測システムにおいて、
前記パラメータ優先度算出部は、
当年度以前の各年度におけるほ場の属性情報の統計情報を比較することにより、前記優先度を算出することを特徴とする収量予測システム。
【請求項８】
請求項１に記載の収量予測システムにおいて、
前記収量推計部は、
当年度の時系列気象データから、前記ほ場の生育ステージが変化する日付を推定する生育ステージ推定部と、を備えることを特徴とする収量予測システム。
【請求項９】
請求項１に記載の収量予測システムにおいて、
前記予測収量を表示する収量推計結果表示部をさらに有することを特徴とする収量予測システム。
【請求項１０】
ほ場の形状情報と位置情報と、当該ほ場が撮影された空撮画像とから、前記ほ場の画像特徴量を出力する画像解析部と、
前記画像特徴量より、前記ほ場の予測収量を出力する収量推計部とを有し、
前記収量推計部は、
予め記憶された気象データの時系列パターンを基に、前記ほ場で栽培される作物の特定の生育ステージごとに、前記作物の生育状況との相関を有するパラメータ群を出力する時系列パターン解析部を備え、
前記画像特徴量と前記パラメータ群を用いて前記ほ場の予測収量を算出すること
を特徴とする収量予測装置。
【請求項１１】
請求項１０に記載の収量予測装置において、
前記予測収量を表示する収量推計結果表示部をさらに有することを特徴とする収量予測装置。
【請求項１２】
請求項１０に記載の収量予測装置において、
前記ほ場の形状情報と位置情報とを記憶する記憶部を備えることを特徴とする収量予測装置。
【請求項１３】
請求項１０に記載の収量予測装置において、
前記空撮画像を受信する受信部を備えることを特徴とする収量予測装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、農作物の収量を推定するシステムや装置に関するものである。
【背景技術】
【０００２】
日本の公的機関や諸外国の損害保険会社などでは、主要農作物に対し農業災害補償制度を取り入れていることが多く、農作地の被害状況に応じて、農家に保険金を支払っている。保険金を支払うためには、農家から被害申告を受けたほ場を、調査員が目視検査、または実測による収量調査を行い、被害の程度から保険金を算出する。
さらに、全ほ場の平均的な収量を得るための被害ほ場以外における実測調査、各ほ場の基準となる単収をあらかじめ見積もるための実測調査なども行われ、調査員によるほ場実測、目視検査は多数回、広範囲におよび、コスト面で大きな負担となり、さらに調査員の高齢化も問題視されている。また人手を介在することによる評価結果の偏りも、農家の大きな不満となっている。
【０００３】
そのため地球表面を広域かつ反復的に撮影できる衛星画像等を用いて収量を一元的に予測するシステムが求められている。
【０００４】
衛星画像を用いて農作物の収量を推定する先行技術として、特許文献１に示すような技術が知られている。
【０００５】
特許文献１では、合成開口レーダ画像を用いて得られる生育期前半における標本圃場での後方散乱強度と、現地調査により得られる茎数等の水稲の生育的特徴との相関に基づき収量予測式を生成する水稲収量予測モデル生成方法、及び水稲収量予測方法が開示されている。
【先行技術文献】
【特許文献】
【０００６】
【特許文献１】特開２０１１−１６７１６３号公報
【発明の概要】
【発明が解決しようとする課題】
【０００７】
衛星画像等のスペクトルを用いて収量推計を行う際、高精度に推計するためには、実測ほ場（学習データ）が大量に必要となる。しかし多くのほ場に対し、現地実測をするのには多くの人手や時間を要する上に高コストとなる。
【０００８】
またスペクトルとしてほぼ同一なほ場でも、被害の有無によって収量推計モデルが異なるため、同一のモデルにより収量推計を行うと精度が下がるという問題がある。
【課題を解決するための手段】
【０００９】
上記課題を解決するための本願発明の一例としては、ほ場を含む空撮画像を受信できる受信部と、当該ほ場を含む複数のほ場の形状情報と、位置情報とをそれぞれ記憶する記憶部と、受信した前記ほ場を含む空撮画像と、当該ほ場の形状情報と位置情報とから、当該ほ場の画像特徴量を算出する画像解析部と、算出された前記画像特徴量から前記ほ場の予測収量を算出する収量推計部とを有し、前記収量推計部は、予め記憶された気象データの時系列パターンを基に、前記ほ場で栽培される作物の特定の生育ステージごとに、前記作物の生育状況と相関を有する第一パラメータ群を出力する時系列パターン解析部を備え、
前記画像特徴量と前記第一パラメータ群を用いて前記ほ場の予測収量を算出する。
【発明の効果】
【００１０】
本発明によると、特定の生育ステージごとの気象データの時系列パターンを用いることにより、被害の有無に関わらず、高精度な収量推計を実現できる。
【図面の簡単な説明】
【００１１】
【図１】本発明の第１の実施形態の収量予測システムの基本構成を表すブロック図である。
【図２】本発明の第１の実施形態の収量予測装置で行われる処理を表すシーケンス図である。
【図３】本発明の第１の実施形態のパラメータ優先度算出部にて行われる処理を表したフローチャート図である。
【図４】本発明の第１と第２の実施形態の画像解析部にて行われる処理を表したフローチャート図である。
【図５】本発明の第１と第２の実施形態の実測ほ場選択部にて行われる処理を表したフローチャート図である。
【図６】本発明の第１と第２の実施形態の収量推計部にて行われる処理を表したフローチャート図である。
【図７】本発明の第１と第２の実施形態のほ場ＧＩＳのデータ構造例である。
【図８】本発明の第１と第２の実施形態のほ場特徴データ、実測ほ場データと非実測ほ場データのデータ構造例である。
【図９】本発明の第１と第２の実施形態の気象ＤＢのデータ構造例である。
【図１０】本発明の第１と第２の実施形態の過去データと、第１の実施形態の統計データと、第１と第２の実施形態の生育ＤＢのデータ構造例である。
【図１１】本発明の第１と第２の実施形態の実測ほ場選択部で行われる処理の説明図である。
【図１２】本発明の第１と第２の実施形態のほ場内画素抽出部で行われる処理の説明図である。
【図１３】本発明の第１と第２の実施形態の時系列パターン解析部で行われる処理の説明図である。
【図１４】本発明の第１と第２の実施形態の実測データ入力部の入力画面である。
【図１５】本発明の第１と第２の実施形態の収量推計結果表示部の出力画面である。
【図１６】本発明の第２の実施形態の収量予測システムである。
【図１７】本発明の第２の実施形態の生育ステージ分類部で行われる処理のフローチャート図である。
【図１８】本発明の第２の実施形態のほ場生育ステージＤＢのデータ構造例である。
【図１９】本発明の第１と第２の実施形態を示し、収量予測装置のハードウェアを示すブロック図である。
【図２０】本発明の第１と第２の実施形態を示し、本発明の運用フローを示すシーケンス図である。
【発明を実施するための形態】
【００１２】
以下、本発明の一実施形態を添付図面に基づいて説明する。
【００１３】
本発明は、ほ場ごとの作物の収量等を推定するシステムに関するものである。ここでいうほ場とは、所定の領域内の農作区域を示すものであり、同一属性情報を持つと想定される空間的に連続した農作区域を含む。
また収量に限らず、窒素、リン酸、カリウム等の土壌成分値や農作物の草丈等の定量的なパラメータを推定する際にも本発明を適用することができる。
【００１４】
本発明に関する詳細については以下の実施形態例において説明するが、本発明は以下の内容を一例として含む。
【００１５】
過去に蓄積された調査地域の画像やほ場の属性情報と、調査対象であるほ場が撮影された空撮画像から、調査すべき実測ほ場を選定するための判断基準となるパラメータ群を決定し、前記パラメータ群ができる限り分散を持つように、実測ほ場を選定し、また調査の負担をできる限り減らすために、実測ほ場の候補ができる限り位置的に集中するように選定する。また気象データの時系列パターンを各生育ステージごとに解析することにより、生育状況に相関のあるパラメータ群を算出し、画像特徴量、ほ場の属性情報と前記パラメータ群を説明変数とした収量推計を実施する。
【００１６】
本発明は、過去に蓄積された調査地域の画像やほ場の属性情報と、調査対象であるほ場が撮影された空撮画像から、調査すべき実測ほ場を選定するための判断基準となるパラメータ群を決定し、前記パラメータ群ができる限り分散を持つように、実測ほ場を選定し、また調査の負担をできる限り減らすために、実測ほ場の候補ができる限り位置的に集中するように選定する。
【００１７】
また空撮画像上ではスペクトルが同様に見えるほ場でも、品種や作付日、被害、気象データの違いにより、収量推計モデルが異なる。そのため、本発明では、気象データの時系列パターンを各生育ステージごとに解析することにより、生育状況に相関のあるパラメータ群を算出し、画像特徴量、ほ場の属性情報と前記パラメータ群を説明変数とした収量推計を実施する。
【００１８】
本発明では、調査員が、本システムによって指定されたほ場の収穫量を調査するだけで、残りのほ場に関しては、衛星や航空機に搭載されたセンサによって得られた調査地域の空撮画像を用いて推定することが可能となることから、調査員の負担や調査にかかる費用を削減し、また収量推定の精度を向上させることができる。
【００１９】
以下、本発明の実施形態を、図面を参照して詳細に説明する。
【００２０】
＜第１の実施形態＞
図１は、収量予測システムの基本構成を表すブロック図である。
【００２１】
本実施形態の収量予測システムは収量予測装置１０を備える。
【００２２】
収量予測装置１０は、過去データ１１から昨年度以前に取得されたデータを、画像ＤＢ１２から収量推計対象であるほ場が撮影された今年度の画像を、ほ場ＧＩＳ１３から収量推計対象であるほ場のデータを受信した後、現地調査員１５へ実測対象として選定されたほ場の情報を送信する機能と、現地調査員１５から現地調査結果を受信した後、気象ＤＢ１４から今年度の時系列気象データを受信した後、収量推計結果１７を出力する機能と、を備える。
【００２３】
上記機能を実現するために、収量予測装置１０は、パラメータ優先度算出部１０１と、画像解析部１０２と、実測ほ場選択部１０３と、ほ場特徴データ１０５と、実測ほ場表示部１０６と、実測データ入力部１０７と、実測ほ場データ１０８と、非実測ほ場データ１０９と、収量推計部１１０と、収量推計結果表示部１１１と、を備える。
【００２４】
パラメータ優先度算出部１０１は、過去データ１１から昨年度以前に取得されたデータを受信し、実測ほ場を選定する判断材料に用いられるパラメータの優先度を、実測ほ場選択部１０３へ送信する機能を備える。上記機能を実現するために、パラメータ優先度算出部１０１は、統計解析部１０１１と、優先度算出部１０１２と、統計データ１０１３と、を備える。各部の具体的な処理内容と、統計データ１０１３に格納されるデータについては後述する。
【００２５】
画像解析部１０２は、画像ＤＢ１２から収量推計対象であるほ場が撮影された今年度の画像を、ほ場ＧＩＳ１３から収量推計対象であるほ場のデータを受信し、ほ場特徴データ１０５に画像から算出されたほ場スペクトル特徴量を格納する。上記機能を実現するために、画像解析部１０２は、ほ場内画素抽出部１０２１と、ほ場内特徴算出部１０２２と、を備える。各部の具体的な処理内容については後述する。
【００２６】
実測ほ場選択部１０３は、ほ場特徴データ１０５からほ場スペクトル特徴量を、パラメータ優先度算出部１０１から実測ほ場を選定する判断材料に用いられるパラメータの優先度を受信し、実測ほ場表示部１０６に、調査地域に存在するほ場の中から選択された実測ほ場の情報を送信する。上記機能を実現するために、実測ほ場選択部１０３は、ほ場選定軸決定部１０３１と、ヒストグラム作成部１０３２と、実測候補選定部１０３３と、実測コスト算出部１０３４と、実測ほ場決定部１０３５と、を備える。各部の具体的な処理内容については後述する。
【００２７】
実測ほ場表示部１０６は、実測ほ場選択部１０３から受信した、対象地域に存在するほ場の中から選択された実測ほ場の情報を現地調査員１５に表示する。具体的な表示画面例については後述する。
【００２８】
実測データ入力部１０７は、現地調査員１５から実測結果を、ほ場特徴データ１０５からほ場スペクトル特徴量を受信し、実測ほ場データ１０８、非実測ほ場データ１０９に格納されたデータを更新する。
【００２９】
収量推計部１１０は、ほ場ＧＩＳ１３から収量推計対象であるほ場のデータを、気象ＤＢ１４から今年度の気象データの時系列パターンを、実測ほ場データ１０８と非実測ほ場データ１０９からほ場のスペクトル特徴量と実測によって得た収量を受信し、非実測ほ場データ１０９に格納されているデータを更新し、収量推計結果表示部１１１へ収量推計の終了フラグを送信する。上記機能を実現するために、収量推計部１１０は、生育ステージ推定部１１０１と、時系列パターン解析部１１０２と、モデル策定部１１０３と、収量推計部１１０４と、生育ＤＢ１１０５とを備える。各部の具体的な処理内容については後述する。
【００３０】
収量推計結果表示部１１１は、収量推計部１１０から収量推計の終了フラグを受信した後、実測ほ場データ１０８から実測によって得た収量を、非実測ほ場データ１０９から推計によって得た収量を受信し、収量推計結果１７を出力する。
【００３１】
過去データ１１、ほ場ＧＩＳ１３、気象ＤＢ１４、ほ場特徴データ１０５、実測ほ場データ１０８、非実測ほ場データ１０９に格納されるデータについては後述する。なお収量予測装置１０は、過去データ１１と、画像ＤＢ１２と、ほ場ＧＩＳ１３と、気象ＤＢ１４とを内部に備えていてもよい。
【００３２】
画像ＤＢ１２は、少なくとも撮影日時、画像位置情報をヘッダ情報に含むような、衛星、または航空機に搭載されたセンサによって得られた画像データを格納する。ここで撮影日時は、当該画像データがいつ撮影されたかを表し、画像位置情報とは、画像の四隅の位置情報を表す。なお、画像位置情報は、例えば、世界測地系で用いられる緯度、経度の情報であってもよい。画像ＤＢ１２に含まれる空撮画像を主に用いることにより、収量推計を実施する。
【００３３】
図２は、収量予測装置１０にて行われる処理例を表すシーケンス図である。以下に具体的な処理について説明する。
【００３４】
まずほ場ＧＩＳ１３は調査地域の収量推計対象ほ場のＩＤと座標情報を画像解析部１０２へ送信する（Ｓ２０１）。この時に、画像解析部１０２は、画像ＤＢ１２から収量推計対象ほ場が撮影された今年度の画像を受信し、収量推計対象ほ場ごとの画像特徴量を算出する。そして収量推計対象ほ場の画像特徴量をほ場特徴データ１０５に格納する（Ｓ２０２）。ほ場特徴データ１０５は、画像解析部１０２から収量推計対象ほ場の画像特徴量を格納されると、前記画像特徴量を実測ほ場選択部１０３へ送信する（Ｓ２０４）。パラメータ優先度算出部１０１は過去データ１１から昨年度以前に取得されたデータを受信し、前記受信データに対して、統計解析を行うことにより、実測ほ場を選定する判断材料に用いられるパラメータ優先度を算出し、実測ほ場選択部１０３へ送信する（Ｓ２０３）。なおＳ２０１→Ｓ２０２→Ｓ２０４の処理と、Ｓ２０３の処理は並列に処理してもよいし、どちらの処理を先に行ってもよい。実測ほ場選択部１０３は、ほ場ＧＩＳ１３から調査地域の収量推計対象ほ場の情報を受信し（Ｓ２０５）、ほ場特徴データ１０５から前記画像特徴量を受信し、パラメータ優先度算出部１０１から前記パラメータ優先度を受信すると、実測ほ場の選定を行い、実測ほ場表示部１０６を通して、現地調査員１５に選定された実測ほ場に関する情報を送信する（Ｓ２０６）。現地調査員１５は、受信した実測ほ場に関する情報に従って現地実測を行い、各ほ場の収量を取得する（Ｓ２０７）。現地調査員１５によって各ほ場の収量を取得した後、実測ほ場データ１０８、非実測ほ場データ１０９は、ほ場特徴データ１０５から各ほ場のＩＤとその画像特徴量を受信し（Ｓ２０９）、また現地調査員１５より各ほ場の収量を、実測データ入力部１０７を通して、受信する（Ｓ２０８）。この際、実測ほ場データ１０８には、現地調査員１５によって調査されたほ場のＩＤと対応する画素特徴量、収量が格納され、非実測ほ場データ１０９には、現地調査員１５によって調査されなかったほ場のＩＤと対応する画素特徴量、収量が格納される。非実測ほ場データ１０９に格納される収量は、Ｓ２０８の段階では情報がないため、ＮＵＬＬを格納する。その後、収量推計部１１０は、実測ほ場データ１０８、非実測ほ場データ１０９から各ほ場のＩＤ、画素特徴量、収量を受信し（Ｓ２１０）、ほ場ＧＩＳ１３から収量推計対象であるほ場のデータを受信し（Ｓ２１１）、気象ＤＢ１４から今年度の気象データの時系列パターンを受信し（Ｓ２１２）、収量推計モデルを策定する。その後、策定された収量推計モデルに従って推定された、各ほ場の収量を実測ほ場データ１０８、非実測ほ場データ１０９に格納する（Ｓ２１３）。その後、収量推計部１１０は、収量推計が終了したフラグを収量推計結果表示部１１１に送信し（Ｓ２１５）、収量推計結果表示部１１１は、実測ほ場データ１０８、非実測ほ場データ１０９より、各ほ場のＩＤ、収量、推計収量等を受信する（Ｓ２１４）。最後に収量推計結果表示部１１１にて、結果を表示して、収量予測装置１０にて行われる処理が終了する。
【００３５】
図３は、パラメータ優先度算出部１０１にて行われる処理例を示すフローチャート図である。パラメータ優先度算出部１０１は、過去データ１１から昨年度以前に取得されたデータを受信し、実測ほ場を選定する判断材料に用いられるパラメータの優先度を、実測ほ場選択部１０３へ送信する。以下に具体的な処理について説明する。
【００３６】
Ｓ３１は、パラメータ優先度算出部１０１の処理開始を表す。
【００３７】
Ｓ３２では、過去データ１１から受信した昨年度以前に取得されたデータを統計解析部１０１１で解析し、統計データ１０１３に格納する。具体的には、昨年度以前に取得されたデータの各パラメータ（１１０３〜１１１１）の統計値を各年ごとに算出する。ここでいう統計値は平均、分散、サンプル数、ヒストグラム、最大値、最小値等のことを指す。ここでいうヒストグラムは、対応するパラメータの階級に応じた度数を表現したものである。なお階級の決め方は任意であり、例えば最小値から最大値の区域を１０分割して決定してもよい。前記パラメータが被害、品種等の離散値である場合、平均、分散、最大値、最小値だけは計算できないために、ＮＵＬＬを格納する。
【００３８】
Ｓ３３では、優先度算出部１０１２が統計データ１０１３から各パラメータの年間の統計値を受信し、各パラメータの優先度を算出し、各パラメータ名と算出された優先度の値のペアを実測ほ場選択部１０３へ送信する。具体的には、各パラメータの年間の分布を比較し、似た分布が多い場合は優先度を高く、年間の分布が大きく異なる場合は、優先度を低めに算出する。例えばパラメータが連続値の場合、合計年数をＩ、年のインデックスをｉ、該当パラメータの平均をMi、全ての年度を合わせた場合の分散をV、サンプル数をNiとすると、優先度は、以下の（数１）としてもよい。
【００３９】
【数１】
2015000049.tif 000003
【００４０】
または分布間距離を用いて表現してもよい。分布間距離の場合は、パラメータが連続値でも離散値でも適用することが可能である。例えば年iのパラメータxの確率分布をpi(x)とし、全ての年度を合わせたパラメータxの確率分布をp(x)とした時に、以下の（数２）としてもよい。
【００４１】
【数２】
2015000049.tif 000004
【００４２】
ここでd(pi(x),pj(x))は確率分布pi(x),pj(x)の距離を測る尺度であり、たとえばカルバック・ライブラー情報量でもよい。
【００４３】
Ｓ３４は、パラメータ優先度算出部１０１の処理終了を表す。
【００４４】
図４は、画像解析部１０２にて行われる処理例を示すフローチャート図である。画像解析部１０２は、画像ＤＢ１２から収量推計対象であるほ場が撮影された今年度の空撮画像を、ほ場ＧＩＳ１３から収量推計対象であるほ場のデータを受信し、ほ場特徴データ１０５に画像から算出されたほ場内の特徴量を格納する。以下に具体的な処理について説明する。
Ｓ４１は、画像解析部１０２の処理開始を表す。
Ｓ４２では、画像ＤＢ１２から受信した空撮画像と、ほ場ＧＩＳ１３から受信した各ほ場の形状情報から、ほ場内に含まれる空撮画像の画素値を抽出する。具体的な抽出方法については、後述する。
Ｓ４３では、Ｓ４２の処理時に抽出された画素数が１以上であればＳ４４へ進み、そうでなければＳ４５へ進む。
【００４５】
Ｓ４４では、ほ場内特徴算出部１０２２が、抽出された画素値の集合から特徴量を計算し、ほ場特徴データ１０５に対応するほ場ＩＤ１０５１を、画素値１０５２にほ場内の画素値を、特徴量１０５３に算出されたほ場内の特徴量を格納する。画素値１０５２に格納する画素値は、ほ場内の全画素値を並べてもよいし、ほ場内の平均画素値を用いてもよい。特徴量１０５３も同様である。ここで特徴量とは、画素値に対し、所定の演算をすることによって得られるパラメータのことであり、例えば、近赤外の画素値から赤の画素値を引いた値と、近赤外の画素値と赤の画素値の和の比である正規化植生指数等を特徴量として用いてもよいし、ほ場内に含まれる画素値の集合に対し、セグメンテーションを行い、各グループに所属する平均画素値に対して所定の演算を行うことにより算出されたパラメータと、各グループに所属する画素数をほ場内の画素数で割ったパラメータを並べたものを特徴量として用いてもよい。
【００４６】
Ｓ４５では、抽出された画素数が０なので、ほ場特徴データ１０５に対応するほ場ＩＤ１０５１を格納し、画素値１０５２と特徴量１０５３にはＮＵＬＬを格納する。
Ｓ４６では、終了判定を行う。ほ場ＧＩＳ１３に含まれる全てのほ場について、Ｓ４４、もしくはＳ４５を行った場合はＳ４７に進み、まだＳ４４、Ｓ４５のいずれも行っていないほ場が存在する場合はＳ４２に進み、次のほ場についてもほ場内画素値を抽出する。
Ｓ４７は、画像解析部１０２の処理終了を表す。
【００４７】
図５は、実測ほ場選択部１０３にて行われる処理例を示すフローチャート図である。実測ほ場選択部１０３は、ほ場特徴データ１０５からほ場スペクトル特徴量を、パラメータ優先度算出部１０１から実測ほ場を選定する判断材料に用いられるパラメータの優先度を受信し、実測ほ場表示部１０６に、調査地域に存在するほ場の中から選択された実測ほ場の情報を送信する。以下に具体的な処理について説明する。
【００４８】
Ｓ５１は、実測ほ場選択部１０３の処理開始を表す。
Ｓ５２では、ほ場特徴データ１０５とほ場ＧＩＳ１３から格納データを受信し、各パラメータに対する今年度の統計値を算出する。前記統計値は統計データ１０１３に格納されているデータと同様である。
【００４９】
Ｓ５３では、実測の際に重要視するパラメータを複数選択する。パラメータ優先度算出部１０１から受信した優先度が一定の閾値以上となったパラメータと、Ｓ５２で算出された各パラメータの統計値の中で分散が一定の閾値以上となったパラメータを選択してもよいし、各パラメータに関し、前記優先度と分散の値を掛けた値を算出し、その値に関して一定の閾値以上となったパラメータを選択してもよい。なおＳ５２で統計値を算出したパラメータが連続値でなく、離散値である場合は、分散の代わりにエントロピー等の統計値を利用してもよい。
Ｓ５２、Ｓ５３の処理は、ほ場選定軸決定部１０３１によって行われる。
Ｓ５４では、Ｓ５３にて選択された複数のパラメータを軸とする多変量ヒストグラムを作成する。
Ｓ５４の処理は、ヒストグラム作成部１０３２によって行われる。
【００５０】
Ｓ５５では、Ｓ５３で作成された多変量ヒストグラムの各軸を複数の区域に分割する。この際の分割方法や分割数は任意で、例えば各パラメータの最小値から最大値までの区域を１０等分してもよいし、各パラメータの、平均、中央値から標準偏差の定数倍引いた値から平均から標準偏差の定数倍足した値までの区域を１０等分してもよい。
Ｓ５６では、各軸について複数の区域に分割された多変量ヒストグラムから、実測ほ場候補の選定を行う。具体的な処理内容については、後述する。
Ｓ５５、Ｓ５６の処理は、実測候補選定部１０３３によって行われる。
【００５１】
Ｓ５７では、Ｓ５６で選ばれた実測ほ場の候補の全組み合わせに対し、実測コストを算出する。実測コストは、選ばれた実測ほ場が空間的にどのような分布を持っているかで決定することができ、基本的に選ばれた実測ほ場同士の距離が小さければ小さいほど、実測コストが小さくなり、前記距離が大きい場合は実測コストが大きくなる。選定ほ場の座標情報（各ほ場に関する座標データ１３０２の平均値）の集合から算出された共分散行列の行列式を実測コストとしてもよい。
Ｓ５７の処理は、実測コスト算出部１０３４によって行われる。
【００５２】
Ｓ５８では、算出した実測コストが最も小さくなるような実測ほ場候補を最終的な実測ほ場として選択し、実測ほ場表示部１０６に送信する。この際に、ほ場ＧＩＳ１３に格納されているデータの内、選択されたほ場ＩＤに対応するほ場リストを送信する。
Ｓ５８の処理は、実測ほ場決定部１０３５によって行われる。
Ｓ５９は、実測ほ場選択部１０３の処理終了を表す。
【００５３】
図６は、収量推計部１１０にて行われる処理例を示すフローチャート図である。ほ場ＧＩＳ１３から収量推計対象であるほ場のデータを、気象ＤＢ１４から今年度の気象データの時系列パターンを、実測ほ場データ１０８と非実測ほ場データ１０９からほ場のスペクトル特徴量と実測によって得た収量を受信し、非実測ほ場データ１０９に格納されているデータを更新し、収量推計結果表示部１１１へ収量推計の終了フラグを送信する。以下に具体的な処理について説明する。
【００５４】
Ｓ６１は、収量推計部１１０の処理開始を表す。
Ｓ６２では、気象ＤＢ１４から時系列気象データを受信し、各ほ場に関して生育ステージの時期を推定する。具体的には各ほ場に対し、田植日から画像撮影日もしくは収穫日までの全ての日について、積算気温、積算日射量、積算降水量等を算出し、生育ＤＢ１１０５の格納データとの比較を行う。これにより、田植日から画像撮影日もしくは収穫日までの全ての日が、どの生育ステージに相当するかを推定する。この際、積算気温、積算日射量、積算降水量等を説明変数、生育ステージが目的変数となる。推定方法としてｋ近傍法やｋ平均法、混合正規分布等のクラスタリング手法等を適用してもよい。なおこの際、用いる気象データは、対応するほ場から最も近い位置で観測される気象データでもよいし、複数地点の気象データを用いて対応するほ場の位置の気象データを補間してもよい。
Ｓ６２の処理は、生育ステージ推定部１１０１によって行われる。
【００５５】
Ｓ６３では、対応するほ場の時系列気象データから収量推計に用いる説明変数を算出する。Ｓ６１で推定された生育ステージの時期ごとに気象データの各パラメータの時系列パターンを解析することによって説明変数を計算する。詳細については後述する。
Ｓ６３の処理は、時系列パターン解析部１１０２によって行われる。
【００５６】
Ｓ６４では、実測ほ場データ１０８から受信した各ほ場のスペクトル特徴量、収量と、またＳ６３で算出された時系列気象データに基づく説明変数を用いることにより、収量推計モデルを策定する。この際、収量推計モデルにおける説明変数は、実測ほ場データ１０８から受信した各ほ場のスペクトル特徴量とＳ６３で算出された時系列気象データに基づく説明変数であり、目的変数は収量となる。なお収量推計モデル策定のために用いる手法は多々存在する。例えば、説明変数のインデックスをi、説明変数をXi、目的変数をYとした時に、(数３)のような線形モデルを立てたり、(数４)や(数５)のような非線形モデルを立てて、各パラメータを求めてもよい。
【００５７】
【数３】
2015000049.tif 000005
【００５８】
【数４】
2015000049.tif 000006
【００５９】
【数５】
2015000049.tif 000007
【００６０】
また品種毎の混合モデルを用いてもよいし、階層ベイズモデルを用いてもよい。
Ｓ６４の処理は、モデル策定部１１０３によって行われる。
【００６１】
Ｓ６５は、実測ほ場データ１０８、非実測ほ場データ１０９から受信した各ほ場のスペクトル特徴量、Ｓ６３で算出された時系列気象データに基づく説明変数を、Ｓ６４で策定された収量推計モデルに代入することにより、各ほ場の収量を推計する。そして各ほ場について推定された収量を、実測ほ場データ１０８、非実測ほ場データ１０９の推定収量１０８５に格納する。その後、収量推計結果表示部１１１に終了フラグを送信する。
Ｓ６５の処理は、収量推計部１１０４によって行われる。
Ｓ６６は、収量推計部１１０の処理終了を表す。
【００６２】
図７は、ほ場ＧＩＳ１３のデータ構造例である。ほ場ＧＩＳ１３は事前に用意されていた各ほ場の座標や形状、属性情報に関するデータ１３０１〜１３０９を格納する。
ほ場ＩＤ１３０１は、収量推計対象であるほ場を識別するためのラベルである。
座標データ１３０２は、各ほ場の形状を多角形とした時の、各頂点の座標を並べたものである。座標データ１３０２に格納される座標は任意の座標系でよく、例えば世界測地系の緯度、経度を格納してもよい。
【００６３】
肥料１３０３は、各ほ場に対する施肥量を表す。肥料１３０３は、例えば栽培前に施された窒素量や燐酸量、カリウム量でもよい。単位は任意であり、例えば１０アール当たりのキログラムでもよい。
品種１３０４は、各ほ場で栽培される作物の品種名が格納される。品種１３０４は、作物名も一緒に格納されてもよい。
被害１３０５は、各ほ場で起きた被害名が格納される。一つのほ場にて複数の被害が起きた際には、起きた被害の中で一番影響が大きいものを格納してもよいし、起きた被害全てを格納してもよい。
【００６４】
基準単収１３０６は、各ほ場で平均的に収穫できると想定される収量が格納される。単位は任意であり、例えば１０アール当たりのキログラムでもよい。
田植期１３０７は、各ほ場における今年度の田植日が格納される。
出穂揃日１３０８は、各ほ場における今年度の出穂揃日が格納される。
収穫日１３０９は、各ほ場における今年度の収穫日が格納される。
【００６５】
図８は、ほ場特徴データ１０５、実測ほ場データ１０８と非実測ほ場データ１０９のデータ構造例である。ほ場特徴データ１０５は、各ほ場の画素特徴量等に関するデータ１０５１〜１０５３を格納する。実測ほ場データ１０８、非実測ほ場データ１０９は、各ほ場の画素特徴量、収量等に関するデータ１０８１〜１０８５を格納する。
【００６６】
ほ場ＩＤ１０５１は、収量推計対象であるほ場を識別するためのラベルであり、ほ場ＧＩＳ１３におけるほ場ＩＤ１３０１と一致する。
【００６７】
画素値１０５２は、調査地域の空撮画像中において対応するほ場内に含まれる画素値である。ほ場内に含まれる画素値全てを格納してもよいし、ほ場内の平均画素値を格納してもよい。
【００６８】
特徴量１０５３は、画素値１０５２から所定の演算をすることによって算出されたパラメータを並べたものである。例えば近赤外の画素値から赤の画素値を引いた値と、近赤外の画素値と赤の画素値の和の比である正規化植生指数等をパラメータとして用いてもよいし、ほ場内に含まれる画素値の集合に対し、セグメンテーションを行い、各グループに所属する平均画素値に対して所定の演算を行うことにより算出されたパラメータと、各グループに所属する画素数をほ場内の画素数で割ったパラメータを並べてもよい。
【００６９】
ほ場ＩＤ１０８１は、収量推計対象であるほ場を識別するためのラベルであり、ほ場ＧＩＳ１３におけるほ場ＩＤ１３０１と一致する。
【００７０】
画素値１０８２、調査地域の空撮画像中において対応するほ場内に含まれる画素値であり、画素値１０５２と一致する。
【００７１】
特徴量１０８３は、画素値１０５２から所定の演算をすることによって算出されたパラメータを並べたものであり、特徴量１０５３と一致する。
【００７２】
収量１０８４は、対応するほ場の収量である。単位は任意であり、例えば１０アール当たりのキログラムでもよい。非実測ほ場データ１０９に格納されるほ場は、調査員によって実測されないほ場なので、非実測ほ場データ１０９の収量１０８４は、ＮＵＬＬが格納される。
【００７３】
推定収量１０８５は、対応するほ場の推定収量である。単位は任意であり、例えば１０アール当たりのキログラムでもよい。推定収量１０８５は収量推計部１１０によって推定された収量が格納される。そのため、Ｓ２０９が行われた直後は、ＮＵＬＬが格納されており、Ｓ２１３が行われた後に推定された収量が格納される。
【００７４】
図９は、気象ＤＢ１４のデータ構造例である。
気象ＤＢ１４は、調査地域にて観測された時系列気象データに関するデータ１４００〜１４１２を格納する。気象ＤＢ１４は複数地点で観測されたデータを格納している。
観測ＩＤ１４００は、観測地点を識別するためのＩＤである。
観測場所１４０１は、対応する観測地点の場所に関する座標情報である。観測場所１４０１に格納される座標は、任意の座標系でよく、例えば世界測地系の緯度、経度を格納してもよい。
【００７５】
日付１４０２は、気象データが観測された日付を表す。
平均気温１４０３は、対応する日付内での観測地点の平均気温を表す。単位は任意であり、例えば℃でもよい。
最高気温１４０４は、対応する日付内での観測地点の最高気温を表す。単位は任意であり、例えば℃でもよい。
最低気温１４０５は、対応する日付内での観測地点の最低気温を表す。単位は任意であり、例えば℃でもよい。
日射量１４０６は、対応する日付の観測地点の平均日射量を表す。なお気温の場合と同じく、最高日射量や最低日射量も同様に格納してもよい。単位は任意であり、例えばＷ／ｍ２でもよい。
【００７６】
湿度１４０７は、対応する日付の観測地点の平均湿度を表す。なお気温の場合と同じく、最高湿度や最低湿度も同様に格納してもよい。単位は任意であり、例えば％でもよい。
【００７７】
気圧１４０８は、対応する日付の観測地点の平均気圧を表す。なお気温の場合と同じく、最高気圧や最低気圧も同様に格納してもよい。単位は任意であり、例えば％でもよい。
【００７８】
降水量１４０９は、対応する日付の観測地点の積算降水量を表す。単位は任意であり、例えばｍｍでもよい。
日照時間１４１０は、対応する日付の観測地点の日照時間を表す。単位は任意であり、例えばｈでもよい。
風向１４１１は、対応する日付の観測地点の風向を表す。
風速１４１２は、対応する日付の観測地点の平均風速を表す。なお瞬間最大風速なども同様に格納してもよい。単位は任意であり、例えばｍ／ｓでもよい。
【００７９】
なお気象ＤＢ１４は、以上の例のように観測地点ごとにそれぞれテーブルを持っていてもよいし、観測ＩＤ１４００と観測場所１４０２をカラムに追加して、全データを１つのテーブルに格納してもよい。
【００８０】
図１０は、過去データ１１と、統計データ１０１３、生育ＤＢ１１０５のデータ構造例である。過去データ１１は、各ほ場について昨年度以前に取得された各ほ場ごとの属性情報や画素特徴量等のデータ１１０１〜１１１２を格納する。基本的にほ場ＧＩＳ１３とほ場特徴データ１０５と一致するが、ほ場ＧＩＳ１３とほ場特徴データ１０５は今年度取得されるものであるという点でのみ異なる。
【００８１】
ほ場ＩＤ１１０１は、収量推計対象である対応する年度のほ場を識別するためのラベルであり、ほ場ＩＤ１３０１と一致する。
座標データ１１０２は、対応する年度の各ほ場の形状を多角形とした時の、各頂点の座標を並べたものであり、座標データ１３０２と一致する。
肥料１１０３は、対応する年度の各ほ場に対する施肥量を表し、肥料１３０３と一致する。
【００８２】
品種１１０４は、対応する年度の各ほ場で栽培される作物の品種名が格納され、品種１３０４と一致する。
被害１１０５は、対応する年度に各ほ場で起きた被害名が格納され、被害１３０５と一致する。
基準単収１１０６は、対応する年度の各ほ場で平均的に収穫できると想定される収量が格納され、基準単収１３０６と一致する。
田植期１１０７は、対応する年度の各ほ場における田植日が格納され、田植期１３０７と一致する。
出穂揃日１１０８は、対応する年度の各ほ場における出穂揃日が格納され、出穂揃日１３０８と一致する。
収穫日１１０９は、対応する年度の各ほ場における収穫日が格納され、収穫日１３０９と一致する。
【００８３】
画素値１１１０は、対応する年度における、調査地域の空撮画像中において対応するほ場内に含まれる画素値であり、画素値１０５２と一致する。
特徴量１１１１は、画素値１１１０から所定の演算をすることによって算出されたパラメータを並べたものであり、特徴量１０５３と一致する。
年１１１２は、対応するデータが取得される年度を表す。
統計データ１０１３は、過去データ１１に関する各年の統計データ１０１３１〜１０１３８を格納する。
年１０１３１は、統計値を算出したデータの取得年を表す。過去データ１１の格納データ１１１２に対応する。
【００８４】
パラメータ１０１３２は、統計値を算出するパラメータの種類を表す。過去データ１１の格納データ１１０２〜１１１１のカラム名に対応する。
平均１０１３３は、対応する年の、対応パラメータに関する平均を表す。なお対応パラメータが品種や被害等の離散値に相当する場合は、平均を計算することができないために、ＮＵＬＬ値を格納する。
【００８５】
分散１０１３４は、対応する年の、対応パラメータに関する分散を表す。同様に、対応パラメータが品種や被害等の離散値に相当する場合は、分散を計算することができないために、ＮＵＬＬ値を格納する。
サンプル数１０１３５は、対応する年の、対応パラメータに関するデータ数を表す。
【００８６】
区域定義１０１３６は、対応する年の、対応パラメータに関するヒストグラムにおける階級の定義域を表す。具体的に説明すると、対応パラメータが連続値である場合、区域定義１０１３６は、各階級の定義域の最大値、最小値のペアを階級の分割数だけ並べたものに相当する。対応パラメータが離散値の場合は、対応パラメータがとりうる値全てを並べたものに相当する。
【００８７】
区域Ａ１０１３７は、対応する年の、対応パラメータに関するデータの中で、一番最初の区域定義に存在するデータ数を格納する。
区域Ｂ１０１３８は、対応する年の、対応パラメータに関するデータの中で、二番目の区域定義に存在するデータ数を格納する。同様に区域Ｃ、区域Ｄ・・・と、区域が階級数の分、存在しうる。例えば階級数が１０だった場合、カラムとして区域Ｉまで存在する。これに対し、例えば対応パラメータが５種類しか存在しない離散値の場合、区域Ｆ〜ＩにはＮＵＬＬ値を格納する。
【００８８】
生育ＤＢ１１０５は、各生育ステージと積算気温、積算日射量、積算降水量、品種との関係をデータベースとして格納している。
生育ステージ１１０５１は、各生育ステージを表す。
積算気温１１０５２には、対応する生育ステージの時の積算気温が格納される。単位は任意であり、例えば℃でもよい。
積算日射量１１０５３には、対応する生育ステージの時の積算日射量が格納される。単位は任意であり、例えばＷ／ｍ２でもよい。
積算降水量１１０５４には、対応する生育ステージの時の積算降水量が格納される。単位は任意であり、例えばｍｍでもよい。
品種１１０５５は、該当データの品種情報を表す。なおここに作物の種類の情報を加えてもよい。
【００８９】
図１１は、実測ほ場選択部１０３のＳ５６で行われる処理の説明図である。以下に具体的な処理内容について説明する。
Ｓ５６では選定されたパラメータの軸を分割し、生成された各グリッドから決められた数のデータを抽出する。図１１は選定されたパラメータが基準単収（軸５６１）と、正規化植生指数（軸５６２）である場合の説明図である。
【００９０】
まず軸５６１、５６２に対し、ほ場特徴データ１０５とほ場ＧＩＳ１３を結合したデータに関して、散布図５６０を作成する。そして軸５６１、５６２を一定数で分割（５６３）する。ここで分割の方法や分割数は任意である。例えばパラメータの最小値から最大値までの区域を１０等分してもよいし、パラメータの、平均、中央値から標準偏差の定数倍引いた値から平均から標準偏差の定数倍足した値までの区域を１０等分してもよい。分割によって複数のグリッド５６４が生成されるため、各グリッドからそれぞれ定数個のデータを抽出（５６５）する。もしグリッド内に存在するデータ数が定数個以下の場合は、存在するデータ数のみ抽出する。そして抽出の全組合せを実測ほ場候補として選定し、Ｓ５７に進む。例えば、抽出するデータ数を３個、全グリッド数が４つ、それぞれのグリッド内に存在するデータ数が１０、１、２、５である場合、一つ目のグリッドでは１０！／（７！×３！）通りの選び方があり、二つ目、三つ目のグリッドでは１通り、四つ目のグリッドでは５！／（２！×３！）通りの選び方があり、合計すると、１２００通りの選び方が存在する。これらのデータの選び方を実測ほ場候補として出力する。ここで、！は階乗を意味する。
【００９１】
図１２は、ほ場内画素抽出部１０２１で行われる処理の説明図である（Ｓ４２）。以下に具体的な処理内容について説明する。
図１２は画像４２１に対し、ほ場４２２内に含まれる画素を抽出する時の例である。画素４２３は、ほ場４２２近傍の画素の集合である。
画像４２１の四隅の位置情報を左上：（ｓｘ１、ｓｙ１）、右上：（ｓｘ２、ｓｙ１）、左下：（ｓｘ１、ｓｙ２）、右下：（ｓｘ２、ｓｙ２）とし、ほ場４２２の各頂点の位置情報を（ｘ１、ｙ１）、（ｘ２、ｙ２）、（ｘ３、ｙ３）、（ｘ４、ｙ４）・・・とすると、画像４２１の位置座標とほ場４２２の位置座標を比較することにより、画像内の各画素がほ場内に含まれるか否かを調べることができる。ほ場内に含まれる画素が、ほ場内画素抽出部１０２１の出力となる。なおほ場４２２をまたぐ画素が存在することが多々ある。この場合、各々の画素の中心がほ場内に含まれる場合にのみ、上記出力に加える。
【００９２】
またほ場ＧＩＳ１３や画像ＤＢ１２の位置精度によっては、この抽出方法だとほ場外のあぜ道等の画素が混在することもある。そのため、ほ場４２２に対し、一定距離のバッファ４２４を設け、新しくほ場４２５を作成し、ほ場内画素を抽出してもよい。
【００９３】
図１３は、時系列パターン解析部１１０２で行われる処理の説明図である（Ｓ６３）。以下に具体的な処理内容について説明する。
【００９４】
Ｓ６３では各気象データの時系列パターンを解析し、作物の生育に相関のあるパラメータを算出し、収量推計に用いる説明変数として出力する。グラフ６３０１は横軸に日付６２０２、縦軸に降水量６３０３としてデータをプロットした図である。日付６３０２に対し、生育ステージを推定しているため、田植期６３０４等の各生育ステージごとに時系列データを分割することができる。また画像取得日６３０５は画像ＤＢ１２から取得される空撮画像の撮影日であり、データのヘッダ情報に付加されている。
【００９５】
時系列パターンを解析するために、まずグラフ６３０１のノイズとして表れている高周波成分を除去する。そのために、例えば移動平均法やローパスフィルタ等により、データを処理してもよい。グラフ６３０６は、高周波成分を除去して得られた時系列パターンである。
【００９６】
グラフ６３０６の各生育ステージに存在する時系列パターンのピーク値６３０７や各ピークの幅６３０８（例えば分散）、ピーク位置６３０９、積算値６３１０等を収量推計に用いる説明変数として出力する。また画像取得日６３０５の前後で上記パラメータを分けて算出し、説明変数として加えてもよい。
【００９７】
この時系列パターン解析による効果を下記に説明する。
グラフ６３１１はほ場ごとの生長パラメータの時系列変化をグラフに表したものである。ここでいう生長パラメータとは、例えば生育の良さと正の相関を持つパラメータで、例えば画像から算出された正規化植生指数等を用いてもよい。
【００９８】
曲線６３１２は作付日は早いが途中で被害を受けてしまっているほ場の生長曲線、曲線６３１３は作付日は遅いが被害を受けずに順調に生長しているほ場の生長曲線である。例えば、曲線６３１２、６３１３の交点６３１４の時点で画像が撮影されると、画像のスペクトルより、収量推計を行うことが難しくなる。そのため、時系列気象データから、画像撮影前後の生長曲線を再現することにより、説明変数を取得することにより、収量推計の精度向上を実現することが可能となる。例えば、出穂期から収穫期にかけて雨が続くと、倒伏の被害が出て減収となったり、幼穂分化期等の時期に高温が続くと、高障害を起こし、生長曲線に影響が出る。そのために、時系列パターンから各生育ステージごとに算出されたパラメータ６３０７〜６３１０を説明変数として収量推計に利用する。
【００９９】
なおこの図１３の処理に入力する気象データとして、天気予報等のような、未来を予測したデータを入力することによって、当年度の収量を作付前に予測することが可能となる。そのため、例えば品種、作付日等、様々な条件ごとのシミュレーションによる収量予測を行い、その結果に基づき、営農支援等に応用することも可能である。
【０１００】
図１４は、実測データ入力部１０７の入力画面の例である。以下に入力画面の詳細について説明する。
【０１０１】
画面１０７１０は、ほ場のポリゴンを表示した図である。ポリゴン表示時に、衛星画像等の空撮画像の上にほ場ポリゴンを重畳表示してもよい。
ほ場１０７０１は現在選択されているほ場であり、ほ場１０７０２は選択されていないほ場である。ほ場１０７０１の情報は表１０７０３に表示される。表１０７０３に表示されるのはほ場ＧＩＳ１３に格納されているデータと、現地調査員によって調査されたほ場の実測収量１０７０４である。実測収量１０７０４に、キーボード等を用いてデータ入力することが可能である。
【０１０２】
矢印１０７０５〜１０７０８を選択することにより、画面１０７１０に表示されている地図の範囲を、矢印の方向に変更することができる。
【０１０３】
ボタン１０７１０〜１０７１７を選択することにより、各処理を行うことが可能である。この際、ボタンはタッチパネルのように画面に触れるだけで選択することができてもよいし、マウスのカーソル１０７０９やキーボード入力等によってボタンを選択することができてもよい。
【０１０４】
ボタン１０７１０を選択することにより、画面１０７００に表示されているポリゴンの表示方法を変更することができる。例えば、表１０７０３の属性値によってポリゴンの色や枠線の太さ等を変更して表示することができたり、実測収量１０７０４を入力済みのほ場と入力していないほ場でポリゴンの色や枠線の太さ等を変えてもよい。また画面１０７００では、ポリゴンの上に対応するほ場ＩＤを重畳表示しているが、重畳表示する属性の種類を変更してもよい。
【０１０５】
ボタン１０７１１を選択することにより、画面１０７００の表示の縮尺を拡大することができる。ボタン１０７１２を選択することにより、画面１０７００の表示の縮尺を縮小することができる。
【０１０６】
ボタン１０７１３を選択することにより、各属性値ごとにほ場ポリゴンを検索することができる。検索後は、画面１０７００が検索されたほ場の近傍を表示する。
【０１０７】
ボタン１０７１４は、ほ場ＧＩＳ１３や実測ほ場データ１０８、非実測ほ場データ１０９にデータを新しく追加することができる。キーボード等で表１０７０３に追記することにより、ほ場が追加される。
【０１０８】
ボタン１０７１５は、ほ場ＧＩＳ１３や実測ほ場データ１０８、非実測ほ場データ１０９にデータを更新することができる。キーボード等で表１０７０３を編集することにより、ほ場が更新される。
【０１０９】
ボタン１０７１６は、ほ場ＧＩＳ１３や実測ほ場データ１０８、非実測ほ場データ１０９にデータを削除することができる。現在選択しているほ場のデータが削除される。
【０１１０】
ボタン１０７１７は、データ入力をキャンセルすることができる。これにより、現地調査員によるデータ入力が終了する。
【０１１１】
図１５は、収量推計結果表示部１１１の出力画面の例である。以下に出力画面の詳細について説明する。
【０１１２】
画面１１１００、ほ場１１１０１、１１１０２、表１１１０３、実測収量１１１０４、矢印１１１０５〜１１１０８、ボタン１１１１０〜１１１１３は、それぞれ画面１０７００、ほ場１０７０１、１０７０２、表１０７０３、実測収量１０７０４、矢印１０７０５〜１０７０８、ボタン１０７１０〜１０７１３と同様である。
【０１１３】
推計収量１１１１４は、実測ほ場データ１０８、非実測ほ場データ１０９に格納されている推定収量１０８５のデータである。
【０１１４】
ボタン１１１１５は、収量推計の精度評価を行う。実測収量１１１０４と推定収量１１１１４がそろっているほ場は、実測ほ場データ１０８に格納されているほ場のみなので、そのほ場に対し、実測収量と推定収量を比較し、統計値を算出することにより、精度評価を行う。例えば、実測収量と推定収量の差の平均値であるＭＡＥ（Ｍｅａｎ  Ａｖｅｒａｇｅ  Ｅｒｒｏｒ）を算出したり、実測収量と推定収量の差が一定の閾値以内となるほ場の割合を算出してもよいし、実測収量と推定収量の相関係数を算出してもよい。
【０１１５】
ボタン１１１１６は、実測収量１１１０４や推定収量１１１１４に関するヒストグラムを生成する。ここでいうヒストグラムとは、対応するパラメータの階級に応じた度数をグラフとして表現したものである。なお階級の決め方は任意であり、例えば最小値から最大値の区域を１０分割して決定してもよい。
【０１１６】
ボタン１１１１７は、表１１１０３に含まれる属性名、属性値をＣＳＶ（Ｃｏｍｍａ−Ｓｅｐａｒａｔｅｄ  Ｖａｌｕｅｓ）ファイルとして出力する。この際に、出力データの属性値の条件を入力することにより、出力するデータをフィルタリングすることができてもよい。
【０１１７】
図１９は、収量予測システムのハードウェア構成例を示すブロック図である。
【０１１８】
収量予測装置１０は、操作部１９０１と、表示部１９０２と、プロセッサ１９０３と、メインメモリ１９０４と、記憶装置１９０５と、を備える計算機である。
【０１１９】
プロセッサ１９０３は、メインメモリ１９０４に格納されたプログラムを実行する。
【０１２０】
メインメモリ１９０４は、例えば半導体メモリであり、プロセッサ１９０３によって実行されるプログラム及びプロセッサ１９０３によって参照されるデータを格納する。具体的には、記憶装置１９０５に格納されたプログラム及びデータの少なくとも一部が、必要に応じてメインメモリ１９０４にコピーされる。
【０１２１】
プロセッサ１９０３は、各機能部のプログラムに従って動作することによって、所定の機能を実現する機能部として動作する。例えば、収量予測プログラムに従って動作することで、プロセッサ１９０３は収量推計部１１０として機能する。他のプログラムについても同様である。収量予測装置１０は、これらの機能部を含む装置である。
【０１２２】
操作部１９０１は、ユーザからの入力操作を受ける。操作部は、例えばキーボードまたはマウス等を含んでもよい。
【０１２３】
表示部１９０２は、ユーザへ情報を出力する。表示部１９０２は、例えば液晶ディスプレイのような画像表示装置であってもよい。
【０１２４】
記憶装置１９０５は、例えばハードディスク装置（ＨＤＤ）又はフラッシュメモリのような不揮発性の記憶装置である。本実施形態の記憶装置１９０５は、少なくともパラメータ優先度算出部１０１、画像解析部１０２、実測ほ場選択部１０３、収量推計部１１０が格納される。なお、各ＤＢや記憶部も、記憶装置１９０５に格納されてもよい。
【０１２５】
また物体認識装置１０の各機能を実現するプログラム、ＤＢ等の情報は、記憶装置１５０５や不揮発性半導体メモリ、ハードディスク装置、ＳＳＤ（Ｓｏｌｉｄ  Ｓｔａｔｅ  Ｄｒｉｖｅ）等の記憶デバイス、または、ＩＣカード、ＳＤカード、ＤＶＤ等の計算機読み取り可能な非一時的データ記憶媒体に格納することができる。
【０１２６】
図２０は、本発明の運用例を示すシーケンス図である。以下に具体的なフローについて説明する。
【０１２７】
まず顧客３０４が収量推計業者３０３に対し、収量推計の発注を行う（３０５）。その際に、顧客３０４が保持しているほ場ＧＩＳ１３のデータを提供する（３０６）。その後、収量推計業者３０３は、空撮画像の撮影業者３０２に画像撮影リクエストを行う（３０７）。撮影業者３０２は、画像撮影リクエストを受けると、衛星サーバ３０１に対し、撮影命令を行う（３０８）。そして画像撮影３０９が開始される。画像が撮影されたタイミングで撮影業者３０２に画像が送信される（３１０）。そして収量推計業者３０３に画像が送信され（３１１）、顧客３０４にも画像が送信される（３１２）。その後、顧客３０４は、撮影された画像で収量推計を行うかどうかを判断する（３１３）。例えば、画像中の雲量が多い場合等は、その画像は使わないと判断し（３１４）、その際、収量推計業者３０３は、画像撮影の再オーダーを行う（３１５）。その後は３０８〜３１３のフローを繰り返す。３１２の後に、画像を使うと判断した場合（３１６）、収量推計業者３０３は、収量推計装置１０により、収量推計を実施する（３１７）。なお３０６は、３１６の後から３１７の前に行われてもよい。その後、顧客３０４に結果を送信する（３１８）。
【０１２８】
＜第２の実施形態＞
以下、本発明の第２の実施形態について図面を参照して説明する。
【０１２９】
本発明の第１の実施形態では、過去に蓄積されたデータのばらつきから実測ほ場の選定を行い、その実測ほ場のリストを現地調査員に送信した。しかし過去に蓄積されたデータが存在しない場合もあるため、第２の実施形態では、収量推計を実施する年のデータのみを用いて、実測ほ場の選定を行う。そのために時系列空撮画像や標高データを用いて、収量推計に用いる空撮画像が撮影されたタイミングの生育ステージを推定し、その生育ステージをほ場を選定する際に用いるパラメータとして用いる。
【０１３０】
第２の実施形態では、パラメータ優先度算出部１０１、過去データ１１以外は第１の実施形態と同じである。具体的には、図４〜９、１１〜１５は第２の実施形態でも適用される。以下、第２の実施形態について、第１の実施形態と異なる点のみを説明する。
【０１３１】
図１６は、収量予測システムの基本構成を表すブロック図である。
【０１３２】
本実施形態の収量予測システムは収量予測装置２０を備える。
【０１３３】
収量予測装置２０は、標高ＤＢ２１から対象地域のＤＥＭ（Ｄｉｇｉｔａｌ  Ｅｌｅｖａｔｉｏｎ  Ｍｏｄｅｌ）データを、時系列画像ＤＢ２２から対象地域の時系列空撮画像を、画像ＤＢ１２から収量推計対象であるほ場が撮影された今年度の画像を、ほ場ＧＩＳ１３から収量推計対象であるほ場のデータを受信した後、現地調査員１５へ実測対象として選定されたほ場の情報を送信する機能と、現地調査員１５から現地調査結果を受信した後、気象ＤＢ１４から今年度の時系列気象データを受信した後、収量推計結果１７を出力する機能と、を備える。
【０１３４】
上記機能を実現するために、収量予測装置１０は、生育ステージ分類部２０１と、画像解析部１０２と、実測ほ場選択部１０３と、ほ場特徴データ１０５と、実測ほ場表示部１０６と、実測データ入力部１０７と、実測ほ場データ１０８と、非実測ほ場データ１０９と、収量推計部１１０と、収量推計結果表示部１１１と、を備える。
【０１３５】
画像解析部１０２と、実測ほ場選択部１０３と、ほ場特徴データ１０５と、実測ほ場表示部１０６と、実測データ入力部１０７と、実測ほ場データ１０８と、非実測ほ場データ１０９と、収量推計部１１０と、収量推計結果表示部１１１の機能は第１の実施形態と同一である。
【０１３６】
生育ステージ分類部２０１は、標高ＤＢ２１から対象地域のＤＥＭ（Ｄｉｇｉｔａｌ  Ｅｌｅｖａｔｉｏｎ  Ｍｏｄｅｌ）データを、時系列画像ＤＢ２２から対象地域の時系列空撮画像を受信し、実測ほ場選択部１０３に実測ほ場を選択する際に用いるパラメータの優先度を送信する機能を備える。上記機能を実現するために、生育ステージ分類部２０１は、気温変動推定部２０１１と、積算気温推定部２０１２と、生育ステージ推定部２０１３と、ほ場生育ステージＤＢ２０１４とを備える。具体的な処理内容については後述する。
【０１３７】
標高ＤＢ２１には対象地域のＤＥＭデータが、時系列画像ＤＢ２２には対象地域の時系列空撮画像が格納される。いずれのデータも画像ＤＢ１２に格納される空撮画像と同一の形式のデータである。但し、ＤＥＭデータの各画素値は、その位置座標における標高値を表している。単位は任意であり、例えばｍでもよい。
【０１３８】
図１７は、生育ステージ分類部２０１にて行われる処理例を表すフローチャート図である。以下に具体的な処理内容について説明する。
【０１３９】
Ｓ２０１０１は、生育ステージ分類部２０１の処理開始を表す。
【０１４０】
Ｓ２０１０２では、各ほ場の位置における気温等の時系列気象データを、標高ＤＢ２１に格納されているＤＥＭデータ、気象ＤＢ１４に格納されている気象データ、時系列画像ＤＢ２２に格納されている時系列空撮画像より推定する。
【０１４１】
気象ＤＢ１４に格納される気象データの観測場所は、ほ場の近傍とは限らず、対応するほ場の最近傍の観測場所に相当する気象データを用いても、実際のほ場の気象データと異なることが多い。そこで、ＤＥＭデータや時系列空撮画像を用いて、対応するほ場位置の田植日から収穫日までの気象データを推定する。例えば、ｊ番目の気象データＴｊの観測場所に対応するＤＥＭデータの画素値Ｄｊ、ｊ番目の気象データＴｊの観測日、観測場所に対応する空撮画像の画素値Ａｊとから、(数６)という線形回帰式を策定してα、β、γを算出し、全ほ場について先述の式を適用することにより、時系列気象データを推定してもよい。
【０１４２】
【数６】
2015000049.tif 000008
【０１４３】
例えばα、β、γのパラメータを求めるために、最小二乗法等を適用してもよい。
【０１４４】
Ｓ２０１０２の処理は、気温変動推定２０１１によって行われる。
Ｓ２０１０３では、各ほ場の作付日や出穂日等から収量推計に用いる空撮画像の撮影日までの気象データ積算値を算出する。
Ｓ２０１０３の処理は、積算気温推定部２０１２によって行われる。
Ｓ２０１０４では、Ｓ２０１０３で算出された気象データ積算値を用いて、各ほ場の生育ステージを推定する。推定方法は、Ｓ６２に記述されている方法でもよい。その後、ほ場生育ステージＤＢ２０１４にほ場ＩＤと推定生育ステージを格納する。
Ｓ２０１０５では、生育ステージに関して優先度を算出する。但し、存在するパラメータは生育ステージのみなので、生育ステージの取りうる値を見て、使うか否かだけを決定すればよい。例えば、ほ場生育ステージＤＢ２０１４を見て、生育ステージが１種類しかない場合は、全く生育ステージにバラつきがないため、優先度を０とし、それ以外の場合は優先度１としてもよい。
Ｓ２０１０４、Ｓ２０１０５の処理は生育ステージ推定部２０１３によって行われる。
Ｓ２０１０６は、生育ステージ分類部２０１の処理終了を表す。
【０１４５】
図１８は、ほ場生育段階データ２０１４のデータ構造例である。
ほ場生育段階データ２０１４は、収量推計用に取得された空撮画像が撮影された時に関して、各ほ場の生育ステージを推定した結果を格納する。
ほ場ＩＤ２０１４１は、収量推計対象であるほ場を識別するためのラベルであり、ほ場ＩＤ１３０１と一致する。
生育ステージ２０１４２は、対応するほ場の生育ステージを表す。
【０１４６】
以上の実施の形態例において説明した内容によれば、調査ほ場数をできる限り抑制しながら、空撮画像と農作物の特定の生育ステージにおける時系列気象データを用いて、収量を高精度に推定することが可能となる。
【符号の説明】
【０１４７】
１０  収量予測装置
１１  過去データ
１２  画像ＤＢ
１３  ほ場ＧＩＳ
１４  気象ＤＢ
１５  現地調査員
１７  収量推計結果
２１  標高ＤＢ
２２  時系列画像ＤＢ
１０１  パラメータ優先度算出部
１０２  画像解析部
１０３  実測ほ場選択部
１０５  ほ場特徴データ
１０６  実測ほ場表示部
１０７  実測データ入力部
１０８  実測ほ場データ
１０９  非実測ほ場データ
１１０  収量推計部
１１１  収量推計結果表示部
２０１  生育ステージ分類部
１０１１  統計解析部
１０１２  優先度算出部
１０１３  統計データ
１０２１  ほ場内画素抽出部
１０２２  ほ場内特徴算出部
１０３１  ほ場選定軸決定部
１０３２  ヒストグラム作成部
１０３３  実測候補選定部
１０３４  実測コスト算出部
１０３５  実測ほ場決定部
１１０１  生育ステージ推定部
１１０２  時系列パターン解析部
１１０３  モデル策定部
１１０４  収量推計部
１１０５  生育ＤＢ
１９０１  操作部
１９０２  表示部
１９０３  プロセッサ
１９０４  メインメモリ
１９０５  記憶装置
２０１１  気温変動推定部
２０１２  積算気温推定部
２０１３  生育ステージ推定部
２０１４  ほ場生育ステージＤＢ
【図１】
2015000049.tif 000009
【図２】
2015000049.tif 000010
【図３】
2015000049.tif 000011
【図４】
2015000049.tif 000012
【図５】
2015000049.tif 000013
【図６】
2015000049.tif 000014
【図７】
2015000049.tif 000015
【図８】
2015000049.tif 000016
【図９】
2015000049.tif 000017
【図１０】
2015000049.tif 000018
【図１１】
2015000049.tif 000019
【図１２】
2015000049.tif 000020
【図１３】
2015000049.tif 000021
【図１４】
2015000049.tif 000022
【図１５】
2015000049.tif 000023
【図１６】
2015000049.tif 000024
【図１７】
2015000049.tif 000025
【図１８】
2015000049.tif 000026
【図１９】
2015000049.tif 000027
【図２０】
2015000049.tif 000028