介绍
在本教程中,您将探索几种用于邻近度分析的技术。 特别是,您将学习如何执行以下操作:
- 测量地图上点之间的距离,以及
- 选择某要素半径内的所有点。
import folium
from folium import Marker, GeoJson
from folium.plugins import HeatMap
import pandas as pd
import geopandas as gpd
您将使用来自美国环境保护署(EPA)的数据集,该数据集跟踪美国宾夕法尼亚州费城的有毒化学物质的释放。
releases = gpd.read_file("../input/geospatial-learn-course-data/toxic_release_pennsylvania/toxic_release_pennsylvania/toxic_release_pennsylvania.shp")
releases.head()
您还将使用包含来自同一城市的空气质量监测站的读数的数据集。
stations = gpd.read_file("../input/geospatial-learn-course-data/PhillyHealth_Air_Monitoring_Stations/PhillyHealth_Air_Monitoring_Stations/PhillyHealth_Air_Monitoring_Stations.shp")
stations.head()
Measuring distance(测量距离)
要测量来自两个不同GeoDataFrame的点之间的距离,我们首先必须确保它们使用相同的坐标参考系统(CRS)。 幸运的是,这里就是这种情况,它们都使用EPSG 2272。
print(stations.crs)
print(releases.crs)
{'init': 'epsg:2272'}
{'init': 'epsg:2272'}
我们还检查CRS,以查看其使用的单位(米,英尺或其他)。 在这种情况下,EPSG 2272具有英尺单位。 (如果您愿意,可以在此处进行检查。)
在GeoPandas中计算距离相对简单。 下面的代码单元计算出last_release
中相对较新的版本事件与station
GeoDataFrame中的每个站点之间的距离(以英尺为单位)。
# Select one release incident in particular
recent_release = releases.iloc[360]
# Measure distance from release to each station
distances = stations.geometry.distance(recent_release.geometry)
distances
out:
0 44778.509761
1 51006.456589
2 77744.509207
3 14672.170878
4 43753.554393
5 4711.658655
6 23197.430858
7 12072.823097
8 79081.825506
9 3780.623591
10 27577.474903
11 19818.381002
dtype: float64
使用计算出的距离,我们可以获得诸如到每个站点的平均距离之类的统计信息。
print('Mean distance to monitoring stations: {} feet'.format(distances.mean()))
Mean distance to monitoring stations: 33516.28487007786 feet
或者,我们可以找到最近的监测站。
print('Closest monitoring station ({} feet):'.format(distances.min()))
print(stations.iloc[distances.idxmin()][["ADDRESS", "LATITUDE", "LONGITUDE"]])
Closest monitoring station (3780.623590556444 feet):
ADDRESS 3100 Penrose Ferry Road
LATITUDE 39.9128
LONGITUDE -75.1854
Name: 9, dtype: object
Creating a buffer(创建一个缓冲区)
如果我们想了解地图上所有距某个点一定半径的点,最简单的方法是创建一个缓冲区。
下面的代码单元创建一个GeoSeries two_mile_buffer
,其中包含12个不同的Polygon对象。 每个多边形是在不同的空气监测站周围2英里(或2 * 5280英尺)的缓冲区。
two_mile_buffer = stations.geometry.buffer(2*5280)
two_mile_buffer.head()
0 POLYGON ((2721944.640797138 257149.3104284704,...
1 POLYGON ((2682494.289907977 271248.9000113755,...
2 POLYGON ((2744886.638220146 280980.2466903776,...
3 POLYGON ((2703638.579968393 233247.1013432145,...
4 POLYGON ((2726959.772827223 251134.9763285518,...
dtype: object
我们使用`folium.GeoJson()`在地图上绘制每个多边形。 请注意,由于filium需要纬度和经度的坐标,因此在绘制之前必须将CRS转换为EPSG 4326。
Create map with release incidents and monitoring stations
m = folium.Map(location=[39.9526,-75.1652], zoom_start=11)
HeatMap(data=releases[['LATITUDE', 'LONGITUDE']], radius=15).add_to(m)
for idx, row in stations.iterrows():
Marker([row['LATITUDE'], row['LONGITUDE']]).add_to(m)
Plot each polygon on the map
GeoJson(two_mile_buffer.to_crs(epsg=4326)).add_to(m)
Show the map
m
现在,要测试在**任何**监测站的2英里范围内是否发生了有毒物质释放,我们可以对每个多边形运行12种不同的测试(以单独检查其是否包含该点)。
但是更有效的方法是先将所有多边形折叠成一个**MultiPolygon**对象。 我们使用`unary_union`属性来执行此操作。
Turn group of polygons into single multipolygon
my_union = two_mile_buffer.geometry.unary_union
print('Type:', type(my_union))
Show the MultiPolygon object
my_union
Type: <class 'shapely.geometry.multipolygon.MultiPolygon'>

我们使用`contains()`方法来检查多边形是否包含一个点。 我们将使用本教程前面部分中的发布事件,我们知道它离最近的监视站大约3781英尺。
The closest station is less than two miles away
my_union.contains(releases.iloc[360].geometry)
`out:True`
但是,并非所有的释放都发生在空气监测站的两英里范围内!
The closest station is more than two miles away
my_union.contains(releases.iloc[358].geometry)
`out:False`
# Your turn
在[最后的练习](https://www.kaggle.com/kernels/fork/5832147)中,您将调查纽约市的医院覆盖率。