中間赤外域バンドを利用した森林の樹種分類
【課題】より高精度の森林植生分類図を作成する方法の提供する。
【解決手段】衛星画像を用いて植生分類図を作成する際に用いられる光および／または電磁波の種類（波長）を適切に選択することによって植生分類図の精度が向上すること、とくに中間赤外波を含む電磁波による衛星画像を用いることを含めて森林の植生分類図を作成する。
【解決手段】衛星画像を用いて植生分類図を作成する際に用いられる光および／または電磁波の種類（波長）を適切に選択することによって植生分類図の精度が向上すること、とくに中間赤外波を含む電磁波による衛星画像を用いることを含めて森林の植生分類図を作成する。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、中間赤外域バンドを利用した森林の樹種分類に関し、とくに衛星データによる森林植生分類図の作成方法および森林植生分類図に関する。
【背景技術】
【０００２】
人工造林地が正常に成林して十分な木材収穫が行われることは、林業経営にとって不可欠である。一方、人工造林地のモニタリングを行うことも、収穫予測の策定するために不可欠であるところ、同モニタリングが正確に行えず、正常に成林しているか否かが不明な林分が多い場合には伐採時期における精度の高い伐採可能量の予測を立てることができない。
【０００３】
また、森林施業計画には伐採可能量の算出が不可欠であるところ、将来的な伐採可能量は森林調査簿のデータを基に算出しているため、同調査簿の精度がきわめて重要である。
しかしながら、森林調査簿の更新は、間伐等の施業が実施された場合には毎年内容が更新されるが、施業間隔は通常１０年程度であるため、人工林が強風による風倒等の何らかの理由により広葉樹が侵入して針広混交林化し減損している場合には、樹種の把握が不正確になっているおそれもある。その結果、伐採可能量が実際の量より過大なものとして算出されてしまう懸念がある。
したがって、上記モニタリングにより伐採可能量の推定を極力正確に行うことが、林業経営等においては極めて重要である。
【０００４】
また、地球温暖化による各種の弊害が全世界的な規模で問題視される一方、森林・樹木による炭酸ガスの同化作用による地球温暖化の防止が重要であることは論を待たない。そのため、森林・樹木による炭酸ガス同化作用の算定の基礎となるバイオマスを増加させるために、森林整備の必要性が叫ばれ、森林整備が盛んに行なわれている。
かかる植林後の樹種やすでに森林を構成している樹種の正確な把握も、バイオマスの推定のために不可欠である。したがって、森林における植生分類のより高い精度による把握は、今後日本国内にとどまらず世界的に重要性が増していくことは確実である。
【０００５】
モニタリングの手法として、理論的には全山踏査して毎木調査を行うことが考えられる。しかしながら、かかる調査には膨大な労力、莫大なコスト・時間を要するため、現実的にはほぼ不可能である。そのため、モニタリングの代替手法として、従来から航空写真による目視判読や衛星画像解析による植生分類が行われている。とくに、衛星画像解析によれば、より広範な地域の情報をより効率的に得ることができる。
【０００６】
たとえば、航空カラー写真や衛星カラー画像をコンピュータに入力し、コンピュータ画像処理により森林地域の樹冠形状画像を求め、この樹冠形状画像と元のカラー画像から求めた樹木の色彩及び輝度等とから森林地域全体の樹木植生の調査・評価、すなわちリモートセンシングが行われている（特許文献１）。
【０００７】
また、特許文献２には、調査対象地域のバンドデータの季節変化を利用しておこなう衛星データによる森林地域の植生分類方法が示されている。
さらに、特許文献３および４には航空機や衛星を利用したリモートセンシングによる画像の利用が示され、非特許文献１および２には、高分解能ＩＫＯＮＯＳ衛星による樹種分類の方法が記載されている。また、非特許文献３には、空中写真を用いた森林地理情報システム（ＧＩＳ）データベースによる森林ゾーニングについて報告されている。
【先行技術文献】
【特許文献】
【０００８】
【特許文献１】特開２００１−３５７３８０号公報
【特許文献２】特開２００６−８５５１７号公報
【特許文献３】再表２００３−０６９３１５号公報
【特許文献４】特開２００８−７９５４９号公報
【非特許文献】
【０００９】
【非特許文献１】「高分解能ＩＫＯＮＯＳ衛星による針広混交林の樹種分類」加藤正人，森林航測第１９８号６〜９（２００２．１１）
【非特許文献２】「トレーニングデータの選定方法がＩＫＯＮＯＳ画像の樹種分類に与える影響」宮地ら，第１１５回日林学術講，１９０頁
【非特許文献３】「朝日の森を対象とするＧＩＳを用いた森林機能評価とゾーニング」呉ら，東京大学農学部演習林報告，１１１，５９−８３（２００４）
【発明の概要】
【発明が解決しようとする課題】
【００１０】
上記のとおり森林における植生分類をより高い精度に行うことは重要であるにも拘らず、従来の手法は、必ずしも高い精度で樹種を分類し、高精度の植生図の作成を可能としたものとはいえない。たとえば、非特許文献２においては、トレーニングデータの選定方法に依存して分類することができる可能性があるスギおよびヒノキの分類精度は、それぞれ７０％および７８％程度に留まっている。
【００１１】
したがって、本発明は、より高精度の森林植生分類図を作成する方法を提供することを目的とする。
【課題を解決するための手段】
【００１２】
従来の衛星画像のセンサの多くは、可視域光（青、緑、赤）に加えて、植生による反射ピークが存在する近赤外の波長帯を用いている。しかしながら、可視光域に対する植生による反射量は小さいため、樹種ごとの差異を検出しがたいと考えられた。また、近赤外波を併せて用いても、得られる植生の分類精度は、前記のとおり必ずしも十分なものではなかった。
そこで、本発明者らは、衛星画像を用いて植生分類図を作成する際に用いられる光および／または電磁波の種類（波長）を適切に選択することによって植生分類図の精度が向上すること、とくに中間赤外波を用いることによって上記課題が解決することを見出し、さらに鋭意研究を進めて本発明を完成した。
すなわち、本発明は、中間赤外波を含む電磁波による衛星画像を用いることを含む、森林の植生分類図を作成する方法に関する。
【００１３】
また、本発明は、中間赤外波の波長が１．３μｍ〜２．０μｍである、記方法に関する。
さらに本発明は、さらに近赤外波および可視光を用いる前記方法に関する。
また、本発明は、衛星画像から、２本以上の単木を含むオブジェクトを作成する工程を含む、前記方法に関する。
【００１４】
さらに、本発明は、分類パラメータを用いてオブジェクトの土地被覆分類を行う工程を含み、該分類パラメータが中間赤外波に対する反射輝度値を用いるパラメータを含む前記方法に関する。
また、本発明は、分類パラメータが正規化水指数ＮＤＷＩ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲおよびＲａｔｉｏ＿ＳＷＩＲからの１種または２種以上である前記方法に関する。
【００１５】
また、本発明は、階層的分類法によるオブジェクトの土地被覆分類を行う工程を含む前記方法に関する。
そして本発明は、前記いずれかの方法を用いて作成された森林の植生分類図に関する。
【発明の効果】
【００１６】
本発明の方法によれば、森林の植生分類図を高精度で作成することができる。
本発明の方法のうち中間赤外波長の波長が１．３μｍ〜２．０μｍであるものによれば、森林の植生分類図をより高精度で作成することができる。
本発明の方法のうちさらに近赤外波および可視光を用いるものによれば、森林の植生分類図をさらに高精度で作成することができる。
本発明の方法のうち近赤外波の波長が０．７５μｍ〜１．０μｍであるものまたは可視光が緑色光および赤色光であるものによれば、森林の植生分類図をより一層高精度で作成することができる。
本発明の方法のうち衛星画像から、２本以上の単木を含むオブジェクトを作成する工程を含むものによれば、高分解能（分解能１ｍ以下）を具備していない衛星による衛星画像によっても高精度の分類が可能となり、植生分類図の作成をより低コストで行うことができる。
【００１７】
本発明の方法のうち分類パラメータを用いてオブジェクトの土地被覆分類を行う工程を含むものによれば、森林の植生分類図をより一層高精度で作成することができる。
本発明の方法のうち分類パラメータが中間赤外波に対する反射輝度値を用いるパラメータを含むもの、とくに分類パラメータが正規化水指数ＮＤＷＩ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲおよびＲａｔｉｏ＿ＳＷＩＲからの１種または２種以上であるもの、または分類パラメータが正規化植生指数ＮＤＶＩ、正規化水指数ＮＤＷＩ、色相Ｈｕｅ、彩度Ｓａｔｕｒａｔｉｏｎ、明度Ｉｎｔｅｎｓｉｔｙ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲ、Ｒａｔｉｏ＿ＳＷＩＲからなるものによれば、森林の植生分類図をさらに高精度で作成することができる。
本発明の方法のうち階層的分類法によるオブジェクトの土地被覆分類を行う工程を含むによれば、森林の植生分類図としてより妥当なものを作成することができる。
本発明の方法のうちスギおよび／またはヒノキの分類精度が８０％以上であるものによれば、スギやヒノキ等の針葉樹についてとくに分類精度が高い森林の植生分類図を作成することができる。
【００１８】
本発明による森林の植生分類図によれば、収穫対象林分と減損林分を正確に把握することができる。
【００１９】
植生による反射ピークは近赤外域に存在するのに対し、中間赤外波は必ずしも植生による反射は大きくなく、少なくとも近赤外波に対する反射より小さい傾向がある。したがって、本発明において中間赤外波を用いてより高精度の森林植生分類図の作成することは、本発明の特徴の一つである。
【図面の簡単な説明】
【００２０】
【図１】ＳＰＯＴ−５による提供データ範囲および切り出し範囲を示す図である。
【図２】スケールパラメータを２５としてオブジェクトを生成した図である。
【図３Ａ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｂ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｃ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｄ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｅ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｆ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｇ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｈ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｉ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｊ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｋ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｌ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｍ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｎ】各メンバーシップ関数および分類項目との関係を示す図である。
【図４】階層的分類法による土地被覆分類の結果を示す図である。
【図５】非階層的分類法による土地被覆分類の結果を示す図である。
【発明を実施するための形態】
【００２１】
（定義）
本明細書において「中間赤外（波）」（「中間赤外域バンド波」）とは、一般的な近赤外波と遠赤外波との中間の波長を有する赤外波を意味する。典型的には、「中間赤外（波）」とは、約１．２μｍ〜４μｍの波長を有する電磁波である。また、青色光、緑色光、赤色光および近赤外波は、それぞれ約０．４μｍ〜約０．５μｍ、約０．５μｍ〜約０．６μｍ、約０．６μｍ〜約０．７μｍおよび０．７μｍ〜１．２μｍの波長を有する。
【００２２】
また、それぞれの波長の境界付近の波長を有する可視光または電磁波は、当該境界を構成する可視光または電磁波のいずれにも分類し得る。たとえば、０．４５μｍ〜０．５２μｍのものは、青色光または緑色光とも称され得る。
【００２３】
本発明の方法のうち、下記（１）〜（４）の各工程を含むものは好ましい。これらの工程をすべて含むものはより好ましい。
（１）衛星画像データおよび地形データの前処理
（２）オブジェクトの作成
（３）グラウンドトゥルース
（４）オブジェクトの土地被覆分類
上記（１）〜（４）の各工程につき、以下にさらに詳細に説明する。
【００２４】
（１）衛星画像データおよび地形データの前処理
解析の対象となる衛星画像データおよび地形データのオリジナルデータを解析により適した形式に変換する工程である。
当該工程は、中間赤外波を含む電磁波での衛星画像を用いるものであれば、他の工程・手法等はとくに限定されない。すなわち、当該衛星画像データは、中間赤外波を含む電磁波での衛星画像を用いたデータであれば、本発明に用いることができる。
【００２５】
本発明の方法における中間赤外波の波長は約１．２μｍ〜４μｍであるところ、当該波長が１．３μｍ〜２．０μｍであるものは好ましい。また、中間赤外波の波長は、より好ましくは１．５μｍ〜１．８μｍであり、最も好ましくは１．５８μｍ〜１．７５μｍである。
【００２６】
本発明の方法は、さらに近赤外波ならびに可視光を用いることができる。近赤外波の波長は０．７μｍ〜１．２μｍであり、好ましくは０．７５μｍ〜１．０μｍであり、より好ましくは０．７８μｍ〜０．８９μｍである。
【００２７】
可視光としては、青色光、緑色光および赤色光からの１種または２種以上が好適であり、緑色光および赤色光を用いたものはより好ましい。
【００２８】
（２）オブジェクトの作成
使用画像の解像度や分類対象物の大きさに応じて適切なスケールパラメータを決定し、画像を解析の単位であるオブジェクト（ピクセルの集合体）に区分する工程である。
すでに衛星画像の分解能は高くなり（１ｍ以下）、単木単位での樹種判読は高い精度で分類が可能となっている。しかしながら、高分解能の衛星画像は高価であるため、広域な森林を調査するには莫大な費用が掛かるという欠点がある。また、ピクセル単位での画像分類は、与えられる情報量は豊富であるが、たとえば樹冠の陰影のような過剰な情報が含まれることがあり、当該目的にはむしろ適さない面もある。
【００２９】
また、ピクセル単位での画像分類においては、解像度が高いことにより単木樹冠が複数のピクセルに分割されることにより、ピクセルベースの分類では複数のクラスに誤分類されるおそれが高く、林分のまとまりを性格に把握することが困難である。
【００３０】
これに対し、オブジェクトベースの分類においては、同様なテクスチャや分光反射特性の領域を単位として認識分類を行うため、林相区分図作成により好適である。
また、ある程度のまとまりをもったピクセルの集合体、すなわちオブジェクトを画像分類の単位とする手法が開発され、新規のアプローチが可能となっている。そのため、本発明においては衛星画像から、２本以上の単木を含むオブジェクトを作成する工程を含むものは好ましいのである。
しかしながら、オブジェクトを作成する工程を含む方法においては、樹種ごとに適切な（分類）パラメータが明らかではない。したがって、分類パラメータの新たな検討・決定が、当該工程を含む方法においては重要である。当該分類パラメータの検討・決定についてはさらに後記する。
【００３１】
（３）グラウンドトゥルース
現地状況を実地に調査する工程である。
グラウンドトゥルースによって得られた現地状況と画像の判読結果からトレーニングエリアが画像から選択され、さらに当該トレーニングエリアは、分類のために用いられる分類用セットと精度評価に用いられる評価用セットとにランダムに振り分けられる。
作成された植生分類図および当該分類の精度評価は、グラウンドトゥルースによって得た現地状況からトレーニングエリアを選択し、当該トレーニングエリアにおける比較対照することによって行うことができる。
したがって、本発明の方法のうち、グラウンドトゥルースデータを得る工程を含むものは好ましい。
【００３２】
（４）オブジェクトの土地被覆分類
分類パラメータを用いる手法等により、作成されたオブジェクトごとに土地被覆分類を決定し、植生分類図を完成する工程である。
本発明のうちオブジェクトを作成する工程を含む方法においては、樹種ごとに適切な分類パラメータの新たな検討・決定を行うことが重要である。したがって、本発明の方法のうち、分類パラメータを用いて土地被覆分類を行う工程を含むものは好ましい。
【００３３】
分類パラメータを用いる手法には、例えば階層的または非階層的分類法を例示することができる。
階層的分類法は、大分類から小分類へと段階的に、各段階でそれぞれ異なったメンバーシップ関数を用いて分類する手法である。階層的分類法は、より具体的には、
・各種分類パラメータの分類クラスごとの平均値、最大値、最小値、標準偏差等のグラフおよび任意の２つの分類パラメータを組み合わせた散布図等を作成し、２以上の土地被覆項目の特性を検討し、さらに、
・上記検討の結果から階層構造およびメンバーシップ関数を設定し、画像全体を分類する方法である。したがって、本明細書においては、階層的分類（法）を「メンバーシップ関数を用いた分類（法）」ということがある。
【００３４】
非階層的分類法は、あらかじめ設定した複数の土地被覆分類項目に、同じ分類パラメータを用いて一度に分類する手法である。
これらの手法のうち、階層的分類法によるオブジェクトの土地被覆分類を行う工程を含む本発明の方法は、分類がより好適に行われるため好ましい。
【００３５】
分類パラメータとしては、正規化植生指数ＮＤＶＩ、正規化水指数ＮＤＷＩ、色相Ｈｕｅ、彩度Ｓａｔｕｒａｔｉｏｎ、明度Ｉｎｔｅｎｓｉｔｙ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲ、Ｒａｔｉｏ＿ＳＷＩＲ等が例示される。
【００３６】
上記各分類パラメータのうち、正規化植生指数ＮＤＶＩ、正規化水指数ＮＤＷＩ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲ、Ｒａｔｉｏ＿ＳＷＩＲは、それぞれ下記のように定義される：
ＮＤＶＩ＝（Ｂａｎｄ３−Ｂａｎｄ２）／（Ｂａｎｄ３＋Ｂａｎｄ２）＊１００＋１００
ＮＤＷＩ＝（Ｂａｎｄ３−Ｂａｎｄ４）／（Ｂａｎｄ３＋Ｂａｎｄ４）＊１００＋１００
ＮＤＩ＿ＳＷＲ＝（Ｂａｎｄ４−Ｂａｎｄ２）／（Ｂａｎｄ４＋Ｂａｎｄ２）＊１００＋１００
Ｒａｔｉｏ＿ＮＩＲ＝Ｂａｎｄ３／（Ｂａｎｄ１＋Ｂａｎｄ２＋Ｂａｎｄ３＋Ｂａｎｄ４）＊１００
Ｒａｔｉｏ＿ＳＷＩＲ＝Ｂａｎｄ４／（Ｂａｎｄ１＋Ｂａｎｄ２＋Ｂａｎｄ３＋Ｂａｎｄ４）＊１００
【００３７】
上記定義中のバンド１〜４（Ｂａｎｄ１〜４）の波長は、それぞれ０．５０μｍ〜０．５９μｍ（青）、０．６１μｍ〜０．６８μｍ（赤）、０．７９μｍ〜０．８９μｍ（近赤外波）、１．５８μｍ〜１．７５μｍ（中赤外波）である
また、色相Ｈｕｅ、彩度Ｓａｔｕｒａｔｉｏｎ、明度ＩｎｔｅｎｓｉｔｙはＴＮＴｍｉｐｓのＲＧＢ−ＨＳＩ変換機能（Ｒａｓｔｅｒ／Ｃｏｍｂｉｎｅ／Ｃｏｎｖｅｒｔ Ｃｏｌｏｒ．．．）を利用して作成される。
上記のとおり、各波長の光（電磁波）に対する反射輝度値を用いる５種類のパラメータのうち４種類に中赤外波（Ｂａｎｄ４）に対する反射輝度値を用い得ることは本発明の特長である。これら４種類のパラメータからの１種または２種以上を用いることは好ましい。さらに、分類パラメータが、正規化植生指数ＮＤＶＩ、正規化水指数ＮＤＷＩ、色相Ｈｕｅ、彩度Ｓａｔｕｒａｔｉｏｎ、明度Ｉｎｔｅｎｓｉｔｙ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲ、Ｒａｔｉｏ＿ＳＷＩＲのみからなる本発明の方法はより好ましい。
【００３８】
なお、ＮＤＶＩは近赤外域と可視光（赤色光）域における反射量の差をみるためのパラメータであるところ、当該パラメータは植物のバイオマス量や活力度と正の相関があることが知られている。また、ＮＤＩは植生指数の一種である。Ｒａｔｉｏは全バンド中において特定のバンドが占める割合であって、一般的な変換方法による反射輝度値を用いた指数である。
【００３９】
本発明の方法は、スギおよび／またはヒノキの他の樹種との分類にとくに好適である。たとえば、スギおよび／またはヒノキの分類精度が８０％以上である本発明の方法は好ましい。
【００４０】
前記のとおり、本発明は上記いずれかに記載の方法を用いて作成された森林の植生分類図にも関するところ、たとえば本植生分類図によれば収穫対象林分と減損林分を正確に把握することができるのである。
【００４１】
本発明につき、実施例に基づきより詳細に以下に説明する。かかる実施例はあくまで本発明を説明するためのものであって、いかなる意味においても本発明を限定するものではない。
【実施例】
【００４２】
ＳＰＯＴ５衛星画像の解析により、四国地域における住友林業株式会社社有林（造林地）の植生分類図を作成し、収穫対象林分と減損林分とにゾーニングした。
なお、衛星画像として用いられたＳＰＯＴ５衛星画像の概要は以下のとおりである：
ＳＰＯＴ５衛星画像「ＳＰＯＴＶｉｅｗ Ｏｒｔｈｏ」
（画像中心座標：３３°４９′３０″・１３３°２３′３０″）
・パンクロ画像（２７×２７ｋｍ、解像度２．５ｍ）１シーン
・カラー画像（２７×２７ｋｍ、解像度１０ｍ）１シーン
【００４３】
また、住友林業株式会社内調査簿データ（２００８年度分；ｅｘｃｅｌファイル）および同社内森林ＧＩＳデータ（ｓｈａｐｅファイル・ｔｉｆｆファイル）も適宜用いた。
【００４４】
Ｉ．解析方法
概要は下記（１）〜（６）のとおりである。
（１）使用ソフトウェア
リモートセンシングデータ解析およびＧＩＳ解析にはＴＮＴｍｉｐｓ ｖｅｒ７．０〜７．４（合衆国ＭｉｃｒｏＩｍａｇｅｓ 社）、オブジェクトベース分類にはｅＣｏｇｎｉｔｉｏｎ ｖｅｒ４．０（ドイツＤｅｆｉｎｉｅｎｃｅ Ｉｍａｇｉｎｇ 社）を利用した。二つのソフト間でのデータ交換は、ＴＮＴｍｉｐｓによりラスタデータはＧｅｏＴｉｆｆ変換、ベクタデータはＳｈａｐｅ ｆｉｌｅ変換を行うことにより行った。
（２）前処理
入手したオリジナルデータをＴＮＴｍｉｐｓ形式に変換し、衛星画像データについては大気補正、切り出し、ＮＤＶＩ等画像作成を、地形データについてはモザイク、切り出し、Ｓｈａｄｉｎｇ計算、斜面位置区分などの処理を行った。
（３）ｅＣｏｇｎｉｔｉｏｎオブジェクト作成
ｅＣｏｇｎｉｔｉｏｎでは、スケールパラメータを変えることにより生成されるオブジェクトの大きさが異なるため、使用画像の解像度や分類対象物の大きさに応じて適切なスケールパラメータを設定する必要がある。本実施例においてはスケールパラメータを１５から５０まで変化させたときのオブジェクトサイズを検討し、今回の対象地に適当なスケールパラメータを決定した。
【００４５】
（４）グラウンドトゥルース
２００８年１１月３日〜６日に、別子山・高藪・大永山の各事業区を中心に現地調査を行い、衛星画像と現地状況の対応を確認した。
（５）非階層的および階層的分類法によるオブジェクトの土地被覆分類
土地被覆分類については、非階層的分類および階層的分類の２種類を検討した。前者は、あらかじめ設定した複数の土地被覆分類項目に、同じ分類パラメータを用いて一度に分類する手法である。後者は、大分類から小分類へと段階的に、各段階でそれぞれ異なったメンバーシップ関数を用いて分類する手法である。後者の手法を実行するためには、分類パラメータおよび閾値を適切に選択する必要があるため、ｅＣｏｇｎｉｔｉｏｎで生成したオブジェクトをＴＮＴｍｉｐｓに取り込み、トレーニングデータを画像中から選択して検討を行った。
（６）分類結果の検討
二つの手法で行った分類結果について、精度評価用サンプルエリアにより分類精度を評価した。
【００４６】
実際に行われた解析の詳細は以下のとおりである。
【００４７】
ＩＩ．衛星データの前処理
（１）ＳＰＯＴ−５データの地理補正とリサンプリング、切り出し
提供されたＳＰＯＴ−５データをｍｉｐｓ形式に変換後、空中写真および林小班界と重ね合わせて表示したところ、地理的な位置ずれはほとんど無視できると判断したため、ＧＣＰを用いた地理補正は省略した。この点については、例えばＪＡＸＡ提供のＡＬＯＳデータと比較した場合の利点であり、作業量を増やすことがない。提供データはＵＴＭ座標系に基づいたものであったため、平面直角座標系（ＩＶ系）、２．５ｍ解像度へとＮｅａｒｅｓｔ Ｎｅｉｇｈｂｏｕｒ（ＮＮ）法によりリサンプリングを行った。さらに、画像サイズを必要最小限にするため、五良津・大永山・別子山・高藪の４事業区を含むように切り出した。ただし、大永山事業区の西端部分は提供データ範囲に含まれていなかったため、切り出し後のデータにも含まれていない。また、石鎚山事業区は提供データ範囲に含まれていなかったため、今回の解析には含んでいない。切り出し範囲の座標値（Ｎ，Ｅ）は以下の通りである（単位ｍ）。左上：（１０２，９７０，−２１，６０２．５）、右上：（１０２，９７０，１，１４７．５）、左下：（７９，４７０，−２１，６０２．５）、右下：（７９，４７０，１，１４７．５）。画像サイズは９，１００ピクセル×９，４００ラインである。図１に提供データ範囲を示す。
（２）大気補正
大気補正は最小値法を利用した。提供データ範囲において、ＳＰＯＴ−５データの各バンド（Ｂａｎｄ１〜４、Ｐａｎ）のヒストグラムを用いて最小値を検索した。各バンドの最小値は表−１に示すとおりである。この表の値から１を引いた値を各バンドから差し引いた。最小値から１を引いたのは、大気補正後の衛星画像データＤＮ値が０になるとＮＤＶＩ計算の際に分母が０となり計算できないことがあるためである。
【００４８】
【表１】
【００４９】
なお、バンド１〜４（Ｂａｎｄ１〜４）の波長は、それぞれ０．５０μｍ〜０．５９μｍ（青）、０．６１μｍ〜０．６８μｍ（赤）、０．７９μｍ〜０．８９μｍ（近赤外波）、１．５８μｍ〜１．７５μｍ（遠赤外波）である。また、Ｐａｎはバンド１〜４のセンサ（ＨＲＧ−Ｘ）とは異なるセンサ（ＨＲＶ−Ｐ）によって検知される０．４８μｍ〜０．７１μｍの波長の可視光である。
【００５０】
（３）各種指数画像の作成
分類に使用する指数として、正規化植生指数ＮＤＶＩ、正規化水指数ＮＤＷＩ、色相Ｈｕｅ、彩度Ｓａｔｕｒａｔｉｏｎ、明度Ｉｎｔｅｎｓｉｔｙ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲ、Ｒａｔｉｏ＿ＳＷＩＲを検討の対象とし、それぞれの画像を作成した。Ｈｕｅ、Ｓａｔｕｒａｔｉｏｎ、Ｉｎｔｅｎｓｉｔｙは（Ｒ，Ｇ，Ｂ）を（Ｂａｎｄ４，Ｂａｎｄ３，Ｂａｎｄ２）に割り当て作成した。以下、各指数の計算式を示す。
ＮＤＶＩ＝（Ｂａｎｄ３−Ｂａｎｄ２）／（Ｂａｎｄ３＋Ｂａｎｄ２）＊１００＋１００
ＮＤＷＩ＝（Ｂａｎｄ３−Ｂａｎｄ４）／（Ｂａｎｄ３＋Ｂａｎｄ４）＊１００＋１００
ＮＤＩ＿ＳＷＲ＝（Ｂａｎｄ４−Ｂａｎｄ２）／（Ｂａｎｄ４＋Ｂａｎｄ２）＊１００＋１００
Ｒａｔｉｏ＿ＮＩＲ＝Ｂａｎｄ３／（Ｂａｎｄ１＋Ｂａｎｄ２＋Ｂａｎｄ３＋Ｂａｎｄ４）＊１００
Ｒａｔｉｏ＿ＳＷＩＲ＝Ｂａｎｄ４／（Ｂａｎｄ１＋Ｂａｎｄ２＋Ｂａｎｄ３＋Ｂａｎｄ４）＊１００
Ｈｕｅ、Ｓａｔｕｒａｔｉｏｎ、ＩｎｔｅｎｓｉｔｙはＴＮＴｍｉｐｓのＲＧＢ−ＨＳＩ変換機能（Ｒａｓｔｅｒ／Ｃｏｍｂｉｎｅ／Ｃｏｎｖｅｒｔ Ｃｏｌｏｒ．．．）を利用して作成した。
【００５１】
ＩＩＩ．標高データの前処理
（１）ＤＥＭデータ
２００８年１０月に国土地理院基盤情報地図ダウンロードサイト（ｈｔｔｐ：／／ｆｇｄ．ｇｓｉ．ｇｏ．ｊｐ／ｄｏｗｎｌｏａｄ／ＧｓｉＤＬＳｅｌＩｔｅｍＳｅｒｖｌｅｔ）に、四国地方の１０ｍＤＥＭデータが掲載されたデータを利用した。１０ｍＤＥＭデータは、国土地理院の独自フォーマット（ＪＰＧＩＳ２．０（ＧＭＬ）形式基盤地図情報ダウンロードデータ）となっており、そのままでは既存のＧＩＳに読み込むことができない。そこで、（株）エコリスの「基盤地図情報標高データ変換ツール」（ｈｔｔｐ：／／ｗｗｗ．ｅｃｒｉｓ．ｃｏ．ｊｐ／ｄｅｍｔｏｏｌ．ｈｔｍｌ）を利用して基盤地図情報ＧＭＬ形式をＧｅｏＴｉｆｆ形式に変換し、ＴＮＴｍｉｐｓに読み込んだ。その後、ＳＰＯＴデータをリファレンスにしてＣｕｂｉｃＣｏｎｖｏｌｕｔｉｏｎにより内挿を行い２．５ｍ解像度のＤＥＭデータを作成した。
【００５２】
（２）Ｓｈａｄｉｎｇ画像作成
ＤＥＭデータより斜面方位Ａｓｐｅｃｔ、斜面傾斜角Ｓｌｏｐｅおよび相対日射量Ｓｈａｄｉｎｇ画像を作成した。Ｓｈａｄｉｎｇ画像作成時の太陽位置は、仰角７１．３°、方位角１１９．６°とした。
画像全体のヒストグラムを計算したところ、太陽仰角が７１°と高いこともあり、Ｓｈａｄｉｎｇ値による解析対象地の絞り込みは行わなかった。
【００５３】
ＩＶ．オブジェクトサイズの検討
ＳＰＯＴデータ（Ｂａｎｄ１〜４およびＰａｎ）により、ｅＣｏｇｎｉｔｉｏｎを利用してオブジェクト生成を行い、スケールパラメータについて検討を行った。オブジェクト生成のテストサイトを、別子山事業区の含まれるデータの南半分とした。
スケールパラメータ以外の設定はＳｈａｐｅ ｆａｃｔｏｒ＝０．２，Ｃｏｍｐａｃｔｎｅｓｓ＝０．５、Ｓｍｏｏｔｈｎｅｓｓ＝０．５とし、スケールパラメータを１５、２０、２５、３５、５０の５通りに設定してオブジェクトを生成した。画像上で判読可能な林分の領域とオブジェクト境界、生成されるオブジェクト数を考慮して、スケールパラメータを２５としてオブジェクト分類を行うのが適当であると結論した（図２）。
【００５４】
Ｖ．グラウンドトゥルース
２００８年１１月３日〜６日に、住友林業（株）新居浜山林事業所および住友林業フォレストサービス（株）新居浜山林事業所の協力の下、別子山・高藪・大永山の各事業区を中心に現地調査を行い、衛星画像と現地状況の対応を確認した。
調査には、ＧＰＳ受信機（Ｇａｒｍｉｎ ＧＰＳｍａｐ６０ＣＳＸ）、ＧＰＳカメラ（Ｒｉｃｏｈ Ｃａｐｌｉｏ ５００ＳＥおよびＢｌｕｅＴｏｏｔｈ ＧＰＳ受信機ＨＡＩＣＯＭ ＨＩ−４０８ＢＴ）および衛星画像を持参し、調査軌跡の記録、撮影写真の位置情報取得を行った。
【００５５】
ＶＩ．トレーニングエリア選択
まずｅＣｏｇｎｉｔｉｏｎでＳＰＯＴデータ（Ｂａｎｄ１〜４、Ｐａｎ）全域をスケールパラメータ２５としてオブジェクトを生成した。オブジェクト数は１１７，４５７であった。次に、ｅＣｏｇｎｉｔｉｏｎで生成されたオブジェクトの境界ベクタをＳｈａｐｅファイルに変換し、ＴＮＴｍｉｐｓにインポートして以下の作業を行った。Ｓｈａｐｅファイル変換の際には、ＳＰＯＴデータの５バンドに加えて５種類の指数画像についても各オブジェクトにおける平均値をＣＳＶファイルおよびＤＢＦファイルとして出力した。
ＴＮＴｍｉｐｓ上でグラウンドトゥルース調査地点を中心に画像判読を行いながら、オブジェクトポリゴンにクラス属性を付与する手法でトレーニングエリアを選択した。分類クラスは、非植生項目として「道路」および「水域」、植生項目として「スギ」、「スギ（緑）」、「ヒノキ」、「スギ（疎）」、「ヒノキ（疎）」、「ヒノキ（幼齢）」、「カラマツ」、「針広混交」、「ザツ（広葉樹）」、「マツ天然」、「草地および裸地」、「ササ地」の１４項目で、合計で４７９オブジェクトとなった。このうち、「スギ（青）」、「カラマツ」および「ヒノキ（幼齢）」の３クラスに関しては、トレーニングエリアとなるオブジェクトが１０個以下しかなかったため、以降の分類作業では分類項目としなかった。
４７９のトレーニングエリアを、分類に使用するための分類用セットと、精度評価に使用するための評価用セットがおよそ６：４になるようにランダムに振り分けた。表−２に各クラスのトレーニングエリア、分類用セット、評価用セットのオブジェクト個数を示す。
【００５６】
【表２】
【００５７】
ＶＩＩ．階層的分類法による土地被覆分類
（１）手法
トレーニングエリアとして指定されたオブジェクト番号をＴＮＴｍｉｐｓからＣＳＶファイルに出力し、ｅＣｏｇｎｉｔｉｏｎからは全オブジェクトについて以下に示す２８種類のパラメータとオブジェクト番号をＣＳＶファイルに出力し、それらをＥｘｃｅｌに読み込み作業を行った：Ｂａｎｄ１平均値、Ｂａｎｄ２平均値、Ｂａｎｄ３平均値、Ｂａｎｄ４平均値、Ｐａｎ平均値、ＮＤＶＩ平均値、ＮＤＷＩ平均値、Ｈｕｅ平均値、Ｓａｔｕｒａｔｉｏｎ平均値、Ｉｎｔｅｎｓｉｔｙ平均値、ＮＤＩ＿ＳＷＲ平均値、Ｒａｔｉｏ＿ＮＩＲ平均値、Ｒａｔｉｏ＿ＳＷＩＲ平均値、標高平均値、Ｂａｎｄ１標準偏差、Ｂａｎｄ２標準偏差、Ｂａｎｄ３標準偏差、Ｂａｎｄ４標準偏差、Ｐａｎ標準偏差、ＮＤＶＩ標準偏差、ＮＤＷＩ標準偏差、Ｈｕｅ標準偏差、Ｓａｔｕｒａｔｉｏｎ標準偏差、Ｉｎｔｅｎｓｉｔｙ標準偏差、ＮＤＩ＿ＳＷＲ標準偏差、Ｒａｔｉｏ＿ＮＩＲ標準偏差、Ｒａｔｉｏ＿ＳＷＩＲ標準偏差、標高標準偏差。なお、平均値、標準偏差とも各オブジェクト内での値である。
各パラメータの分類クラス毎の平均値・最大値・最小値・標準偏差のグラフ、任意の２
パラメータを組み合わせた散布図を作成し、１４種類の土地被覆項目の特性を検討した。
【００５８】
（２）メンバーシップ関数
検討の結果、表−３に示すような階層構造およびメンバーシップ関数を設定した。それぞれの関数および分類項目との関係を図−６（図３Ａ〜図３Ｎ）に示す。図３Ｂ〜図３Ｎにおいて図示された直線は、各メンバーシップ関数を表すものである。すなわち、当該直線によって、各分類項目の２群が属する座標平面上の領域が示されている。
なお、この中で「マツ天然」に関しては、「スギ」と「マツ天然」が衛星画像データのみでの判別が不可能であったため、現地調査での知見をもとに、「マツ天然」クラスは尾根（稜線）部のみしか出現しないという仮定を設けた。ＤＥＭデータから以下のアルゴリズム（ｈｔｔｐ：／／ｗｗｗ．ｏｐｅｎｇｉｓ．ｃｏ．ｊｐ／ｈｔｍ／ｋａｋｏ＿ｍａｉｌ／ｍａｉｌ＿ｍａｇ＿１３７．ｈｔｍ）を用いて尾根部を判別し、メンバーシップ関数に組み込んだ。
（ｉ）ＴＮＴｍｉｐｓのＲａｓｔｅｒ／Ｅｌｅｖａｔｉｏｎ／Ｗａｔｅｒｓｈｅｄ．．．により、谷線を表すＳＴＤＦＬＯＷＰＡＴＨベクタと流域界を表すＳＴＤＢＡＳＩＮベクタを作成する。
（ｉｉ）ＴＮＴｍｉｐｓのＧｅｏｍｅｔｒｉｃ／Ｃｏｍｐｕｔｅ／Ｄｉｓｔａｎｃｅ Ｒａｓｔｅｒ．．．により、（ｉ）で作成した二つのベクタからの距離ラスタを発生させ、それぞれをＤＩＳＴ＿ＢＡＳＩＮとＤＩＳＴ＿ＦＬＯＷとする。
（ｉｉｉ）以下のジオフォーミュラにより尾根部、斜面部、谷部を分ける。
ｄｉｓｔ＿ｐａｒａ＝２；
ｉｆ（ｓｑｒｔ（（ＤＩＳＴ＿ＢＡＳＩＮ−ＤＩＳＴ＿ＦＬＯＷ）＾２）＜（（（ＤＩＳＴ＿ＢＡＳＩＮ＋ＤＩＳＴ＿ＦＬＯＷ）／２）／ｄｉｓｔ＿ｐａｒａ））｛ＳｌｏｐｅＺｏｎｅ＝２；｝ｅｌｓｅ｛ｉｆ（ＤＩＳＴ＿ＢＡＳＩＮ＜ＤＩＳＴ＿ＦＬＯＷ）｛ＳｌｏｐｅＺｏｎｅ＝１；｝ｅｌｓｅ｛ＳｌｏｐｅＺｏｎｅ＝３；｝｝
このアルゴリズムで、ＳｌｏｐｅＺｏｎｅが１が尾根部を表す。
【００５９】
また、分類の際に「針広混交」クラスは、Ｓａｔｕｒａｔｉｏｎの値により３段階に区分を行った。
Ｓａｔｕｒａｔｉｏｎが０〜５０の場合は針葉樹混交率が高い「ＮＬ＿ｍｉｘｅｄ Ｈ」、Ｓａｔｕｒａｔｉｏｎが５０〜１５０の場合は針葉樹混交率が中庸の「ＮＬ＿ｍｉｘｅｄ Ｍ」、Ｓａｔｕｒａｔｉｏｎが１５０〜１９０の場合は針葉樹混交率の低い「ＮＬ＿ｍｉｘｅｄ Ｌ」とした。現地調査やトレーニングエリアのクラスでは混交率の差について考慮していないため、ここで設けた混交率による３クラスはあくまでもＳａｔｕｒａｔｉｏｎの大小による分類クラスである。なぜ、このように混交率クラスを設けたかについては、グラフ（図−６（ｅ）（図３Ｆ））上ではＳａｔｕｒａｔｉｏｎが同じ「針広混交」クラスでも大きく異なっていたこと、Ｓａｔｕｒａｔｉｏｎ値が針葉樹クラスでは小さく「ザツ（広葉樹）」クラスでは大きな値を持っており、以下の図−６（ｅ）（図３Ｆ）に示すように針葉樹クラスから「針広混交」クラス、「ザツ（広葉樹）」クラスまでＳａｔｕｒａｔｉｏｎ値が連続的であるように見えたためである。
【００６０】
【表３】
【００６１】
（３）分類結果と精度評価
表−３に示したメンバーシップ関数を用いて画像全体を分類した結果を図−８（図４）に示す。
トレーニングエリアの評価用セットを用いたエラーマトリクスを表−４に示す。この中で、「針広混交」クラスについては、前述の通り混交率によるトレーニングエリアがないので、混交率Ｍ（ＮＬ＿ｍｉｘｅｄ Ｍ）を「針広混交」クラスとした形で分類精度の計算を行った。また、「Ｎｏｎ＿ｖｅｇｅｔａｔｉｏｎ」：非植生は、「道路」および「水域」をまとめたクラスである。なお、数値はピクセル数である。
【００６２】
【表４】
【００６３】
全体の分類精度が８１．９％、カッパ係数も７８．２％であり、評価用セットにおいては分類は妥当に行われたと判断できる。
針葉樹種（密）については分類精度がヒノキ８９％、スギ９８％、図化精度がヒノキ９１％、スギ１００％と非常に高い。ヒノキと分類されたオブジェクトには、「針広混交Ｈ」および「針広混交Ｍ」が若干含まれているが、スギと分類されたオブジェクト誤分類は全くなく、スギについての信頼度がかなり高いことが推測できる。針葉樹種（疎）については、ヒノキは「草地および裸地」「針広混交Ｍ」「針広混交Ｈ」への誤分類および未分類があり分類精度５１％であるが、スギは誤分類が全くない。分類結果では他のトレーニングエリアが針葉樹種（疎）と分類されたものが全くないため、図化精度はどちらも１００％となっている。
【００６４】
分類項目をまとめて、より大分類にすることにより精度がどのように変化するかを検討した。ここでは、「Ｎｏｎ＿ｓｔｏｃｋｅｄ」：未立木地（「非植生」、「草地および裸地」および「ササ」）、「ザツ（広葉樹）」、「針広混交」（「針広混交」および「マツ天然」）、「スギ」（「スギ」および「スギ（疎）」）、「ヒノキ」（「ヒノキ」および「ヒノキ（疎）」）の５項目にまとめた。分類精度を表−５に示す。
全体精度は９１％、「ザツ」を除く分類項目の分類精度は約９０％以上、針葉樹種項目の図化精度も９０％以上となった。
【００６５】
【表５】
【００６６】
ＶＩＩＩ．非階層的分類法による土地被覆分類
トレーニングエリアの分類用セットを教師データとし、階層的分類法で使用した１０パラメータ（地形条件を除く）、Ｂａｎｄ１、Ｂａｎｄ３、Ｂａｎｄ４、Ｐａｎ、ＮＤＶＩ、ＮＤＷＩ、Ｈｕｅ、Ｓａｔｕｒａｔｉｏｎ、Ｉｎｔｅｎｓｉｔｙ、Ｒａｔｉｏ＿ＳＷＲを用いて、ｅＣｏｇｎｉｔｉｏｎのＮｅａｒｅｓｔ Ｎｅｉｇｈｂｏｒ（ＮＮ）法により土地被覆分類を行った。全域の分類結果を図−１０（図５）に示す。
また、トレーニングエリアの評価用セットを用いたエラーマトリクスを表−６に示す。数値はピクセル数である。
【００６７】
【表６】
【００６８】
全体の分類精度は８６％、カッパ係数も８３％であり、分類精度に関してはメンバーシップ関数による分類結果よりもよい結果となった。
針葉樹種（密）に関しては、分類精度がヒノキ８４％、スギ９８％、図化精度も同様の値である。「ヒノキ」のトレーニングエリアは「スギ」、「ヒノキ（疎）」および「針広混交」に誤分類が行われているが、「スギ」のトレーニングエリアは「スギ（疎）」への誤分類のみである。一方、「ヒノキ」と分類されたオブジェクトには「ザツ（広葉樹）」、「マツ天然」および「針広混交」が若干含まれており、「スギ」と分類されたオブジェクトには「ヒノキ」および「マツ天然」が若干含まれている。メンバーシップ関数による分類と比較すると、メンバーシップ関数による分類の方がクラス間の混同が少なかった。
【００６９】
前章と同様に、分類項目をまとめてより大分類にすることにより精度がどのように変化するかを検討した。ここでは、「Ｎｏｎ＿ｓｔｏｃｋｅｄ」：未立木地（「非植生」、「草地および裸地」および「ササ」）、「ザツ（広葉樹）」、「針広混交」（「針広混交」および「マツ天然」）、「スギ」（「スギ」および「スギ（疎）」）、「ヒノキ」（「ヒノキ」および「ヒノキ（疎）」）の５項目にまとめた。分類精度を表−７に示す。
全体の分類精度はやや向上したが、メンバーシップ関数による分類精度（表−５）の方が、「針広混交」の分類精度および「ヒノキ」の図化精度がより高かった。
【表７】
【００７０】
ＩＸ．考察
分類結果の画像とＳＰＯＴ画像、現地の状況を比較検討すると、今回の２手法ではメンバーシップ関数による分類結果の方が、全体的により妥当であると判断された。
今回の精度評価は、選択したトレーニングエリアをランダムに２つに分け、その片方の評価用セットで行った。トレーニングエリアは、現地調査を参考に画像を判読しながら行うが、比較的判別しやすいところを選んだ。
【産業上の利用可能性】
【００７１】
本発明によれば、森林の植生分類図をより高精度で作成することができる。したがって、本発明は、林業および環境関連産業等の発展に寄与するところ大である。
【技術分野】
【０００１】
本発明は、中間赤外域バンドを利用した森林の樹種分類に関し、とくに衛星データによる森林植生分類図の作成方法および森林植生分類図に関する。
【背景技術】
【０００２】
人工造林地が正常に成林して十分な木材収穫が行われることは、林業経営にとって不可欠である。一方、人工造林地のモニタリングを行うことも、収穫予測の策定するために不可欠であるところ、同モニタリングが正確に行えず、正常に成林しているか否かが不明な林分が多い場合には伐採時期における精度の高い伐採可能量の予測を立てることができない。
【０００３】
また、森林施業計画には伐採可能量の算出が不可欠であるところ、将来的な伐採可能量は森林調査簿のデータを基に算出しているため、同調査簿の精度がきわめて重要である。
しかしながら、森林調査簿の更新は、間伐等の施業が実施された場合には毎年内容が更新されるが、施業間隔は通常１０年程度であるため、人工林が強風による風倒等の何らかの理由により広葉樹が侵入して針広混交林化し減損している場合には、樹種の把握が不正確になっているおそれもある。その結果、伐採可能量が実際の量より過大なものとして算出されてしまう懸念がある。
したがって、上記モニタリングにより伐採可能量の推定を極力正確に行うことが、林業経営等においては極めて重要である。
【０００４】
また、地球温暖化による各種の弊害が全世界的な規模で問題視される一方、森林・樹木による炭酸ガスの同化作用による地球温暖化の防止が重要であることは論を待たない。そのため、森林・樹木による炭酸ガス同化作用の算定の基礎となるバイオマスを増加させるために、森林整備の必要性が叫ばれ、森林整備が盛んに行なわれている。
かかる植林後の樹種やすでに森林を構成している樹種の正確な把握も、バイオマスの推定のために不可欠である。したがって、森林における植生分類のより高い精度による把握は、今後日本国内にとどまらず世界的に重要性が増していくことは確実である。
【０００５】
モニタリングの手法として、理論的には全山踏査して毎木調査を行うことが考えられる。しかしながら、かかる調査には膨大な労力、莫大なコスト・時間を要するため、現実的にはほぼ不可能である。そのため、モニタリングの代替手法として、従来から航空写真による目視判読や衛星画像解析による植生分類が行われている。とくに、衛星画像解析によれば、より広範な地域の情報をより効率的に得ることができる。
【０００６】
たとえば、航空カラー写真や衛星カラー画像をコンピュータに入力し、コンピュータ画像処理により森林地域の樹冠形状画像を求め、この樹冠形状画像と元のカラー画像から求めた樹木の色彩及び輝度等とから森林地域全体の樹木植生の調査・評価、すなわちリモートセンシングが行われている（特許文献１）。
【０００７】
また、特許文献２には、調査対象地域のバンドデータの季節変化を利用しておこなう衛星データによる森林地域の植生分類方法が示されている。
さらに、特許文献３および４には航空機や衛星を利用したリモートセンシングによる画像の利用が示され、非特許文献１および２には、高分解能ＩＫＯＮＯＳ衛星による樹種分類の方法が記載されている。また、非特許文献３には、空中写真を用いた森林地理情報システム（ＧＩＳ）データベースによる森林ゾーニングについて報告されている。
【先行技術文献】
【特許文献】
【０００８】
【特許文献１】特開２００１−３５７３８０号公報
【特許文献２】特開２００６−８５５１７号公報
【特許文献３】再表２００３−０６９３１５号公報
【特許文献４】特開２００８−７９５４９号公報
【非特許文献】
【０００９】
【非特許文献１】「高分解能ＩＫＯＮＯＳ衛星による針広混交林の樹種分類」加藤正人，森林航測第１９８号６〜９（２００２．１１）
【非特許文献２】「トレーニングデータの選定方法がＩＫＯＮＯＳ画像の樹種分類に与える影響」宮地ら，第１１５回日林学術講，１９０頁
【非特許文献３】「朝日の森を対象とするＧＩＳを用いた森林機能評価とゾーニング」呉ら，東京大学農学部演習林報告，１１１，５９−８３（２００４）
【発明の概要】
【発明が解決しようとする課題】
【００１０】
上記のとおり森林における植生分類をより高い精度に行うことは重要であるにも拘らず、従来の手法は、必ずしも高い精度で樹種を分類し、高精度の植生図の作成を可能としたものとはいえない。たとえば、非特許文献２においては、トレーニングデータの選定方法に依存して分類することができる可能性があるスギおよびヒノキの分類精度は、それぞれ７０％および７８％程度に留まっている。
【００１１】
したがって、本発明は、より高精度の森林植生分類図を作成する方法を提供することを目的とする。
【課題を解決するための手段】
【００１２】
従来の衛星画像のセンサの多くは、可視域光（青、緑、赤）に加えて、植生による反射ピークが存在する近赤外の波長帯を用いている。しかしながら、可視光域に対する植生による反射量は小さいため、樹種ごとの差異を検出しがたいと考えられた。また、近赤外波を併せて用いても、得られる植生の分類精度は、前記のとおり必ずしも十分なものではなかった。
そこで、本発明者らは、衛星画像を用いて植生分類図を作成する際に用いられる光および／または電磁波の種類（波長）を適切に選択することによって植生分類図の精度が向上すること、とくに中間赤外波を用いることによって上記課題が解決することを見出し、さらに鋭意研究を進めて本発明を完成した。
すなわち、本発明は、中間赤外波を含む電磁波による衛星画像を用いることを含む、森林の植生分類図を作成する方法に関する。
【００１３】
また、本発明は、中間赤外波の波長が１．３μｍ〜２．０μｍである、記方法に関する。
さらに本発明は、さらに近赤外波および可視光を用いる前記方法に関する。
また、本発明は、衛星画像から、２本以上の単木を含むオブジェクトを作成する工程を含む、前記方法に関する。
【００１４】
さらに、本発明は、分類パラメータを用いてオブジェクトの土地被覆分類を行う工程を含み、該分類パラメータが中間赤外波に対する反射輝度値を用いるパラメータを含む前記方法に関する。
また、本発明は、分類パラメータが正規化水指数ＮＤＷＩ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲおよびＲａｔｉｏ＿ＳＷＩＲからの１種または２種以上である前記方法に関する。
【００１５】
また、本発明は、階層的分類法によるオブジェクトの土地被覆分類を行う工程を含む前記方法に関する。
そして本発明は、前記いずれかの方法を用いて作成された森林の植生分類図に関する。
【発明の効果】
【００１６】
本発明の方法によれば、森林の植生分類図を高精度で作成することができる。
本発明の方法のうち中間赤外波長の波長が１．３μｍ〜２．０μｍであるものによれば、森林の植生分類図をより高精度で作成することができる。
本発明の方法のうちさらに近赤外波および可視光を用いるものによれば、森林の植生分類図をさらに高精度で作成することができる。
本発明の方法のうち近赤外波の波長が０．７５μｍ〜１．０μｍであるものまたは可視光が緑色光および赤色光であるものによれば、森林の植生分類図をより一層高精度で作成することができる。
本発明の方法のうち衛星画像から、２本以上の単木を含むオブジェクトを作成する工程を含むものによれば、高分解能（分解能１ｍ以下）を具備していない衛星による衛星画像によっても高精度の分類が可能となり、植生分類図の作成をより低コストで行うことができる。
【００１７】
本発明の方法のうち分類パラメータを用いてオブジェクトの土地被覆分類を行う工程を含むものによれば、森林の植生分類図をより一層高精度で作成することができる。
本発明の方法のうち分類パラメータが中間赤外波に対する反射輝度値を用いるパラメータを含むもの、とくに分類パラメータが正規化水指数ＮＤＷＩ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲおよびＲａｔｉｏ＿ＳＷＩＲからの１種または２種以上であるもの、または分類パラメータが正規化植生指数ＮＤＶＩ、正規化水指数ＮＤＷＩ、色相Ｈｕｅ、彩度Ｓａｔｕｒａｔｉｏｎ、明度Ｉｎｔｅｎｓｉｔｙ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲ、Ｒａｔｉｏ＿ＳＷＩＲからなるものによれば、森林の植生分類図をさらに高精度で作成することができる。
本発明の方法のうち階層的分類法によるオブジェクトの土地被覆分類を行う工程を含むによれば、森林の植生分類図としてより妥当なものを作成することができる。
本発明の方法のうちスギおよび／またはヒノキの分類精度が８０％以上であるものによれば、スギやヒノキ等の針葉樹についてとくに分類精度が高い森林の植生分類図を作成することができる。
【００１８】
本発明による森林の植生分類図によれば、収穫対象林分と減損林分を正確に把握することができる。
【００１９】
植生による反射ピークは近赤外域に存在するのに対し、中間赤外波は必ずしも植生による反射は大きくなく、少なくとも近赤外波に対する反射より小さい傾向がある。したがって、本発明において中間赤外波を用いてより高精度の森林植生分類図の作成することは、本発明の特徴の一つである。
【図面の簡単な説明】
【００２０】
【図１】ＳＰＯＴ−５による提供データ範囲および切り出し範囲を示す図である。
【図２】スケールパラメータを２５としてオブジェクトを生成した図である。
【図３Ａ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｂ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｃ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｄ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｅ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｆ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｇ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｈ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｉ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｊ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｋ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｌ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｍ】各メンバーシップ関数および分類項目との関係を示す図である。
【図３Ｎ】各メンバーシップ関数および分類項目との関係を示す図である。
【図４】階層的分類法による土地被覆分類の結果を示す図である。
【図５】非階層的分類法による土地被覆分類の結果を示す図である。
【発明を実施するための形態】
【００２１】
（定義）
本明細書において「中間赤外（波）」（「中間赤外域バンド波」）とは、一般的な近赤外波と遠赤外波との中間の波長を有する赤外波を意味する。典型的には、「中間赤外（波）」とは、約１．２μｍ〜４μｍの波長を有する電磁波である。また、青色光、緑色光、赤色光および近赤外波は、それぞれ約０．４μｍ〜約０．５μｍ、約０．５μｍ〜約０．６μｍ、約０．６μｍ〜約０．７μｍおよび０．７μｍ〜１．２μｍの波長を有する。
【００２２】
また、それぞれの波長の境界付近の波長を有する可視光または電磁波は、当該境界を構成する可視光または電磁波のいずれにも分類し得る。たとえば、０．４５μｍ〜０．５２μｍのものは、青色光または緑色光とも称され得る。
【００２３】
本発明の方法のうち、下記（１）〜（４）の各工程を含むものは好ましい。これらの工程をすべて含むものはより好ましい。
（１）衛星画像データおよび地形データの前処理
（２）オブジェクトの作成
（３）グラウンドトゥルース
（４）オブジェクトの土地被覆分類
上記（１）〜（４）の各工程につき、以下にさらに詳細に説明する。
【００２４】
（１）衛星画像データおよび地形データの前処理
解析の対象となる衛星画像データおよび地形データのオリジナルデータを解析により適した形式に変換する工程である。
当該工程は、中間赤外波を含む電磁波での衛星画像を用いるものであれば、他の工程・手法等はとくに限定されない。すなわち、当該衛星画像データは、中間赤外波を含む電磁波での衛星画像を用いたデータであれば、本発明に用いることができる。
【００２５】
本発明の方法における中間赤外波の波長は約１．２μｍ〜４μｍであるところ、当該波長が１．３μｍ〜２．０μｍであるものは好ましい。また、中間赤外波の波長は、より好ましくは１．５μｍ〜１．８μｍであり、最も好ましくは１．５８μｍ〜１．７５μｍである。
【００２６】
本発明の方法は、さらに近赤外波ならびに可視光を用いることができる。近赤外波の波長は０．７μｍ〜１．２μｍであり、好ましくは０．７５μｍ〜１．０μｍであり、より好ましくは０．７８μｍ〜０．８９μｍである。
【００２７】
可視光としては、青色光、緑色光および赤色光からの１種または２種以上が好適であり、緑色光および赤色光を用いたものはより好ましい。
【００２８】
（２）オブジェクトの作成
使用画像の解像度や分類対象物の大きさに応じて適切なスケールパラメータを決定し、画像を解析の単位であるオブジェクト（ピクセルの集合体）に区分する工程である。
すでに衛星画像の分解能は高くなり（１ｍ以下）、単木単位での樹種判読は高い精度で分類が可能となっている。しかしながら、高分解能の衛星画像は高価であるため、広域な森林を調査するには莫大な費用が掛かるという欠点がある。また、ピクセル単位での画像分類は、与えられる情報量は豊富であるが、たとえば樹冠の陰影のような過剰な情報が含まれることがあり、当該目的にはむしろ適さない面もある。
【００２９】
また、ピクセル単位での画像分類においては、解像度が高いことにより単木樹冠が複数のピクセルに分割されることにより、ピクセルベースの分類では複数のクラスに誤分類されるおそれが高く、林分のまとまりを性格に把握することが困難である。
【００３０】
これに対し、オブジェクトベースの分類においては、同様なテクスチャや分光反射特性の領域を単位として認識分類を行うため、林相区分図作成により好適である。
また、ある程度のまとまりをもったピクセルの集合体、すなわちオブジェクトを画像分類の単位とする手法が開発され、新規のアプローチが可能となっている。そのため、本発明においては衛星画像から、２本以上の単木を含むオブジェクトを作成する工程を含むものは好ましいのである。
しかしながら、オブジェクトを作成する工程を含む方法においては、樹種ごとに適切な（分類）パラメータが明らかではない。したがって、分類パラメータの新たな検討・決定が、当該工程を含む方法においては重要である。当該分類パラメータの検討・決定についてはさらに後記する。
【００３１】
（３）グラウンドトゥルース
現地状況を実地に調査する工程である。
グラウンドトゥルースによって得られた現地状況と画像の判読結果からトレーニングエリアが画像から選択され、さらに当該トレーニングエリアは、分類のために用いられる分類用セットと精度評価に用いられる評価用セットとにランダムに振り分けられる。
作成された植生分類図および当該分類の精度評価は、グラウンドトゥルースによって得た現地状況からトレーニングエリアを選択し、当該トレーニングエリアにおける比較対照することによって行うことができる。
したがって、本発明の方法のうち、グラウンドトゥルースデータを得る工程を含むものは好ましい。
【００３２】
（４）オブジェクトの土地被覆分類
分類パラメータを用いる手法等により、作成されたオブジェクトごとに土地被覆分類を決定し、植生分類図を完成する工程である。
本発明のうちオブジェクトを作成する工程を含む方法においては、樹種ごとに適切な分類パラメータの新たな検討・決定を行うことが重要である。したがって、本発明の方法のうち、分類パラメータを用いて土地被覆分類を行う工程を含むものは好ましい。
【００３３】
分類パラメータを用いる手法には、例えば階層的または非階層的分類法を例示することができる。
階層的分類法は、大分類から小分類へと段階的に、各段階でそれぞれ異なったメンバーシップ関数を用いて分類する手法である。階層的分類法は、より具体的には、
・各種分類パラメータの分類クラスごとの平均値、最大値、最小値、標準偏差等のグラフおよび任意の２つの分類パラメータを組み合わせた散布図等を作成し、２以上の土地被覆項目の特性を検討し、さらに、
・上記検討の結果から階層構造およびメンバーシップ関数を設定し、画像全体を分類する方法である。したがって、本明細書においては、階層的分類（法）を「メンバーシップ関数を用いた分類（法）」ということがある。
【００３４】
非階層的分類法は、あらかじめ設定した複数の土地被覆分類項目に、同じ分類パラメータを用いて一度に分類する手法である。
これらの手法のうち、階層的分類法によるオブジェクトの土地被覆分類を行う工程を含む本発明の方法は、分類がより好適に行われるため好ましい。
【００３５】
分類パラメータとしては、正規化植生指数ＮＤＶＩ、正規化水指数ＮＤＷＩ、色相Ｈｕｅ、彩度Ｓａｔｕｒａｔｉｏｎ、明度Ｉｎｔｅｎｓｉｔｙ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲ、Ｒａｔｉｏ＿ＳＷＩＲ等が例示される。
【００３６】
上記各分類パラメータのうち、正規化植生指数ＮＤＶＩ、正規化水指数ＮＤＷＩ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲ、Ｒａｔｉｏ＿ＳＷＩＲは、それぞれ下記のように定義される：
ＮＤＶＩ＝（Ｂａｎｄ３−Ｂａｎｄ２）／（Ｂａｎｄ３＋Ｂａｎｄ２）＊１００＋１００
ＮＤＷＩ＝（Ｂａｎｄ３−Ｂａｎｄ４）／（Ｂａｎｄ３＋Ｂａｎｄ４）＊１００＋１００
ＮＤＩ＿ＳＷＲ＝（Ｂａｎｄ４−Ｂａｎｄ２）／（Ｂａｎｄ４＋Ｂａｎｄ２）＊１００＋１００
Ｒａｔｉｏ＿ＮＩＲ＝Ｂａｎｄ３／（Ｂａｎｄ１＋Ｂａｎｄ２＋Ｂａｎｄ３＋Ｂａｎｄ４）＊１００
Ｒａｔｉｏ＿ＳＷＩＲ＝Ｂａｎｄ４／（Ｂａｎｄ１＋Ｂａｎｄ２＋Ｂａｎｄ３＋Ｂａｎｄ４）＊１００
【００３７】
上記定義中のバンド１〜４（Ｂａｎｄ１〜４）の波長は、それぞれ０．５０μｍ〜０．５９μｍ（青）、０．６１μｍ〜０．６８μｍ（赤）、０．７９μｍ〜０．８９μｍ（近赤外波）、１．５８μｍ〜１．７５μｍ（中赤外波）である
また、色相Ｈｕｅ、彩度Ｓａｔｕｒａｔｉｏｎ、明度ＩｎｔｅｎｓｉｔｙはＴＮＴｍｉｐｓのＲＧＢ−ＨＳＩ変換機能（Ｒａｓｔｅｒ／Ｃｏｍｂｉｎｅ／Ｃｏｎｖｅｒｔ Ｃｏｌｏｒ．．．）を利用して作成される。
上記のとおり、各波長の光（電磁波）に対する反射輝度値を用いる５種類のパラメータのうち４種類に中赤外波（Ｂａｎｄ４）に対する反射輝度値を用い得ることは本発明の特長である。これら４種類のパラメータからの１種または２種以上を用いることは好ましい。さらに、分類パラメータが、正規化植生指数ＮＤＶＩ、正規化水指数ＮＤＷＩ、色相Ｈｕｅ、彩度Ｓａｔｕｒａｔｉｏｎ、明度Ｉｎｔｅｎｓｉｔｙ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲ、Ｒａｔｉｏ＿ＳＷＩＲのみからなる本発明の方法はより好ましい。
【００３８】
なお、ＮＤＶＩは近赤外域と可視光（赤色光）域における反射量の差をみるためのパラメータであるところ、当該パラメータは植物のバイオマス量や活力度と正の相関があることが知られている。また、ＮＤＩは植生指数の一種である。Ｒａｔｉｏは全バンド中において特定のバンドが占める割合であって、一般的な変換方法による反射輝度値を用いた指数である。
【００３９】
本発明の方法は、スギおよび／またはヒノキの他の樹種との分類にとくに好適である。たとえば、スギおよび／またはヒノキの分類精度が８０％以上である本発明の方法は好ましい。
【００４０】
前記のとおり、本発明は上記いずれかに記載の方法を用いて作成された森林の植生分類図にも関するところ、たとえば本植生分類図によれば収穫対象林分と減損林分を正確に把握することができるのである。
【００４１】
本発明につき、実施例に基づきより詳細に以下に説明する。かかる実施例はあくまで本発明を説明するためのものであって、いかなる意味においても本発明を限定するものではない。
【実施例】
【００４２】
ＳＰＯＴ５衛星画像の解析により、四国地域における住友林業株式会社社有林（造林地）の植生分類図を作成し、収穫対象林分と減損林分とにゾーニングした。
なお、衛星画像として用いられたＳＰＯＴ５衛星画像の概要は以下のとおりである：
ＳＰＯＴ５衛星画像「ＳＰＯＴＶｉｅｗ Ｏｒｔｈｏ」
（画像中心座標：３３°４９′３０″・１３３°２３′３０″）
・パンクロ画像（２７×２７ｋｍ、解像度２．５ｍ）１シーン
・カラー画像（２７×２７ｋｍ、解像度１０ｍ）１シーン
【００４３】
また、住友林業株式会社内調査簿データ（２００８年度分；ｅｘｃｅｌファイル）および同社内森林ＧＩＳデータ（ｓｈａｐｅファイル・ｔｉｆｆファイル）も適宜用いた。
【００４４】
Ｉ．解析方法
概要は下記（１）〜（６）のとおりである。
（１）使用ソフトウェア
リモートセンシングデータ解析およびＧＩＳ解析にはＴＮＴｍｉｐｓ ｖｅｒ７．０〜７．４（合衆国ＭｉｃｒｏＩｍａｇｅｓ 社）、オブジェクトベース分類にはｅＣｏｇｎｉｔｉｏｎ ｖｅｒ４．０（ドイツＤｅｆｉｎｉｅｎｃｅ Ｉｍａｇｉｎｇ 社）を利用した。二つのソフト間でのデータ交換は、ＴＮＴｍｉｐｓによりラスタデータはＧｅｏＴｉｆｆ変換、ベクタデータはＳｈａｐｅ ｆｉｌｅ変換を行うことにより行った。
（２）前処理
入手したオリジナルデータをＴＮＴｍｉｐｓ形式に変換し、衛星画像データについては大気補正、切り出し、ＮＤＶＩ等画像作成を、地形データについてはモザイク、切り出し、Ｓｈａｄｉｎｇ計算、斜面位置区分などの処理を行った。
（３）ｅＣｏｇｎｉｔｉｏｎオブジェクト作成
ｅＣｏｇｎｉｔｉｏｎでは、スケールパラメータを変えることにより生成されるオブジェクトの大きさが異なるため、使用画像の解像度や分類対象物の大きさに応じて適切なスケールパラメータを設定する必要がある。本実施例においてはスケールパラメータを１５から５０まで変化させたときのオブジェクトサイズを検討し、今回の対象地に適当なスケールパラメータを決定した。
【００４５】
（４）グラウンドトゥルース
２００８年１１月３日〜６日に、別子山・高藪・大永山の各事業区を中心に現地調査を行い、衛星画像と現地状況の対応を確認した。
（５）非階層的および階層的分類法によるオブジェクトの土地被覆分類
土地被覆分類については、非階層的分類および階層的分類の２種類を検討した。前者は、あらかじめ設定した複数の土地被覆分類項目に、同じ分類パラメータを用いて一度に分類する手法である。後者は、大分類から小分類へと段階的に、各段階でそれぞれ異なったメンバーシップ関数を用いて分類する手法である。後者の手法を実行するためには、分類パラメータおよび閾値を適切に選択する必要があるため、ｅＣｏｇｎｉｔｉｏｎで生成したオブジェクトをＴＮＴｍｉｐｓに取り込み、トレーニングデータを画像中から選択して検討を行った。
（６）分類結果の検討
二つの手法で行った分類結果について、精度評価用サンプルエリアにより分類精度を評価した。
【００４６】
実際に行われた解析の詳細は以下のとおりである。
【００４７】
ＩＩ．衛星データの前処理
（１）ＳＰＯＴ−５データの地理補正とリサンプリング、切り出し
提供されたＳＰＯＴ−５データをｍｉｐｓ形式に変換後、空中写真および林小班界と重ね合わせて表示したところ、地理的な位置ずれはほとんど無視できると判断したため、ＧＣＰを用いた地理補正は省略した。この点については、例えばＪＡＸＡ提供のＡＬＯＳデータと比較した場合の利点であり、作業量を増やすことがない。提供データはＵＴＭ座標系に基づいたものであったため、平面直角座標系（ＩＶ系）、２．５ｍ解像度へとＮｅａｒｅｓｔ Ｎｅｉｇｈｂｏｕｒ（ＮＮ）法によりリサンプリングを行った。さらに、画像サイズを必要最小限にするため、五良津・大永山・別子山・高藪の４事業区を含むように切り出した。ただし、大永山事業区の西端部分は提供データ範囲に含まれていなかったため、切り出し後のデータにも含まれていない。また、石鎚山事業区は提供データ範囲に含まれていなかったため、今回の解析には含んでいない。切り出し範囲の座標値（Ｎ，Ｅ）は以下の通りである（単位ｍ）。左上：（１０２，９７０，−２１，６０２．５）、右上：（１０２，９７０，１，１４７．５）、左下：（７９，４７０，−２１，６０２．５）、右下：（７９，４７０，１，１４７．５）。画像サイズは９，１００ピクセル×９，４００ラインである。図１に提供データ範囲を示す。
（２）大気補正
大気補正は最小値法を利用した。提供データ範囲において、ＳＰＯＴ−５データの各バンド（Ｂａｎｄ１〜４、Ｐａｎ）のヒストグラムを用いて最小値を検索した。各バンドの最小値は表−１に示すとおりである。この表の値から１を引いた値を各バンドから差し引いた。最小値から１を引いたのは、大気補正後の衛星画像データＤＮ値が０になるとＮＤＶＩ計算の際に分母が０となり計算できないことがあるためである。
【００４８】
【表１】
【００４９】
なお、バンド１〜４（Ｂａｎｄ１〜４）の波長は、それぞれ０．５０μｍ〜０．５９μｍ（青）、０．６１μｍ〜０．６８μｍ（赤）、０．７９μｍ〜０．８９μｍ（近赤外波）、１．５８μｍ〜１．７５μｍ（遠赤外波）である。また、Ｐａｎはバンド１〜４のセンサ（ＨＲＧ−Ｘ）とは異なるセンサ（ＨＲＶ−Ｐ）によって検知される０．４８μｍ〜０．７１μｍの波長の可視光である。
【００５０】
（３）各種指数画像の作成
分類に使用する指数として、正規化植生指数ＮＤＶＩ、正規化水指数ＮＤＷＩ、色相Ｈｕｅ、彩度Ｓａｔｕｒａｔｉｏｎ、明度Ｉｎｔｅｎｓｉｔｙ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲ、Ｒａｔｉｏ＿ＳＷＩＲを検討の対象とし、それぞれの画像を作成した。Ｈｕｅ、Ｓａｔｕｒａｔｉｏｎ、Ｉｎｔｅｎｓｉｔｙは（Ｒ，Ｇ，Ｂ）を（Ｂａｎｄ４，Ｂａｎｄ３，Ｂａｎｄ２）に割り当て作成した。以下、各指数の計算式を示す。
ＮＤＶＩ＝（Ｂａｎｄ３−Ｂａｎｄ２）／（Ｂａｎｄ３＋Ｂａｎｄ２）＊１００＋１００
ＮＤＷＩ＝（Ｂａｎｄ３−Ｂａｎｄ４）／（Ｂａｎｄ３＋Ｂａｎｄ４）＊１００＋１００
ＮＤＩ＿ＳＷＲ＝（Ｂａｎｄ４−Ｂａｎｄ２）／（Ｂａｎｄ４＋Ｂａｎｄ２）＊１００＋１００
Ｒａｔｉｏ＿ＮＩＲ＝Ｂａｎｄ３／（Ｂａｎｄ１＋Ｂａｎｄ２＋Ｂａｎｄ３＋Ｂａｎｄ４）＊１００
Ｒａｔｉｏ＿ＳＷＩＲ＝Ｂａｎｄ４／（Ｂａｎｄ１＋Ｂａｎｄ２＋Ｂａｎｄ３＋Ｂａｎｄ４）＊１００
Ｈｕｅ、Ｓａｔｕｒａｔｉｏｎ、ＩｎｔｅｎｓｉｔｙはＴＮＴｍｉｐｓのＲＧＢ−ＨＳＩ変換機能（Ｒａｓｔｅｒ／Ｃｏｍｂｉｎｅ／Ｃｏｎｖｅｒｔ Ｃｏｌｏｒ．．．）を利用して作成した。
【００５１】
ＩＩＩ．標高データの前処理
（１）ＤＥＭデータ
２００８年１０月に国土地理院基盤情報地図ダウンロードサイト（ｈｔｔｐ：／／ｆｇｄ．ｇｓｉ．ｇｏ．ｊｐ／ｄｏｗｎｌｏａｄ／ＧｓｉＤＬＳｅｌＩｔｅｍＳｅｒｖｌｅｔ）に、四国地方の１０ｍＤＥＭデータが掲載されたデータを利用した。１０ｍＤＥＭデータは、国土地理院の独自フォーマット（ＪＰＧＩＳ２．０（ＧＭＬ）形式基盤地図情報ダウンロードデータ）となっており、そのままでは既存のＧＩＳに読み込むことができない。そこで、（株）エコリスの「基盤地図情報標高データ変換ツール」（ｈｔｔｐ：／／ｗｗｗ．ｅｃｒｉｓ．ｃｏ．ｊｐ／ｄｅｍｔｏｏｌ．ｈｔｍｌ）を利用して基盤地図情報ＧＭＬ形式をＧｅｏＴｉｆｆ形式に変換し、ＴＮＴｍｉｐｓに読み込んだ。その後、ＳＰＯＴデータをリファレンスにしてＣｕｂｉｃＣｏｎｖｏｌｕｔｉｏｎにより内挿を行い２．５ｍ解像度のＤＥＭデータを作成した。
【００５２】
（２）Ｓｈａｄｉｎｇ画像作成
ＤＥＭデータより斜面方位Ａｓｐｅｃｔ、斜面傾斜角Ｓｌｏｐｅおよび相対日射量Ｓｈａｄｉｎｇ画像を作成した。Ｓｈａｄｉｎｇ画像作成時の太陽位置は、仰角７１．３°、方位角１１９．６°とした。
画像全体のヒストグラムを計算したところ、太陽仰角が７１°と高いこともあり、Ｓｈａｄｉｎｇ値による解析対象地の絞り込みは行わなかった。
【００５３】
ＩＶ．オブジェクトサイズの検討
ＳＰＯＴデータ（Ｂａｎｄ１〜４およびＰａｎ）により、ｅＣｏｇｎｉｔｉｏｎを利用してオブジェクト生成を行い、スケールパラメータについて検討を行った。オブジェクト生成のテストサイトを、別子山事業区の含まれるデータの南半分とした。
スケールパラメータ以外の設定はＳｈａｐｅ ｆａｃｔｏｒ＝０．２，Ｃｏｍｐａｃｔｎｅｓｓ＝０．５、Ｓｍｏｏｔｈｎｅｓｓ＝０．５とし、スケールパラメータを１５、２０、２５、３５、５０の５通りに設定してオブジェクトを生成した。画像上で判読可能な林分の領域とオブジェクト境界、生成されるオブジェクト数を考慮して、スケールパラメータを２５としてオブジェクト分類を行うのが適当であると結論した（図２）。
【００５４】
Ｖ．グラウンドトゥルース
２００８年１１月３日〜６日に、住友林業（株）新居浜山林事業所および住友林業フォレストサービス（株）新居浜山林事業所の協力の下、別子山・高藪・大永山の各事業区を中心に現地調査を行い、衛星画像と現地状況の対応を確認した。
調査には、ＧＰＳ受信機（Ｇａｒｍｉｎ ＧＰＳｍａｐ６０ＣＳＸ）、ＧＰＳカメラ（Ｒｉｃｏｈ Ｃａｐｌｉｏ ５００ＳＥおよびＢｌｕｅＴｏｏｔｈ ＧＰＳ受信機ＨＡＩＣＯＭ ＨＩ−４０８ＢＴ）および衛星画像を持参し、調査軌跡の記録、撮影写真の位置情報取得を行った。
【００５５】
ＶＩ．トレーニングエリア選択
まずｅＣｏｇｎｉｔｉｏｎでＳＰＯＴデータ（Ｂａｎｄ１〜４、Ｐａｎ）全域をスケールパラメータ２５としてオブジェクトを生成した。オブジェクト数は１１７，４５７であった。次に、ｅＣｏｇｎｉｔｉｏｎで生成されたオブジェクトの境界ベクタをＳｈａｐｅファイルに変換し、ＴＮＴｍｉｐｓにインポートして以下の作業を行った。Ｓｈａｐｅファイル変換の際には、ＳＰＯＴデータの５バンドに加えて５種類の指数画像についても各オブジェクトにおける平均値をＣＳＶファイルおよびＤＢＦファイルとして出力した。
ＴＮＴｍｉｐｓ上でグラウンドトゥルース調査地点を中心に画像判読を行いながら、オブジェクトポリゴンにクラス属性を付与する手法でトレーニングエリアを選択した。分類クラスは、非植生項目として「道路」および「水域」、植生項目として「スギ」、「スギ（緑）」、「ヒノキ」、「スギ（疎）」、「ヒノキ（疎）」、「ヒノキ（幼齢）」、「カラマツ」、「針広混交」、「ザツ（広葉樹）」、「マツ天然」、「草地および裸地」、「ササ地」の１４項目で、合計で４７９オブジェクトとなった。このうち、「スギ（青）」、「カラマツ」および「ヒノキ（幼齢）」の３クラスに関しては、トレーニングエリアとなるオブジェクトが１０個以下しかなかったため、以降の分類作業では分類項目としなかった。
４７９のトレーニングエリアを、分類に使用するための分類用セットと、精度評価に使用するための評価用セットがおよそ６：４になるようにランダムに振り分けた。表−２に各クラスのトレーニングエリア、分類用セット、評価用セットのオブジェクト個数を示す。
【００５６】
【表２】
【００５７】
ＶＩＩ．階層的分類法による土地被覆分類
（１）手法
トレーニングエリアとして指定されたオブジェクト番号をＴＮＴｍｉｐｓからＣＳＶファイルに出力し、ｅＣｏｇｎｉｔｉｏｎからは全オブジェクトについて以下に示す２８種類のパラメータとオブジェクト番号をＣＳＶファイルに出力し、それらをＥｘｃｅｌに読み込み作業を行った：Ｂａｎｄ１平均値、Ｂａｎｄ２平均値、Ｂａｎｄ３平均値、Ｂａｎｄ４平均値、Ｐａｎ平均値、ＮＤＶＩ平均値、ＮＤＷＩ平均値、Ｈｕｅ平均値、Ｓａｔｕｒａｔｉｏｎ平均値、Ｉｎｔｅｎｓｉｔｙ平均値、ＮＤＩ＿ＳＷＲ平均値、Ｒａｔｉｏ＿ＮＩＲ平均値、Ｒａｔｉｏ＿ＳＷＩＲ平均値、標高平均値、Ｂａｎｄ１標準偏差、Ｂａｎｄ２標準偏差、Ｂａｎｄ３標準偏差、Ｂａｎｄ４標準偏差、Ｐａｎ標準偏差、ＮＤＶＩ標準偏差、ＮＤＷＩ標準偏差、Ｈｕｅ標準偏差、Ｓａｔｕｒａｔｉｏｎ標準偏差、Ｉｎｔｅｎｓｉｔｙ標準偏差、ＮＤＩ＿ＳＷＲ標準偏差、Ｒａｔｉｏ＿ＮＩＲ標準偏差、Ｒａｔｉｏ＿ＳＷＩＲ標準偏差、標高標準偏差。なお、平均値、標準偏差とも各オブジェクト内での値である。
各パラメータの分類クラス毎の平均値・最大値・最小値・標準偏差のグラフ、任意の２
パラメータを組み合わせた散布図を作成し、１４種類の土地被覆項目の特性を検討した。
【００５８】
（２）メンバーシップ関数
検討の結果、表−３に示すような階層構造およびメンバーシップ関数を設定した。それぞれの関数および分類項目との関係を図−６（図３Ａ〜図３Ｎ）に示す。図３Ｂ〜図３Ｎにおいて図示された直線は、各メンバーシップ関数を表すものである。すなわち、当該直線によって、各分類項目の２群が属する座標平面上の領域が示されている。
なお、この中で「マツ天然」に関しては、「スギ」と「マツ天然」が衛星画像データのみでの判別が不可能であったため、現地調査での知見をもとに、「マツ天然」クラスは尾根（稜線）部のみしか出現しないという仮定を設けた。ＤＥＭデータから以下のアルゴリズム（ｈｔｔｐ：／／ｗｗｗ．ｏｐｅｎｇｉｓ．ｃｏ．ｊｐ／ｈｔｍ／ｋａｋｏ＿ｍａｉｌ／ｍａｉｌ＿ｍａｇ＿１３７．ｈｔｍ）を用いて尾根部を判別し、メンバーシップ関数に組み込んだ。
（ｉ）ＴＮＴｍｉｐｓのＲａｓｔｅｒ／Ｅｌｅｖａｔｉｏｎ／Ｗａｔｅｒｓｈｅｄ．．．により、谷線を表すＳＴＤＦＬＯＷＰＡＴＨベクタと流域界を表すＳＴＤＢＡＳＩＮベクタを作成する。
（ｉｉ）ＴＮＴｍｉｐｓのＧｅｏｍｅｔｒｉｃ／Ｃｏｍｐｕｔｅ／Ｄｉｓｔａｎｃｅ Ｒａｓｔｅｒ．．．により、（ｉ）で作成した二つのベクタからの距離ラスタを発生させ、それぞれをＤＩＳＴ＿ＢＡＳＩＮとＤＩＳＴ＿ＦＬＯＷとする。
（ｉｉｉ）以下のジオフォーミュラにより尾根部、斜面部、谷部を分ける。
ｄｉｓｔ＿ｐａｒａ＝２；
ｉｆ（ｓｑｒｔ（（ＤＩＳＴ＿ＢＡＳＩＮ−ＤＩＳＴ＿ＦＬＯＷ）＾２）＜（（（ＤＩＳＴ＿ＢＡＳＩＮ＋ＤＩＳＴ＿ＦＬＯＷ）／２）／ｄｉｓｔ＿ｐａｒａ））｛ＳｌｏｐｅＺｏｎｅ＝２；｝ｅｌｓｅ｛ｉｆ（ＤＩＳＴ＿ＢＡＳＩＮ＜ＤＩＳＴ＿ＦＬＯＷ）｛ＳｌｏｐｅＺｏｎｅ＝１；｝ｅｌｓｅ｛ＳｌｏｐｅＺｏｎｅ＝３；｝｝
このアルゴリズムで、ＳｌｏｐｅＺｏｎｅが１が尾根部を表す。
【００５９】
また、分類の際に「針広混交」クラスは、Ｓａｔｕｒａｔｉｏｎの値により３段階に区分を行った。
Ｓａｔｕｒａｔｉｏｎが０〜５０の場合は針葉樹混交率が高い「ＮＬ＿ｍｉｘｅｄ Ｈ」、Ｓａｔｕｒａｔｉｏｎが５０〜１５０の場合は針葉樹混交率が中庸の「ＮＬ＿ｍｉｘｅｄ Ｍ」、Ｓａｔｕｒａｔｉｏｎが１５０〜１９０の場合は針葉樹混交率の低い「ＮＬ＿ｍｉｘｅｄ Ｌ」とした。現地調査やトレーニングエリアのクラスでは混交率の差について考慮していないため、ここで設けた混交率による３クラスはあくまでもＳａｔｕｒａｔｉｏｎの大小による分類クラスである。なぜ、このように混交率クラスを設けたかについては、グラフ（図−６（ｅ）（図３Ｆ））上ではＳａｔｕｒａｔｉｏｎが同じ「針広混交」クラスでも大きく異なっていたこと、Ｓａｔｕｒａｔｉｏｎ値が針葉樹クラスでは小さく「ザツ（広葉樹）」クラスでは大きな値を持っており、以下の図−６（ｅ）（図３Ｆ）に示すように針葉樹クラスから「針広混交」クラス、「ザツ（広葉樹）」クラスまでＳａｔｕｒａｔｉｏｎ値が連続的であるように見えたためである。
【００６０】
【表３】
【００６１】
（３）分類結果と精度評価
表−３に示したメンバーシップ関数を用いて画像全体を分類した結果を図−８（図４）に示す。
トレーニングエリアの評価用セットを用いたエラーマトリクスを表−４に示す。この中で、「針広混交」クラスについては、前述の通り混交率によるトレーニングエリアがないので、混交率Ｍ（ＮＬ＿ｍｉｘｅｄ Ｍ）を「針広混交」クラスとした形で分類精度の計算を行った。また、「Ｎｏｎ＿ｖｅｇｅｔａｔｉｏｎ」：非植生は、「道路」および「水域」をまとめたクラスである。なお、数値はピクセル数である。
【００６２】
【表４】
【００６３】
全体の分類精度が８１．９％、カッパ係数も７８．２％であり、評価用セットにおいては分類は妥当に行われたと判断できる。
針葉樹種（密）については分類精度がヒノキ８９％、スギ９８％、図化精度がヒノキ９１％、スギ１００％と非常に高い。ヒノキと分類されたオブジェクトには、「針広混交Ｈ」および「針広混交Ｍ」が若干含まれているが、スギと分類されたオブジェクト誤分類は全くなく、スギについての信頼度がかなり高いことが推測できる。針葉樹種（疎）については、ヒノキは「草地および裸地」「針広混交Ｍ」「針広混交Ｈ」への誤分類および未分類があり分類精度５１％であるが、スギは誤分類が全くない。分類結果では他のトレーニングエリアが針葉樹種（疎）と分類されたものが全くないため、図化精度はどちらも１００％となっている。
【００６４】
分類項目をまとめて、より大分類にすることにより精度がどのように変化するかを検討した。ここでは、「Ｎｏｎ＿ｓｔｏｃｋｅｄ」：未立木地（「非植生」、「草地および裸地」および「ササ」）、「ザツ（広葉樹）」、「針広混交」（「針広混交」および「マツ天然」）、「スギ」（「スギ」および「スギ（疎）」）、「ヒノキ」（「ヒノキ」および「ヒノキ（疎）」）の５項目にまとめた。分類精度を表−５に示す。
全体精度は９１％、「ザツ」を除く分類項目の分類精度は約９０％以上、針葉樹種項目の図化精度も９０％以上となった。
【００６５】
【表５】
【００６６】
ＶＩＩＩ．非階層的分類法による土地被覆分類
トレーニングエリアの分類用セットを教師データとし、階層的分類法で使用した１０パラメータ（地形条件を除く）、Ｂａｎｄ１、Ｂａｎｄ３、Ｂａｎｄ４、Ｐａｎ、ＮＤＶＩ、ＮＤＷＩ、Ｈｕｅ、Ｓａｔｕｒａｔｉｏｎ、Ｉｎｔｅｎｓｉｔｙ、Ｒａｔｉｏ＿ＳＷＲを用いて、ｅＣｏｇｎｉｔｉｏｎのＮｅａｒｅｓｔ Ｎｅｉｇｈｂｏｒ（ＮＮ）法により土地被覆分類を行った。全域の分類結果を図−１０（図５）に示す。
また、トレーニングエリアの評価用セットを用いたエラーマトリクスを表−６に示す。数値はピクセル数である。
【００６７】
【表６】
【００６８】
全体の分類精度は８６％、カッパ係数も８３％であり、分類精度に関してはメンバーシップ関数による分類結果よりもよい結果となった。
針葉樹種（密）に関しては、分類精度がヒノキ８４％、スギ９８％、図化精度も同様の値である。「ヒノキ」のトレーニングエリアは「スギ」、「ヒノキ（疎）」および「針広混交」に誤分類が行われているが、「スギ」のトレーニングエリアは「スギ（疎）」への誤分類のみである。一方、「ヒノキ」と分類されたオブジェクトには「ザツ（広葉樹）」、「マツ天然」および「針広混交」が若干含まれており、「スギ」と分類されたオブジェクトには「ヒノキ」および「マツ天然」が若干含まれている。メンバーシップ関数による分類と比較すると、メンバーシップ関数による分類の方がクラス間の混同が少なかった。
【００６９】
前章と同様に、分類項目をまとめてより大分類にすることにより精度がどのように変化するかを検討した。ここでは、「Ｎｏｎ＿ｓｔｏｃｋｅｄ」：未立木地（「非植生」、「草地および裸地」および「ササ」）、「ザツ（広葉樹）」、「針広混交」（「針広混交」および「マツ天然」）、「スギ」（「スギ」および「スギ（疎）」）、「ヒノキ」（「ヒノキ」および「ヒノキ（疎）」）の５項目にまとめた。分類精度を表−７に示す。
全体の分類精度はやや向上したが、メンバーシップ関数による分類精度（表−５）の方が、「針広混交」の分類精度および「ヒノキ」の図化精度がより高かった。
【表７】
【００７０】
ＩＸ．考察
分類結果の画像とＳＰＯＴ画像、現地の状況を比較検討すると、今回の２手法ではメンバーシップ関数による分類結果の方が、全体的により妥当であると判断された。
今回の精度評価は、選択したトレーニングエリアをランダムに２つに分け、その片方の評価用セットで行った。トレーニングエリアは、現地調査を参考に画像を判読しながら行うが、比較的判別しやすいところを選んだ。
【産業上の利用可能性】
【００７１】
本発明によれば、森林の植生分類図をより高精度で作成することができる。したがって、本発明は、林業および環境関連産業等の発展に寄与するところ大である。
【特許請求の範囲】
【請求項１】
中間赤外波を含む電磁波による衛星画像を用いることを含む、森林の植生分類図を作成する方法。
【請求項２】
中間赤外波の波長が１．３μｍ〜２．０μｍである、請求項１に記載の方法。
【請求項３】
電磁波がさらに近赤外波および可視光を用いる、請求項１または２に記載の方法。
【請求項４】
衛星画像から、２本以上の単木を含むオブジェクトを作成する工程を含む、請求項１〜３のいずれかに記載の方法。
【請求項５】
分類パラメータを用いてオブジェクトの土地被覆分類を行う工程を含み、該分類パラメータが中間赤外波に対する反射輝度値を用いるパラメータを含む、請求項１〜４のいずれかに記載の方法。
【請求項６】
分類パラメータが正規化水指数ＮＤＷＩ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲおよびＲａｔｉｏ＿ＳＷＩＲからの１種または２種以上である、請求項５に記載の方法。
【請求項７】
階層的分類法によるオブジェクトの土地被覆分類を行う工程を含む、請求項５または６に記載の方法。
【請求項８】
請求項１〜７のいずれかに記載の方法を用いて作成された森林の植生分類図。
【請求項１】
中間赤外波を含む電磁波による衛星画像を用いることを含む、森林の植生分類図を作成する方法。
【請求項２】
中間赤外波の波長が１．３μｍ〜２．０μｍである、請求項１に記載の方法。
【請求項３】
電磁波がさらに近赤外波および可視光を用いる、請求項１または２に記載の方法。
【請求項４】
衛星画像から、２本以上の単木を含むオブジェクトを作成する工程を含む、請求項１〜３のいずれかに記載の方法。
【請求項５】
分類パラメータを用いてオブジェクトの土地被覆分類を行う工程を含み、該分類パラメータが中間赤外波に対する反射輝度値を用いるパラメータを含む、請求項１〜４のいずれかに記載の方法。
【請求項６】
分類パラメータが正規化水指数ＮＤＷＩ、ＮＤＩ＿ＳＷＲ、Ｒａｔｉｏ＿ＮＩＲおよびＲａｔｉｏ＿ＳＷＩＲからの１種または２種以上である、請求項５に記載の方法。
【請求項７】
階層的分類法によるオブジェクトの土地被覆分類を行う工程を含む、請求項５または６に記載の方法。
【請求項８】
請求項１〜７のいずれかに記載の方法を用いて作成された森林の植生分類図。
【図１】
【図２】
【図３Ａ】
【図３Ｂ】
【図３Ｃ】
【図３Ｄ】
【図３Ｅ】
【図３Ｆ】
【図３Ｇ】
【図３Ｈ】
【図３Ｉ】
【図３Ｊ】
【図３Ｋ】
【図３Ｌ】
【図３Ｍ】
【図３Ｎ】
【図４】
【図５】
【図２】
【図３Ａ】
【図３Ｂ】
【図３Ｃ】
【図３Ｄ】
【図３Ｅ】
【図３Ｆ】
【図３Ｇ】
【図３Ｈ】
【図３Ｉ】
【図３Ｊ】
【図３Ｋ】
【図３Ｌ】
【図３Ｍ】
【図３Ｎ】
【図４】
【図５】
【公開番号】特開２０１１−２４４７１（Ｐ２０１１−２４４７１Ａ）
【公開日】平成２３年２月１０日（２０１１．２．１０）
【国際特許分類】
【出願番号】特願２００９−１７２７３７（Ｐ２００９−１７２７３７）
【出願日】平成２１年７月２４日（２００９．７．２４）
【新規性喪失の例外の表示】特許法第３０条第１項適用申請有り 掲載年月日：２００９年３月１７日 掲載アドレス：ｈｔｔｐ：／／ｗｗｗ．ｊｓｔａｇｅ．ｊｓｔ．ｇｏ．ｊｐ／ａｒｔｉｃｌｅ／ｊｆｓｃ／１２０／０／５２０／＿ｐｄｆ／−ｃｈａｒ／ｊａ／
【出願人】（０００１８３４２８）住友林業株式会社 (540)
【出願人】（５０９２０９９９３）
【Ｆターム（参考）】
【公開日】平成２３年２月１０日（２０１１．２．１０）
【国際特許分類】
【出願日】平成２１年７月２４日（２００９．７．２４）
【新規性喪失の例外の表示】特許法第３０条第１項適用申請有り 掲載年月日：２００９年３月１７日 掲載アドレス：ｈｔｔｐ：／／ｗｗｗ．ｊｓｔａｇｅ．ｊｓｔ．ｇｏ．ｊｐ／ａｒｔｉｃｌｅ／ｊｆｓｃ／１２０／０／５２０／＿ｐｄｆ／−ｃｈａｒ／ｊａ／
【出願人】（０００１８３４２８）住友林業株式会社 (540)
【出願人】（５０９２０９９９３）
【Ｆターム（参考）】
[ Back to top ]