Rで解析:地図の描写が面白い「tmap」パッケージ


ggplot2パッケージと同様レイヤーの概念を持ち、数量データを示す地図を作成するのに便利なパッケージの紹介です。データをもとに分割して地図を作成することが可能です。

付与するデータ次第で面白い表現が可能と考えます。

パッケージバージョンは1.10。windows 10のR version 3.3.3で動作を確認しています。


シェープファイルを準備する

先日紹介した「Rで解析:「leaflet」パッケージで使える日本のシェープファイルの作成例」で作成したシェープファイルに気象庁より観測所一覧を取得してシェープファイルに付与します。

・Rで解析:「leaflet」パッケージで使える日本のシェープファイルの作成例
 https://www.karada-good.net/analyticsr/r-594

#パッケージの読み込み
library("tcltk")

#作成したシェープファイルを読み込む
#参考:https://www.karada-good.net/analyticsr/r-594
load(paste0(as.character(tkgetOpenFile(title = "ファイルを選択",
                                       filetypes = '{"ファイル" {".*"}}',
                                       initialfile = c("*.*")))))

###気象庁より観測所一覧を取得#####
#Tmpフォルダに保存し解凍後に読み込む
temp <- tempfile()
download.file(&quot;http://www.jma.go.jp/jma/kishou/know/amedas/ame_master.zip&quot;, temp, mode = &quot;wb&quot;)
Kansokujyo <- read.csv(unzip(temp))
unlink(temp)
#必要データを日本地図シェープファイルに付与して保存
MasterKansokujyo <- Kansokujyo&#91;, c(1, 4, 7:10)&#93;
#北海道地方を北海道に置換
MasterKansokujyo&#91;, 1&#93; <- c(rep(&quot;北海道&quot;, 227), as.character(MasterKansokujyo&#91;228:nrow(MasterKansokujyo), 1&#93;))
#緯度経度を変換
LatData <- type.convert(paste0(MasterKansokujyo&#91;, 3&#93; + MasterKansokujyo&#91;, 4&#93;/60))
LonData <- type.convert(paste0(MasterKansokujyo&#91;, 5&#93; + MasterKansokujyo&#91;, 6&#93;/60))
#シェープファイルに付与
library(&quot;stringi&quot;)
attributes(JPNPref)$Kansokujyo <- data.frame(&quot;都道府県&quot; = stri_encode(MasterKansokujyo&#91;, 1&#93;, &quot;&quot;, &quot;utf-8&quot;),
                                             &quot;観測所名&quot; = stri_encode(MasterKansokujyo&#91;, 2&#93;, &quot;&quot;, &quot;utf-8&quot;),
                                             &quot;経度&quot; = LatData, &quot;緯度&quot; = LonData)
#シェープファイルを保存
save(JPNPref, file = &quot;JPNPref.shp&quot;)
&#91;/code&#93;
<hr />

<h2>気象庁より最高気温最新を取得する</h2>
作成した日本地図シェープファイルを読み込み後、以下のコードを実行します。例では日本地図シェープファイルのデータ名は「JPNPref」です。


###気象庁より最高気温最新を取得#####
#http://www.data.jma.go.jp/obd/stats/data/mdrr/docs/csv_dl_readme.html
MaxTemp <- read.csv(&quot;http://www.data.jma.go.jp/obd/stats/data/mdrr/tem_rct/alltable/mxtemsadext00_rct.csv&quot;)

#各都道府県の最高気温を抽出
#データ保管用変数
NewMaxTemp <- data.frame()
#処理
for(n in 1:47){
#都道府県を抽出
GetPrefData <- MaxTemp&#91;which(MaxTemp&#91;, 2&#93; %in% grep(JPNPref@data$gns_name&#91;n&#93;, MaxTemp&#91;, 2&#93;, value = TRUE)),&#93;
#最高気温を降順で並び替え
GetPrefData <- GetPrefData&#91;order(GetPrefData&#91;, 10&#93;, decreasing = TRUE),&#93;
#観測所の対象都道府県
GetKansokujyo <- JPNPref@Kansokujyo&#91;which(JPNPref@Kansokujyo&#91;, 1&#93; %in%
                                            grep(substr(JPNPref@data$gns_name&#91;n&#93;, 1, 2),
                                                 as.character(JPNPref@Kansokujyo&#91;, 1&#93;), value = TRUE)),&#93;
#観測所の経度緯度を取得
LatLonData <- GetKansokujyo&#91;which(GetKansokujyo&#91;, 2&#93; %in% GetPrefData&#91;1, 3&#93;),&#93;&#91;1, 4:3&#93;
#1位の地点名と最高気温,経度緯度を結合
NewMaxTemp <- rbind(NewMaxTemp, cbind(JPNPref@data$gns_name&#91;n&#93;, GetPrefData&#91;1, c(3, 10)&#93;, LatLonData))}
#1列目と2列目を結合後に1列目を置き換え
NewMaxTemp&#91;, 1&#93; <- paste0(NewMaxTemp&#91;, 1&#93;, &quot;:&quot;, NewMaxTemp&#91;, 2&#93;)
#2列目を削除
NewMaxTemp <- NewMaxTemp&#91;, -2&#93;
#データ名を変更
colnames(NewMaxTemp) <- c(&quot;都道府県:地点名&quot;, &quot;気温&quot;, &quot;lat&quot;, &quot;lon&quot;)
#JPNPref@dataへ追加
JPNPref@data <- cbind(JPNPref@data, NewMaxTemp&#91;, 1:2&#93;)
#SpatialPointsDataFrameを作成
MaxSpDF <- SpatialPointsDataFrame(SpatialPoints(NewMaxTemp&#91;, 3:4&#93;,
                                                proj4string = CRS(&quot;+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs +towgs84=0,0,0&quot;)),
                                  data = data.frame(id = as.factor(NewMaxTemp&#91;, 1&#93;)))

#MaxSpDF,NewMaxTemp以外を削除
rm(GetPrefData, MaxTemp, n)
&#91;/code&#93;
<hr />

<h2>パッケージのインストール</h2>
下記コマンドを実行してください。


#パッケージのインストール
install.packages(&quot;tmap&quot;)

コマンドの紹介

詳細はコマンド、パッケージのヘルプを確認してください。

#パッケージの読み込み
library("tmap")

#プロット出力を選択する:tmap_modeコマンド
#プロット作成前に実行します
#ggplot2でプロット
#オプションをplotにする
tmap_mode("plot")
#leafletでプロット
#オプションをviewにする
tmap_mode("view")

#簡易プロット:qtmコマンド
qtm(JPNPref, fill = "気温", text = "都道府県:地点名", style = "gray",
text.root = 3, legend.show = TRUE)

#tm_shapeで指定した情報をもとにtm_XXXXコマンドレイヤーを重ねる
tm_shape(JPNPref) +
tm_polygons("気温", id = "都道府県:地点名", legend.show = TRUE) +
tm_shape(MaxSpDF) +
tm_bubbles(size = .01, "red") +
tm_style_gray()

#例えば気温で分割してプロット
#気温を区分分け:「fancycutp」パッケージを利用
#https://www.karada-good.net/analyticsr/r-544
#install.packages("fancycut")
library("fancycut")
#分類
CutLabel <- fancycut(x = JPNPref@data[, 5], intervals = c("[10, 20)", "[20, 25)", "[25, 30)", "[30, 35]"), buckets = c("10-20", "20-25", "25-30", "30-35"), unmatched.bucket = "範囲外") #結合 JPNPref@data <- cbind(JPNPref@data, CutLabel) #プロット tm_shape(JPNPref) + tm_polygons("気温", id = "都道府県:地点名", legend.show = TRUE) + tm_facets(by = "CutLabel") + tm_shape(MaxSpDF) + tm_bubbles(size = .01, "red") + tm_style_gray() [/code]


出力例

2017/05/13 13:00頃に取得したデータです。
・簡易プロット:qtmコマンド

・例えば気温で分割してプロット

・例えば気温で分割してインタラクティブにプロット
いろいろいじってみてください。
 ・別ウィンドウで開く
 https://www.karada-good.net/wp/wp-content/uploads/2017/05/tmapfacet.html


少しでも、あなたの解析が楽になりますように!!けいおんは素晴らしい。

スポンサードリンク

関連コンテンツ


スポンサードリンク