UAV 搭載型 LiDAR およびマルチスペクトルセンサを用いた地すべり活動の把握

北海道大学 笠井美青
北海道大学 早川裕弌
東京大学 飯塚浩太郎

令和 5 年 9 月

目次
1．はじめに ................................................................ 1
2.

3.

対象地すべりの概要 ..................................................... 2
2.1

オタリマップ地すべり ................................................ 3

2.2

研究林地すべり ...................................................... 4

計測 ................................................................... 5
3.1

LiDAR 計測 .......................................................... 5

3.2

マルチスペクトル計測 ................................................ 7

3.2.1. マルチスペクトル画像解析 ........................................ 7
3.2.2. オタリマップ地すべり ............................................ 8
3.2.3. 研究林地すべり .................................................. 9
3.3

VLF 計測 ........................................................... 103.3.1

地中センサー ................................................... 10

3.3.2

機材 ........................................................... 10

3.3.3

データ観測・処理方法 ........................................... 10

3.3.4

データの可視化 ................................................. 11

3.3.5

VLF 計測結果 ................................................... 114.

地すべり動態の把握 ................................................... 12
4.1

オタリマップ地すべり ............................................... 12

4.1.1. オタリマップにおける点群解析 ................................... 124.1.2勾配の変化 ..................................................... 124.1.3

PIV 解析........................................................ 134.2
5.

研究林地すべり ..................................................... 13

結果 .................................................................. 14
5.1

オタリマップ地すべり ............................................... 14

5.1.1

地すべりの活動箇所の分布 ....................................... 14

5.1.2

地すべりの活動の誘因事象 ....................................... 19

5.1.3

NDVI の空間及び季節における分布の特徴 .......................... 23

5.1.4

地すべりの変動と各樹木の NDVI の変化 ............................ 29

5.2

6.

研究林地すべり ..................................................... 34

5.2.1

地すべり動態 ................................................... 34

5.2.2

地すべりと NDVI の関係 .......................................... 34

まとめと考察 .......................................................... 35
6.1

地すべり活動の把握 ................................................. 35

6.2

地すべり地における NDVI ............................................ 36

参考文献 .................................................................. 37

1．はじめに地すべり多発地域において、現在も活動を続ける地すべりを適切に抽出し、その動態を定期的にモニタリングすることは、施設の建設や維持管理を安全かつ効率よく実施する上で重要な事項である。地すべりは、地下の岩盤に形成されたすべり面に沿って、斜面の土塊が下方に移動する現象である。すべり面は、物理変成や化学変成を受けて微細な粒子成分からなる難透水性の層に形成されやすく、そのような脆弱な地質が分布する一帯では広域に地すべりが多発する。地すべりは、豪雨、急激な融雪、地震により、すべり面上部の土塊中の間隙水圧が高まり、土塊のせん断強度が低下した際に滑動する。その速度は、急勾配の斜面で表層の土層が数秒で滑る表層崩壊よりも、かなり遅い場合が多く、また誘因となるイベントの発生により断続的に移動する傾向がある。しかし上述の誘因によって急
速に地すべり土塊内の水分が増加すると，厚い土塊が一度に長距離滑ることがある。そのような際には、人命や財産の損失、及びインフラの破壊に加えて，地域社会や経済に長期に及ぶ多大な損害を与える場合がある．従って地すべり災害が発生しやすい場所を予め的確に把握し，効率よく被害が最小限に抑える対策を施すことは，地域社会の安全を維持していく上で不可欠な事項である。
地すべりは土塊内の間隙水圧が高まる度にすべり面上をすべる現象であり、その移動期
間は長期間に及ぶ。その為、その移動履歴は頭部の滑落崖や，中腹部の緩斜面、末端部の急斜面や土塊の境界部に発達する沢などの、所謂「地すべり地形」として現れる。従って、まずは地形図等を用いて地形判読を実施し、地すべりの存在を認識することが、災害を軽減する上での不可欠な手段となる。しかし地形判読からは地すべりの存在は分かっても、それらの地すべりが現時点で活発に移動しているとはいえない。地質が脆弱な地すべり多発地域において，道路やダムなどの建設や維持管理を効率よく進める為には，「地すべり地形」が過去に終了した，もしくは現在進行中の滑動の結果生じているのか識別し、地すべりの活動状態を把握することが、地すべり対策の優先順位や手法を決定する上で重要である。これまで現在の活動状況を示す指標としては、地すべり土塊の移動に伴って形成される亀裂や段差、小丘等の数メートルスケールの地形（微地形）の分布を現在の活動状況を示す指標とみなした研究例もある。しかし、一旦形成された微地形は、風化し、消滅するまでに時間を要する為、それらの微地形が示唆する活動の期間には、必ずしも現在が含まれているとは限らない。移動を続ける地すべりの表面では、亀裂などの微細な地形が出現・消失しつつ変化していく為、微地形を現在の地
すべりの活動度の指標とするには、植生下にあっても高精細な地形データが得られる航空レーザー(LIDAR)計測を定期的に実施して、解析をする必要がある。しかし広域に渡る地すべり多発地域を対象として、LIDAR 計測を何度も実施す
ることは経済的な面での負担が大きく、広域を対象にして安価に効率よく活動中の地すべりを抽出できる手段の登場が望まれる。そこで本研
究では，地すべり土塊上に生育する植生を、地す
図 2.1 対象地すべりの位置

1

べり活動状況を推し量るための指標として活用することを提案する。
緩慢に滑動し続ける地すべり上には樹木も生育し，土塊の移動に合わせて根曲
がりや幹折れ等の生育異常が見られることが多い．地すべり地に生育する樹木の樹幹から円盤を切り出し，「あて」が形成された年代を年輪から判読し，滑動履歴を調べる研究はこれまでも数十年に渡り実施されてきた．ただしこの手法は、個別の地すべりを対象にしているため、広域にて現在活動中の地すべりを抽出する
目的には適していない。一方で，衛星デ図 2.1.1 オタリマップ地すべりータを用いて，滑動中の地すべり地では
植生の活性度が低くなることを示唆した研究もある。例えば，後藤（1985）は，
Landsat 衛星データを用いて植生の活性度 で あ る NDVI （ Normalized Difference
Vegetation Index：正規化植生指標))を求め、活動中の地すべり上の植生活性は周囲と比べて低いこと、また、植生活性度が弱
まる秋が、NDVI を用いた地すべりの活動度の推定に、最も適した季節であるこ

図 2.1.2 地すべり末端に露出する凝灰質泥岩と も 示 し た 。 近 年 で は Jacquemart andTiampo（2020）が、InSAR データと Sentinel-2 衛星データを組み合わせることにより、活動的な地すべり地における NDVI は常に周囲と比較して低いこと、また地すべり活動がよ
り活発になると、地すべり上の植生が示す NDVI は低下することを示した。しかし、衛星
データの解像度や地上データの欠落から，NDVI を用いて活動的な地すべりを抽出できるまでには現在至っていない。衛星データの入手と解析が一般にも容易になってきている近年，ここで NDVI と地すべり活動との関係を明確に示すことにより、広域を対象にして、安価に効率よく滑動中の地すべりを抽出することに貢献できる可能性がある．そこで本研
究ではこの関係を明確に示すことを目的に、現在活動中の地すべり地にて UAV 搭載型
LiDAR およびマルチスペクトルセンサによる定期計測を実施し、NDVI と地すべり地形変化との関係を求める。また、植生により NDVI が異なること、また地形変化は地すべり移動以外の原因でも発生することから、現地にて同時に植生調査と地形調査も行った。地形変化については、土塊内の間隙水圧が大きく関与する地すべりによることを示す為に、
UAV 搭載型 VLF 計測も同時に実施した。

2. 対象地すべりの概要
本研究は、北海道沙流郡平取町内オタリマップ地すべりと、北海道天塩郡幌内町問寒別
2地区内北大研究林地すべりにて実施した(図 2.1)。以下に各地すべりについて概要を示す。
2.1

オタリマップ地すべり

オタリマップ地すべり（面積 3 ha）は、沙流川水系オタリマップ川左岸に位置する（図2.1.1）。周辺には神居古潭変成帯の新第三期層の滝ノ上層の泥岩堆積物と，中生代後期白亜紀に形成された超苦鉄質岩類が分布し、地すべり地形が多く存在する（防災科学技術研究
所：地すべり地形分布図）。オタリマップ川河岸には暗灰色を呈した凝灰質の泥岩が露出する（図 2.1.2）。周辺のカラマツとトドマツの植林は 1951 年頃から始まり、造林地外は広葉樹類で構成される（日高北部森林管理署国有林野データ, 2018）。2016 年と 2017 年に北海道森林管理局日高北部森林管理署が実施したオタリマップ地すべりの調査報告書によると，地すべりは風化岩盤すべりに分類され、土塊の移動が年 1 m 以上観測されている。地すべ
り地内の 4 点のボーリング孔の調査結果からは、すべり面深度は滑落崖直下では約 4.7m、末端付近では約 8m 位置にあると推定される。地すべり土塊の頂部には明瞭な滑落崖があ図 2.1.3 オタリマップ地すべりの詳細

3

図 2.2.1 研究林地すべり

図 2.2.2 地すべり末端に露出する蛇紋岩り、その上部にも滑落崖と並行に亀裂の形成が見られる（図 2.1.3a）。滑落崖の直下や土塊
上には湿地や湧水が見られ（図 2.1.3b）、地下水位の高さを示唆している。また、段差や亀裂も広く分布するとともに、植生の傾倒や幹曲がりも多く見られる（図 2.1.3c）。地すべり
末端部はオタリマップ川に押し出し、その水みちを湾曲させている（図 2.1.3d）。
2.2

研究林地すべり

研究林地すべり（面積 9.5 ha）(図 2.2.1)は、北海道幌延町問寒別地区の北海道大学天塩研究林内のパンケルペシペ川の右岸に位置し、その末端は川に張り出ている(図 2.2.1a)。頭部の滑落崖や，滑落崖直下の緩斜面、緩斜面上の線状凹地の発達など、典型的な地すべり地形を示す。また、その周囲にも多くの地すべり地形が存在する。周辺の地質は、神居古
潭変成帯の中生代白亜紀の苦鉄質超塩基性岩である蛇紋岩が分布し、河岸には青白色を呈した蛇紋岩が露出する(図 2.2.2)。この地域は明治末期から数回にわたり山火事の被害を受4けており 、ササ 地が 広 がってい る (図2.2.3)。またその他にも、70 から 100 年
生のアカエゾマツの植林木(図 2.2.3)や、
ダケカンバを優占種とした広葉樹も見られる。研究林地すべりは小ブロックに分化しており、内部には副次滑落崖(図2.2.1b)や、ブロック境界を形成する小沢(図 2.2.1c)が複数見られる。また葦が生育する湿地もあり(図 2.2.1d)、現在も地下水面が高いことを示す。しかし 1947
図 2.2.3 地すべり内のササ地とアカエゾマツ林

年より以降に取得された空中写真からは、地すべり地形に顕著な変化は見られない(図 2.2.4)。

3. 計測
オタリマップ地すべりでの UAV 計測
は、2021 年 7 月（マルチスペクトル）
2021 年 10 月 （ マ ル チ ス ペ ク ト ル ,LiDAR）, 2022 年 6 月（マルチスペクト
ル, LiDAR, VLF）, 2022 年 10 月（マルチスペクトル, LiDAR, VLF）, 2023 年 6 月（マルチスペクトル, LiDAR, VLF）に実施した。一方で、研究林地すべりでの計
測は、2021 年 7 月（マルチスペクトル）と 2021 年 10 月（マルチスペクトル,
LiDAR）に実施した。
3.1

LiDAR 計測

2021 年 10 月及び 2022 年 6 月に実施したレーザー測量では、DJI 社製 Matrice
600 Pro ド ロ ー ン 及 び 、 GreenValley
International 社製レーザー測量機 LiAir50を使用した。M600 Pro は産業用途に設計された高度な多回転式ドローンで、6 枚のブレードや、積載量 6kg までのペイロード、6 つのバッテリーによる先進的な飛行制御システムとプロペラ設計によって、優れた飛行性能と高度な安全機能を提供している。LiAir50 は M600 Pro 用
図 2.2.4 研究林地すべりの変遷に搭載可能なシステムになっており、
5

Velodyne VLP-16 レーザースキャナー、
KVH 1750 IMU（慣性計測装置）、および 2 つのアンテナを備えた NovAtel
GNSS（全球測位システム）レシーバーが搭載されている。また別枠の基地局を用いて RTK 通信を行い観測精度を高めている。本システムを用いて、対象地域を対地高度 60m 及び、飛行速度
を 4m/s (14.4km/h)に設定し、シングル
グリッドフライトを用いて（図 3.1.1）

図 3.1.1 ドローン飛行経路

点群データを取得した。またサイドラ

ップは 50%となるように設定した。取得された点群データに対して Progressive TIN
densification filtering algorithm を用い、すべての点群情報から純粋な地形だけの点群を分類し（図 3.1.2）、後に数値地形データ(Digital Terrain Model: DTM)を作成した。DTM の作成に
は、地形の点群データを用いて Natural Neighbor 法による内挿処理を施し、解像度 0.3m のラスターデータを作成した。なお解像度はより高くできるが、データ処理や容量の観点か
ら 0.3m に留めている。過去の地形データ及び、未来の地形データを用いて、本対象エリアの地形変化を解析する。
2022 年 10 月・2023 年 6 月の計測にて使用した UAS-Lidar は，無人航空機の機体が DJI
社 Matrice 300 RTK，搭載された Lidar センサが DJI Zenmuze L1 であり，これにより得られる地表・地物を含む 3 次元点群の元データを，フィルタリング処理を行うことで地表と地物に分離した。フィルタリング処理には CloudCompare ソフトウェアに備わる Cloth
Simulation Filtering アルゴリズム（Zhang et al., 2016）のプラグイン機能，およびオンライ
ン点群処理ツールである ScanX を，2022 年 10 月・2023 年 6 月のそれぞれのデータに適用した。

図 3.1.2 オリジナルの LiDAR 点群データ及び、地形情報抽出後の点群

6

3.2

マルチスペクトル計測

3.2.1. マルチスペクトル画像解析
地すべり移動体および周辺域における植生の活動度を把握するためには，近赤外域のバンドを含むマルチスペクトルカメラによるオルソ補正画像の生成が有効である（Huang etal. 2021; Reilly et al. 2021）。このような背景から，本研究では地すべり移動体とその周辺地
域における植生活動の詳細なモニタリングを行うため，低空からの撮影が可能なマルチス
ペクトルセンサを搭載した UAS（Unmanned Aircraft System）として，DJI P4 Multispectral（以下，P4M）を使用する。
P4M は，通常のカメラが可視光域の波長を 1 つの RGB センサで捉えるのに対し，RGB
センサに加えて赤色，青色，緑色，レッドエッジ，近赤外の波長帯をそれぞれ個別に捉えることができるマルチスペクトルセンサを搭載している。これにより，より豊富なスペク
トル情報を収集し，植生活動の微細な変化を捉えることができる。さらに，P4M には太陽による陰陽の補正を行う照度計が機体上部に搭載されており，太陽光の強弱による反射率の変化を補正するための情報を付与する機能をもつ。これにより，異なる光環境下での正確なデータ収集が可能である。また，撮影前後の反射率較正を行うために，一定の反射率を備えた白板を使用することもできる。これにより，P4M による反射率データをより正確に校正し，信頼性の高い情報を得ることができる。
さらに，位置情報の精度向上のため，飛行中にリアルタイムで取得される GNSS（Global
Navigation Satellite System）位置補正情報を用いて，RTK（Real-time Kinematic）補正を行う
ことで，cm 単位の精度をもった値を得ることができる。機体の中での GNSS アンテナ位置とカメラセンサの位置との相対関係を加味し，撮影される画像にもそれぞれ RTK 補正された位置情報が付与される。高精度な RTK 補正を行うために，飛行地点付近に一時的な基準局として GNSS 受信アンテナ（DJI D-RTK2）を設置し，機体のコントローラを介した無線通信で位置補正情報を伝送する。また，一時的な基準局の地理座標については，D-RTK2
により得られたアンテナ位置における GNSS 受信情報を，近隣の国土地理院電子基準点による受信情報を用いて後処理補正することで，cm〜mm 単位での精度をもった位置情報として得る。なお，GNSS により得られる高さ情報は基本的に楕円体高であるが，調査対象地におけるジオイド高はおよそ 29.15 m であり，これを楕円体高に加算することで調査対象地付近における標高が得られる。
こ こ で ， 植 生 の 活 性 度 を 示 す 代 表 的 な 指 標 と し て ， 正 規 化 植 生 活 性 指 数 （ NDVI:Normalized Diffierential Vegetation Index）を用いる（Eq. 3-1）。
NDVI = (NIR - R) / (NIR + R) (Eq. 3-1)
ここで NIR：近赤外バンド（Near-Infrared）における反射率，R: 赤色バンド（Red）における反射率をそれぞれ示す。NDVI は，式 3-1 で示されるように近赤外および赤色の波長帯の差を正規化した指標であり，生態学や農林学などの研究分野で広く用いられている。一般に，地球観測衛星や UAS などを用いて収集されたリモートセンシングデータから計算され，植生の活性度を評価するために利用される。NDVI の値は通常，-1 から+1 の範囲となる。この範囲内での値の解釈が重要で，+1 に近づくほど植生の活性度が高く，健康的な7

成長や生育状態を表す。一方，-1 に近づくほど，植生の活性が低く，病気や枯死などの問題が発生していることも示唆する。これは植生活性が高いと赤色光の反射率が低く，近赤外光の反射率が高くなるという植物の分光反射特性を利用している。すなわち，健康な植物は，赤色光を吸収し，近赤外光を反射するため，NIR と R の反射率の差が大きくなり，NDVI の値が高くなる。一方，病気や水不足などのストレスを受けた植物は，この反射パター
ンが変化し，NDVI の値が低下する。
この特性を活用し，NIR と R との差分を一定の基準で標準化して評価することで，広範囲の植生の健康状態を効率的に監視することができる。
Fig. 3.2.2.1. Distribution of NDVI values in the
Otarimappu landslide. (a) July 17, 2021. (b) October
5, 2021. Green indicates higher NDVI values (up to+1), while red indicates lower values (down to -1).

NDVI の応用は多岐にわたる。たとえば，農業分野では，畑や作物の成長
モニタリングに NDVI が使用され，水や肥料の効率的な使用を支持する。環境モニタリングでは，森林や生態系の変化を追跡し，環境保護活動に役立てられる。都市計画では，都市緑化の進捗状況を評価するのにも NDVI が適用される。さらに，気候変動の研究や緑地管理，自然災害の予測にも NDVI は欠かせない指標となっている。このよ
うに，NDVI は地球観測と環境モニタ
リングの分野で幅広く活用され，植生の健康や生育状態の評価に貢献している。UAS によりこの指標が高精細かつ広範囲に得られるようになったことで，地すべりの活動度に対する植生の活性度の面的な評価が可能となる。
3.2.2. オタリマップ地すべり

Fig. 3.2.2.2. Distribution of 25-m cells with lowNDVI values in the Otarimappu landslide. (a)July 17, 2021 (NDVI values <= 0.826). (b)October 5, 2021 (NDVI values <=0.404).

8

地すべりの活動と季節による NDVI
値の変動との関係をみるため，P4M による観測を 2021 年の夏季と秋季の 2

Fig. 3.2.3.1. Distribution of NDVI values in
the Teshio landslide. (a) June 29, 2021. (b)October 15, 2021. Green indicates higher
NDVI values (up to +1), while red indicateslower values (down to -1).

Fig. 3.2.3.2. Distribution of 25-m cells with
low-NDVI values in the Teshio landslide.
(a) June 29, 2021 (NDVI values <= 0.695). (b)October 15, 2021 (NDVI values <=0.429).

時期に実施した。近隣の国土地理院電子基準点には平取（950133 平取，世界測地系（測地成果 2011）における北緯 42°43′38″.0658，東経 142°17′44″.8919，標高 92.444 m，楕円体高121.16 m）を用いた。夏季は 2021 年 7 月 17 日，2022 年 6 月 25 日、2023 年 6 月 6 日、秋
季は 2021 年 10 月 5 日、2022 年 10 月 13 日の晴天日に計測した（図 3.2.2.1）。解像度はそれぞれおよそ 10 cm である。取得された RGB, RE, NIR を含む反射率のデータから NDVI
値を算出する。また，25 ｍ 四方の方形内における NDVI の低い値の密度を求めた。解析
範囲内の全 NDVI 値の 0.5％に相当する値は，2021 年 7 月 17 日は 0.826，2021 年 10 月 5 日
は 0.404 であった（図 3.2.2.2）。なおカラマツ分布範囲に限定した場合には、NDVI の値の
下位 0.5％に当たる値は、2021 年 7 月計測では 0.90，10 月は 0.55 であった。ここでカラマツ林の抽出は，国有林林班図（林野庁日高北部森林管理署）と UAS による可視光（RGB）のオルソ補正画像を用いて行った。植生活性が低い箇所は、植生被覆が乏しいブロック境
界を示す小沢や、地すべり滑落崖や末端の崩壊下の崖錐斜面に出現する。
3.2.3. 研究林地すべり
地すべりの活動と季節による NDVI 値の変動との関係をみるため，P4M による観測を
2021 年の夏季と秋季の 2 時期に実施した。夏季は 2021 年 6 月 29 日，秋季は 2021 年 10 月
15 日の晴天日に計測した（図. 3.2.3.1）。解像度はそれぞれおよそ 10 cm である。取得された RGB, RE, NIR を含む反射率のデータから NDVI 値を算出し、25 ｍ 四方の方形内にお9

ける NDVI の低い値の密度を求めた。解析範囲内の全 NDVI 値の 0.5％に相当する値は，夏季は 0.695，秋季は 0.429 であった（図 3.2.3.2）。またササとアカエゾマツの混在域ではそれぞれ 0.698 と 0.398 であった。低 NDVI 値は、地すべり内の滑落崖や，線状凹地とその中に形成された湧水の周辺等と空間的に対応することが示された。
3.3 VLF 計測
VLF 計測は、2022 年 6 月 9 日、2022 年 10 月 13 日、2023 年 6 月 6 日に、オタリマップ
地すべりにて実施した。
3.3.1

地中センサー

地すべり発生地域における地中情報の観測のために、国立台湾大学の朱子豪教授が開発している超低周波を用いた地中観測センサーを用いて対象地域のモニタリングを行った。
本センサーは、地中内部から反射もしくは透過される 0.3Hz～3kHz の電磁波を観測し、データを比抵抗値び GeoStress に変換するための新規センサー及び解析システムである。仕様上の性能としては深度方向 800m までを 1m 解像度で解析、または深度 6 ㎞までは 10m
解像度で解析などを行うことが可能である（unpublished）。地下構造物を非破壊的な方法で解析可能な本手法は、地滑り発生地域の内部構造を把握するのに最適であるため、実験的ではあるが、導入をした。
3.3.2

機材

DJI 社製 Matrice 210 ドローンと、その脚部にはめ込んで固定されるようの土台を 3D プ
リンターを用いて作成した。土台に基盤を搭載し、落ちないように固定をした。M210 のソ
ケットに取り付け可能なジンバルを用意し、センサー部分をジンバルに搭載した。飛行中は直下を常に向くように設定した。本システムはスマートフォン（以下スマホ）のアプリを用いたログ記録を行い、さらに緯度経度情報などもスマホの情報を利用されるため、より座標精度を高めるために別途、Emlid 社製の GNSS 機である Reach をドローン上部に設置した（図 3.3.2.1）。
3.3.3

データ観測・処理方法

ドローンの飛行ルートは、M600 の LiDAR 取得のための飛行と同じようにシングルグリ
ッドフライトを行った。一方で飛行時間の制限から飛行速度を 5 m/s (18 ㎞/h)に設定し、飛図 3.3.2.1 地中構造解析システムを搭載した状態

10図 3.3.5.1 地中観測センサーによる超低周波の観測結果。データは Geostress 値に変換されている。

行経路のライン間隔は 25 m に設定した。本センサーは高さに依存しないため、飛行高度は割愛する。地中観測センサーのログ頻度は 1 秒に設定し、飛行経路上を直下方向に観測している。また、同時に Reach GNSS も 1 Hz（1 秒）の観測間隔にし、後に地中観測センサ
ーのタイムスタンプログを用いて数値を入れ替えることで座標精度を高める。Reach GNSSの観測衛星には GPS、GLONASS、GALILEO、BEIDOU を用いている。観測後、取得された Reach の GNSS データは、近くの電子基準点から取得できる衛星信号データを用いてPost Processing Kinematics (PPK)処理を施し、より精度の高い座標値を計算している。
3.3.4

データの可視化

解析処理を施した地中観測センサーの生データを可視化するのに Golden Software 社製の Voxler3D を用いる。本ソフトウェアは様々な 3 次元データの可視化を可能とするものである。今回は取得されたデータの GeoStress 値を可視化し、地中の状態や構造が 2022 年 6月以降のモニタリング結果と比較してどうなっているかを確認していく。
3.3.5

VLF 計測結果

図 3.3.5.1 に地中観測センサーから取得したデータを可視化したものを示す。2022 年 6 月
9 日、2022 年 10 月 13 日、そして 2023 年 6 月 6 日の 3 機関での Geostress の変化を見たところ。特に注目したいのは 2022 年 10 月に極端に Geostress 値があがっていることである。
2022 年 6 月及び 2023 年 6 月の数値ではほとんどの数値が 0.2～1 の間に分布しているが、2022 円 10 月では 10~60 に分布しており、数十倍に数値があがっている。平面的に全体を見ても Geostress の値が相対的に高いところが増えており、このエリアにて何かしらの影響が出ているように見られる。直前に多量の雨量が観測されているため、それらが影響して数値として現れている可能性はある。
図 3.3.5.2 は 2023 年 6 月時点で観測している地滑りの著しかった箇所のプロファイルを取って、東側から断面を観察したものである。日高北部森林管理署の報告書で示されているボーリング調査によるすべり面の位置と、高 geostress 値と低 geostress 値との範囲の境界がおおよそ符号しているように見える。
11

図 3.3.5.2 該当プロファイルを東側から見た図。6 月の 2 つの数値スケールは一緒だが、10 月のスケール
は絶対値が違うことに注意（相対的なスケールを似せている）

4.
地すべり動態の把握
4.1 オタリマップ地すべり
4.1.1. オタリマップにおける点群解析
無人航空機（UAS: Unmanned Aircraft System）に搭載されたレーザ測量センサ（Lidar: Lightdetection and ranging）による反射データから得られる地表・地物の 3 次元点群データを用
いて，2022 年 10 月 13 日および 2023 年 6 月 6 日に取得した UAS-Lidar の点群データの位
置関係を比較し，オタリマップ地すべりの約半年間における地表変動を解析した。地表を示す点群同士の位置関係の比較，すなわち地表変化の抽出には，M3C2 距離アルゴリズム（Lague et al. 2013）を用いた。これは，2 時期の点群のうち，片方の点群を基準として，その表面の法線方向におけるもう片方の点群までの距離を一定の範囲を基準として示すもの
であり，地すべりのような移動する物体や，侵食・堆積などの 3 次元的な変化抽出を行うために用いられる。
4.1.2

勾配の変化

2009 年と 2015 年に北海道開発局によって取得された LiDAR データから 1 m セル毎に勾配を求め、各計測の間に発生した地形変化を、勾配の変化量として求めた。あるセルにおける勾配は、周囲の 8 セルも含んで計算されるために、微地形の形状を表す指標として適12

していると考えられる。本研究では、微地形の形状の変化を確実に把握する目的で、一般に用いられる標高ではなく勾配の変化量を採用した。この変化量を、
LiDAR データから作成した地形図と比較したところ、±6 度以上の勾配の変化が表れた場所にて、地すべりによる地形変化が顕著に見られることが分かった。微地形の分布を把握するにあたり、一辺が25 m サイズのウィンドウが適していたため、このウィンドウにて±6 度以上の勾図 4.2.1 年輪サンプリング

配の変化が表れたセルの密度を求め、地
すべり動態の指標とした。

4.1.3

PIV 解析

2010 年と 2015 年に北海道開発局によって有人飛行にて取得された LiDAR データ、及び2021 年 10 月と 2023 年 6 月に UAV によって得られた LiDAR データからそれぞれ 1 m セル
サイズの DEM を作成し、PIV 解析を実施した。まず DEM の値を、近傍 5 m を対象に平滑化し、この平滑化した値にラプラシアンフィルターをかけて、地表のエッジを強調したラプラシアン図を作成した。ここでラプラシアンフィルターとは、対象セルの縦方向と横方
向の 2 次微分を足したフィルターのことであり、得られたラプラシアン図では、地すべりに出現する崖や線状凹地が強調される。この図について、PIV.lab（Thielicke, 2021）を用いて Matlab 上で PIV 解析を実施した。この PIV 解析では、5 m 四方を単位としたある時期の地形について、次の時期において認識された類似の地形までの距離と方向を算出し、変位ベクトルとして表すことができる。
4.2

研究林地すべり

研究林地すべりは、典型的な地すべり地形を呈してはいるものの、土塊の移動は顕著には見られない。ササが密生していることから、レーザー光も地上に届きにくく、LIDAR データを用いた数値解析を実施し、地形変化を捉えることは困難であった。そこで 1947 年から始まる一連の空中写真及び LIDAR データから作成した地形図を用いて亀裂や線状凹地の伸長を判読し、長期的な地形変化を捉えることを試みた。
加えて、地すべりの活動に伴い，地すべり地上に生育する樹木の傾倒や，樹木の根のせん断による成長の劣化が発生し、その履歴が年輪に反映すると考えられることから，地すべり内のアカエゾマツの年輪判読も実施した。ここで年輪のサンプリングは、成長錘を用
いて，斜面下方側から材を採取した（図 4.2.1）。そして顕微鏡とノギスを用いて、年輪幅の計測を実施した。

13

5. 結果5.1 オタリマップ地すべり
5.1.1 地すべりの活動箇所の分布
各 LIDAR から求められた勾配図を、図 5.1.1.1 に示す。等高線は、1m 間隔で示している。赤線で囲った範囲は、2023 年 6 月 6 日時点の勾配図と等高線より判読されたオタリマ
ップ地すべりの範囲である。2009 年 11 月の時点から比較すると、末端である a 行の部分
にて地すべりがオタリマップ川に向かって徐々に押し出ていっていることが分かる。なお、2021 年 10 月から 2022 年 10 月までは、例えば b-3 に見られるように局所的な地形の変化のみが確認されており、地すべり土塊が大きく移動したのは 2022 年 10 月から 2023 年 6 月であった。

図 5.1.1.1 オタリマップ地すべり周辺の勾配図14図 5.1.1.2 上段：25 m 四方において、期間中の勾配の変化が絶対値 6 度以上のセルの出現数、下段：
PIV 解析結果

図 5.1.1.2 には、上段に 25 m2 の範囲にて、絶対値としての勾配の変化が 6 度以上であった 1 m セルの数を、下段には PIV 解析の結果を示す。2021 年と 2023 年については、UAV計測を始めた時点及び最近の計測時の LIDAR データを用いた解析結果を示している。どちらの図においても、図 5.1.1.1 のみからは判読しにくかった地形の変化が明確に示されている。これらの図からは、2009 年から 2023 年にかけて、地すべり土塊の移動とともに、地形が激しく変化した場所も斜面下方側にずれていっていることが分かる。勾配を用いた解析結果からは、そのような場所が認識できるのみである一方で、PIV 解析においては、
地すべり土塊上の各地形の変動距離と方向が示される。図 5.1.1.3 には変動方向の変遷を示す。ここで矢印の大きさは、変動距離の大きさを指す。ただし、各期間において矢印の大きさを示す距離の範囲が異なることに留意する。黒枠と白線、および行と列の記号と番号15図 5.1.1.3 PIV 解析による地すべり土塊の変動方向の変遷

は、図 5.1.1.1 及び 5.1.1.2 と同じである。各図の背景には、それぞれ 2015 年、2021 年、2023 年の勾配図を用いている。2009 年から 2015 年までは、地すべり内の滑落崖付近（図5.1.1.2 中の c-1 と c-2 の範囲における黒線、図 5.1.1.3 中の青矢印で示した付近）にて最も変動が大きい。図 5.1.1.1 における c-2 をみても、滑落崖下の地形が下方に移動しているこが示唆される。おそらくこの土塊の移動に伴い、滑落崖上方の斜面も下方に移動したと考えられる。この期間においては、地すべり頭部から続く右側の崖からも（図 5.1.1.3 中の黒枠矢印で示した付近）土塊の北西方向への移動が確認できる。滑落崖付近の変動は、2015
年から 2021 年までの期間も認識されるものの、2021 年から 2023 年までの間は顕著ではなかった。一方で、c-2 付近の頭部付近の右側の崖付近においては、地すべり土塊の北東方向への変位が 2009 年以降の期間を通じて示唆される。また同期間を通じては、b-2 と b-3 にかかる範囲にて、北西方向への変位が顕著であった。この変位に伴い、地すべり土塊上に黄矢印で示される沢も西方向へと位置を変えている。末端部（a 行）においても、地形の変位方向は、2009 年から期間を通じて大きな変化は発生していない。この範囲では、赤線で示した軸を境に、土塊の移動方向が分裂しており、図 5.1.1.2 を見ると、二つのブロックとして分裂していることが示唆されている（白矢印）。図 5.1.1.1 では、この軸を境に末端部の崩壊の向きが異なり、かつ亀裂もその上方の土塊上に形成されている（赤矢印）。また緑矢印で示されている箇所付近に今後副次滑落崖が発達していくと考えられる。なお図5.1.1.4 には、22 年 10 月に UAV で取得されたオルソ写真から無作為に抽出した 56 本の樹
木の位置を示しており、星印が 23 年 6 月のオルソ写真でも位置が変わらなかった樹木、丸印が移動した樹木になる。図中の紫色の直線は、樹木の移動を示している。最も移動距離が長かった b-2 の左上側にある樹木以外は、概ね移動方向は PIV 解析で得た移動方向と一致していた。この移動距離が長かった樹木については、PIV 解析による移動方向と直行
した方向に移動している。土塊の移動方向の変換点に位置していることと何等かの関係が16

図 5.1.1.4 2021 年 10 月から 2023 年 6 月までの PIV 解析による地すべり土塊の変動方向と、2022 年 10 月から 2023 年 6 月までの立木の移動

図 5.1.1.5 b-3 付近の地すべり地表面の変遷あるとも考えられるが、現時点ではその理由は不明である。なお b-3 にて地すべり内に残る塩ビパイプ付近で撮影されている写真（図 5.1.1.5）では、時間の経過とともに周辺の樹木が傾倒し、パイプも浮き上がっていることが分かる。
また 2022 年 10 月と 2023 年 6 月の計測で取得された点群については、M3C2 距離アルゴ
リズム（Lague et al. 2013）を用いることで得られた各点における変化方向と距離の情報を，
17

Fig. 5.1.1.6. Three-dimensional movement of the surface of the Otarimappu landslide from October 13,
2022 to June 06, 2023, using the ground surface data extracted from the UAS-Lidar point cloud.

QGIS を用いて矢印として地図上に示している（図 5.1.1.6）。いくつかの点はデータや M3C2
距離計算の誤差により過大な長さをもつ矢印が散在するが，これらを除くと，地すべり移
動体の変化の状況をかなり詳細に把握することができる。たとえば，地すべり移動体の末端北西面では，その移動方向は概ね北向きを示すのに対し（図 5.1.1.6 a, b），末端北東面ではその移動方向は東向きとなっている（図 5.1.1.6 d,e,f）。一方，地すべり全体の移動方向は北東向きであるが，移動体内部の各所においては，それぞれより局所的な変動を示し，より複雑な変動を示している（図 5.1.1.6 g, h, i, j, k）。滑落崖の付近においても，基本的には
崖面の急傾斜方向に沿った変動を示すが，凹地の形成される場所についてはその変動は一
18

Fig. 5.1.1.6. Three-dimensional movement of the surface of the Otarimappu landslide from October 13,
2022 to June 06, 2023, using the ground surface data extracted from the UAS-Lidar point cloud.

様ではなく，また末端部に比べ変動幅は小さい（図 5.1.1.6 l, m, n）。全体的には移動体が引き伸ばされ，その表面が複雑に起伏を再構成するような変動を示している。
5.1.2

地すべりの活動の誘因事象

オタリマップ地すべりから 2.9ｋｍ北西にある平取町内の気象庁旭気象観測地点にて記
録されている月降水量のグラフを図 5.1.2.1 に示す。日高北部森林管理署による報告書では、2016 年 6 月から 10 月の約 4 か月間（緑枠）では、地すべり地内に 53 本設置された移動杭が平均で約 2.5 m 移動したことが記載されている。この期間中に発生した 2016 年 8 月
19Fig. 5.1.1.6. Three-dimensional movement of the surface of the Otarimappu landslide from October 13,
2022 to June 06, 2023, using the ground surface data extracted from the UAS-Lidar point cloud.

の台風 7 号により（最大日雨量：120 mm）では、地内のボーリング孔に設置したひずみ計が破断している。また同月には台風 9 号も発生しており、地すべり土塊の移動に寄与したと考えられる。また 2018 年 9 月 6 日に発生した北海道胆振東部地震でも、周辺は震度 6 弱の揺れが記録されていることから、この地震が土塊を移動させた可能性も高い。
一方、UAV を用いた観測を開始した 2021 年 7 月以降については、2022 年 8 月に月降水量が 400 m（最大日雨量：81 mm）を超えた。7 月の降水量は平年並みであったものの、6月の降水量についても、2016 年時と同様に多かった。しかし上述したように、地すべり土塊の変動が最も顕著に確認されたのは、夏の多雨期ではなく、2022 年 10 月から 2023 年 620

Fig. 5.1.1.6. Three-dimensional movement of the surface of the Otarimappu landslide from October 13,
2022 to June 06, 2023, using the ground surface data extracted from the UAS-Lidar point cloud.

月にかけての期間であった。その誘因として、2023 年 3 月に急速に進行した融雪が考えられた（図 5.1.2.2）。図 5.1.2.2 では、オタリマップ地すべりから 19 km 南東にある新冠町新和に設置された気象庁の観測地点における 2023 年 2 月 1 日から 3 月 31 日までの一日あた
りの最低、平均、最高気温と、旭気象観測地点で記録されている積雪深を示している。気温については、5 日間の移動平均も点線にて表している。気温の観測地点は地すべりから距離が離れているものの、大体の変化の傾向は同様であったと考えられる。グラフからは、
2 月の下旬から気温が急激に上昇し、積雪深も減少したことが分かる。特に 3 月 8 日以降
は最高気温が 10 度に届く日も増え、最低気温が 0 度を下回らない日も出てくる。2 月 20
21

Fig. 5.1.1.6. Three-dimensional movement of the surface of the Otarimappu landslide from October 13,
2022 to June 06, 2023, using the ground surface data extracted from the UAS-Lidar point cloud.

日は 60 cm あった積雪深が、3 月 15 日には 0 cm となる。この観測地点とオタリマップ地
すべりにおける積雪深は全く同じではないものの、雪の密度を 0.5 g/cm3 とした場合には、
この融雪量は 300 mm の降水量に相当することになる。
なお 2009 年 11 月から 2015 年 6 月にかけての期間については、地すべり土塊の移動の誘因は不明であるものの、2010 年夏の多雨である可能性がある。また 2012 年 4 月には北海道内にて急激な融雪が発生し、各所で斜面崩壊や地すべりが記録されている。上述の、2022 年 10 月から 2023 年 6 月までの期間における地すべり土塊の移動が融雪であると考えられることから、この時期においても同様の現象が発生した可能性が非常に高い。
22

Fig. 5.1.1.6. Three-dimensional movement of the surface of the Otarimappu landslide from October 13,
2022 to June 06, 2023, using the ground surface data extracted from the UAS-Lidar point cloud.

5.1.3

NDVI の空間及び季節における分布の特徴

NDVI の値には樹種による影響があると考えられることから、ここではカラマツ林に限定し、値の分布を求めた。樹種の判定には，国有林林班図（林野庁日高北部森林管理署）
及び 2021 年 10 月 5 日と 30 日に取得した UAV によるオルソ画像（m セル）を用いた。計測した NDVI は、1 m セルごとに表した。
図 5.1.3.1 に、各計測日の NDVI 値の分布を示す。赤線で 2015 年以降に地すべりの滑動が明確に示された範囲を囲っており、図 5.1.1.2 で示した範囲よりも狭いことに留意する。
23

Fig. 5.1.1.6. Three-dimensional movement of the surface of the Otarimappu landslide from October 13,
2022 to June 06, 2023, using the ground surface data extracted from the UAS-Lidar point cloud.

図 5.1.2.1 気象庁旭気象観測地点における月降水量

また図に示した青（地すべり 1）、黒（地すべり 2）、ピンク（地すべり 3）の線で囲まれた範囲はそれぞれ地すべり地形を有している。地すべり 1 と 2 については顕著な土塊の移動が近年は確認されていない一方で、地すべり 3 については、オタリマップ地すべりより移動距離は少ないものの、頭部における亀裂の拡大が、2022 年 10 月と 2023 年６月に取得されたオルソ写真より確認されている。図 5.1.3.1 においては、6 月または 7 月と比較して気温が低い 10 月の計測は、値も全体的に低くなることが分かる。特に 2022 年 10 月の計測
では、2021 年 10 月よりも植生の活性がかなり低い。一方で、2023 年 6 月の計測でも、前年および前々年と比較して NDVI 値が全体的に低い。新和気象観測地点でのマルチスペク
24トル計測日における 13 時の気温は、2021 年 7 月 17 日は 30.2 度（12 時から
13 時までの日照時間 0.6 時間）、2021
年 10 月 5 日は 17.3 度（12 時から 13 時までの日照時間 1 時間）、2022 年 6 月
25 日は 23.9 度（12 時から 13 時までの日照時間 0.7 時間）、2022 年 10 月 13 日
は 17.6 度（12 時から 13 時までの日照時間 1 時間）、2023 年 6 月 6 日は 21.8
度（12 時から 13 時までの日照時間 0.9
図 5.1.2.2 2023 年 2 月及び 3 月の気象庁旭気象観測地点における積雪深と、新和観測地点における日平均気温時間）であった。観測地点が離れているため、これらの値を直接対象

図 5.1.3.1 各計測日の NDVI 値の分布25

図 5.1.3.2 新和気象観測地点における各年の日平均気温. 上：5 月から 7 月、下：9 月から 10 月

地に適用することはできないが、気温が同様に変化する傾向にあったとするならば、6 月
または 7 月の計測日については、21 年が最も気温が高く、23 年が最も低かったことになる。また 10 月の計測日については、気温に大差はなかった。一方で、それぞれの計測日周
辺の気温の変化を見ると（図 5.1.3.2）、2023 年 6 月計測日の以前は、その他の年の 6 月もしくは 7 月の計測日比較して気温がかなり低く推移していたことが分かる。従って 2023 年
6 月計測日時点においては、冬季からの植生の回復が他の計測日よりも進んでいなかった可能性がある。また 10 月の計測日については、前述したように 13 時時点の気温に差はな
かったものの、2022 年計測日の前には 2021 年計測日よりも気温が低く推移しており、植生の活性が前年より落ちていたことがうかがえる。なお図 5.1.3.1 では、オタリマップ地すべり内にて極端に植生の活性の低い範囲が各所に見受けられるものの、各地すべり内やその他のエリアにおける NDVI 値の明確な違いは確認できない。
図 5.1.3.3 には、図 5.1.1.2 と同じ範囲におけるオルソ写真及び NDVI 値の分布を示す。
図 5.1.3.1 にて局所的に NDVI 値が低かった範囲は、c-１や c-2 の滑落崖直下や、滑落崖下から発達する沢沿い（c-2, c-3, b3）、a-3 内にある植生に乏しい箇所であることが分かる。秋
季は落葉により地上が露出し、NDVI が低い範囲が増える。また a-3 と a-4 の境界にある、
図 5.1.1.2 で示唆された末端ブロックの境界部分についても、秋季には低い NDVI 値が出現する。
NDVI 値について、図 5.1.3.4 では、連続する計測日における差を図 5.1.3.1 と同じ範囲にて示す。オタリマップ地すべり内の NDVI 値は、図にて NDVI 値が低かった場所にて顕著
26

図 5.1.3.3 オタリマップ地すべり周辺の各計測日の NDVI 値に変動する。夏季から秋季にかけての植生活性の変化の傾向には、上述した計測日による活性度の違いはあるものの、オタリマップ地すべりと、地すべり外や他の地すべりとの間に特徴的な違いは見られない。その一方で、秋季から夏季にかけての変化を比較すると、地すべりの移動が生じた 2022 年から 2023 年にかけては、植生が活性化した箇所が増加したように見受けられる。しかし地すべり 1, 2, 3 およびオタリマップ地すべり内の NDVI 値の分布について、平均値及び標準偏差を示した図 5.1.3.5 からは、その傾向を読み取ることは困難である。グラフからは、オタリマップ地すべりには植生が乏しい箇所が存在するた
めに、値の平均値が低くなるとともに、そのばらつきも大きくなることが分かる。特に 2022
年 10 月の計測では、その特徴が明白に表れる。これらの傾向は、地すべり 2 と 3 についても見られるものの、オタリマップ地すべりにおける NDVI 値の期間を通じての変動は、他の地すべりと比較して最も大きい。なお図 5.1.3.6 では、2022 年 10 月から 2023 年 6 月の期間を対象にした PIV 解析において、各地すべり内での平均移動速度と、2023 年 6 月時点
での平均 NDVI 値との関係を示している。このグラフでは、地すべりの移動距離が大きいほど NDVI の値が小さくなることが示唆されているが、NDVI 値の範囲が狭く、かつ対象となる地すべりの数が限られているために、明確な関係とは言い難い。ここで地すべり 3については解析範囲外の為、グラフには示されていない。
図 5.1.3.7 には、NDVI 値について、夏季及び秋季を対象にした計測年度間の差を示す。
秋季については地すべり間や地すべり内外の違いが明白ではないが、2021 年と 2022 年の夏季の計測については、図 5.1.3.3 において植生が乏しく活性が低かった場所が明確に示される。一方で、2022 年 6 月と 2023 年同月の NDVI の差は、オタリマップ地すべり内で顕著に表れる。特に植生活性が増加している範囲も多く、図 5.1.3.1 にて示された 2023 年 6月における NDVI が 2022 年 6 月よりも低い傾向にあった結果と反する。そこで同じ箇所
27

図 5.1.3.4 連続する計測日間の NDVI 値の差における 2022 年 6 月と 2023 年同月の NDVI の差と、同計測日におけるオルソ画像を比べたところ、地すべりの変動に伴い出現した亀裂や裸地にて NDVI 値は低下したことが分かった（図 5.1.3.8）。また、樹木が移動したことで出現したギャップでは、下層の植生の活性が高くなった一方で、移動した樹木自体については活性が低下した例も多く見られた。ただしこの樹木の活性の低下については、前述したように計測日以前の気温による影響も大きいと考えられる。地すべりの変動による樹木の活性への影響については、次章で検討する。
28

図 5.1.3.5 各地すべり内における計測時の平均 NDVI 値及びその標準偏差

図 5.1.3.6 各地すべり内における 2022 年 10 月から 2023 年までの地表の平均移動速度及び 2023 年 6月計測時の平均 NDVI 値

地すべり活動が与えた NDVI の空間分布への影響を強調するために、図 5.1.3.7 で示した
1 m セルごとの NDVI の差を、フォーカル統計にて 10 m 四方の平均値として求めた図を図5.1.3.9 に示す。この図では、図 5.1.1.2 で示された地すべりの変動が大きかった場所にて、2022 年と 2023 年における夏季の値の分布が、2021 年と 2022 年の夏季の分布の傾向と異なっていたことが、鮮明に示されている。
5.1.4

地すべりの変動と各樹木の NDVI の変化本研究では、UAV を用いることにより、マルチスペクトル計測による詳細な NDVI の分
布と LIDAR による地形の変化、及びオルソ画像を取得できた。これらのデータを用いて、樹木ごとに位置と NDVI の変化を、地形の変動とともに追跡することが可能となった。前述したとおり、2023 年 2 月から 3 月の融雪時に発生したと考えられる地すべり土塊の移動
方向は、オルソ画像から無作為に抽出した樹木の移動方向と概ね一致していた（図 5.1.1.4）。
これらの樹木について、林冠の中心部における NDVI の値を 2022 年および 2023 年の 6 月に取得された NDVI 画像から取得し、その変化を調べた。解析の対象となった樹木及び29

図 5.1.3.7 各年の夏季および秋季の NDVI 値の差

NDVI 値を図 5.1.4.1 に示す。上段の図の背景には、2022 年 6 月及び 2023 年 6 月計測時のオルソ写真をそれぞれ用いている。また樹木ごとに NDVI 値を記している。下段は、背景が 2022 年 6 月と 2023 年 6 月の差であり、正の値が 2023 年にて植生活性が高くなったことを表す。赤と青の矢印は、それぞれ図 5.1.3.8 に示した画像の位置である。紫色の直線が樹木の移動方向であり、移動先の樹木は、2023 年 6 月と 2022 年 6 月における NDVI 値の比として表されている。１を超えると植生活性が前年より高くなったことを示し、暖色で示している。図中の星印が、2022 年 10 月と 2023 年の 6 月にて位置が変わらなかった樹
木、丸印が変わった樹木になる。地すべり１と２、及び滑落崖上部の樹木については位置の変化は見られなかった。選択した樹木のうち最も長く移動したのは b-3 に位置し、その30

図 5.1.3.8 オタリマップ地すべりにおける 2022 年 6 月と 2023 年 6 月の NDVI 値の差と、同箇所のオルソ写真

距離は約 14 m であった。2022 年と 2023 年の計測日周辺の気象条件の影響もあり、各樹木の活性は前年と比較して概ね低いが、オタリマップ地すべり外には活性が高い樹木もある。
また地すべり内においても、a-3 や b-1 に位置する地すべり末端の樹木の活性は高い。オルソ写真を見ると、この 2 本については、位置の移動に伴い近接した樹木がなくなっている
為、日当たりが以前より良好な状況におかれていたことも考えられる。
図 5.1.4.2 には、2022 年 10 月から 2023 年 6 月までに位置の移動があった樹木、および
移動しなかった樹木の NDVI について示している。ここで示した全ての樹木については、
2021 年 10 月から 2022 年 10 月までは、ほぼ同位置にあったことがオルソ写真より確認されている。また移動しなかった樹木は、全てオタリマップ地すべりの活動範囲外にある。
このグラフより、オタリマップ地すべり内にある樹木の方が、植物活性の季節変動が激しいことが示唆され、特に夏季における植生の活性が高い傾向にある。これは裸地が多いた
めに植生の密度が疎であることに加え、地すべり地は地下水面が高いために、生育に十分な水分も供給されているためとも考えられる。地すべり移動に伴う樹木の変形や傾倒で、活性が弱くなることも当初予想されたが、2023 年 6 月の時点でも地すべり外の植生と比較すると、活性は高い傾向にある。ただし 2021 年や 2022 年と比較すると 2023 年では地す
べり内外の樹木の活性は類似している。一方で秋季については、2021 年 10 月の計測時では NDVI の値に差異は見られなかったものの、2022 年については活性が低い樹木が多くなる。この理由としては、植生密度が疎であることや地下水位が高いという地すべり地特有
31

図 5.1.3.9 各年の夏季および秋季の NDVI 値の差についての 10 m 平均値の条件と関係があるとは考えられるが、現時点では不明であり、更なる検討が必要である。
図 5.1.4.3 には、各樹木における 2021 年と 2022 年夏季、及び 2022 年と 2023 年夏季のNDVI の比の関係を示す。地すべり土塊の滑動に伴い位置が移動したケースのいくつかの樹木については、地すべり外にて移動が発生しなかった樹木と同様の関係にある一方で、2023 年の値がこの関係より低い値をとる樹木も多い。すなわち、冬季からの活性の回復が遅れていると解釈できる。図 5.1.4.4 には、2022 年 10 月と 2023 年 6 月に取得されたオルソ画像から判読された移動距離と、2022 年と 2023 年夏季の NDVI の比の関係を示す。移動距離が 3 m を超えると、位置が移動しなかった樹木と比較して前年度の夏季の NDVI に対する今年度の値が低い傾向にある。
図 5.1.4.3 及び 5.1.4.4 において、2022 年と 2023 年夏季の NDVI 比が最も低かった樹木に32

図 5.1.4.1 上：2022 年と 2023 年の夏季における立木の NDVI、下： 各立木の 2022 年と 2023 年の夏季における NDVI の比ついて、2022 年 6 月と 2023 年 6 月に取得されたオルソ画像を図 5.1.4.5 に示す。この樹木の位置を図 5.1.4.1 の下図において紫色の矢印で示す（a-3）。この樹木は、この期間に約 6m 位置が移動した。黄色の矢印で移動方向を示す。図 5.1.4.5 では、樹冠の右斜め上方向に幹が伸びている様子が見られ、この樹木が地すべり移動にともない傾倒したことがうかがえる。またこの樹木と右斜め上の樹木との間のギャップが大きくなっていることから、こ33

図 5.1.4.2 2022 年 10 月から 2023 年 6 月の期間において移動が発生した樹木と移動しなかった樹木における各計測時の NDVI 値のスペースには亀裂が新たに生成された可能性も高い。実際に現場では地すべりの土塊の移動に伴う亀裂や階段状の地形の形成のために、根が切れている樹木も多く確認されている（図 5.1.4.6）。この樹木についても、同様の事象が発生したことにより植生の活性が低下したと考えられる。
5.2

図 5.1.4.3 2022 年 10 月から 2023 年 6 月の期間において移動が発生した樹木と移動しなかった樹木における、各年の夏季の NDVI の比研究林地すべり5.2.1

地すべり動態研究林地すべりは、典型的な地すべり地形を示し、地すべり内には複数のブロックが判読できるものの、ブロック境界を形成する小沢の若干の伸長の他は、土塊としての顕著な動きを確認することは1947 年以降の空中写真からは困難であった。採取した全てのアカエゾマツの年輪
サンプルからも、地すべり土塊が移動した際に見られる年輪幅の急激な変化は見
図 5.1.4.4 2022 年 10 月から 2023 年 6 月の期間における各樹木の移動距離と 2022 年と 2023 年の夏季の NDVI の比られなかった。
5.2.2

地すべりと NDVI の関係

研究林地すべりの内部には複数の小ブロックの他にも微細な亀裂などの線状凹地が存在する。また湿地もあることから、地下水面が高いことも分かる。これらの地形的特徴はいずれも地すべり活動に伴って現れるが、研究林地すべりにて土塊の顕著な移動があったとすれば、それは最も古い空中写真（1947 年取得）以前であったと考えられる。一連の空中写真や LIDAR 画像、年輪判読からは、現在土塊が移動中であるとしても、かなり緩慢な34

図 5.1.4.5

図 5.1.4.4 にて 2022 年と 2023 年の夏季の NDVI の比が最も低かった樹木の移動例

速度であることがうかがえる。この地すべりが植生活性に与える影響は、小ブロックを形成する小沢や、崖下など、地すべり活動の結果副次的に出現する地形の周辺に限られており、これらの地点では、水流や崖錐堆積物の蓄積により、植生の生育が阻害されていたことが現地調査にて確認された。
図 5.1.4.6 地すべり土塊の移動に伴う階段状の地形の形成と樹木の根切れ（2023 年 6 月撮影）

6. まとめと考察
6.1 地すべり活動の把握
本研究では、地すべり活動の把握に際し、航空写真や有人航空機搭載型のレーザー計測に加え、UAV 搭載型のレーザー計測を採用した。レーザー計測では、樹木が生育していても地表の微地形を詳細に捉えることが出来るため、地すべりのように線状凹地や湿地の出没と消滅を繰り返しながら移動していく現象を計測するに適した手段である。しかしササのように地表をほとんど被覆するような植生が繁茂する場所では、レーザー光が地表に到達しないため、地表の地形の計測は困難な状況である。その為、研究林地すべりでは空中写真判や年輪の判読が地すべり動態の把握に用いられた。
一方、オタリマップ地すべりでは、土塊の移動に伴う地形変化を把握するにあたり、
LIDAR 計測が強力な手段となった。測量データから得られた地形のイメージ図（本研究では勾配図）からでは移動箇所の判読による移動箇所の把握は若干困難ではあるものの、移動箇所の把握については、時期の異なる測量データから得られた勾配の差分を用いることで可能であった。なおこの点については DEM の差分を用いる手法も一般的ではあるものの、今回は標高値の誤差を考慮し、形状の変化を把握することに重点をおいたため勾配の35差分を採用している。LIDAR 計測による DEM を用いた PIV 解析からは、少なくとも 2009年以降は地すべり土塊の移動の中心は、斜面下方へ移ってきていることが明らかになった。
すなわち、地すべり土塊が下方へゆっくりと流れていっている状況を示しているともいえる。ただし VLF 計測の結果からは、オタリマップ地すべりでは、頭部滑落崖のその情報（図 5.1.1.2 の d-1）においても地下水面がまだ高いことが確認される。C-2 でも滑落崖の直上に亀裂の発生も見られることから、いずれ 2009 年から 2015 年までに発生した状況と同様の規模の地すべり移動が起こる可能性もある。
UAV によって取得した LIDAR データでは詳細な地形が表現されているために、計測期間において最も土塊が移動した現象を的確に把握することができたと考える。M3C2 距離
アルゴリズム（Lague et al. 2013）を用いた解析では、土塊上の詳細な変動を捉えることが可能であった。また PIV 解析とともに、b-3 にて地すべりの末端部が二つのブロックに分裂し、オタリマップ川に向かって別方向に移動していく状況を明らかにすることができた。
この末端部における二つのブロックの移動方向は 2009 年以降変わっていないことから、
オタリマップ地すべりは今後も引き続き、川に向かって扇状に広がっていくと考えられる。
なお本研究は LIDAR データを地形解析に用いて地すべりの動態を検討したが、航空写
真では植生の位置の変化が分かるため、植生の移動を追跡することで、地すべりの移動方
向と距離の把握が可能になることも分かった。航空写真では、土塊の移動に伴い発生した亀裂などが樹冠間のギャップの形成からも推測できることから、LIDAR データが存在しない地域でも、航空写真は地すべりの動態の把握の為に活用できることが改めて示された。
6.2

地すべり地における NDVI

後藤（1985）は，活動中の地すべり上の植生活性は周囲と比べて低いこと、また、植生
活性度が弱まる秋が、NDVI を用いた地すべりの活動度の推定に最も適した季節であることを、Landsat 衛星データを解析して示した。また Jacquemart and Tiampo（2020）も、Sentinel2 衛星データから、活動的な地すべり地における NDVI 値は常に周囲と比較して低いこと、
また地すべり活動がより活発になると、NDVI 値は更に低下することを示している。本研究でも、当初は、地すべり土塊の移動に伴い、地表の樹木が変形したり傾倒したりすることで、地すべり地では植生活性が低くなると想定していた。しかし結果からは、地すべり
内では地すべり外と比較して NDVI が低くなることが示されたものの、その主な理由は、土塊の移動に伴って形成された沢や崖、また崖下の崖錐堆積地などの植生が育ちにくい場所が存在することにあった。特に活動的な地すべりでは、そのような場所も頻繁に生成されるために、地すべり全体を対象にした場合の NDVI の平均値は低下することになる。衛
星データを使用している場合、その解像度から植生のみが NDVI の値に直接反映しているとは言い難く、低い値はそれらの地すべり特有の地形による可能性が高い。本研究で UAVを用いて、各樹木についての NDVI 値が判明できる解像度のマルチスペクトル計測を行ったところ、活動的な地すべり地における樹木の植生活性は、地すべり外と比較し夏季は寧ろ高い結果となった。また秋季においては、地すべり内では活性がかなり下がる樹木もあるものの、地すべり外と同様の活性を示す樹木も存在する。活発な地すべり内には裸地に加え活性の高い植生も存在するために、緩慢に滑動する地すべりよりも NDVI の平均値は低くはなるが、その値のばらつきは大きくなる。
36

その一方で、地すべりの移動に伴い位置が移動した地すべりでは、夏季になっても植生の活性が回復しがたい例が見られた。本研究では無作為に選択した樹木に対して、現地にてその傾倒の程度や根切れについて計測していない為、この例について言及することは出来ない。しかし図 5.1.4.3 や 5.1.4.4 で示された結果については、現地における樹木の状態も観測することで、より明確に地すべり移動による植生活性の影響を示すことが出来ると考える。
図 5.1.3.1 からは、NDVI を利用して活動的な地すべりを抽出する場合には、斜面中腹に出現する値が低い範囲に注目すべきことが分かる。地すべりの滑落崖と直下であれば、植生が生育しないために NDVI の値は低いはずであり、それらは衛星画像からでも容易に抽出できるサイズである。しかし研究林地すべりのように、活動がかなり緩慢であっても、滑落崖や湿地、線状凹地や地すべりブロック境界に発達する沢などの地すべり特有の地形が出現する場所において、植生が生育しにくい例もある。本研究の結果からは、最近移動が発生した地すべりの抽出にあたっては、夏季における前年との NDVI の差を求めることが有効な手段であることが示唆された。夏季においては、地すべり内の樹木の活性が高くなる傾向があるが、移動した地すべりについては、移動がなかった時に獲得できるはずの活性を得ていない可能性がある。また、樹木が地すべり土塊の上で形状や根を維持したま
ま移動した場合は、移動した場において、前年度より植生活性が上昇する場合もありえる。
移動で生じた裸地においては、草本が生育する例もあるが、周辺からの水土砂の供給が頻繁にあるならば、裸地がそのまま維持される場もある。そのような地すべり地において、夏季における前年との NDVI の差は、各箇所によってかなりばらつきを持つことになる。
このばらつきが大きい範囲が、最近滑動した地すべりである可能性がある。この範囲は、
10m 四方における平均値を、フォーカル統計を用いて求めることで、より明瞭に把握することができる。
なお、地すべり土塊の顕著な移動が発生した後について、秋季における前年との NDVIの差については不明なままである。地すべり移動によりダメージを受けた樹木の活性が夏季に回復しないままであれば、秋季の NDVI 値はより低下することも考えられる。従って
秋季の NDVI 値の差分についても、夏季と同様に最近土塊の移動があった地すべりの場所を示唆できる可能性がある。

参考文献
後藤ら (1996) 人工衛星データによる植物活性を導入した地すべりの発生予測法, 地すべり, 第 33 巻, 第 2 号, p. 25-34, https://doi.org/10.3313/jls1964.33.25
Huang, S., Tang, L., Hupy, J.P. et al. (2021) A commentary review on the use of normalized
difference vegetation index (NDVI) in the era of popular remote sensing. J. For. Res. 32, 1–6(2021). https://doi.org/10.1007/s11676-020-01155-1
Jacquemart, M. and Tiampo, K. (2021) Leveraging time series analysis of radar coherence andnormalized difference vegetation index ratios to characterize pre-failure activity of the Mud
Creek landslide, California Natural Hazards and Earth System Sciences, 21, 629–642.
Lague D., Brodu N., and Leroux J. (2013) Accurate 3D comparison of complex topography with37

terrestrial laser scanner: Application to the Rangitikei canyon (N-Z). ISPRS Journal of
Photogrammetry and Remote Sensing, 82, 10-26.
https://doi.org/10.1016/j.isprsjprs.2013.04.009.
Reilly, S., Clark, M.L., Bentley, L.P., et al.. (2021) The Potential of Multispectral Imagery and 3D
Point Clouds from Unoccupied Aerial Systems (UAS) for Monitoring Forest Structure and the
Impacts of Wildfire in Mediterranean-Climate Forests, Remote Sensing 13, no. 19: 3810.
https://doi.org/10.3390/rs13193810
Thielicke,W.

(2022),

PIVlabparticle

image

velocimetry

(PIV)

tool

with

GUI,

https://jp.mathworks.com/matlabcentral/fileexchange/27659
Zhang W, Qi J, Wan P, Wang H, et al. An Easy-to-Use Airborne LiDAR Data Filtering Method Basedon Cloth Simulation. Remote Sensing. 2016; 8(6):501.

38