BigQueryにShpeFileデータを取り込む

shapefileを読み込んだ後のsfクラス内から該当データを直接取り出して加工します。geometryをキャラクターに変換するのはst_as_text関数を使います。

### ShpeFileからBigQuery用のWKTを作る
# 目的:Google datapotalでの行政界ポリゴンがひどすぎる(超きたない)から自分で作る
# ESRI市区町村界ShapeFileを使っているが、この方法では規定違反になるので取り扱い注意
#
library(sf)
x <- st_read("japan_ver821/japan_ver821.shp")
x2 <- x$geometry
x3 <- st_as_text(x2)
z2 <- x$JCODE
z3 <- as.data.frame(z2)
df1 <- cbind(z3,x3)
df2 <- na.omit(df1)
write.table(df2,"japan_ver821.csv",row.names = F,col.names = F, sep = ",")

GoogelCloud上での後処理

GoogleCloudではGoogleCloudStrageとBigQueryを使います。 生成したcsvをGoogleCloudStrageにアップロードしてから、BigQueryでテーブルを作成します。 csvが10M以下であればBigQueryに直接アップロードできますが、GCS>BQの手順で覚えたほうが良いです。 テーブル作成時にはスキーマを以下のように手動で作ってください。 (スキーマ作成を自動でやると、JCODEを数値と認識してしまい、またWKT変換エラーが発生します。)

# スキーマ
# JCODE	STRING	NULLABLE	
# WKT_TEXT	STRING	NULLABLE

WKT_TEXTはSTRINGで読み込まれていますので、別途クエリでWKTに変換することでBigQuery用のWKTになります。スキーマ指定で直接WKTを生成するとエラーになるためにこの方法を使います。

# WKT変換用クエリ
# SELECT *,ST_GEOGFROMTEXT(WKT_TEXT,make_valid => TRUE) AS WKT FROM `データセット名.テーブル名` 

BigQueryのST_GEOGFROMTEXT関数を使います。make_valid => TRUE の部分でポリゴンがループしている部分のエラーを修正します。直接WKTを生成するとこの機能がないためエラーになります。

※ESRIのシェープファイルはMULTIPOLYGONになっています。JCODEとマルチポリゴンが一対です。

国土地理院の行政区域データを小さくする(ローポリ化)

なぜ小さくするのかというと、国土地理院の行政区域データは非常に細かく出来ていますが、利用したいデータポータルでWKTを表示する場合ポリゴンのポイント総数が10000ポイント以内という制限があるためローポリ化して表示させる必要があります。 表示系に制限がなければ小さくする必要はありません。

ローポリ化の方法はマップシェーパーを使います。

https://mapshaper.org/

やりかたは、以下のページを参照してください。非常に簡単にローポリ化できます。

https://humo-life.net/memo/doku.php?id=%E3%82%B5%E3%83%BC%E3%83%90%E3%82%BD%E3%83%95%E3%83%88%E3%82%A6%E3%82%A7%E3%82%A2:mapshaper:web%E7%89%88%E3%81%AE%E4%BD%BF%E3%81%84%E6%96%B9

ローポリ後のshp,shx,dbfを使ってWKTに変換する。(ここでは沖縄の行政区域データ)

library(sf)
xx <- st_read("N03-20_47_200101_2/N03-20_47_200101.shp")
xx2 <- xx$geometry
xx3 <- st_as_text(xx2)
zz2 <- xx$N03_007   #JCODEの場所
zz3 <- as.data.frame(zz2)
ddf1 <- cbind(zz3,xx3)
ddf2 <- na.omit(ddf1)
write.table(ddf2,"okinawa.csv",row.names = F,col.names = F, sep = ",")

※国土地理院のシェープファイルはPOLYGONになっています。そのため同じJCODEが複数あります。

rshapemapperでローポリ化も一発でやる方法

CRANにはrshapemapperパッケージがあるので、それを利用します。

library(sf)
library(rmapshaper)
yy <- st_read("N03-20_47_200101_2/N03-20_47_200101.shp",options = "ENCODING=CP932") #エンコード指定しないと変換時にエラーになる
yy1 <- ms_simplify(yy,keep = 0.05 )  #ここでローポリ化の指数を入れる。0.05は元データの5%に縮小される。
yy3 <- yy1$geometry
yy4 <- st_as_text(yy3)
zz2 <- yy1$N03_007
zz3 <- as.data.frame(zz2)
dddf1 <- cbind(zz3,yy4)
dddf2 <- na.omit(dddf1)
write.table(dddf2,"okinawa.csv",row.names = F,col.names = F, sep = ",")

コメント欄


2164


トップ   新規 一覧 検索 最終更新   ヘルプ   最終更新のRSS