Access via Dask read_parquet

[mikeskaug]

This is a wonderful project. Thanks for your work on this.

I'm primarily a python user, so I wanted to explore accessing the data using dask.dataframe.read_parquet or pandas.read_parquet. Given the size of the dataset, I decided to start with Dask because it allows loading select row groups based on filters.

Loading the places data using Dask's read_parquet kind of works, but I've encountered a couple challenges. This is where I am at:

df = dd.read_parquet(
    's3://overturemaps-us-west-2/release/2023-07-26-alpha.0/theme=places/type=*/*',
    filters=[
        ('bbox.minx', '>', -105.296),
        ('bbox.minx', '<', -105.251),
        ('bbox.miny', '>', 40.000),
        ('bbox.miny', '<', 40.033)
    ],
    columns=['names', 'bbox.minx', 'bbox.miny'],
    engine='fastparquet',
    storage_options={"anon": True}
)

Two main challenges:

1. The filters argument doesn't do anything. df.compute() returns a dataframe with ~59 million rows, which I believe is the full dataset? I'm not a parquet expert, but I think the filters will only be applied if there are column stats available for each partition.
2. I haven't figured out how to load the "geometry" column. I understand that the geometry column is WKB, so I've tried loading it as bytes, but dask always tries to decode it as UTF-8, which fails.

Any advice would be appreciated. I also wanted to put this here in case it helps anyone else trying the same approach.

[TomAugspurger]

The filters argument doesn't do anything

I can't tell from the documentation at https://arrow.apache.org/docs/python/generated/pyarrow.parquet.ParquetDataset.html whether referencing struct members like 'bbox.min' is supported. Maybe @jorisvandenbossche knows off hand? Perhaps you need to make some kind of pyarrow.compute.Expression, but I couldn't sort it out.

I haven't figured out how to load the "geometry" column. I understand that the geometry column is WKB, so I've tried loading it as bytes, but dask always tries to decode it as UTF-8, which fails.

If you load the data as binary, you can use geopandas.Geoseries.from_wkb to parse it into polygons:

import dask.dataframe as dd
import geopandas
import dask_geopandas

df = dd.read_parquet(
    'az://release/2023-07-26-alpha.0/theme=places/type=place/*',
    columns=['bbox', 'geometry'],
    engine='pyarrow',
    dtype_backend="pyarrow",
    storage_options={"anon": True, "account_name": "overturemapswestus2"},
    parquet_file_extensions=False,
)

geometry = df["geometry"].map_partitions(geopandas.GeoSeries.from_wkb, meta=geopandas.GeoSeries(name="geometry"))
gdf = dask_geopandas.from_dask_dataframe(df, geometry=geometry)
gdf.head()

outputs

                                                bbox                     geometry
0  {'minx': -4.175979, 'maxx': -4.175979, 'miny':...    POINT (-4.17598 53.21174)
1  {'minx': -58.53305, 'maxx': -58.53305, 'miny':...  POINT (-58.53305 -34.58090)
2  {'minx': -46.5936956, 'maxx': -46.5936956, 'mi...  POINT (-46.59370 -23.51118)
3  {'minx': 99.555783, 'maxx': 99.555783, 'miny':...    POINT (99.55578 11.55832)
4  {'minx': 98.7590515, 'maxx': 98.7590515, 'miny...     POINT (98.75905 8.41062)

[mikeskaug]

Thank you for the pointers. That definitely solves parsing the geometry. I didn't know about the dask_geopandas project, so thanks for pointing that out!

Here's where I'm at now. The result is filtered by bbox, but I assume this is getting applied only after reading full partitions, so not benefiting from knowing column metadata.

import dask.dataframe as dd
import geopandas
import dask_geopandas

df = dd.read_parquet(
    's3://overturemaps-us-west-2/release/2023-07-26-alpha.0/theme=places/type=*/*',
    columns=['bbox', 'geometry'],
    engine='pyarrow',
    dtype_backend="pyarrow",
    storage_options={"anon": True},
)

def get_field(s, field):
    return s.apply(lambda v: v[field])
    
df['bbox.minx'] = df['bbox'].map_partitions(get_field, 'minx', meta=pd.Series(dtype=float))
df['bbox.miny'] = df['bbox'].map_partitions(get_field, 'miny', meta=pd.Series(dtype=float))

geometry = df["geometry"].map_partitions(geopandas.GeoSeries.from_wkb, meta=geopandas.GeoSeries(name="geometry"))
gdf = dask_geopandas.from_dask_dataframe(df, geometry=geometry)

gdf = gdf[
    (gdf['bbox.minx'] > -105.296)
    & (gdf['bbox.minx'] < -105.251)
    & (gdf['bbox.miny'] > 40.000)
    & (gdf['bbox.miny'] < 40.033)
]

gdf.head()

output:

bbox | geometry | bbox.minx | bbox.miny
{'minx': -105.28189, 'maxx': -105.28189, 'miny... | POINT (-105.28189 40.01656) | -105.281890 | 40.016560
{'minx': -105.252388, 'maxx': -105.252388, 'mi... | POINT (-105.25239 40.01725) | -105.252388 | 40.017250
{'minx': -105.276607, 'maxx': -105.276607, 'mi... | POINT (-105.27661 40.01818) | -105.276607 | 40.018185
{'minx': -105.2618899, 'maxx': -105.2618899, '... | POINT (-105.26189 40.02301) | -105.261890 | 40.023010
{'minx': -105.26339, 'maxx': -105.26339, 'miny... | POINT (-105.26339 40.02377) | -105.263390 | 40.023772

[ksmithNau]

I was also able to filter down using a shapely Polygon and the GeoDataframe.geometry.within() method, same idea.

import dask.dataframe as dd
import dask_geopandas as dgpd
import geopandas as gpd
from shapely.geometry import Polygon

minx, miny, maxx, maxy = -152.37741936472457, 60.27885301659569, -147.8292574609619, 61.99921515608003
bbox = Polygon([
    (minx, miny),
    (minx, maxy),
    (maxx, maxy),
    (maxx, miny),
    (minx, miny)
])

df = dd.read_parquet('s3://overturemaps-us-west-2/release/2023-07-26-alpha.0/theme=places/type=*/*',
    engine='fastparquet',
    columns=['geometry','id','names','sources','categories.main','categories.alternate','websites','socials', 'brand.names','brand.wikidata'],
    storage_options={"anon": True}
    )

geometry = df["geometry"].map_partitions(gpd.GeoSeries.from_wkb, meta=gpd.GeoSeries(name="geometry"))
gdf = dgpd.from_dask_dataframe(df, geometry=geometry)

anchorage_gdf = gdf[gdf.geometry.within(bbox)]
print(anchorage_gdf.head())

[mikeskaug]

I like that. Using a Polygon is more flexible as well.