Rでラスターの解像度を変更する

Rのstarsパッケージはtidyverseとの相性もよく、またsfパッケージと同じ感覚で操作できます。ggplotとの相性が良い点でも好ましいパッケージです。しかしながら、巨大なラスタを操作しようとするとRstudioが落ちてしまいます。今回は、terraパッケージを使ってラスタ解像度を調整しました。

library(sf)
library(stars)
library(tidyverse)

基盤地図情報の10mメッシュと国土数値情報の行政界、北海道教育委員会の包蔵地GISデータを読み込みます。

# データ読込
dem <-
  read_stars("merge3_utm54")
# 行政界
hym <-
  st_read("Kokudo_gyousei_H28_utm54.shp") %>% 
  filter(
    N03_002 == "渡島総合振興局" |  # 渡島振興局と檜山振興局の領域を選択
      N03_002 == "檜山振興局"
  ) 
# 北の遺跡案内
site <-
  st_read("site_JGD2000_UTMzone54.gpkg") 

starsパッケージのオリジナルの描画はシンプルで好ましいです。

plot(dem)

遺跡データはこんな感じです。

plot(site)

DEMを渡島檜山の領域でクリップします。starsパッケージはsfパッケージのような操作感で記述できるので使いやすいです。

# 渡島檜山の領域でDEMをクロップ  
dem <-
  st_crop(dem, hym)
plot(dem)

包蔵地も渡島檜山の領域でクリップします。

# 渡島檜山の領域で遺跡をクロップ
site <-
  st_intersection(site, hym) 
plot(site)

ggplotで描画する

せっかくstarsパッケージを利用するのでさくっとggplotで描画したいです。シンプルな記述に感動します。

ggplot() +
  geom_stars(data = dem) +
  geom_sf(data = site)

止まりました・・・
どうやらメモリを消費しすぎているみたいです。道南全体の10mDEMが重すぎるのだと思います。なので解像度を変更したいと思います。

ラスターの解像度を変える

試行錯誤したのですが、starsパッケージで読み込んだDEMの解像度を変えようとすると、結局巨大なファイルサイズがネックとなり、途中でRstudioが止まってしまいます。そこで、terraパッケージで読み込んで、解像度を変更した後、stars形式に戻します。

terraパッケージは高速で動作する上に、なぜかメモリの負荷も小さいみたいです。

library(terra)
# terra経由で読み込み
dem_ag <-
  rast("../../../GIS_Data/kiban10mDEM_dounan_utm/merge3_utm54") %>% # terra形式で読み込み
  aggregate(fact = 20) %>% # terraパッケージのaggregate関数で解像度を20分1にする
  st_as_stars() %>% # stars形式に変換
  st_crop(hym) # 渡島檜山の領域でクリップ
plot(dem_ag)

先ほどと同じように描画されました。ただし、解像度は20倍したので200mメッシュになっているはずです。

いよいよ、ggplotで描画します。

# ggplotに描画
ggplot() +
  geom_stars(data = dem_ag) +
  geom_sf(data = site)

今度は無事に描画されました。

見栄えを整えるため、cptctyパッケージとggthemesパッケージを使います。cptcityパッケージの使い方とカラーパレットの実例はこちら（地図データの可視化に役立つcptcityのカラーパレットをためす）にまとめてあります。

library(cptcity)
library(ggthemes)
# ggplotに描画
ggplot() +
  geom_stars(data = dem_ag) +
  geom_sf(data = site,
          size = 0.8,
          alpha = 0.7) +
  scale_fill_gradientn(
    colours = cptcity::cpt(
      pal = "esri_hypsometry_eu_europe_5",
      n = 10)
  ) +
  theme_map() +
  theme(
    legend.position = "none",
    panel.background = element_rect(fill = "white")　# NA値を白塗り
  )

きれいに描画されました。

まとめ

RのGISは重いDEMの操作が苦手な傾向にあります。今回のように単に背景地図として利用したい場合は必要な解像度まで落として利用することになりますが、あまりにもデータが大きい場合は解像度を落とすことさえできないことがあります。

今回利用したterraパッケージは高速に動作することで知られていますが、Rasterパッケージやstarsパッケージでは処理が落ちてしまうような場合にも利用することができました。Rのラスタデータ処理は決定版的なパッケージはなく、当面は得意不得意に応じてパッケージを使い分けていくことになるのだと思います。