【GEE】CARTで土地被覆分類を行う【衛星データ解析】

SHARE

これまで、Google Earth Engine (以下GEE)を利用した衛星データの主な利用方法について紹介してきましたが、さらに発展させて機械学習を利用した解析について紹介します。

機械学習を利用することで衛星データから土地被覆分類や時系列予測などができるようになります。

土地被覆分類に利用される手法を中心に、前回はランダムフォレストについて紹介しました。今回は別の手法のCART法をGEEで利用する方法について紹介します。

動作検証済み環境

Python 3.11.4, earthengine-api 0.1.377
Google Earth Engineを用いた衛星データ解析入門 ~タイムラプスを作成し地球の環境変動を観察~

Pythonを使った衛星データ解析の初心者向け技術書を販売中！

¥2,500 → ¥800 今なら約70%OFF!!

Google Earth Engineを使えば簡単に
自宅にいながら地球環境の変動を観察できます

技術書ページへ

土地被覆分類とは？

衛星画像や航空写真などのリモートセンシングデータを分析し、地表のさまざまなカバータイプ（例えば、森林、農地、水域、都市地域など）を識別するプロセスです。土地被覆分類は地球環境のモニタリング、土地利用計画、生態系の管理、気候変動の研究など、多くの応用分野で利用されます。

以前の記事で利用した環境省生物多様性センターの植生図をはじめ、国土地理院の地理院地図、産総研のシームレス地質図などが作成され公開されています。

GEEには土地被覆分類を行うためのオブジェクトとして**ee.Classifier**というものが用意されています。これを利用することで機械学習を利用した土地被覆分類を行うことができます。

今回は、CARTを利用して、ローカルな土地被覆分類を行ってみます。

CARTとは？

CART（Classification And Regression Trees）は、教師あり学習アルゴリズムの一種で、分類問題と回帰問題の両方に使用することができます。このアルゴリズムは、データを最適に分割する決定木を作成することで、複雑なデータセットのパターンを捉え、予測を行います。

以下がCARTの主な特徴です。

  1. 二分木構造： CARTは、データを二分岐することにより決定木を成長させていきます。各ノード（分岐点）で、最も有効な分割を見つけ出すことにより、データをより均一なサブセットに分割します。このプロセスは、予め設定された停止条件を満たすまで繰り返されます。
  2. 分類と回帰の対応： CARTはその名が示す通り、分類問題だけでなく回帰問題にも適用可能です。分類木は、カテゴリカルな目的変数を予測するために使用され、回帰木は数値的な目的変数の予測に使用されます。
  3. 過学習の制御： CARTでは、木の成長後に枝刈りを行うことで過学習を防ぎます。枝刈りは、木の複雑さとトレーニングデータへの適合度とのバランスを取ることにより、モデルの一般化能力を向上させます。
  4. 特徴の選択： 分割の際、CARTは使用可能な特徴の中から最適な分割を選択するため、特徴選択のプロセスを自然に組み込んでいます。これにより、予測に不要な特徴がモデルから自動的に除外されることがあります。
  5. 解釈の容易さ： 決定木は視覚的に表現することができ、生成された木の構造を分析することで、どの特徴が予済みに重要な役割を果たしているかを理解しやすいです。これは、モデルの意思決定過程を説明する場合に特に有用です。

CARTはその柔軟性と解釈の容易さから、多くの分野でのデータ分析や予測モデリングにおいて広く利用されています。

CARTとランダムフォレストの違い

  • 複雑さと精度：ランダムフォレストは、複数の決定木を組み合わせることでCARTよりも高い予測精度を提供しますが、モデルの複雑さと計算コストが増加します。
  • 過学習の耐性：ランダムフォレストは過学習に対してより強い耐性を持ちますが、CARTは適切な枝刈りを行うことで過学習を制御できます。
  • 解釈性：CARTの単純な二分木構造は解釈が容易ですが、ランダムフォレストの複数の決定木を組み合わせる方法は解釈を難しくします。

CARTは解釈性が重要な場合や単純な問題に適しており、ランダムフォレストは精度を最大化したい複雑な問題や過学習に対する耐性が必要な場合に適しています。

今回はランダムフォレストとの比較はせず、GEEでCARTを利用する単純な方法について紹介します。

CARTによる土地被覆分類の実装

今回もGoogleのサンプルコードを参考に以下のように実装しました。今回はgeemapというライブラリを使う関係でJupyterLab上で実行します。

ESAのWorldCoverというデータセットを教師データとして、函館市付近の土地被覆分類を行ってみます。

import ee
import geemap

ee.Initialize()

# Sentinel-2の表面反射率イメージの読み込み
img = ee.Image('COPERNICUS/S2_SR_HARMONIZED/20200417T012649_20200417T013323_T54TVM') \
    .select('B.*')

# ESAのWorldCoverデータセットを使用
lc = ee.Image('ESA/WorldCover/v100/2020')

# 土地被覆クラスのリマッピング
class_values = [10, 20, 30, 40, 50, 60, 70, 80, 90, 95, 100]
remap_values = ee.List.sequence(0, 10)
label = 'lc'
lc = lc.remap(class_values, remap_values).rename(label).toByte()

# データセットの結合
sample_img = img.addBands(lc)

# トレーニングデータのサンプリング
roi = ee.Geometry.Rectangle([140.625, 41.721, 140.868, 41.862])
sample = sample_img.stratifiedSample(
    numPoints=5000,
    classBand=label,
    region=roi,
    scale=10,
    geometries=True
)

# データセットをトレーニングセットとバリデーションセットに分割
sample = sample.randomColumn()
training_sample = sample.filter('random <= 0.8')
validation_sample = sample.filter('random > 0.8')

# 分類器のトレーニング (CART)
trained_classifier_cart = ee.Classifier.smileCart(10).train(
    features=training_sample,
    classProperty=label,
    inputProperties=img.bandNames()
)

# トレーニングサンプルの精度評価 (CART)
train_accuracy_cart = trained_classifier_cart.confusionMatrix()
print('CART Training accuracy:', train_accuracy_cart.accuracy().getInfo())

# バリデーションサンプルの分類と精度評価 (CART)
validation_sample_cart = validation_sample.classify(trained_classifier_cart)
validation_accuracy_cart = validation_sample_cart.errorMatrix(label, 'classification')
print('CART Validation accuracy:', validation_accuracy_cart.accuracy().getInfo())

# イメージの分
classified_img_cart = img.classify(trained_classifier_cart)

# レイヤーの色パレットの設定
class_vis = {
    'min': 0,
    'max': 10,
    'palette': [
        '006400',
        'ffbb22',
        'ffff4c',
        'f096ff',
        'fa0000',
        'b4b4b4',
        'f0f0f0',
        '0064c8',
        '0096a0',
        '00cf75',
        'fae6a0',
    ],
}

# 地図の作成とレイヤーの追加
m = geemap.Map()
m.set_center(140.75, 41.79, 12)
m.addLayer(classified_img_cart, class_vis, 'CART Classified')
m

ref: ee.Classifier.smileCart

実行結果

まず、Accuracy（精度）を確認します。今回は以下のようになりました。

これらの値については以下の通りです。

  1. トレーニング精度（Training accuracy）: 約58.52%
    • これは、トレーニングデータセットに対する分類器の性能を示しており、分類器がトレーニング中に使用したデータをどれだけうまく分類できたかを意味します。
    • 一般に、トレーニング精度が高いことは良いことですが、過剰適合（オーバーフィッティング）の可能性も考慮する必要があります。これは、モデルがトレーニングデータに対して過度に最適化されており、新しい未知のデータに対してはうまく機能しない状況を指します。
  2. バリデーション精度（Validation accuracy）: 約58.48%
    • バリデーション精度は、トレーニングに使用されていないデータ（バリデーションデータセット）に対する分類器の性能を示します。
    • バリデーション精度が低い場合、モデルの一般化能力に問題があるか、データセットの特徴が十分でない可能性があります。

Accuracy（精度）の確認が完了し、問題ない場合は作成した土地被覆分類図を確認します。

今回は以下のような地図が表示されていれば成功です。

緑色が植生、赤色が市街地、灰色が建造物など、黄色が裸地や砂地、ピンクが農地、青が水域に分類されています。

前回のランダムフォレストと比較してどちらのスコアも下がっています。これは単純比較でこの処理方法における土地被覆分類はランダムフォレストのほうが適していると言えます。ただし、これはあくまで訓練に用いるデータ（今回はESAのデータセット）や決定木の数によっても変わってくるので一概には言えません。

参考文献のShetty, S. (2019)ではGEEを利用した土地被覆分類ではランダムフォレストが高精度になることから推奨されています。

コードの解説

上に書いたソースコードの解説をしていきます。

ライブラリのインポートと初期化:

import ee
import geemap

ee.Initialize()
  • ee と geemap をインポートします。ee は Earth Engine API、geemap は地図可視化用のライブラリです。
  • ee.Initialize() で Earth Engine のセッションを初期化します。

2.　データの取得と前処理:

# Sentinel-2の表面反射率イメージの読み込み
img = ee.Image('COPERNICUS/S2_SR_HARMONIZED/20200417T012649_20200417T013323_T54TVM') \
    .select('B.*')

# ESAのWorldCoverデータセットを使用
lc = ee.Image('ESA/WorldCover/v100/2020')

# 土地被覆クラスのリマッピング
class_values = [10, 20, 30, 40, 50, 60, 70, 80, 90, 95, 100]
remap_values = ee.List.sequence(0, 10)
label = 'lc'
lc = lc.remap(class_values, remap_values).rename(label).toByte()

# データセットの結合
sample_img = img.addBands(lc)
  • Sentinel-2の表面反射率イメージを特定の日付で読み込み、B.* バンド（全てのバンド）を選択します。
  • ESAのWorldCoverデータセットを読み込み、土地被覆クラスを0から始まる連続したシリーズにリマッピングします。
  • class_values は、WorldCoverデータセットに含まれる元の土地被覆クラスを定義しています（例えば、10は樹木、20は草地など）。
  • remap_values は、これらのクラス値を0から始まる連続した整数にリマップするために使用される新しい値のシリーズです。
  • lc.remap(class_values, remap_values).rename(label).toByte() は、元のクラス値を新しい連続した値にリマップし、それを新しいバンド名**'lc'**（土地被覆）でイメージに追加して、バイト型に変換しています。
  • sample_img = img.addBands(lc) は、Sentinel-2のイメージ（img）にWorldCoverの土地被覆イメージ（lc）をバンドとして追加しています。これにより、反射率データと土地被覆データが組み合わされた新しいイメージが作成されます。

3.　トレーニングデータの準備:

# トレーニングデータのサンプリング
roi = ee.Geometry.Rectangle([140.625, 41.721, 140.868, 41.862]) 
sample = sample_img.stratifiedSample(
    numPoints=5000,
    classBand=label,
    region=roi,
    scale=10,
    geometries=True
)

# データセットをトレーニングセットとバリデーションセットに分割
sample = sample.randomColumn()
training_sample = sample.filter('random <= 0.8')
validation_sample = sample.filter('random > 0.8')
  • 特定の地理的領域（roi）内のデータから5000点のサンプルを抽出します。
  • 抽出されたサンプルにランダムな列を追加し、80%をトレーニング用、20%をバリデーション用に分割します。

4.　ランダムフォレスト分類器のトレーニング:

# 分類器のトレーニング (CART)
trained_classifier_cart = ee.Classifier.smileCart(10).train(
    features=training_sample,
    classProperty=label,
    inputProperties=img.bandNames()
)
  • ee.Classifier.smileCart(10)で決定木を10本に設定し、分類器をトレーニングサンプルからトレーニングします。

5.　トレーニングとバリデーション精度の計算:

# トレーニングサンプルの精度評価 (CART)
train_accuracy_cart = trained_classifier_cart.confusionMatrix()
print('CART Training accuracy:', train_accuracy_cart.accuracy().getInfo())

# バリデーションサンプルの分類と精度評価 (CART)
validation_sample_cart = validation_sample.classify(trained_classifier_cart)
validation_accuracy_cart = validation_sample_cart.errorMatrix(label, 'classification')
print('CART Validation accuracy:', validation_accuracy_cart.accuracy().getInfo())
  • トレーニングサンプルに基づく混同行列を作成し、その精度を計算します。
  • バリデーションサンプルを分類し、その精度を計算します。

6.　分類されたイメージの生成と表示:

# イメージの分
classified_img_cart = img.classify(trained_classifier_cart)

# レイヤーの色パレットの設定
class_vis = {
    'min': 0,
    'max': 10,
    'palette': [
        '006400',
        'ffbb22',
        'ffff4c',
        'f096ff',
        'fa0000',
        'b4b4b4',
        'f0f0f0',
        '0064c8',
        '0096a0',
        '00cf75',
        'fae6a0',
    ],
}

# 地図の作成とレイヤーの追加
m = geemap.Map()
m.set_center(140.75, 41.79, 12)
m.addLayer(classified_img_cart, class_vis, 'CART Classified')
m
  • トレーニングされた分類器を使用して、元のイメージを分類します。
  • 分類されたイメージを地図に追加し、特定の色パレットで表示します。
    • '006400'（濃い緑色）: 植生や森林。
    • 'ffbb22'（オレンジ色）: 農地や草地。
    • 'ffff4c'（黄色）: 裸地や砂地など、植生が少ない地域。
    • 'f096ff'（ピンク色）: 農地、湿地や水域。
    • 'fa0000'（赤色）: 都市化された地域や市街地。
    • 'b4b4b4'（灰色）: 都市や建造物。
    • 'f0f0f0'（薄灰色）: 雪や氷。
    • '0064c8'（青色）: 水域。
    • '0096a0'（シアン）: 河川や湖。
    • '00cf75'（薄緑色）: 植生や農地。
    • 'fae6a0'（ベージュ色）: 砂漠や裸地。
  • 最後にgeemapで地図を生成します。

Accuracyが低い場合は、まずは3.のnumPointsや4.のsmileCart()の数値を調整してみてください。

最後に

Google Earth Engine（GEE）のAPIを活用することで、事前にチューニングされた機械学習モデルを容易に利用できるようになります。しかし、高品質な土地被覆分類を達成するためには、適切な訓練データの選択やパラメータの設定が重要です。

リモートセンシングの分野では、機械学習を用いたデータ分析が近年、活発化しています。GEEの利用は入門者にとっても敷居が低く、非常に便利です。

この機会にぜひ初めて見てはどうでしょうか。

参考文献

Google Earth Engineを用いた衛星データ解析入門 ~タイムラプスを作成し地球の環境変動を観察~

Pythonを使った衛星データ解析の初心者向け技術書を販売中！

¥2,500 → ¥800 今なら約70%OFF!!

Google Earth Engineを使えば簡単に
自宅にいながら地球環境の変動を観察できます

技術書ページへ

コメントを残す

メールアドレスが公開されることはありません。 ※ が付いている欄は必須項目です

前の記事

プベルル酸の毒性予測と治療薬への可能性の模索【in silic…

次の記事

【PySCF】Pythonで始める量子化学計算【一点計算】
【PySCF】Pythonで始める量子化学計算【一点計算】