Go to Qiita Advent Calendar 2024 Top
LoginSignup
27
23

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

More than 1 year has passed since last update.

@nokonoko_1203(のこのこ)inMIERUNE

GeoPandasをやるならFlatGeobufより10倍早いGeoParquetを使おう！

Last updated at Posted at 2022-06-06

GeoParquetとは

GIS界隈では、AIやビッグデータ処理などデータサイエンスの分野で利用されるようになっている高速かつ効率の良いデータ処理手法などを全て取り込めてる訳ではなく、若干遅れをとっている状態と言われています。

これは地理空間情報というもの自体が複雑な仕様をもつためかとおもいますが、近年はGISでもGB/TB級のデータを取り扱う機会が増えており、Hadoopに代表されるような表形式データの分散処理基盤などを有効活用できない点において、「どげんかせんといかんなぁ」と考える方々も増えています。

そんなモチベーションからここ最近活発に開発されている形式のデータがGeoParquetです。

GeoParquet: https://github.com/opengeospatial/geoparquet

Apache Parquetという、こちらも近年活発に開発されている形式のデータを地理空間情報に応用しよう！として開発されているもので、シンプルに言えば「高効率に圧縮された表形式のバイナリデータ」ですが、ここに(またもや活発に開発されているApache Arrow形式の地理空間活用である)Geo Arrowを取り込み、列形式の地理空間情報をより上手に処理しようという試みのようです。

Apache Parquet: https://parquet.apache.org/docs/overview/motivation/

Geo Arrow: https://github.com/geopandas/geo-arrow-spec

現状、使えそうな感じ？

2022/06/07時点ではv0.4と実験的なファイル形式ですが、数ヶ月以内のメジャーアップデートを目標としており、現時点でも触ってみた感じかなり効率よくデータを処理することができました。

  • 地理座標系/投影座標系を処理できる
  • geometry列を複数持てる(！！！)
  • 高圧縮
  • チャンク/パーティションを区切って効率よく並列で計算・処理
  • 空間インデックスが有効

FlatGeobufと比較して、どうよ

Web上でのクラウドネイティブなベクトルデータとしての評判は、以下のように見えます。

  • GeoParquetはまだ実験的な形式であるため、FlatGeobufを使うのが良いだろう
  • ただし、将来的にはApache Parquetのデータ処理基盤を利用できる点、優れた圧縮技術やパーティションごとのストリーミングが可能な点などで優位に働くだろう

他のファイル形式と比較してみよう

ということで早速適当にデータを作成して、いろいろな形式でread/wtierしてみましょう。

1000万レコードのポイントデータを作成する

まずはデータを作成していきましょう。

jupyterを利用していきますが、環境構築がなされていない方はこちらを参考にどうぞ。

以下のようにして、GeoPandas・numpy・shapelyを用いてランダムでポイントデータを作成していきます。

%%time
import geopandas as gpd
import numpy as np
from shapely.geometry import box

%%time
box = box(-180, -90, 180, 90)
data = {"geometry": [box.wkt]}

poly = gpd.GeoSeries.from_wkt(data["geometry"]).unary_union

x_min, y_min, x_max, y_max = poly.bounds
# CPU times: user 1.28 ms, sys: 1.09 ms, total: 2.37 ms
# Wall time: 1.48 ms

%%time
n = 10_000_000

x = np.random.uniform(x_min, x_max, n)
y = np.random.uniform(y_min, y_max, n)
# CPU times: user 150 ms, sys: 30.1 ms, total: 180 ms
# Wall time: 179 ms

%%time
points = gpd.GeoSeries(gpd.points_from_xy(x, y))
points = points[points.within(poly)]
# CPU times: user 3.39 s, sys: 387 ms, total: 3.78 s
# Wall time: 3.78 s

%%time
gdf = gpd.GeoDataFrame(geometry=points)
# CPU times: user 406 ms, sys: 108 ms, total: 513 ms
# Wall time: 512 ms

4~5秒程度で1000万レコードを作成し終えました。早いですね。

作成したGeoDataFrameをいろんな形式で出力する

  • GeoParquet
%%time
# GeoParquet
gdf.to_parquet(
    "export.parquet",
    index=False,
    compression="brotli",
)
# CPU times: user 28.8 s, sys: 2.2 s, total: 31 s
# Wall time: 31.5 s
  • GeoPackage
%%time
# GeoPackage
gdf.to_file("export.gpkg", index=False, driver="GPKG", layer="layer")
# CPU times: user 11min 7s, sys: 2min 28s, total: 13min 36s
# Wall time: 13min 37s
  • FlatGeobuf
%%time
# FlatGeobuf
gdf.to_file("export.fgb", index=False, driver="FlatGeobuf", spatial_index="NO")
# CPU times: user 7min 55s, sys: 5.25 s, total: 8min
# Wall time: 8min 6s
  • GeoJSON
%%time
# GeoJSON
gdf.to_file("export.geojson", index=False, driver="GeoJSON")
# CPU times: user 9min 47s, sys: 8.96 s, total: 9min 56s
# Wall time: 10min
  • Shapefile
%%time
# Shapefile
gdf.to_file("export.shp", index=False, driver="ESRI Shapefile")
# CPU times: user 8min 29s, sys: 1min, total: 9min 29s
# Wall time: 9min 30s

こんな関数を作ってファイルサイズを人間にわかりやすい形式でファイルサイズを出力してみましょう。

参考: https://pystyle.info/python-data-size-conversion/

import math

def convert_size(size):
    units = ("B", "KB", "MB", "GB", "TB", "PB", "EB", "ZB")
    i = math.floor(math.log(size, 1024)) if size > 0 else 0
    size = round(size / 1024**i, 2)

    return f"{size} {units[i]}"

出力結果はこのようになりました。

import os

geo_parquet = os.path.getsize("export.parquet")
geo_package = os.path.getsize("export.gpkg")
flat_geobuf = os.path.getsize("export.fgb")
geojson = os.path.getsize("export.geojson")
shapefile = (
    os.path.getsize("export.shp")
    + os.path.getsize("export.shx")
    + os.path.getsize("export.cpg")
    + os.path.getsize("export.dbf")
)

print(f"geo_parquet={convert_size(geo_parquet)}")
print(f"geo_package={convert_size(geo_package)}")
print(f"flat_geobuf={convert_size(flat_geobuf)}")
print(f"geojson={convert_size(geojson)}")
print(f"shapefile={convert_size(shapefile)}")

# geo_parquet=150.0 MB
# geo_package=903.49 MB
# flat_geobuf=610.35 MB
# geojson=1.26 GB
# shapefile=457.76 MB

出力したファイルを読み込んでみる

出力したファイルをGeoDataFrameに読み込ませて時間を計測してみましょう。

  • GeoParquet
%%time
# GeoParquet
geo_parquet_gdf = gpd.read_parquet("export.parquet", columns=["geometry"])
geo_parquet_gdf.head()
# CPU times: user 4.66 s, sys: 960 ms, total: 5.62 s
# Wall time: 5.67 s
# geometry
# 0	POINT (84.21447 -60.10291)
# 1	POINT (-72.92868 53.36820)
# 2	POINT (-64.87842 44.19289)
# 3	POINT (-179.93497 -72.37761)
# 4	POINT (26.41531 3.05427)
  • GeoPackage
%%time
# GeoPackage
geo_package_gdf = gpd.read_file("export.gpkg")
geo_package_gdf.head()
# CPU times: user 3min 8s, sys: 3.86 s, total: 3min 12s
# Wall time: 3min 12s
# geometry
# 0	POINT (84.21447 -60.10291)
# 1	POINT (-72.92868 53.36820)
# 2	POINT (-64.87842 44.19289)
# 3	POINT (-179.93497 -72.37761)
# 4	POINT (26.41531 3.05427)
  • FlatGeobuf
%%time
# FlatGeobuf
flat_geobuf_gdf = gpd.read_file("export.fgb")
flat_geobuf_gdf.head()
# CPU times: user 2min 54s, sys: 2.9 s, total: 2min 57s
# Wall time: 2min 57s
# geometry
# 0	POINT (84.21447 -60.10291)
# 1	POINT (-72.92868 53.36820)
# 2	POINT (-64.87842 44.19289)
# 3	POINT (-179.93497 -72.37761)
# 4	POINT (26.41531 3.05427)
  • GeoJSON
%%time
# GeoJSON
geojson_gdf = gpd.read_file("export.geojson")
geojson_gdf.head()
# CPU times: user 3min 57s, sys: 3.74 s, total: 4min 1s
# Wall time: 4min 1s
# geometry
# 0	POINT (84.21447 -60.10291)
# 1	POINT (-72.92868 53.36820)
# 2	POINT (-64.87842 44.19289)
# 3	POINT (-179.93497 -72.37761)
# 4	POINT (26.41531 3.05427)
  • Shapefile
%%time
# Shapefile
shapefile_gdf = gpd.read_file("export.shp")
shapefile_gdf.head()
# CPU times: user 3min 40s, sys: 4.54 s, total: 3min 45s
# Wall time: 3min 46s
# FID	geometry
# 0	0	POINT (84.21447 -60.10291)
# 1	1	POINT (-72.92868 53.36820)
# 2	2	POINT (-64.87842 44.19289)
# 3	3	POINT (-179.93497 -72.37761)
# 4	4	POINT (26.41531 3.05427)

作成したGeoDataFrameの最初の5行を取得してみましたが、すべて同じジオメトリが出力されていますね。
レコード数を確認しましたが、特にかけている様子もありませんでした。

%%time
len(geo_parquet_gdf)
# CPU times: user 11 µs, sys: 2 µs, total: 13 µs
# Wall time: 16.2 µs
# 10000000

%%time
len(geo_package_gdf)
# CPU times: user 13 µs, sys: 1e+03 ns, total: 14 µs
# Wall time: 17.2 µs
# 10000000

%%time
len(flat_geobuf_gdf)
# CPU times: user 10 µs, sys: 1e+03 ns, total: 11 µs
# Wall time: 13.8 µs
# 10000000

%%time
len(geojson_gdf)
# CPU times: user 11 µs, sys: 1e+03 ns, total: 12 µs
# Wall time: 15.3 µs
# 10000000

%%time
len(shapefile_gdf)
# CPU times: user 12 µs, sys: 1e+03 ns, total: 13 µs
# Wall time: 32.9 µs
# 10000000

結果

read/writeの結果、以下のようになりました。

FlatGeobufのread/writeが思ったより遅いかったりと、気になる点はありますが、どちらにしろGeoParquetが圧倒的なパフォーマンスで出力できることがわかったかと思います。

今回使ったnotebookはこちらに置いておきます。

write

  1. GeoParquet: 31.5 s
  2. FlatGeobuf: 8min 6s
  3. Shapefile: 9min 30s
  4. GeoJSON: 10min
  5. GeoPackage: 13min 37s

read

  1. GeoParquet: 5.67 s
  2. FlatGeobuf: 2min 57s
  3. GeoPackage: 3min 12s
  4. Shapefile: 3min 46s
  5. GeoJSON: 4min 1s

と、いうことでみなさんも積極的にGeoParquetを利用してみてはいかがでしょうか？

余談

ちなみに、ogr2ogrコマンドで利用するには2022/05上旬にリリースされた最新のGDALv3.5から利用できるようです。
https://gdal.org/drivers/vector/parquet.html?highlight=parquet

v3.4では利用できませんでした。

$ ogr2ogr --version
GDAL 3.4.3, released 2022/04/22

$ ogr2ogr -f "Parquet" data.Parquet data.geojson
ERROR 1: Unable to find driver `Parquet'.
27
23
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
27
23

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?