JP2019144607A

JP2019144607A - 衛星画像を用いた樹種推定方法、および、樹種推定された樹木の健全度判定方法

Info

Publication number: JP2019144607A
Application number: JP2018025376A
Authority: JP
Inventors: 靖博松田; Yasuhiro Matsuda; 謙治橋爪; Kenji Hashizume; 詳悟林; Shogo Hayashi
Original assignee: West Nippon Expressway Engineering Shikoku Co Ltd; West Nippon Expressway Co Ltd
Current assignee: West Nippon Expressway Engineering Shikoku Co Ltd; West Nippon Expressway Co Ltd
Priority date: 2018-02-15
Filing date: 2018-02-15
Publication date: 2019-08-29
Anticipated expiration: 2038-02-15
Also published as: JP7046432B2

Abstract

【課題】低コストで、簡単かつ確実に、樹種を推定できるようにする。【解決手段】衛星を介し、解析対象の地表を撮影し、樹木を含むパンクロマチック画像を生成し、生成されたパンクロマチック画像から、画像特徴から樹木を自動抽出し、円で樹木を囲み、抽出円の内側領域の同時生起行列を算出し、円中心のマルチスペクトル正規化データを取得し、正規化処理を実行し、教師データと抽出された樹木を付き合わせ、多変量解析モデルにより、樹種推定モデルを作成し、解析モデルにより、解析対象の地表のパンクロマチック画像から抽出された樹木の樹種を推定する。【選択図】図１

Description

衛星画像の解析により、樹木の種類を推定する方法、および、樹種推定された樹木の健全度判定方法に関する。

山林や森林を構成する樹木の種類や、樹木の活性度を調べる場合、一般的に、現地に出向いて作業員が確認する必要がある。

しかしながら、作業員が容易に近寄ることができない区域や、調査対象区域が広大な場合、コストと時間が膨大にかさみ、現実的ではない。
このため、低コストで、簡単かつ確実に、樹種を推定し、樹木の健全度を判定する方法が望まれている。

特開2014-100099

本発明は、上記の問題を解決するためになされたものであって、低コストで、簡単かつ確実に、樹種を推定し、樹木の活性度を判定できるようにすることを目的とする。

本発明にかかる衛星画像を用いた樹種推定方法、および、樹種推定された樹木の健全度判定方法は、
衛星画像を用いた樹種推定方法であって、
衛星を介し、解析対象の地表を撮影し、樹木を含むパンクロマチック画像を生成する撮影ステップと、
生成されたパンクロマチック画像から、画像特徴から樹木を自動抽出し、円で樹木を囲む樹木抽出ステップと、
抽出円の内側領域の同時生起行列を算出し、円中心のマルチスペクトル正規化データを取得し、正規化処理を実行する解析データ整理ステップと、
教師データと抽出された樹木を付き合わせ、多変量解析モデルにより、樹種推定モデルを作成する解析モデル作成ステップと、
解析モデルにより、解析対象の地表のパンクロマチック画像から抽出された樹木の樹種を推定する衛星画像解析ステップと、
を含むことを特徴とする。

本発明にかかる衛星画像を用いた樹種推定方法、および、樹種推定された樹木の健全度判定方法は、前記の樹種推定方法において、
教師データが、
衛星画像中に含まれる樹木に対応する位置に植生する樹木の種類を現地で確認し、
確認された樹木の種類と、確認された樹木に対応する、衛星画像中に含まれる樹木のマルチスペクトル正規化データとを関連付け、
樹種に対応したマルチスペクトル正規化データをデータベース化することにより作成される
ことを特徴とする。

本発明にかかる衛星画像を用いた樹種推定方法、および、樹種推定された樹木の健全度判定方法は、前記の樹種推定方法において、
樹木抽出ステップにおいて、画像中に、ブロッコリー状の画像特徴が認められたときに、これを樹木として自動抽出することを特徴とする。

本発明にかかる衛星画像を用いた樹種推定方法、および、樹種推定された樹木の健全度判定方法は、前記の樹種推定方法において、
解析データ整理ステップが、数値データに対して、撮影時の衛星位置情報、太陽と地球との距離の情報に基づく補正を行い、マルチスペクトル画像の各ＢＡＮＤの反射率を求めるものであって、
解析データ整理ステップにおける正規化処理が、
天頂角θを求めるステップと、
Earth-Sun Distanceを求めるステップと、
地表面からの放射輝度Ｌを求めるステップと、
反射率を求めるステップと、
反射率を正規化処理するステップと
を含み、
各ＢＡＮＤの反射率の計算式は、

であり、
マルチスペクトル正規化データが、マルチスペクトルの反射率の正規化処理であり、反射輝度を、最小値が０、最大値が１になるように正規化することを示すことを特徴とする。

本発明では、衛星画像の解析により、衛星画像中の樹木の種類を簡単かつ確実に推定でき、さらには、樹種が推定された樹木の活性度を判定できるようになる。

図１は、本発明にかかる衛星画像を用いた樹種推定方法の一実施例の流れを示すフロー図である。図２は、図１の樹種推定方法における解析に用いた情報の一覧である。図３は、図１の樹木の抽出に用いられる画像特徴量を抽出した項目の一覧である。図４は、図１の解析データの整理において実行されるマルチスペクトル画像の取得の様子を示す概念図である。図５は、図１の解析データの整理において必要な天頂角の模式図である。図６は、図１の解析データの整理において必要な画像の撮影日時を含む情報の一例である。図７は、図１の解析データの整理において必要な衛星の各ＢＡＮＤの定数の一覧である。図８は、図１の解析データの整理において必要な反射率の正規化の一例を示すグラフである。図９は、図１の解析モデルの作成の概念図である。図１０は、図１に示した樹種推定方法による教師木の樹種推定結果を示す表である。図１１は、図１０の推定結果を画像化した一例の画像である。図１２は、本発明の衛星画像を用いて樹種推定された樹木の健全度判定方法の流れを示すフロー図である。

以下、図面に基づき、本発明を詳細に説明する。
図１は、本発明にかかる衛星画像を用いた樹種推定方法の流れを示すフロー図、図２は、図１の樹種推定方法における解析に用いた情報の一覧、図３は、図１の樹木の抽出に用いられる画像特徴量を抽出した項目の一覧、図４は、図１の解析データの整理において実行されるマルチスペクトル画像の取得の様子を示す概念図、図５は、図１の解析データの整理において必要な天頂角の模式図、図６は、図１の解析データの整理において必要な画像の撮影日時を含む情報の一例、図７は、図１の解析データの整理において必要な衛星の各ＢＡＮＤの定数の一覧、図８は、図１の解析データの整理において必要な反射率の正規化の一例を示すグラフ、図９は、図１の解析モデルの作成の概念図、図１０は、図１に示した樹種推定方法による教師木の樹種推定結果を示す表、図１１は、図１０の推定結果を画像化した一例の画像、図１２は、本発明の衛星画像を用いて樹種推定された樹木の健全度判定方法の流れを示すフロー図である。

以下、図面に基づき、本発明を実施例を用いて説明する。
まず、本実施例では、衛星を介してパンクロマチック画像を生成する（撮影ステップＳ０１）。

本実施例では、既に、衛星軌道上で運用されている衛星によって撮影された衛星画像を使用した。
具体的には、撮影衛星としては、WorldView-2およびWorldView-3によって撮影された撮影データであって、パンクロマチック画像を使用した。

このパンクロマチック画像を確認すると、色味の情報だけではなく、樹木の形状によって、見え方、即ち、画像の輝度分布が異なる。

そこで、本発明では、高分解能（５０ｃｍ／画素程度）衛星パンクロマチック画像（４５０〜８００ｎｍ）の輝度分布から、同時生起行列（ＧＬＣＭ）を取得し、これを画像特徴量として用いることにより、樹種推定を行った。
解析に用いた情報の一覧は、図２に示したとおりである。

次に、パンクロマチック画像から、ブロッコリー状の画像特徴を検出し、ブロッコリー状の画像特徴を樹木として自動抽出し、自動抽出された樹木を円で囲む（樹木抽出ステップＳ０２）。
このとき、誤検出、未検出、過検出については、手動で修正することが望ましい。

この樹種抽出ステップＳ０２では、樹木周りの影の領域で樹木が分割されることに注目して、この影を基準に樹木の抽出を行った。
抽出では、画像の解像度が高いほど検出結果が優れる傾向が認められた。
なお、樹木自体が影に埋もれている場合、抽出が困難であり、また、樹木以外の部分も過検出されることもあった。

次に、抽出円の内側領域の同時生起行列（ＧＬＣＭ）を算出し、円中心のマルチスペクトル正規化データを取得し、正規化処理を実行する（解析データ整理ステップＳ０３）。
同時生起行列（ＧＬＣＭ：Gray-Level Co-occurrence Matrix）とは、画像のテクスチャ解析手法の一つであり、離れた２つの画素対における濃度の配置具合を調べる方法である。

本実施例では、囲み円内の輝度統計量として５項目、ＧＬＣＭから１４項目、合計１９項目の画像特徴量を抽出した。
この特徴量抽出の一覧は図３に示したとおりである。
なお、ＧＬＣＭの算出には、「新編画像解析ハンドブック」１２６０頁から１２６４頁の手法を用いた。

また、この解析データ整理ステップＳ０３では、マルチスペクトルデータからの反射率と、正規化反射率の抽出も実行する。
具体的には、図４に示されるように、可視画像上で、抽出された円から、円の中心座標位置を取得し、マルチスペクトル画像の同一座標から、マルチスペクトル画像の各ＢＡＮＤの１１ｂｉｔの数値データを取得した。

そして、解析データ整理ステップＳ０３では、数値データに対して、撮影時の衛星位置情報、太陽と地球との距離の情報に基づく補正を行い、各ＢＡＮＤの反射率を求めた。
この反射率の計算手順は、以下のとおりである。

１）天頂角度θを求める
２）Earth-Sun Distance（dES）を求める。
３）地表面からの放射輝度Ｌを求める。
４）反射率を求める。
５）反射率を正規化処理する。

また、各ＢＡＮＤの反射率の計算式には以下の数式（１）を用いた。

以下、上記反射率の計算の手順を具体的に説明する。
まず、天頂角θ_Ｓは、撮影日時から求める。
この天頂角の概念は、図５に示すとおりである。

次に、Earth-Sun Distance（ｄ_ＥＳ）の求め方を説明する。
衛星情報には、撮影日時のデータが記録されている。
この情報から、地球と太陽との距離を求める。

例えば、図６に示される情報、具体的には、衛星情報のfirstLineTimeによると、データは２０１６年３月２１日の２時１２分２．３１０２２２秒に撮影されていることがわかる。
即ち、firstLineTime = 2016-03-21T02:12:02.310222Zである。

この情報を元にｄ_ＥＳを求める。
具体的には、上記の場合は、以下の数式（２）となる。

次に、地表面からの放射輝度Ｌの求め方を説明する。
衛星の各ＢＡＮＤの定数（GAIN、absFac、eBw、Offset）と太陽光の成分（Esun）から、画像の１１ｂｉｔデータ（ＤＮ）を変換し、放射輝度Ｌを求める。
各ＢＡＮＤの定数は、図７に示したとおりである。
放射輝度Ｌは、以下の数式（３）によって求められる。

次に、各ＢＡＮＤの反射率について説明する。
各バンドの放射輝度Ｌに、天頂角（θ_Ｓ）と地球と太陽の距離（ｄ_ＥＳ）から求められる反射率係数を乗じて、反射率を得る。

次に、マルチスペクトル反射率の正規化処理について説明する。
日陰の影響など、日射の当たり方の影響をできるだけ除去するため、日射の当たり方の影響を受ける反射輝度を、最小値が０、最大値が１になるように正規化した。
この正規化の一例は、図８に示したとおりである。

上記の解析データ整理ステップＳ０３続き、解析モデルを作成する（解析モデル作成ステップＳ０４）。
この解析モデル作成ステップＳ０４は、教師データと抽出された樹木を付き合わせ、多変量解析モデルにより、樹種推定モデルを作成することによって実行される。

具体的には、解析は、７８０本の樹木を教師データとして解析を実行した。
解析手法として、ニューラルネットワークモデル、決定木分析など、数種類の多変量解析を行い、その結果、決定木分析が最も相関が高かったため、決定木分析を分析手法として採用した。
従属変数である樹種分類については、作業員が確認し易い地域で行い、検出を実行した。
この解析モデルの作成の概念は、図９に示したとおりである。

上記の解析モデル作成ステップＳ０４の後、その解析モデルを用い、解析対象の地表の衛星画像（パンクロマチック画像）から抽出された樹木の樹種を推定する（衛星画像ステップＳ０５）。

この生成画像ステップＳ０５を実際に行った試験結果を以下に示す。
ここでは、教師データとした全７８０本について、精度検証を実施した。
即ち、教師木となった樹木に対し、上記の各ステップを実行し、樹種の推定を行った。

この教師木の樹種の推定結果は、図１０に示したとおりである。
図１０の表のうち、アスタリスクがマークされたテーブルが正解であり、その本数は、７８０本中、７１６本であった。

この結果から、樹種の判定精度は、９２％である。
この推定結果を画像化した画像は、図１１に示したとおりである。
この画像では、樹種の推定結果を、円の囲みの近傍にローマ字表示した。

上記の試験の結果から、上記の本発明に係る衛星画像を用いた樹種推定方法は、実用に耐えうる程度の高精度に、衛星画像から樹種を簡単かつ確実に推定できる。

次に、上記の衛星画像を用いた樹種推定方法によって推定された樹種と、ＮＤＶＩ(Normalized Difference Vegetation Index：正規化植生指標)とを組み合わせた樹木の健全度判定方法について説明する。

ＮＤＶＩとは、植物による光の反射特性に基づき、衛星画像を用いて植生の状況を判断する指標の１種であり、植物の活力を表す従来公知の指標である。
ここでは、ＮＤＶＩの詳細な説明は省略する。

まず、樹種に応じた健全度について考察する。
まず、クズのような樹木を被圧する植物に覆われた樹木については、ＮＤＶＩ指標に基づく活性度の推定が困難であるため、樹木の健全度の判定が困難である、と考えられる。
このため、樹種がクズ、または、これに類する樹木を被圧する植物の場合、その樹木は、「現状確認」と判定し、現場での目視、または、経過観察とする。
樹木の健全度と、ＮＤＶＩの活性度の相関性であるが、ＮＤＶＩが０．６０超の樹木は、健全な樹木であることが確認されたため、ＮＤＶＩが０．０６超の樹木は「健全」と判定する。

ＮＤＶＩが０．５５超の樹木は、健全に近い状態の樹木であることが確認されたため、ＮＤＶＩが０．５５超の樹木は、「観察」と判定する。
ＮＤＶＩが０．４５以下の樹木は、健全度が著しく損なわれているため、「要注意」と判定する。

ＮＤＶＩが０．４５超であって、その樹木の樹種が、例えばアラカシのように、ナラ枯れ樹種の場合、その樹木は「要注意」と判定する。
これは、ナラ枯れは、進行性倒木のリスクが大きいためである。

ナラ枯れ樹種ではないと判定された樹木のうち、その樹木が、マツ枯れ樹種と判定された樹木を「注意」と判定する。
これは、マツ枯れは、倒木のリスクが比較的低いためである。

マツ枯れ樹種ではないと判定されなかった樹木のうち、その樹木が、例えばアカメガシワやサクラのような侵入木と判定された樹木を「注意」と判定する。
これは、侵入木は、育成が早い樹種であり、サクラなどのように幹が柔らかく、虫害を受け易い樹種が多いためである。

侵入木と判定されなかった樹木は、その他の樹木とし、「観察」と判定する。

上記の樹木の健全度の判定方法を整理すると、図１２に示したようになる。
ここでは、樹木の健全度を示す項目を「健全」、「現状確認」、「観察」、「注意」および「要注意」とする。

まず、樹種がクズのような樹木を被圧するような植物（単にクズという）のとき、クズに被圧された樹木は「観察」と判定する。
次に、樹種がクズでなく、ＮＤＶＩが０．６０超の場合、「健全」と判定する。

次に、ＮＤＶＩが０．６０以下の樹木のうち、ＮＤＶＩが０．５５超の場合、「観察」と判定する。

次に、ＮＤＶＩが０．５５以下の樹木のうち、ＮＤＶＩが０．４５以下の場合、その樹木は全て「要注意」と判定する。

次に、ＮＤＶＩが０．４５超であって、ナラ枯れ樹種の場合、「要注意」と判定する。

次に、ナラ枯れ樹種とは判定されなかった樹木のうち、マツ枯れ樹種と判定された樹木を「注意」と判定する。

次に、マツ枯れ樹種と判定されなかった樹木のうち、侵入木と判定された樹木を「注意」と判定する。

次に、侵入木と判定されなかった樹木を、その他の樹種として「観察」と判定する。

ここで、「要注意」と判定された場合は、例えば、直ちに伐採除去等の措置をとる。

また、「注意」と判定された場合は、例えば、時間の経過に伴い、「要注意」となる可能性がある樹木であることを、データベース等に登録しておく。

また、「観察」と判定された場合、例えば、「要注意」に変化する可能性がある樹木であることを、データベース等に登録しておく。

また、「健全」と判定された場合、例えば、特段の措置を講じる必要はない。

なお、上記の実施例では、解析対象は、衛星によって撮影された衛星画像を用いたが、この衛星画像は、上空から空撮された画像であればどのようなものであってもよく、例えば、衛星画像ではなく、例えば、無人航空機を用いて撮影された航空写真であってもよい。

本発明は、教師データさえ作成すれば、衛星から撮影し得る区域であれば、樹木の種類や、樹木の健全度を、作業員の現地調査なしに、簡単かつ確実に把握できるようになり、森林や山林の管理における負担を軽減させることができる。

Claims

衛星画像を用いた樹種推定方法、および、樹種推定された樹木の健全度判定方法であって、
衛星を介し、解析対象の地表を撮影し、樹木を含むパンクロマチック画像を生成する撮影ステップと、
生成されたパンクロマチック画像から、画像特徴から樹木を自動抽出し、円で樹木を囲む樹木抽出ステップと、
抽出円の内側領域の同時生起行列を算出し、円中心のマルチスペクトル正規化データを取得し、正規化処理を実行する解析データ整理ステップと、
教師データと抽出された樹木を付き合わせ、多変量解析モデルにより、樹種推定モデルを作成する解析モデル作成ステップと、
解析モデルにより、解析対象の地表のパンクロマチック画像から抽出された樹木の樹種を推定する衛星画像解析ステップと、
を含む上記の衛星画像を用いた樹種推定方法、および、樹種推定された樹木の健全度判定方法。
教師データが、
衛星画像中に含まれる樹木に対応する位置に植生する樹木の種類を現地で確認し、
確認された樹木の種類と、確認された樹木に対応する、衛星画像中に含まれる樹木のマルチスペクトル正規化データとを関連付け、
樹種に対応したマルチスペクトル正規化データをデータベース化することにより作成される
請求項１に記載の衛星画像を用いた樹種推定方法、および、樹種推定された樹木の健全度判定方法。
樹木抽出ステップにおいて、画像中に、ブロッコリー状の画像特徴が認められたときに、これを樹木として自動抽出する請求項１または２に記載の衛星画像を用いた樹種推定方法、および、樹種推定された樹木の健全度判定方法。
解析データ整理ステップが、数値データに対して、撮影時の衛星位置情報、太陽と地球との距離の情報に基づく補正を行い、マルチスペクトル画像の各ＢＡＮＤの反射率を求めるものであって、
解析データ整理ステップにおける正規化処理が、
天頂角θを求めるステップと、
Earth-Sun Distanceを求めるステップと、
地表面からの放射輝度Ｌを求めるステップと、
反射率を求めるステップと、
反射率を正規化処理するステップと
を含み、
各ＢＡＮＤの反射率の計算式は
であり、
マルチスペクトル正規化データが、マルチスペクトルの反射率の正規化処理であり、反射輝度を、最小値が０、最大値が１になるように正規化することを示す
請求項１から３のいずれかに記載の衛星画像を用いた樹種推定方法、および、樹種推定された樹木の健全度判定方法。
衛星画像の解析により推定された樹種と、ＮＤＶＩ（正規化植生指標）を組み合わせ、樹木の健全度を樹種に応じて、樹木の健全度を判定する衛星画像を用いた樹種推定方法、および、樹種推定された樹木の健全度判定方法。
樹木の健全度を示す項目を「健全」、「現状確認」、「観察」、「注意」および「要注意」とし、
樹種がクズのとき、クズに被圧された樹木は「観察」と判定し、
樹種がクズでなく、ＮＤＶＩが０．６０超の場合、「健全」と判定し、
ＮＤＶＩが０．６０以下の樹木のうち、ＮＤＶＩが０．５５超の場合、「観察」と判定し、
ＮＤＶＩが０．５５以下の樹木のうち、ＮＤＶＩが０．４５以下の場合、その樹木は全て「要注意」と判定し、
ＮＤＶＩが０．４５超であって、ナラ枯れ樹種の場合、「要注意」と判定し、
ナラ枯れ樹種とは判定されなかった樹木のうち、マツ枯れ樹種と判定された樹木を「注意」と判定し、
マツ枯れ樹種と判定されなかった樹木のうち、侵入木と判定された樹木を「注意」と判定し、
侵入木と判定されなかった樹木を、その他の樹種として「観察」と判定する
請求項５に記載の衛星画像を用いた樹種推定方法、および、樹種推定された樹木の健全度判定方法。
衛星画像が、無人航空機を用いて撮影された撮影画像である請求項１から６のいずれかに記載の衛星画像を用いた樹種推定方法、および、樹種推定された樹木の健全度判定方法。