寒い日が続いていますが、もうしばらくすると暖かい春がやってきます。春といえば、花粉の時期ですね。私は花粉症がひどいので、春の訪れとともにマスク必須の生活が始まります。

そこで、本記事では花粉の少ない街を衛星画像を使って探すことにチャレンジします！

※本記事は宙畑メンバーが気になったヒト・モノ・コトを衛星画像から探す不定期連載「宇宙データ使ってみた-Space Data Utilization-」の第19弾です。まだまだ修行中の身のため至らない点があるかと思いますがご容赦・アドバイスいただけますと幸いです！

（１）解析の全体像

花粉の少ない街を探す！という今回のチャレンジですが、まずは何のデータで花粉を測定するのが良いか考えてみます。

調べてみると、環境省が花粉観測システム「はなこさん」というサイトを公開していました。そこで紹介されていたのが、上図のように直接花粉を採取する方法です。

他にも調べて見ましたが、ほとんどが直接花粉を採取する方法で、なかなか広域で花粉の状態を観測することはできなさそうです。

そこで、我らが衛星データの出番です！衛星であれば、広域で様子を観測することができます。

衛星データの種類は様々ありますが、空気中の粒子の様子をみるとすると「ライダー（LIDER）」と呼ばれるセンサが最適です。

ライダーでの花粉観測は、地上での観測で研究が行われている様ですが、まだ発展途上の様です。

花粉がどこを飛んでいるのかを観測するレーザー技術とは？（夢ナビ）

そこで、今回は花粉そのものではなく、花粉を飛ばす大元であるスギに注目。どのあたりにスギが多いのかを調査し、スギの場所以外の要素（コンクリートが多いところだと被害が大きいなどよく聞きますね）と組み合わせて、直接花粉が測定できなくても、花粉の多さを推定できないかを考えます。

まとめると、以下の3つのステップに分けて解析を進めて行きたいと思います。

①衛星データからスギ（木の種類）が分かるか
②衛星データからスギが多い街・少ない街を調べる
③他のデータと合わせて、花粉が少ない街を見つける

今回の記事ではステップ①「衛星データからスギ（木の種類）が分かるか」について、調査していきます。

（２）モノの色が違う理由（原理）

これから衛星画像からスギの木を見つける方法を考えてみたいと思いますが、その前に人間の目でモノが見える仕組みについてご紹介したいと思います。

光の波長と色の関係～目に見えるモノが全てじゃない！～

光は電波の一種なので様々な波長があり、波長が610nm -750nmなら赤、波長が500nm-560nmなら緑といった具合に、光の波長に応じて私たちの目に見える色が変わります。私たちの目に見える範囲の光のことを「可視光（かしこう）」と言います。

そして、表で示すように私たちの目には見えない範囲にも光は存在しています。
赤外線や紫外線などと呼ばれる領域です。

光学センサを使うと、可視光以外の光も捉えることができます。

モノの色が違うのは「反射特性」の違い

私たちの目にモノが見えるのはモノに光が反射し、その光が目に入ることでモノが見えるためです。

モノに光が当たったとき、モノは全ての光を反射するのではなく、一部の光を吸収し一部の光を反射します。例えば、リンゴなら赤色の光を反射し、緑色や青色の光を吸収します。そのため、リンゴは赤く見えます。

モノ(リンゴなど)によって、いろいろな波長の光を反射/吸収しておりその特性（反射特性と呼びます）は異なるため、モノの色は異なって見えます。

そしてこれは、可視光の領域だけでなく、赤外線や紫外線の領域でも同じことが言えます。

樹種によって反射特性が異なる

上図は、様々な樹木に光を当てた時の樹木の光の反射特性を表しています。

例えば、ヒバは800nm~900nm間の波長をよく反射し、750nm以下の波長はあまり反射しないことが図から読み取れます。

また、それぞれの樹木でおおよその特性は同じですが、樹木の種類によって反射特性が異なることが分かります。

750nm以上の波長は人の目では見えないですが、専用のセンサを使うと見ることができ、それを使うと樹木の違いを区別できるということです。

衛星から木の反射特性を調べる

衛星に搭載されているセンサも、赤、青、緑といった可視光だけではなく、様々な波長帯の光を観測しています。それぞれの波長帯のことを「バンド」呼び、各バンドに対応した衛星画像を撮っています。

つまり、Sentinel-2衛星の場合、1場面の観測で12枚の衛星画像が生成されるということです。

例えば、ヨーロッパの地球観測衛星Sentinel-2の場合、観測バンドB1では上表の通り443nm付近の波長を観測しています。

先ほどの、樹木の反射特性で違いが顕著だった800～900nmはSentinel-2ではB8、B8a、B9に該当するということになります。

（３）衛星データからスギを見つける

さていよいよ、衛星データを使ってスギを識別してみたいと思います。

解析の流れは大きく以下の5ステップです。

① スギ林探し

まずは、学習と検証を行うためのスギ林を調査します。
今回は以下の2箇所を選定しました。

エリアⅠ（学習用）：大雄山のスギ

エリアⅡ（検証用）：高尾山のスギ

②衛星画像の入手

スギは図4より、800nm~900nmの波長で反射値が大きくなっています。従って、表2より観測バンドB8、B8a、B9で撮影した衛星写真を使うとスギが発見できそうです。

今回はSentinel-2という衛星の画像を使うことにしました。Sentinel-2は分解能が高い（=解像度の高い）画像を幅広いバンドで撮影することができ、無償で使えるためです。
反射特性は時期によっても変化するため、撮影時期も考慮する必要があります。今回は、スギの活動が最も活発と考えられる夏のデータとしました。

衛星画像は、ESAが民間に公開している衛星画像取得のプラットフォーム『Copernicus』から取得しました。

ここで簡単に、『Copernicus』から画像を取得する手順を紹介します。

尚、衛星データを利用する際の注意点として、衛星データの「処理レベル」が挙げられます。反射特性を得るためには反射率変換処理後のデータが適しています。

今回使用するCopernicusで得られる衛星データの処理レベルはL1Cであり、大気上端での反射率（Top of Atmosphere (TOA) Reflectance）のデータです（ユーザーガイド）。

処理レベルについて詳しく知りたい方は以下の記事を参照ください。

こちらが大雄山のスギ周辺の衛星写真です。

③反射特性の調査（学習）

エリアⅠ：大雄山のスギ周辺の衛星画像を用いて、スギの反射特性を学習します。
反射特性は、画像の各ピクセル値（画素値）を調べることで分かります。今回は、画像のそれぞの画素値をQGISというソフトを使って確認します。

QGISの地物情報を確認する機能を使うと画像に格納されている画素値を確認することができます。画素値とは、画素に含まれる各波長のデジタル値を表しています。

オレンジ点付近の各バンドにおける画素値は、先ほどダウンロードした画像ファイルをQGIS上にドラッグ&ドロップし、図6の手順に従うと確認できます。

オレンジ点付近で、数点画素値を取得した結果、各バンドにおける画素値は下記の様になります。

表3　大雄山のスギの反射特性

表3の様に画素値に幅があるのは、図7の様にスギが生えている箇所の色見が微妙に異なる為です。

④検証用データへの適用

表3の画素値の条件を満たす場所を、検証用のエリアⅡ：高尾山の衛星画像に適用してみます。
エリアⅡのB8、B8a、B9の画像をそれぞれ表3の条件に当てはまる部分だけ白く抜き出す処理を施します。この処理はQGISの「ラスタ」＞「ラスタ計算機」を使います。

⑤検証

④で抽出した場所に本当にスギがあるのかを確認します。

画像左上オレンジで〇をつけてある箇所が、高尾山のスギがあると確認できている箇所です。

残念ながら、オレンジ丸の周辺は白くなっておらず、高尾山のスギは捉える事が出来なかった事がわかります。なぜなのでしょうか…

⑥考察

理由仮説１：スギの反射特性は他の種と区別しづらい

図4の反射特性のグラフをもう一度良く見てみると、スギは他の樹種と反射特性が似ています。これでは、反射値の幅で抽出を試みても、なかなかスギだけを捉えるのは難しそうです。

理由仮説２：反射特性のばらつきを拾い切れていない

もう一つ理由として考えられるのが、反射特性の拾い方です。図6のように、画像をある点をクリックして拾っていく手法では、データ点の数に限りがあります。

点ではなく、面で値を取得する必要があると考えられます。

（４）衛星データからアカエゾマツを見つける！

理由仮説１の正しさを確認するために、他の樹木と比較し反射値が小さく、検出が比較的簡単そうなアカエゾマツを捉えてみることにしてみました。

①アカエゾマツ林探し

アカエゾマツでも同様に、生息エリアの調査から始めます。

アカエゾマツの植林地を探してみると北海道に多くあることが分かりました。下図の赤い部分がアカエゾマツが生息している部分になります。

アカエゾマツの場合、エリアを2つ見つけることができなかったため、上図の中のエリアを学習用と検証用に分け、解析を進めていくことにします。

②衛星画像の入手

スギの時と同様に、Copernicusより衛星画像を入手します。
今回は時期に依る特性も確認するため、雪が降り表面が見えなくなってしまう冬を除いた、春・夏・秋の画像を使用します。

表4　アカエゾマツ反射特性

③反射特性の調査（学習）

理由仮説２より、点ではなくアカエゾマツであると分かっているエリア全体の反射値を取得し、反射値のばらつきまで考慮した条件の設定を行います。

まずは、アカエゾマツと分かっているエリアの形のデータ（geojsonファイルと言います。）を作ります。geojsonファイルは、こちらのサイトを使って作成しました。

作成したgeojsonファイルは以下のように地図上に投影できます。

このgeojsonファイルを使って、対象の画像ファイルを切り出します。また画像ファイルの切り出しはバンド毎に行わなければなりません。
切り出した後の画像は以下のようになります

切り出した後の画像の反射値は、QGIS上で、対象のデータの[プロパティ]から得ることができます。結果は以下の通りです。

表5　アカエゾマツの反射特性

④検証用データへの適用

表5の条件で、検証用に残しておいたエリアⅡを確認してみます。

左側は私たち人間の目で見た色に近い赤青緑のバンドで合成した「True Color画像」で、右側の白黒画像の白い点の部分が、表5にの条件に適合する部分になります。

こうしてみると、True Color画像では同じように緑に見える部分でも、赤外線の領域の反射特性で抽出をしてみると、差があることが分かります。

そして、その差は夏の画像で最も大きく、秋になると抽出が難しくなるようです。

⑤検証

④の結果を先ほどの保護林位置を比較してみます。

図12の保護林位置図のアカエゾマツエリアを図中オレンジ色の線で記載しています。

季節で比較すると、最もうまく抽出できているのが夏の画像でした。

春・秋の画像では右上の畑の裸地と思われる部分も検出してしまっており、植物の中の乾燥等が影響を与えてしまっているのではないかと考えられます。

また、南側中央のエリアで検出できていないのは、True Color画像を見てみるとそもそもあまり緑が深くない様子です。

正解データとして用いた保護林のエリアについて、保護林以外にもアカエゾマツがある可能性と逆に保護林の中に他の樹木がある可能性があり、正解データ自体の精査をする必要があると思われます。

（５）まとめ

今回の記事では、衛星データを使ってスギを捉える事に挑戦しました。

今回分かったことは、

・スギが反射しやすい波長は他の樹木も反射する波長である為、1時期の近赤外線領域の衛星画像だけスギを見つけることは難しい。

・アカエゾマツなど反射特性が特徴的な樹木は、衛星画像を使って植物を抽出することができる。

です。

一方、今回の課題として考えられるのは、

・正解データの正確な位置情報の把握：

今回調べてみて分かったのは正確な位置情報まで含めた杉林の位置のオープンデータが意外とないということでした。
正確な教師データがもう少しあれば、条件の精度を高めることができるかもしれません。

・時期による反射特性の変化まで考慮した樹種判別：

後半のアカエゾマツでは、春夏秋の3時期の画像で解析を行いました。植物は時期によって反射特性が変化するため、どの時期の画像が最も識別しやすいのか、もしくは時期を組み合わせることにより詳細な識別が可能になるのか、検討の余地があると思います。

・他の波長でも検証してみる：

今回は、近赤外の領域に絞って調査を行いましたが、他の波長で特徴的な値となるかは検証してみる必要があります。特に、反射値の立ち上がりの赤の波長（Red Edgeと呼ばれています）などは試してみる価値があると思われます。

さらに、現時点で無料で手に入る衛星データではありませんが、ハイパースペクトルカメラと呼ばれる、100を超えるバンド（観測波長）を持つセンサHISUIが2019年12月に宇宙ステーションに設置されています（宇宙実証用ハイパースペクトルセンサ「HISUI」が「きぼう」船外実験プラットフォームへ設置されました/JAXA）。HISUIの情報が公開されれば、さらに詳細に調査が可能になると考えられます。

尚、同じ内容をTellus内に格納されているLandat-8という衛星データで実施することが可能です。Tellusを使えばもっと広範囲な解析をおこなうことができます。こちらのコードも近々別記事で公開予定です。

冒頭でお伝えした通り、本記事は3つのステップから構成されています。

①衛星データからスギ（木の種類）が分かるか
②衛星データからスギが多い街・少ない街を調べる
③他のデータと合わせて、花粉が少ない街を見つける

次回、衛星データでスギを見つけつつ、他のデータの組み合わせまで行けるのか、ご期待ください！（一緒にやってくれる方も絶賛募集中です！）