从 python 中的地理数据框中的点创建地图边界

Question

我有一个名为 map 的地理数据框，其中包含点列表和一列，Closest_TrainStation_name 其中包含离该点最近的火车站的名称。

	geometry	Closest_TrainStation_name
1	POINT(1,1)	Station_1
2	POINT(10,10)	Station_2
...	...	...

是否可以像下图那样创建包含每个组的边界多边形？每个多边形都具有原始文件中最近的火车站的名称。

这些点是使用最近邻算法生成的，因此它们不会相互交叉。

我还有一个名为 boundary 的整个国家/地区的边界地理数据框，我认为可能需要它来定义此地图的外部边界。还有火车站的档案。

我能找到的所有方法都是关于从边界上的点创建边界的，例如convex_hull.

Answer 1

已经生成了一些示例数据
使用 https://geopandas.org/en/stable/docs/user_guide/aggregation_with_dissolve.html 的组合生成 MultiPoint shapely 对象然后 convex_hull 会给你一个边界关联点的多边形

import io
import pandas as pd
import geopandas as gpd
import shapely.wkt

df = pd.read_csv(
    io.StringIO(
        """geometry,stationname
POINT (6.84365570015329 53.29316918283323),station 0
POINT (6.84397077786025 53.34219781242178),station 0
POINT (6.86411903336562 53.31761406936343),station 0
POINT (6.818510924831408 53.31265038241855),station 0
POINT (6.56484740742038 53.29280081541883),station 1
POINT (6.626702808494024 53.28117175507543),station 1
POINT (6.593563692631294 53.3071306760684),station 1
POINT (6.578341784305515 53.33493210605594),station 1
POINT (6.622799032233983 53.28010824623466),station 1
POINT (6.613683981274606 53.29073387043917),station 1
POINT (7.132754969933528 53.04549586270261),station 2
POINT (7.198050554432563 53.1025440484535),station 2
POINT (7.129751539288898 53.05572095828281),station 2
POINT (7.176645494114577 53.11826239684777),station 2
POINT (7.197200682695362 53.07260800358519),station 2
POINT (7.162805175737367 53.10861252890756),station 2
POINT (7.109351059746762 53.13349222863155),station 2
POINT (7.06391508890864 53.11204313452662),station 2
POINT (7.184832998096279 53.05647939687099),station 2
POINT (7.159339780758176 53.06749816133893),station 2
POINT (7.107083795756558 53.08269279330896),station 2
POINT (7.023547595125693 53.07391566675645),station 2
POINT (7.150033380903857 53.1339122880752),station 2
POINT (7.087567734179277 53.13629412269268),station 2
POINT (7.082082576050039 53.09928180316115),station 2
POINT (7.163037435775678 53.10643603535302),station 2
POINT (7.188994818425656 53.11808478608855),station 2
POINT (7.096403367400944 53.09224186363327),station 2
POINT (7.166715752727275 53.13534692250609),station 2
POINT (7.153873895859038 53.15187182174655),station 2"""
    )
)

gdf = gpd.GeoDataFrame(df, geometry=df["geometry"].apply(shapely.wkt.loads)).dissolve(
    "stationname"
)["geometry"].apply(lambda mp: mp.convex_hull)

gdf.plot()

stationname	geometry
station 0	POLYGON ((6.84365570015329 53.29316918283323, 6.818510924831408 53.31265038241855, 6.84397077786025 53.34219781242178, 6.86411903336562 53.31761406936343, 6.84365570015329 53.29316918283323))
station 1	POLYGON ((6.622799032233983 53.28010824623466, 6.56484740742038 53.29280081541883, 6.578341784305515 53.33493210605594, 6.626702808494024 53.28117175507543, 6.622799032233983 53.28010824623466))
station 2	POLYGON ((7.132754969933528 53.04549586270261, 7.023547595125693 53.07391566675645, 7.087567734179277 53.13629412269268, 7.153873895859038 53.15187182174655, 7.188994818425656 53.11808478608855, 7.198050554432563 53.1025440484535, 7.197200682695362 53.07260800358519, 7.184832998096279 53.05647939687099, 7.132754969933528 53.04549586270261))

Answer 2

使用英格兰的医院作为数据源（有效的车站和医院，可互换）
具有作为英格兰边界的几何图形（剪切 Voronoi 多边形）
然后生成代表所有医院的多边形变得简单
使用sjoin()将多边形与承载所有关联属性的点相关联
已使用 plotly 来可视化和演示多边形具有与之关联的医院属性

import shapely.ops

# points for hospitals
gdf = gpd.GeoDataFrame(
    geometry=dfhos.loc[:, ["Longitude", "Latitude"]].apply(
        shapely.geometry.Point, axis=1
    )
)

# generate voroni polygon for each hospital
gdfv = gpd.GeoDataFrame(
    geometry=[
        p.intersection(uk)
        for p in shapely.ops.voronoi_diagram(
            shapely.geometry.MultiPoint(gdf["geometry"].values)
        ).geoms
    ]
)

# spatial join polygons to points to pick up full details of hospital
gdf3 = gpd.sjoin(gdfv, gdf, how="left").merge(
    dfhos, left_on="index_right", right_index=True
)
gdf3["Color"] = pd.factorize(gdf3["Postcode"], sort=True)[0]

# and visualize
fig = (
    px.choropleth_mapbox(
        gdf3,
        geojson=gdf3.__geo_interface__,
        locations=gdf3.index,
        hover_data=["OrganisationCode","OrganisationName","Postcode"],
        color="Color",
        color_continuous_scale="phase",
    )
    .update_layout(
        mapbox={
            "style": "carto-positron",
            "center": {
                "lon": sum(gdf3.total_bounds[[0, 2]]) / 2,
                "lat": sum(gdf3.total_bounds[[1, 3]]) / 2,
            },
            "zoom": 5,
        },
        margin={"l": 0, "r": 0, "t": 0, "b": 0},
        coloraxis={"showscale":False}
    )
)

fig

获取英格兰医院的位置和英格兰边界的多边形

import geopandas as gpd
import shapely.geometry
import numpy as np
import plotly.express as px
import requests, io
from pathlib import Path
from zipfile import ZipFile
import urllib
import pandas as pd

# fmt: off
# uk geometry
url = "http://geoportal1-ons.opendata.arcgis.com/datasets/687f346f5023410ba86615655ff33ca9_1.zip"
f = Path.cwd().joinpath(urllib.parse.urlparse(url).path.split("/")[-1])

if not f.exists():
    r = requests.get(url, stream=True, headers={"User-Agent": "XY"})
    with open(f, "wb") as fd:
        for chunk in r.iter_content(chunk_size=128):
            fd.write(chunk)
    zfile = ZipFile(f)
    zfile.extractall(f.stem)

f2 = Path.cwd().joinpath("uk.geojson")
if not f2.exists():
    gdf2 = gpd.read_file(list(f.parent.joinpath(f.stem).glob("*.shp"))[0])
    gdf2 = gdf2.loc[gdf2["ctyua16cd"].str[0] == "E"]
    uk = gpd.GeoDataFrame(geometry=[p for p in shapely.ops.unary_union(gdf2.to_crs(gdf2.estimate_utm_crs())["geometry"].values).simplify(5000).geoms]).set_crs(gdf2.estimate_utm_crs()).to_crs("EPSG:4326")

    uk.to_file(Path.cwd().joinpath("uk.geojson"), driver='GeoJSON')
uk = gpd.read_file(f2)
uk = shapely.geometry.MultiPolygon(uk["geometry"].values)
# fmt: on

# get hospitals in UK
dfhos = pd.read_csv(io.StringIO(requests.get("https://assets.nhs.uk/data/foi/Hospital.csv").text),sep="Č",engine="python",)
dfhos = dfhos.loc[lambda d: d["Sector"].eq("NHS Sector") & d["SubType"].eq("Hospital")].groupby("ParentODSCode").first()

从 python 中的地理数据框中的点创建地图边界

Create map boundaries from points within a geodataframe in python

python

mapping

gis

geopandas

获取英格兰医院的位置和英格兰边界的多边形