[数据分析与可视化] Python绘制数据地图4-MovingPandas入门指北

MovingPandas是一个基于Python和GeoPandas的开源地理时空数据处理库,用于处理移动物体的轨迹数据。它提供了一组强大的工具,可以轻松地加载、分析和可视化移动物体的轨迹。通过使用MovingPandas,用户可以轻松地处理和分析移动对象数据,并从中提取有关行为、模式和趋势的见解。无论是处理交通流量数据、物流轨迹数据还是动物迁徙数据,MovingPandas都是一个强大的地理可视化工具。

MovingPandas库具有以下主要特性:

  • 轨迹对象:MovingPandas引入了Trajectory对象,用于表示和管理轨迹数据。该对象包含了时间、位置和属性信息,方便对轨迹进行操作和分析。
  • 时空切片:MovingPandas支持时空切片操作,可以按照时间和空间条件对轨迹数据进行筛选和分割。
  • 轨迹分析:MovingPandas提供了丰富的轨迹分析功能,包括计算轨迹长度、速度、方向、加速度等指标,以及轨迹聚类、交互和相似性分析。
  • 轨迹可视化:MovingPandas可以方便地将轨迹数据可视化,支持在地图上绘制轨迹线、点、热力图等,帮助用户更直观地理解移动物体的行为。
  • 数据格式兼容:MovingPandas支持多种常见的地理数据格式,包括GeoJSON、Shapefile等,方便用户加载和保存轨迹数据。

MovingPandas官方仓库地址为:movingpandas。MovingPandas官方示例代码仓库地址为:movingpandas-examples。本文所有实验数据来自于:movingpandas-examples-data。

本文所有代码见:Python-Study-Notes

MovingPandas作者推荐在Python 3.8及以上环境下安装MovingPandas,并建议使用conda进行安装。可以使用以下命令来安装MovingPandas:

conda install -c conda-forge movingpandas

由于MovingPandas的依赖环境较为复杂,所以不推荐使用pip进行安装。如果坚持使用pip进行安装,可以按照以下命令来安装MovingPandas:

pip install movingpandas
# 本文必安装第三方库
pip install contextily
# 以下第三方库可选
pip install hvplot
pip install cartopy
pip install geoviews

下面的代码展示了MovingPandas的版本信息,本文所用Python版本为Python3.10。

# jupyter notebook环境去除warning
import warnings
warnings.filterwarnings("ignore")# 查看movingpandas版本及依赖库版本信息
import movingpandas as mpd
mpd.show_versions()
MovingPandas 0.16.1SYSTEM INFO
-----------
python     : 3.10.10 (main, Mar 21 2023, 18:45:11) [GCC 11.2.0]
executable : /opt/conda/envs/python35-paddle120-env/bin/python
machine    : Linux-5.4.0-104-generic-x86_64-with-glibc2.31GEOS, GDAL, PROJ INFO
---------------------
GEOS       : None
GEOS lib   : None
GDAL       : 3.6.4
GDAL data dir: /opt/conda/envs/python35-paddle120-env/lib/python3.10/site-packages/fiona/gdal_data
PROJ       : 9.2.1
PROJ data dir: /opt/conda/envs/python35-paddle120-env/lib/python3.10/site-packages/pyproj/proj_dir/share/projPYTHON DEPENDENCIES
-------------------
geopandas  : 0.13.2
pandas     : 2.0.2
fiona      : 1.9.4.post1
numpy      : 1.24.3
shapely    : 2.0.1
rtree      : 1.0.1
pyproj     : 3.6.0
matplotlib : 3.7.1
mapclassify: None
geopy      : 2.3.0
holoviews  : None
hvplot     : None
geoviews   : None
stonesoup  : None

文章目录

  • 1 使用入门
    • 1.1 基础结构
      • 1.1.1 Trajectory类
      • 1.1.2 TrajectoryCollection类
    • 1.2 基础函数
      • 1.2.1 移动信息添加函数
      • 1.2.2 轨迹位置提取函数
      • 1.2.3 与GeoPandas交互
    • 1.3 轨迹处理
      • 1.3.1 轨迹提取
      • 1.3.2 轨迹压缩与平滑
      • 1.3.3 轨迹停留检测
      • 1.3.4 轨迹聚类
      • 1.3.5 轨迹距离计算
  • 2 轨迹绘图实例之海鸥迁徙轨迹分析
  • 3 参考

1 使用入门

1.1 基础结构

# 加载第三方库
import pandas as pd
import geopandas as gpd
import movingpandas as mpd
import shapely as shp
# hvPlot是建立在数据可视化库holoviews的一个高级封装库
# import hvplot.pandas from geopandas import GeoDataFrame, read_file
from shapely.geometry import Point, LineString, Polygon
from datetime import datetime, timedelta
# hvplot通常与HoloViews集成使用
# from holoviews import opts# 设置绘图界面
# opts.defaults(opts.Overlay(active_tools=['wheel_zoom'], frame_width=500, frame_height=400))

1.1.1 Trajectory类

MovingPandas提供Trajectory对象用于表示单个轨迹的对象。它由一系列坐标点组成。Trajectory对象初始化参数如下:

class Trajectory:def __init__(self,df,traj_id,obj_id=None,t=None,x=None,y=None,crs="epsg:4326",parent=None,):

其中各个参数解释如下:

  • df:具有GeoPandas的geometry坐标列和时间戳索引的GeoDataFrame,或Pandas的DataFrame。必填参数。
  • traj_id:任意类型,表示轨迹的唯一标识符。必填参数。
  • obj_id:任意类型,表示移动物体的唯一标识符。默认为 None。
  • t:表示包含时间戳的列名,默认为 None。
  • x:表示包含x坐标的列名,使用Pandas的DataFrame需指定。默认为 None。
  • y:表示包含y坐标的列名,使用Pandas的DataFrame需指定。默认为 None。
  • crs:表示 x/y 坐标的坐标参考系统。默认为 “epsg:4326”,即 WGS84。
  • parent:一个Trajectory 对象,表示父轨迹。默认为 None。

由Trajectory轨迹对象的初始化参数可知,MovingPandas用相关点位置信息和其时间信息来表示轨迹。下面用一个实例代码来展示Trajectory对象的使用。

首先创建一个用于表示位置和其对应时间信息的GeoDataFrame对象。

import pandas as pd
from shapely.geometry import Point
from geopandas import GeoDataFrame
import movingpandas as mpd
from datetime import datetime# 创建DataFrame对象
df = pd.DataFrame([{'geometry': Point(0, 0), 't': datetime(2023, 7, 1, 12, 0, 0)},{'geometry': Point(6, 0), 't': datetime(2023, 7, 1, 12, 6, 0)},{'geometry': Point(6, 6), 't': datetime(2023, 7, 1, 12, 10, 0)},{'geometry': Point(9, 9), 't': datetime(2023, 7, 1, 12, 15, 0)}])
# 设置索引轴
df = df.set_index('t')# 创建GeoDataFrame对象
# EPSG:3857是在线地图常用的投影坐标系,以米为单位
gdf = GeoDataFrame(df, crs='EPSG:3857')
gdf
geometry
t
2023-07-01 12:00:00POINT (0.000 0.000)
2023-07-01 12:06:00POINT (6.000 0.000)
2023-07-01 12:10:00POINT (6.000 6.000)
2023-07-01 12:15:00POINT (9.000 9.000)

然后基于GeoDataFrame对象创建轨迹对象。

# 创建Trajectory对象
# 1表示轨迹标识符
toy_traj = mpd.Trajectory(gdf, 1)
toy_traj
Trajectory 1 (2023-07-01 12:00:00 to 2023-07-01 12:15:00) | Size: 4 | Length: 16.2m
Bounds: (0.0, 0.0, 9.0, 9.0)
LINESTRING (0 0, 6 0, 6 6, 9 9)

接下来可以可视化轨迹路径,横纵轴表示点的xy坐标。

# 可视化路径
toy_traj.plot()
<Axes: >

png

# 交互式展示
# toy_traj.hvplot()

此外也可以直接调用toy_traj的GeoPandas对象进行可视化展示。

toy_traj.df.plot()
<Axes: >

png

1.1.2 TrajectoryCollection类

TrajectoryCollection类是一个集合,用于存储多个Trajectory 对象。它提供了对整个轨迹集合进行操作和分析的功能。可以使用 TrajectoryCollection来合并多个轨迹为一个对象、筛选轨迹、进行空间查询、可视化等。TrajectoryCollection对象初始化参数如下:

class TrajectoryCollection:def __init__(self,data,traj_id_col=None,obj_id_col=None,t=None,x=None,y=None,crs="epsg:4326",min_length=0,min_duration=None,)

其中大部分参数与Trajectory对象的参数一致,不同参数解释如下:

  • min_length: 轨迹的最小长度,低于该长度的轨迹将被丢弃。
  • min_duration: 轨迹的期望最小持续时间,低于该续时间的轨迹将被丢弃。

直接创建TrajectoryCollection对象

# 创建两条轨迹
df = pd.DataFrame({"t": pd.date_range("2023-07-01", periods=5, freq="min"),"trajectory_id": [1, 1, 2, 2, 2],"geometry": [Point(0, 0), Point(0, 1), Point(1, 2), Point(1, 3), Point(2, 4)],}
)
gdf = gpd.GeoDataFrame(df, crs='EPSG:3857')
gdf
ttrajectory_idgeometry
02023-07-01 00:00:001POINT (0.000 0.000)
12023-07-01 00:01:001POINT (0.000 1.000)
22023-07-01 00:02:002POINT (1.000 2.000)
32023-07-01 00:03:002POINT (1.000 3.000)
42023-07-01 00:04:002POINT (2.000 4.000)
# 创建TrajectoryCollection对象,用trajectory_id作为单个轨迹的唯一标识
tc = mpd.TrajectoryCollection(gdf, traj_id_col='trajectory_id', t='t')
tc
TrajectoryCollection with 2 trajectories
# 绘制轨迹, column指定轨迹对象,legend展示轨迹持续时间
import matplotlib.cm as cm
tc.plot(column='trajectory_id', legend=True)
<Axes: >

png

# 交互式展示轨迹
# tc.hvplot()

此外我们还可以从TrajectoryCollection提取单个轨迹,或者筛选所需要的轨迹。

# 从TrajectoryCollection提取单个轨迹
my_traj = tc.trajectories[1]
print(my_traj)
Trajectory 2 (2023-07-01 00:02:00 to 2023-07-01 00:04:00) | Size: 3 | Length: 2.4m
Bounds: (1.0, 2.0, 2.0, 4.0)
LINESTRING (1 2, 1 3, 2 4)
# 筛选路径持续时间大于一分钟的
min_duration = timedelta(minutes=1)
tc.trajectories = [traj for traj in tc if traj.get_duration() > min_duration]
print(tc.trajectories)
[Trajectory 2 (2023-07-01 00:02:00 to 2023-07-01 00:04:00) | Size: 3 | Length: 2.4m
Bounds: (1.0, 2.0, 2.0, 4.0)
LINESTRING (1 2, 1 3, 2 4)]

从文件数据创建TrajectoryCollection对象

TrajectoryCollection可以通过GeoPandas函数从各种地理空间数据存储中读取数据来创建,也可以利用Pandas从各类文件类型中读取数据进行创建。

下面展示了直接通过csv文件创建TrajectoryCollection对象,数据下载地址为:movingpandas-examples-data

# 读取数据
df = pd.read_csv('data/geolife_small.csv', delimiter=';')
df.head()
XYfididsequencetrajectory_idtrackert
0116.39130539.8985731111192008-12-11 04:42:14+00
1116.39131739.8986172221192008-12-11 04:42:16+00
2116.39092839.8986133331192008-12-11 04:43:26+00
3116.39083339.8986354441192008-12-11 04:43:32+00
4116.38941039.8987235551192008-12-11 04:43:47+00
# 转换为traj_collection对象
traj_collection = mpd.TrajectoryCollection(df, 'trajectory_id', t='t', x='X', y='Y')
print(traj_collection)
TrajectoryCollection with 5 trajectories
# 绘制轨迹
traj_collection.plot(column='trajectory_id', legend=True, figsize=(9,5))
<Axes: >

png

1.2 基础函数

MovingPandas中的函数可以同时适用于单个轨迹(Trajectory)对象和轨迹集合(TrajectoryCollection)对象,以实现对单条和多条轨迹的分析。

1.2.1 移动信息添加函数

MovingPandas提供add_speed、add_distance、add_direction、add_timedelta、add_timedelta和add_acceleration等移动信息添加函数来处理移动数据。这些函数的接口如下:

add_speed(overwrite=False, name=SPEED_COL_NAME, units=UNITS())
add_distance(overwrite=False, name=DISTANCE_COL_NAME, units=None)
add_direction(overwrite=False, name=DIRECTION_COL_NAME)
add_angular_difference(overwrite=False,name=ANGULAR_DIFFERENCE_COL_NAM):
add_timedelta(overwrite=False, name=TIMEDELTA_COL_NAME)
add_acceleration(overwrite=False, name=ACCELERATION_COL_NAME, units=UNITS())

所有函数均使用系统默认值或全局变量,具体函数介绍如下:

  • add_speed: 用于计算移动对象的速度。根据两个时间点之间的距离和时间差来计算速度。参数overwrite表示是否覆盖现有的速度列,name表示速度列的名称,units表示速度的单位,名字和单位不指定则使用默认值。
  • add_distance: 用于计算移动对象的距离。根据两个时间点之间的直线距离异来计算距离。参数overwrite表示是否覆盖现有的距离列,name表示距离列的名称,units表示距离的单位。
  • add_direction: 用于计算移动对象的方向。根据两个时间点之间的坐标差异来计算方向。参数overwrite表示是否覆盖现有的方向列,name表示方向列的名称。方向值以度为单位,从北开始顺时针旋转,范围值为[0,360)。
  • add_angular_difference: 用于计算两个对象之间的角度差。根据移动对象的方向信息来计算两个对象之间的角度差。参数overwrite表示是否覆盖现有的方向列,name表示方向列的名称。角度范围值为[0,180]。
  • add_timedelta: 用于计算移动对象的时间差。根据时间戳之间的差异来计算时间差。参数overwrite表示是否覆盖现有的时间差列,name表示时间差列的名称。时间差为当前行与前一行之差。
  • add_acceleration: 用于计算移动对象的加速度。根据速度和时间差来计算加速度。参数overwrite表示是否覆盖现有的加速度列,name表示加速度列的名称,units表示加速度的单位。

以下代码详细介绍了这些函数的使用。

# 创建轨迹
df = pd.DataFrame([{'geometry':Point(0,0), 't':datetime(2023,7,1,12,0,0)},{'geometry':Point(6,0), 't':datetime(2023,7,1,12,0,6)},{'geometry':Point(6,6), 't':datetime(2023,7,1,12,0,11)},{'geometry':Point(9,9), 't':datetime(2023,7,1,12,0,14)}
]).set_index('t')
gdf = GeoDataFrame(df, crs='EPSG:3857')
toy_traj = mpd.Trajectory(gdf, 1)
toy_traj
Trajectory 1 (2023-07-01 12:00:00 to 2023-07-01 12:00:14) | Size: 4 | Length: 16.2m
Bounds: (0.0, 0.0, 9.0, 9.0)
LINESTRING (0 0, 6 0, 6 6, 9 9)
# 查看轨迹的坐标系,以米为单位
toy_traj.crs
<Projected CRS: EPSG:3857>
Name: WGS 84 / Pseudo-Mercator
Axis Info [cartesian]:
- X[east]: Easting (metre)
- Y[north]: Northing (metre)
Area of Use:
- name: World between 85.06°S and 85.06°N.
- bounds: (-180.0, -85.06, 180.0, 85.06)
Coordinate Operation:
- name: Popular Visualisation Pseudo-Mercator
- method: Popular Visualisation Pseudo Mercator
Datum: World Geodetic System 1984 ensemble
- Ellipsoid: WGS 84
- Prime Meridian: Greenwich
# 添加速度列,不指定单位则以米为单位
toy_traj.add_speed(overwrite=True)
toy_traj.df
geometryspeed
t
2023-07-01 12:00:00POINT (0.000 0.000)1.000000
2023-07-01 12:00:06POINT (6.000 0.000)1.000000
2023-07-01 12:00:11POINT (6.000 6.000)1.200000
2023-07-01 12:00:14POINT (9.000 9.000)1.414214
# 添加自定义定单位的速度列
toy_traj.add_speed(overwrite=True, name="speed (ft/min)", units=("ft", "min"))
toy_traj.df
geometryspeedspeed (ft/min)
t
2023-07-01 12:00:00POINT (0.000 0.000)1.000000196.850394
2023-07-01 12:00:06POINT (6.000 0.000)1.000000196.850394
2023-07-01 12:00:11POINT (6.000 6.000)1.200000236.220472
2023-07-01 12:00:14POINT (9.000 9.000)1.414214278.388497
# 其他列的添加
toy_traj.add_distance(overwrite=True, name="distance (m)", units="m")
toy_traj.add_direction(overwrite=True)
toy_traj.add_timedelta(overwrite=True)
toy_traj.add_angular_difference(overwrite=True)
toy_traj.add_acceleration(overwrite=True, name="acceleration (mph/s)", units=("mi", "h", "s"))
toy_traj.df
geometryspeedspeed (ft/min)distance (m)directiontimedeltaangular_differenceacceleration (mph/s)
t
2023-07-01 12:00:00POINT (0.000 0.000)1.000000196.8503940.00000090.0NaT0.00.000000
2023-07-01 12:00:06POINT (6.000 0.000)1.000000196.8503946.00000090.00 days 00:00:060.00.000000
2023-07-01 12:00:11POINT (6.000 6.000)1.200000236.2204726.0000000.00 days 00:00:0590.00.089477
2023-07-01 12:00:14POINT (9.000 9.000)1.414214278.3884974.24264145.00 days 00:00:0345.00.159727
# 可视化某列结果绘制结果
toy_traj.plot(column='speed', linewidth=5, capstyle='round', figsize=(9,3), legend=True, vmax=20)
<Axes: >

png

1.2.2 轨迹位置提取函数

MovingPandas提供add_speed、add_distance、add_direction、add_timedelta、add_timedelta和add_acceleration等轨迹位置提取函数用于提取轨迹位置相关信息数据。具体函数介绍如下:

  • get_start_location(): 该函数用于获取移动对象的起始位置。它返回移动对象的起始位置坐标。
  • get_end_location(): 该函数用于获取移动对象的结束位置。它返回移动对象的结束位置坐标。
  • get_position_at(t, method="interpolated"): 该函数用于在给定的时间点(t)处获取移动对象的位置。参数t表示所需位置的时间点,而可选的参数method指定了用于获取位置的方法,默认为"interpolated",表示使用插值方法获取位置。
  • get_segment_between(t1, t2): 这个函数用于获取两个给定时间点(t1和t2)之间的移动对象的片段。它返回包含这段时间内移动对象的子集。
  • clip(polygon, point_based=False): 该函数用于剪裁移动对象,使其仅保留在指定多边形区域内的部分。参数polygon表示用于剪裁的多边形区域,而可选的参数point_based指定剪裁是否基于点而不是线段。默认情况下,剪裁是基于线段的。此外,MovingPandas也提供了intersection函数来算轨迹数据与空间几何图形之间的交集,效果类似clip函数,具体使用可以参考官方文档。

以下代码详细介绍了这些函数的使用。

# 创建轨迹
df = pd.DataFrame([{'geometry':Point(0,0), 't':datetime(2023,7,1,12,0,0)},{'geometry':Point(6,0), 't':datetime(2023,7,1,12,6,6)},{'geometry':Point(6,6), 't':datetime(2023,7,1,12,10,11)},{'geometry':Point(9,9), 't':datetime(2023,7,1,12,19,14)}
]).set_index('t')
gdf = GeoDataFrame(df, crs='EPSG:3857')
toy_traj = mpd.Trajectory(gdf, 1)
toy_traj
Trajectory 1 (2023-07-01 12:00:00 to 2023-07-01 12:19:14) | Size: 4 | Length: 16.2m
Bounds: (0.0, 0.0, 9.0, 9.0)
LINESTRING (0 0, 6 0, 6 6, 9 9)

起始位置和结束位置

ax = toy_traj.plot()
# 提取轨迹起始点和结束点
gpd.GeoSeries(toy_traj.get_start_location()).plot(ax=ax, color='blue')
gpd.GeoSeries(toy_traj.get_end_location()).plot(ax=ax, color='red')
<Axes: >

png

中间位置

# 提取轨迹中的一个时间点的位置
t = datetime(2023,7,1,12,6,0)
# 改时间点不同位置信息获取方式
print(toy_traj.get_position_at(t, method="nearest"))
print(toy_traj.get_position_at(t, method="interpolated"))
print(toy_traj.get_position_at(t, method="ffill"))
print(toy_traj.get_position_at(t, method="bfill"))
POINT (6 0)
POINT (5.9016393442622945 0)
POINT (0 0)
POINT (6 0)
point = toy_traj.get_position_at(t, method="interpolated")
ax = toy_traj.plot()
gpd.GeoSeries(point).plot(ax=ax, color='red', markersize=100)
<Axes: >

png

轨迹片段

# 提取轨迹中的一个片段
segment = toy_traj.get_segment_between(datetime(2023,7,1,12,6,0), datetime(2023,7,1,12,12,0))
print(segment)
Trajectory 1_2023-07-01 12:06:00 (2023-07-01 12:06:06 to 2023-07-01 12:10:11) | Size: 2 | Length: 6.0m
Bounds: (6.0, 0.0, 6.0, 6.0)
LINESTRING (6 0, 6 6)
ax = toy_traj.plot()
segment.plot(ax=ax, color='red', linewidth=5)
<Axes: >

png

轨迹区域

xmin, xmax, ymin, ymax = 2, 8, -10, 5
# 设置区域
polygon = Polygon([(xmin, ymin), (xmin, ymax), (xmax, ymax), (xmax, ymin), (xmin, ymin)])
# 提取位于该区域的轨迹
intersections = toy_traj.clip(polygon)
intersections
TrajectoryCollection with 1 trajectories
# 绘制轨迹区域
ax = toy_traj.plot()
gpd.GeoSeries(polygon).plot(ax=ax, color='lightgray')
intersections.plot(ax=ax, color='red', linewidth=5, capstyle='round')
<Axes: >

png

# 单独绘制轨迹区域
intersections = toy_traj.clip(polygon)
intersections.plot()
<Axes: >

png

1.2.3 与GeoPandas交互

MovingPandas提供to_point_gdf,to_line_gdf和to_traj_gdf函数以按照不同规则将轨迹类Trajectories和TrajectoryCollections转换回为GeoPandas的GeoDataFrames结构。这些函数解释如下:

  • to_point_gdf(return_orig_tz=False),将轨迹对象转换为与每个点相关的GeoDataFrame。return_orig_tz用于指定是否返回与原始时区匹配的时间戳。默认为False,时间戳会根据本地时区进行转换。
  • to_line_gdf(),将轨迹对象转换为每条线相关的GeoDataFrame。
  • to_traj_gdf(wkt=False, agg=False),返回一个GeoDataFrame,其每一行表示一条轨迹。wkt用于指定是否将轨迹几何表示为Well-Known Text (WKT)格式。默认为False,返回几何对象。agg用于指定是否对轨迹进行聚合。默认为False,不进行聚合,agg可选参数见代码。
# 读取数据
df = pd.read_csv('data/geolife_small.csv', delimiter=';')
df.head()
# 转换为traj_collection对象
traj_collection = mpd.TrajectoryCollection(df, 'trajectory_id', t='t', x='X', y='Y')
print(traj_collection)
TrajectoryCollection with 5 trajectories
# 转换为点GeoDataFrame
traj_collection.to_point_gdf().head()
fididsequencetrajectory_idtrackergeometry
t
2008-12-11 04:42:14111119POINT (116.39131 39.89857)
2008-12-11 04:42:16222119POINT (116.39132 39.89862)
2008-12-11 04:43:26333119POINT (116.39093 39.89861)
2008-12-11 04:43:32444119POINT (116.39083 39.89863)
2008-12-11 04:43:47555119POINT (116.38941 39.89872)
# 转换为线GeoDataFrame
traj_collection.to_line_gdf().head()
fididsequencetrajectory_idtrackertprev_tgeometry
02221192008-12-11 04:42:162008-12-11 04:42:14LINESTRING (116.39131 39.89857, 116.39132 39.8...
13331192008-12-11 04:43:262008-12-11 04:42:16LINESTRING (116.39132 39.89862, 116.39093 39.8...
24441192008-12-11 04:43:322008-12-11 04:43:26LINESTRING (116.39093 39.89861, 116.39083 39.8...
35551192008-12-11 04:43:472008-12-11 04:43:32LINESTRING (116.39083 39.89863, 116.38941 39.8...
46661192008-12-11 04:43:502008-12-11 04:43:47LINESTRING (116.38941 39.89872, 116.39052 39.8...
# 每条轨迹单独聚合GeoDataFrame
traj_collection.to_traj_gdf()
traj_idstart_tend_tgeometrylengthdirection
012008-12-11 04:42:142008-12-11 05:15:46LINESTRING (116.39131 39.89857, 116.39132 39.8...6207.020261186.681376
122009-06-29 07:02:252009-06-29 11:13:12LINESTRING (116.59096 40.07196, 116.59091 40.0...38764.575483250.585295
232009-02-04 04:32:532009-02-04 11:20:12LINESTRING (116.38569 39.89977, 116.38565 39.8...12745.157506304.115160
342009-03-10 10:36:452009-03-10 12:01:07LINESTRING (116.38805 39.90342, 116.38804 39.9...14363.780551300.732843
452009-02-25 09:47:032009-02-25 14:31:24LINESTRING (116.38526 39.90027, 116.38525 39.9...39259.779560305.200501

1.3 轨迹处理

1.3.1 轨迹提取

MovingPandas提供了TemporalSplitter、ObservationGapSplitter,StopSplitter,SpeedSplitter类来根据不同规则从轨迹中提取指定轨迹片段。具体如下:

  • TemporalSplitter:使用规则的时间间隔将轨迹拆分为子区间。
  • ObservationGapSplitter:根据观测时间之间的间隔将轨迹拆分为子区间。如果两个连续观测之间的间隔超过指定的阈值,则认为是该轨迹需要拆分。
  • StopSplitter:根据停留点的定义将轨迹拆分为子区间。停留点是指轨迹在相对较小的区域内停留了一段时间。
  • SpeedSplitter:据速度阈值将轨迹拆分为子区间。如果轨迹在某段速度低于指定阈值,则需要在此拆分轨迹。

这些类都只需待处理轨迹进行初始化,然后调用split函数进行轨迹提取。具体使用见下列代码。

import matplotlib.pyplot as plt# 读取数据
gdf = read_file('data/geolife_small.gpkg').to_crs('EPSG:3857')
tc = mpd.TrajectoryCollection(gdf, 'trajectory_id', t='t')
# 共五条轨迹
tc.to_traj_gdf()
traj_idstart_tend_tgeometrylengthdirection
012008-12-11 04:42:142008-12-11 05:15:46LINESTRING (12956620.805 4851214.113, 12956622...8101.428690186.679744
122009-06-29 07:02:252009-06-29 11:13:12LINESTRING (12978845.964 4876404.973, 12978840...50621.731208250.500509
232009-02-04 04:32:532009-02-04 11:20:12LINESTRING (12955995.635 4851388.236, 12955991...16626.383723304.099365
342009-03-10 10:36:452009-03-10 12:01:07LINESTRING (12956258.794 4851917.156, 12956257...18739.337154300.716597
452009-02-25 09:47:032009-02-25 14:31:24LINESTRING (12955947.434 4851460.354, 12955946...51327.869585305.185128

TemporalSplitter

# 将轨迹切分为以每小时为单位的轨迹片段,min_length表示删除轨迹长度小于该值的轨迹片段,min_lwngth取决于轨迹的单位
split = mpd.TemporalSplitter(tc).split(mode="hour", min_length=20000)
split.to_traj_gdf()
traj_idstart_tend_tgeometrylengthdirection
02_2009-06-29 07:00:002009-06-29 07:02:252009-06-29 07:59:55LINESTRING (12978845.964 4876404.973, 12978840...40106.105641245.123482
15_2009-02-25 10:00:002009-02-25 10:00:042009-02-25 10:54:47LINESTRING (12955251.019 4851210.485, 12955248...25455.673254337.129155
# 绘制轨迹
fig, axes = plt.subplots(nrows=1, ncols=len(split), figsize=(6,4))
for i, traj in enumerate(split):traj.plot(ax=axes[i], linewidth=3.0, capstyle='round', column='speed', vmax=20)

png

ObservationGapSplitter

# 如果两个连续观测超过间隔gap,如30分钟,则认为该轨迹需要拆分
split = mpd.ObservationGapSplitter(tc).split(gap=timedelta(minutes=30),min_length=0)
split.to_traj_gdf()
traj_idstart_tend_tgeometrylengthdirection
01_02008-12-11 04:42:142008-12-11 05:15:46LINESTRING (12956620.805 4851214.113, 12956622...8101.428690186.679744
12_02009-06-29 07:02:252009-06-29 08:20:15LINESTRING (12978845.964 4876404.973, 12978840...47469.321958252.952783
22_12009-06-29 10:57:172009-06-29 11:13:12LINESTRING (12948649.439 4867034.108, 12948650...3040.348707139.615947
33_02009-02-04 04:32:532009-02-04 04:35:03LINESTRING (12955995.635 4851388.236, 12955991...311.72923142.937430
43_12009-02-04 10:03:212009-02-04 11:20:12LINESTRING (12956011.999 4851497.646, 12956043...16228.264596303.229612
54_02009-03-10 10:36:452009-03-10 12:01:07LINESTRING (12956258.794 4851917.156, 12956257...18739.337154300.716597
65_02009-02-25 09:47:032009-02-25 10:54:47LINESTRING (12955947.434 4851460.354, 12955946...27149.500896335.377580
75_12009-02-25 13:30:222009-02-25 14:31:24LINESTRING (12945965.972 4873487.115, 12945952...24074.321167165.727187

StopSplitter

# 如果轨迹在某一点为圆心,直径为10范围内停留60s,则认为轨迹需要在该段分割
split = mpd.StopSplitter(tc).split(max_diameter=10, min_duration=timedelta(seconds=60))
split.to_traj_gdf()
traj_idstart_tend_tgeometrylengthdirection
01_2008-12-11 04:42:142008-12-11 04:42:142008-12-11 05:15:46LINESTRING (12956620.805 4851214.113, 12956622...8101.428690186.679744
12_2009-06-29 07:02:252009-06-29 07:02:252009-06-29 07:05:30LINESTRING (12978845.964 4876404.973, 12978840...608.23301629.527683
22_2009-06-29 07:06:552009-06-29 07:06:552009-06-29 08:02:40LINESTRING (12979026.970 4876730.251, 12979022...41655.491556246.215181
32_2009-06-29 08:03:402009-06-29 08:03:402009-06-29 11:13:12LINESTRING (12949605.674 4863764.794, 12949579...8333.942283357.660458
43_2009-02-04 04:32:532009-02-04 04:32:532009-02-04 11:20:12LINESTRING (12955995.635 4851388.236, 12955991...16626.383723304.099365
54_2009-03-10 10:36:452009-03-10 10:36:452009-03-10 12:01:07LINESTRING (12956258.794 4851917.156, 12956257...18739.337154300.716597
65_2009-02-25 09:47:032009-02-25 09:47:032009-02-25 14:31:24LINESTRING (12955947.434 4851460.354, 12955946...51327.869585305.185128

SpeedSplitter

# 把超过30分钟速度低于1地理系统单位/秒的轨迹分开
split = mpd.SpeedSplitter(tc).split(speed=1, duration=timedelta(minutes=30))
split.to_traj_gdf()
traj_idstart_tend_tgeometrylengthdirection
01_02008-12-11 04:42:142008-12-11 05:15:46LINESTRING (12956620.805 4851214.113, 12956622...8048.160604186.679744
12_02009-06-29 07:02:252009-06-29 08:20:15LINESTRING (12978845.964 4876404.973, 12978840...47336.977010252.952783
22_12009-06-29 10:57:222009-06-29 11:10:07LINESTRING (12948650.441 4867044.718, 12948642...2915.988294138.780873
33_02009-02-04 04:32:532009-02-04 04:35:03LINESTRING (12955995.635 4851388.236, 12955991...310.44078042.937430
43_12009-02-04 10:03:232009-02-04 11:20:12LINESTRING (12956043.836 4851524.490, 12956025...15962.930350302.882421
54_02009-03-10 10:36:452009-03-10 12:01:07LINESTRING (12956258.794 4851917.156, 12956257...18349.431950300.716597
65_02009-02-25 09:47:032009-02-25 10:54:47LINESTRING (12955947.434 4851460.354, 12955946...27081.554127335.377580
75_12009-02-25 13:30:232009-02-25 14:31:24LINESTRING (12945952.057 4873516.928, 12945956...24006.683028165.708568

1.3.2 轨迹压缩与平滑

MovingPandas提供了各种轨迹压缩类和轨迹平滑类,以减少轨迹对象的大小(点的数量),从而加快处理速度。

压缩轨迹类都只需待处理轨迹进行初始化,然后调用generalize函数进行轨迹压缩。用于压缩轨迹的类有:

  • MinDistanceGeneralizer(最小距离):通过删除原始数据中的一些点来压缩轨迹,被删除的点与相邻点之间的距离必须小于指定的最小距离。
  • DouglasPeuckerGeneralizer(道格拉斯-普克):道格拉斯-普克算法根据指定的阈值,逐渐删除原始数据中的部分点,从而生成一个近似的简化线串或轨迹。
  • MaxDistanceGeneralizer(最大距离):删除原始数据中与相邻点之间距离超过指定最大距离的点,从而实现轨迹压缩。
  • MinTimeDeltaGeneralizer(最小时间间隔):通过删除两个连续时间戳之间时间间隔小于指定最小时间间隔的数据点来实现轨迹压缩。
  • TopDownTimeRatioGeneralizer(自顶向下时间比率):根据预先设定的时间比率,逐渐删除原始数据中的数据点。

平滑轨迹类目前只提供了KalmanSmootherCV类进行轨迹平滑,KalmanSmootherCV需要额外安装第三库且处理流程麻烦,所以一般都是效果相近的压缩轨迹类。

# 读取数据
gdf = read_file('data/geolife_small.gpkg').to_crs('EPSG:3857')
tc = mpd.TrajectoryCollection(gdf, 'trajectory_id', t='t')
# 提取单条轨迹进行操作
original_traj = tc.trajectories[1]
# 可以看到轨迹包括897个点
print(original_traj)
Trajectory 2 (2009-06-29 07:02:25 to 2009-06-29 11:13:12) | Size: 897 | Length: 50621.7m
Bounds: (12948595.449314836, 4861831.088215791, 12979030.643375682, 4877940.244020148)
LINESTRING (12978845.964340456 4876404.972802613, 12978840.175726935 4876411.664456934, 12978837.281
# tolerance设置连续点之间的最小距离
mpd.MinDistanceGeneralizer(original_traj).generalize(tolerance=1.0)
Trajectory 2 (2009-06-29 07:02:25 to 2009-06-29 11:13:12) | Size: 818 | Length: 50611.6m
Bounds: (12948595.449314836, 4861831.088215791, 12979030.643375682, 4877940.244020148)
LINESTRING (12978845.964340456 4876404.972802613, 12978840.175726935 4876411.664456934, 12978837.281
# tolerance表示距离公差,具体使用见算法介绍
mpd.DouglasPeuckerGeneralizer(original_traj).generalize(tolerance=1.0)
Trajectory 2 (2009-06-29 07:02:25 to 2009-06-29 11:13:12) | Size: 564 | Length: 50606.8m
Bounds: (12948595.449314836, 4861831.088215791, 12979030.643375682, 4877940.244020148)
LINESTRING (12978845.964340456 4876404.972802613, 12978837.281420175 4876414.573873267, 12978852.086
# tolerance设置连续点之间的最大距离
mpd.MaxDistanceGeneralizer(original_traj).generalize(tolerance=1.0)
Trajectory 2 (2009-06-29 07:02:25 to 2009-06-29 11:13:12) | Size: 324 | Length: 50166.8m
Bounds: (12948595.449314836, 4861831.088215791, 12979029.752819754, 4877940.244020148)
LINESTRING (12978845.964340456 4876404.972802613, 12978837.281420175 4876414.573873267, 12978852.086
# tolerance连续轨迹最小的时间距离
mpd.MinTimeDeltaGeneralizer(original_traj).generalize(tolerance=timedelta(minutes=1))
Trajectory 2 (2009-06-29 07:02:25 to 2009-06-29 11:13:12) | Size: 76 | Length: 47565.6m
Bounds: (12948636.414887449, 4862053.18880025, 12979029.30754179, 4877912.7458354365)
LINESTRING (12978845.964340456 4876404.972802613, 12978815.79675845 4876446.868451926, 12978780.3971
# tolerance距离公差,该类处理速度很慢
mpd.TopDownTimeRatioGeneralizer(original_traj).generalize(tolerance=1.0)
Trajectory 2 (2009-06-29 07:02:25 to 2009-06-29 11:13:12) | Size: 756 | Length: 50616.4m
Bounds: (12948595.449314836, 4861831.088215791, 12979030.643375682, 4877940.244020148)
LINESTRING (12978845.964340456 4876404.972802613, 12978840.175726935 4876411.664456934, 12978837.281

1.3.3 轨迹停留检测

MovingPandas提供了轨迹停留检测的功能,事实上轨迹停留并不意味着轨迹运动速度为零,事实上由于跟踪不准确,物体移动速度很少达到0,例如GPS轨迹往往会在物体的停止位置周围不断移动。所以MovingPandas检测轨迹停留的方式为如果物体停留在指定大小的区域内超过一段时间,则判定为轨迹停留。具体使用如下:

# 读取数据
gdf = read_file('data/geolife_small.gpkg').to_crs('EPSG:3857')
tc = mpd.TrajectoryCollection(gdf, 'trajectory_id', t='t')
# 使用单个轨迹进行运算,也可以轨迹集合
my_traj = tc.trajectories[0]
my_traj
Trajectory 1 (2008-12-11 04:42:14 to 2008-12-11 05:15:46) | Size: 466 | Length: 8101.4m
Bounds: (12955985.950308602, 4845963.532957837, 12956871.051579898, 4851235.877839145)
LINESTRING (12956620.805364596 4851214.112520242, 12956622.141198486 4851220.49700885, 12956578.8379
# 创建停留检测器
detector = mpd.TrajectoryStopDetector(tc)
# 获取在某个直径为10的区域停留60s以上的轨迹时间范围
stop_time_ranges = detector.get_stop_time_ranges(min_duration=timedelta(seconds=60), max_diameter=10)
for x in stop_time_ranges: print(x)
Traj 2: 2009-06-29 07:05:30 - 2009-06-29 07:06:55 (duration: 0 days 00:01:25)
Traj 2: 2009-06-29 08:02:40 - 2009-06-29 08:03:40 (duration: 0 days 00:01:00)
# 获取在某个直径为10的区域停留60s以上的轨迹点
stop_points = detector.get_stop_points(min_duration=timedelta(seconds=60), max_diameter=10)
stop_points
geometrystart_timeend_timetraj_idduration_s
stop_id
2_2009-06-29 07:05:30POINT (12979027.415 4876729.960)2009-06-29 07:05:302009-06-29 07:06:55285.0
2_2009-06-29 08:02:40POINT (12949605.785 4863764.939)2009-06-29 08:02:402009-06-29 08:03:40260.0
# 获取在某个直径为10的区域停留60s以上的轨迹片段
stops = detector.get_stop_segments(min_duration=timedelta(seconds=60), max_diameter=10)
stops
TrajectoryCollection with 2 trajectories

1.3.4 轨迹聚类

Tmovingpandas提供TrajectoryCollectionAggregator类,用于将多个轨迹集合(TrajectoryCollection)合并为一个更大的轨迹集合。具体使用如下:

# 读取数据
gdf = read_file('data/geolife_small.gpkg').to_crs('EPSG:3857')
tc = mpd.TrajectoryCollection(gdf, 'trajectory_id', t='t')
# 简化轨迹
generalized = mpd.MinDistanceGeneralizer(tc).generalize(tolerance=100)
# 聚类
aggregator = mpd.TrajectoryCollectionAggregator(generalized, max_distance=1000, min_distance=100, min_stop_duration=timedelta(minutes=5))
# 从合并的轨迹集合中提取显著点,并返回一个Geopandas的GeoDataFrame对象
pts = aggregator.get_significant_points_gdf()
pts.head()
geometry
0POINT (12956620.805 4851214.113)
1POINT (12956054.412 4846377.879)
2POINT (12956409.855 4851235.878)
3POINT (12956533.642 4851120.522)
4POINT (12956436.794 4851148.817)
# 合并的轨迹集合视为一个整体,并将其划分为不同的簇
# 返回每个簇的质心结果
clusters = aggregator.get_clusters_gdf()
clusters.head()
geometryn
0POINT (12955779.665 4851264.376)30
1POINT (12956533.123 4846314.334)6
2POINT (12956511.072 4849943.391)4
3POINT (12956768.526 4848886.514)2
4POINT (12956668.895 4847819.306)2
# 将合并后的轨迹数据转换为一个GeoDataFrame对象,其中每一条记录代表着一次移动(即一段轨迹)
flows = aggregator.get_flows_gdf()
flows.head()
geometryweight
0LINESTRING (12955779.665 4851264.376, 12956511...1
1LINESTRING (12956511.072 4849943.391, 12956768...1
2LINESTRING (12956768.526 4848886.514, 12956668...1
3LINESTRING (12956668.895 4847819.306, 12956533...1
4LINESTRING (12978848.347 4876830.605, 12978353...1

1.3.5 轨迹距离计算

MovingPandas提供了distance函数和hausdorff_distance函数来计算两个几何对象的距离。distance函数基于shapely-objectdistance计算距离,而hausdorff_distance返回豪斯多夫距离。使用示例如下:

# 定义轨迹1
df = pd.DataFrame([{'geometry':Point(0,0), 't':datetime(2023,7,1,12,0,0)},{'geometry':Point(6,0), 't':datetime(2023,7,1,12,6,0)},{'geometry':Point(6,6), 't':datetime(2023,7,1,12,10,0)},{'geometry':Point(9,9), 't':datetime(2023,7,1,12,15,0)}
]).set_index('t')
# 单位为米
geo_df = GeoDataFrame(df, crs='EPSG:3857')
toy_traj = mpd.Trajectory(geo_df, 1)# 定义轨迹2
df = pd.DataFrame([{'geometry':Point(3,3), 't':datetime(2023,7,1,12,0,0)},{'geometry':Point(3,9), 't':datetime(2023,7,1,12,6,0)},{'geometry':Point(2,9), 't':datetime(2023,7,1,12,10,0)},{'geometry':Point(0,7), 't':datetime(2023,7,1,12,15,0)}
]).set_index('t')
geo_df = GeoDataFrame(df, crs='EPSG:3857')
toy_traj2 = mpd.Trajectory(geo_df, 1)
# 展示轨迹
ax = toy_traj.plot()
toy_traj2.plot(ax=ax, color='red')
<Axes: >

png

# 计算距离,默认单位由坐标系单位决定,该坐标系单位为米
print(toy_traj.distance(toy_traj2))
# 计算hausdorff距离,默认单位由坐标系单位决定。
print(toy_traj.hausdorff_distance(toy_traj2))
3.0
6.082762530298219
# 计算距离,units自定义单位为厘米
print(toy_traj.distance(toy_traj2, units="cm"))
# 计算hausdorff距离,units自定义单位为英寸
print(toy_traj.hausdorff_distance(toy_traj2, units="ft"))
300.0
19.95656998129337

2 轨迹绘图实例之海鸥迁徙轨迹分析

基于海鸥迁徙数据来分析其移动轨迹。数据下载地址:movingpandas-examples-data。

step1 加载数据

df = read_file('data/gulls.gpkg')
# 展示数据
df.head()
event-idvisibletimestamplocation-longlocation-latsensor-typeindividual-taxon-canonical-nametag-local-identifierindividual-local-identifierstudy-namegeometry
01082620685true2009-05-27 14:00:00.00024.5861761.24783gpsLarus fuscus9173291732ANavigation experiments in lesser black-backed ...POINT (24.58617 61.24783)
11082620686true2009-05-27 20:00:00.00024.5821761.23267gpsLarus fuscus9173291732ANavigation experiments in lesser black-backed ...POINT (24.58217 61.23267)
21082620687true2009-05-28 05:00:00.00024.5313361.18833gpsLarus fuscus9173291732ANavigation experiments in lesser black-backed ...POINT (24.53133 61.18833)
31082620688true2009-05-28 08:00:00.00024.5820061.23283gpsLarus fuscus9173291732ANavigation experiments in lesser black-backed ...POINT (24.58200 61.23283)
41082620689true2009-05-28 14:00:00.00024.5825061.23267gpsLarus fuscus9173291732ANavigation experiments in lesser black-backed ...POINT (24.58250 61.23267)
# 查看数据维度
df.shape
(89867, 11)
# 展示数据坐标点
df.plot()
<Axes: >

png

# 海鸥唯一id
df['individual-local-identifier'].unique()
array(['91732A', '91733A', '91734A', '91735A', '91737A', '91738A','91739A', '91740A', '91741A', '91742A', '91743A', '91744A','91745A', '91746A', '91747A', '91748A', '91749A', '91750A','91751A', '91752A', '91754A', '91755A', '91756A', '91758A','91759A', '91761A', '91762A', '91763A', '91764A', '91765A','91766A', '91767A', '91769A', '91771A', '91774A', '91775A','91776A', '91777A', '91778A', '91779A', '91780A', '91781A','91782A', '91783A', '91785A', '91786A', '91787A', '91788A','91789A', '91794A', '91795A', '91797A', '91798A', '91799A','91800A', '91802A', '91803A', '91807A', '91809A', '91810A','91811A', '91812A', '91813A', '91814A', '91815A', '91816A','91819A', '91821A', '91823A', '91824A', '91825A', '91826A','91827A', '91828A', '91829A', '91830A', '91831A', '91832A','91835A', '91836A', '91837A', '91838A', '91839A', '91843A','91845A', '91846A', '91848A', '91849A', '91852A', '91854A','91858A', '91861A', '91862A', '91864A', '91865A', '91866A','91870A', '91871A', '91872A', '91875A', '91876A', '91877A','91878A', '91880A', '91881A', '91884A', '91885A', '91893A','91894A', '91897A', '91900A', '91901A', '91903A', '91907A','91908A', '91910A', '91911A', '91913A', '91916A', '91918A','91919A', '91920A', '91921A', '91924A', '91929A', '91930A'],dtype=object)
# 海鸥数量
len(df['individual-local-identifier'].unique())
126
# 查看海鸥个体轨迹记录数,可以看到呈长尾分布
df['individual-local-identifier'].value_counts().plot(kind='bar', figsize=(16,4))
<Axes: xlabel='individual-local-identifier'>

png

# 创建轨迹集合,以海鸥编号为轨迹编号
tc = mpd.TrajectoryCollection(df, traj_id_col='individual-local-identifier', t='timestamp', min_length=100, crs='EPSG:3857')    
tc
TrajectoryCollection with 125 trajectories
# 压缩轨迹
tc = mpd.MinTimeDeltaGeneralizer(tc).generalize(tolerance=timedelta(days=1))

step2 个体分析

# 挑选个体编号91916A的海鸥轨迹,该海鸥轨迹记录数最多
filtered = tc.filter('individual-local-identifier', '91916A')
my_traj = filtered.trajectories[0].copy()
my_traj
Trajectory 91916A (2009-08-15 15:00:00 to 2015-08-27 09:00:00) | Size: 2134 | Length: 77841407.0m
Bounds: (7.915, 20.70533, 39.51317, 61.54817)
LINESTRING (7.915 54.18533, 9.44183 54.87233, 9.4425 54.87283, 9.41167 54.8555, 9.39583 54.88, 9.396
# 展示轨迹信息
# my_traj.hvplot(title='Movement', line_width=2) 
my_traj.plot()
<Axes: >

png

由上图可以看到该海鸥一直呈现季节性来回移动。所以按照年份,将该海鸥的轨迹信息进行拆分。

trips_by_year = mpd.TemporalSplitter(filtered).split(mode='year')
trips_by_year.to_traj_gdf()
traj_idstart_tend_tgeometrylengthdirection
091916A_2009-12-31 00:00:002009-08-15 15:00:002009-12-31 19:00:00LINESTRING (7.91500 54.18533, 9.44183 54.87233...5.352375e+06131.939002
191916A_2010-12-31 00:00:002010-01-01 19:00:002010-12-31 07:00:00LINESTRING (39.18417 21.17583, 39.18833 21.170...1.232428e+07281.047091
291916A_2011-12-31 00:00:002011-01-01 07:00:002011-12-31 04:00:00LINESTRING (39.17000 21.18017, 39.16883 21.156...1.441793e+07189.238229
391916A_2012-12-31 00:00:002012-01-01 04:00:002012-12-31 19:00:00LINESTRING (39.16933 21.16667, 39.17567 21.142...1.324811e+0762.132640
491916A_2013-12-31 00:00:002013-01-01 19:00:002013-12-31 13:00:00LINESTRING (39.18167 21.17333, 39.18100 21.173...1.293261e+07273.165851
591916A_2014-12-31 00:00:002014-01-01 13:00:002014-12-31 19:00:00LINESTRING (39.17083 21.15400, 39.17100 21.157...1.300973e+0733.742967
691916A_2015-12-31 00:00:002015-01-01 19:00:002015-08-27 09:00:00LINESTRING (39.18167 21.16967, 39.18233 21.169...6.551740e+06343.887905
# 按照年份显示轨迹
from matplotlib.cm import tab20
ax = None
for index,tmp in enumerate(trips_by_year):if ax is None:ax = tmp.plot(color=tab20(index))else:ax= tmp.plot(ax=ax,color=tab20(index))

png

我们可以单独提取某一年的数据进行可视化,如下展示某一年速度:

one_year = trips_by_year.get_trajectory('91916A_2010-12-31 00:00:00')
one_year.add_speed(overwrite=True,units=("km", "h"))
# 可视化当年的轨迹
one_year.plot(column='speed', cmap='tab20', legend=True)
<Axes: >

png

当然也可以可是轨迹中某一天的结果,如下所示:


fig, ax = plt.subplots(figsize=(5, 5))
one_year.plot(ax=ax, linewidth=2.0, color='black')
GeoDataFrame([one_year.get_row_at(datetime(2010,9,1))]).plot(ax=ax, color='red',zorder=2)
GeoDataFrame([one_year.get_row_at(datetime(2010,10,1))]).plot(ax=ax, color='red',zorder=2)
GeoDataFrame([one_year.get_row_at(datetime(2010,11,1))]).plot(ax=ax, color='red',zorder=2)
<Axes: >

png

此外也可以选择一个区域,比如栖息地。以确定某一轨迹到达该区域的次数及进出时间。

area_of_interest = Polygon([(30, 25), (50, 25), (50, 15), (30, 15), (30, 25)])
plotted_area_of_interest = GeoDataFrame(pd.DataFrame([{'geometry': area_of_interest, 'id': 1}]), crs='EPSG:3857')
# 输入的是单个轨迹得到的单个轨迹位于该区域的片段
arrivals = [traj for traj in my_traj.clip(area_of_interest)]
# 位于该区域的轨迹片段数
print("Found {} arrivals".format({len(arrivals)}))
# 提取每条轨迹到达该区域的时间
for traj in arrivals:name = traj.df['individual-local-identifier'].iloc[0]start_time = traj.get_start_time()end_time = traj.get_end_time()print("Individual {} arrived at {} left at {}".format(name, start_time, end_time))
Found {6} arrivals
Individual 91916A arrived at 2009-12-23 02:58:36.946571 left at 2010-04-17 14:31:42.170632
Individual 91916A arrived at 2010-10-30 20:55:36.697832 left at 2011-04-03 11:40:12.803103
Individual 91916A arrived at 2011-11-09 20:25:19.550486 left at 2012-04-04 09:12:15.852988
Individual 91916A arrived at 2012-10-14 11:25:28.063070 left at 2013-03-30 04:22:18.125679
Individual 91916A arrived at 2013-10-08 04:17:33.524488 left at 2014-03-29 18:26:54.311106
Individual 91916A arrived at 2014-10-28 19:05:32.941608 left at 2015-04-07 22:59:45.284499

step3 群体分析

step2中对单个轨迹Trajectory类进行分析,同样的函数也可以应用到轨迹集合TrajectoryCollection类中以对多条轨迹进行分析。

year_of_interest = 2010
# 输入的是多个轨迹得到的位于该区域的单个轨迹
trajs_in_aoi = tc.clip(area_of_interest)
# 提取2010到达该区域的轨迹
relevant = [ traj for traj in trajs_in_aoi if traj.get_start_time().year <= year_of_interest and traj.get_end_time().year >= year_of_interest]
# 查看有多少只海鸥,也就是多少条单个轨迹到达该区域
print("Found {} arrivals".format(len(relevant)))
Found 16 arrivals
# 查看哪些海鸥到达该区域以及停留时间
for traj in relevant:print("Individual '{}' arrived at {} (duration: {})".format(traj.df['individual-local-identifier'].iloc[0], traj.get_start_time().date(), traj.get_end_time()-traj.get_start_time()))
Individual '91732A' arrived at 2010-04-10 (duration: 5 days, 21:27:11.670985)
Individual '91737A' arrived at 2009-12-08 (duration: 140 days, 11:11:29.473206)
Individual '91761A' arrived at 2010-04-11 (duration: 12 days, 6:10:03.857850)
Individual '91761A' arrived at 2010-10-04 (duration: 6 days, 10:42:00.340093)
Individual '91762A' arrived at 2010-04-19 (duration: 42 days, 1:28:22.699066)
Individual '91771A' arrived at 2009-12-23 (duration: 66 days, 8:05:31.053782)
Individual '91789A' arrived at 2009-11-10 (duration: 550 days, 20:39:18.769409)
Individual '91824A' arrived at 2010-05-06 (duration: 12:53:53.409236)
Individual '91832A' arrived at 2010-04-21 (duration: 3 days, 5:48:46.276706)
Individual '91832A' arrived at 2010-09-23 (duration: 1 day, 1:40:25.175880)
Individual '91837A' arrived at 2010-05-04 (duration: 1 day, 18:38:46.170554)
Individual '91846A' arrived at 2010-05-15 (duration: 10 days, 2:50:28.505732)
Individual '91862A' arrived at 2010-01-06 (duration: 248 days, 6:10:16.962620)
Individual '91910A' arrived at 2010-09-28 (duration: 2 days, 20:34:31.117736)
Individual '91916A' arrived at 2009-12-23 (duration: 115 days, 11:33:05.224061)
Individual '91916A' arrived at 2010-10-30 (duration: 154 days, 14:44:36.105271)

对于’91761A’这只海鸥,可以看到进入该区域两次。我们可以展示这只海鸥在当年的飞行轨迹并添加背景地图以更好展示这只海鸥的行进情况。

my_traj = tc.get_trajectory('91761A')
# 提取当年飞行结果
my_traj = my_traj.get_segment_between(datetime(year_of_interest,1,1), datetime(year_of_interest,12,31))
# 压缩轨迹,tolerance越大压缩比例越大
my_traj = mpd.DouglasPeuckerGeneralizer(my_traj).generalize(tolerance=2.0)
my_traj.plot()
<Axes: >

png

通过如下代码绘制该海鸥的飞行轨迹,可以看到其在去途和归途都穿过了该地区,其中添加的背景地图具体使用方法见GeoPandas叠加背景地图。

import contextily as cx
fig, ax = plt.subplots(figsize=(3, 5))
# 添加飞行轨迹
ax = my_traj.plot(ax=ax,color='blue', capstyle='round')
# 添加停留区域
ax = gpd.GeoSeries(area_of_interest).plot(ax=ax, color='lightgray')
# 绘制1月到11月的位置,12月该鸟飞回原地
for i in range(1,12):# 红点表示上半年,绿点表示下半年color = 'red' if i<6 else 'green'GeoDataFrame([my_traj.get_row_at(datetime(year_of_interest,i,1))]).plot(ax=ax, color=color,zorder=2)# 添加背景地图,zoom越大越精细,这里使用自适应zoom。
cx.add_basemap(ax = ax , source=cx.providers.Stamen.TonerLite,crs=my_traj.df.crs, zoom='auto',attribution="")
# 保存图片
fig.savefig("res.jpg",dpi=300)

png

3 参考

  • movingpandas
  • movingpandas-examples
  • movingpandas-examples-data
  • shapely-objectdistance
  • 豪斯多夫距离
  • GeoPandas叠加背景地图

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.mzph.cn/news/19828.shtml

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈email:809451989@qq.com,一经查实,立即删除!

相关文章

[JAVA基础]自动拆装箱NPE问题

1.自动拆装箱场景 自动装箱 当把字面量转换成包装类的时候会自动装箱 比如&#xff1a; Integer a 1; Integer b 1; 自动拆箱 当你对包装类的对象进行运算&#xff08;如加法、减法等&#xff09;时&#xff0c;Java会自动进行拆箱操作。拆箱是将包装类型的对象转换为相应的基…

登月再进一步:Apollo自动驾驶的里程碑

前言 「作者主页」&#xff1a;雪碧有白泡泡 「个人网站」&#xff1a;雪碧的个人网站 「推荐专栏」&#xff1a; ★java一站式服务 ★ ★前端炫酷代码分享 ★ ★ uniapp-从构建到提升★ ★ 从0到英雄&#xff0c;vue成神之路★ ★ 解决算法&#xff0c;一个专栏就够了★ ★ 架…

数据结构 | 递归

目录 一、何谓递归 1.1 计算一列数之和 1.2 递归三原则 1.3 将整数转换成任意进制的字符串 二、栈帧&#xff1a;实现递归 三、递归可视化 四、谢尔平斯基三角形 五、复杂的递归问题 六、动态规划 一、何谓递归 递归是解决问题的一种办法&#xff0c;它将问题不断地分…

Google Colab 部署 Stable Diffusion Webui

1. 下载 vae 访问: https://huggingface.co/stabilityai/sd-vae-ft-mse-original/blob/main/vae-ft-mse-840000-ema-pruned.ckpt 根据下图的红圈部分&#xff0c;点击下载。 将其上传至GoogleDrive根目录 点击新建后&#xff0c;然后点击上传文件。 2. 访问 colab 首先第一步…

数据结构——单链表

无头单向非循环链表的建立 前言——什么链表链表形象图链表分类 一、Single_linked_list.h头文件的建立二、Single_linked_list.c功能函数的定义Single_linked_list_test.c主函数的定义四、代码运行测试五、Single_linked_list完整代码演示&#xff1a;总结 前言——什么链表 链…

一条自由游动的鲸鱼

先看效果&#xff1a; 再看代码&#xff1a; <!DOCTYPE html> <html lang"en"> <head><meta charset"UTF-8"><title>鲸鱼</title><style>#canvas-container {width: 100%;height: 100vh;overflow: hidden;}&l…

HTML之表单标签

目录 表单标签 Form表单 定义&#xff1a; 基本语法结构&#xff1a; form属性&#xff1a; enctyoe属性 fieldeset标签 fieldeset属性 legend标签 label标签 优势 label属性 input标签 input属性 input标签中的type属性 text text输入框有以下配套属性 searc bu…

太魔幻了!计算机就业哀鸿遍野,而高考生疯狂涌向计算机专业

在张雪峰推荐的几大专业里&#xff0c;计算机专业是其中之一。近几年&#xff0c;计算机专业报考热度不减&#xff0c;但就业前景却令人堪忧&#xff0c;互联网裁员接二连三&#xff0c;许多码农找不到工作。 一位网友感叹&#xff1a;一边是计算机就业哀鸿遍野&#xff0c;一…

IDEA中怎么使用git下载项目到本地,通过URL克隆项目(gitee\github)

点击 新建>来自版本控制的项目 点击后会弹出这样一个窗口 通过URL拉取项目代码 打开你要下载的项目仓库 克隆>复制 gitee github也是一样的 返回IDEA 将刚刚复制的URL粘贴进去选择合适的位置点击克隆 下载完成

Mac电脑目录

System&#xff08;系统&#xff09;Applications&#xff08;应用程序&#xff09;应用程序目录&#xff0c;默认所有的GUI应用程序都安装在这里User&#xff08;用户&#xff09;存放用户的个人资料和配置。每个用户有自己的单独目录Library&#xff08;资料库&#xff09;系…

Android的Handler消息通信详解

目录 背景 1. Handler基本使用 2. Handler的Looper源码分析 3. Handler的Message以及消息池、MessageQueue 4. Handler的Native实现 4.1 MessageQueue 4.2 Native结构体和类 4.2.1 Message结构体 4.2.2 消息处理类 4.2.3 回调类 4.2.5 ALooper类 5. 总结&…

基于fpga_EP4CE6F17C8_秒表计数器

文章目录 前言实验手册一、实验目的二、实验原理1&#xff0e;理论原理2&#xff0e;硬件原理 三、系统架构设计四、模块说明1&#xff0e;模块端口信号列表dig_driver(数码管驱动模块)key(按键消抖模块)top(顶层模块) 2&#xff0e;状态转移图3&#xff0e;时序图五、仿真波形…

目标检测与跟踪 (2)- YOLO V8配置与测试

系列文章目录 第一章 目标检测与跟踪 &#xff08;1&#xff09;- 机器人视觉与YOLO V8 目标检测与跟踪 &#xff08;1&#xff09;- 机器人视觉与YOLO V8_Techblog of HaoWANG的博客-CSDN博客3D物体实时检测、三维目标识别、6D位姿估计一直是机器人视觉领域的核心研究课题&a…

7.物联网操作系统互斥信号量

1.使用互斥信号量解决信号量导致的优先级反转&#xff0c; 2.使用递归互斥信号量解决互斥信号量导致的死锁。 3.高优先级主函数中多次使用同一信号量的使用&#xff0c;使用递归互斥信号量&#xff0c;但要注意每个信号量的使用要对应一个释放 优先级翻转问题 优先级翻转功能需…

gateway过滤器没生效,特殊原因

看这边文章的前提&#xff0c;你要会gateway&#xff0c;知道过滤器怎么配置&#xff1f; 直接来看过滤器&#xff0c;局部过滤器 再来看配置 请求路径 http://127.0.0.1:8080/appframework/services/catalog/catalogSpecials.json?pageindex1&pagesize10&pkidd98…

网络安全知识点整理(作业2)

目录 一、js函数声明->function 第一种 第二种 第三种 二、this关键字 this使用场合 1.全局环境 2.构造函数 3.对象的方法 避免多层this 三、js的同步与异步 定时器 setTimeout和setInterval 同步与异步的例子 四、宏任务与微任务 分辨宏任务与微任务 一、js…

深入解析项目管理中的用户流程图

介绍用户流程图 用户流程图的定义 用户流程图(User Flow Diagram)是一种可视化工具&#xff0c;它描绘了用户在应用或网站上完成任务的过程。这些任务可以是购物、注册账户、查找信息等&#xff0c;任何需要用户交互的动作都可以在用户流程图中找到。 用户流程图的重要性 用…

我对排序算法的理解

排序算法一直是一个很困惑我的问题&#xff0c;早在刚开始接触 数据结构的时候&#xff0c;这个地方就很让我不解。就是那种&#xff0c;总是感觉少了些什么的感觉。一开始&#xff0c;重新来过&#xff0c;认真来学习这一部分&#xff0c;也总是学着学着就把概念记住了。过了一…

案例|会展大数据服务系统 智能服务体系建设实践

根据相关的市场调查发现在国内进行的2000多个会展项目中&#xff0c;仅有15%的项目能够提供相关主题网站基本业务和服务&#xff0c;且提供的服务并不能够满足会展各方的应用&#xff0c;剩余的大部门会展项目信息系统构建的不够完善、技术性差、应用复杂等&#xff0c;甚至提供…

三、JVM-如何判断对象已死问题

内存模型以及如何判定对象已死问题 体验与验证 2.4.5.1 使用visualvm visualgc插件下载链接 &#xff1a;https://visualvm.github.io/pluginscenters.html 选择对应JDK版本链接—>Tools—>Visual GC 若上述链接找不到合适的&#xff0c;大家也可以自己在网上下载对应…