交通流合成数据生成原理及实现代码

移动数据是设备的地理位置，通过正常活动被动产生。 它具有从交通规划到迁移预测的重要应用。 由于移动数据稀有且难以收集，研究人员已开始探索综合生成移动数据的解决方案。

在本文中，我将讨论一种用于合成移动数据的简单解决方案。 该合成数据可用于研究目的和训练/微调算法。 例如，可以综合生成标记的移动数据，并训练模型来预测城市交通拥堵。 然后，训练好的模型可以应用于现实生活中的数据。

代码可以在这里找到，你可以使用这个 colab 笔记本自己尝试一下。

 在线工具推荐： Three.js AI纹理开发包 - YOLO合成数据生成器 - GLTF/GLB在线编辑 - 3D模型格式在线转换 - 可编程3D场景编辑器

1、合成数据结构

合成数据将代表从手机设备收集的位置数据记录。 通常，此类数据包含以下属性：

• phone_id — 手机的唯一标识符
• phone_type — 手机操作系统（iOS / Android）
• 时间戳（以纪元时间为单位）
• 纬度
• 经度
• 精度（以米为单位）

phone_id,phone_type,timestamp,latitude,longitude,accuracy
cfcd208495d565ef66e7dff9f98764da,iOS,2022-01-01 00:00:10,30.34639898170681,-97.73158031325373,81.12640643460236
cfcd208495d565ef66e7dff9f98764da,iOS,2022-01-01 00:09:13,30.34642963771432,-97.73140260170493,25.094041385659516
cfcd208495d565ef66e7dff9f98764da,iOS,2022-01-01 00:18:55,30.346380588102306,-97.73148816652471,122.98442053366325
cfcd208495d565ef66e7dff9f98764da,iOS,2022-01-01 00:28:16,30.34643576891581,-97.7315210760708,185.08383036465293
cfcd208495d565ef66e7dff9f98764da,iOS,2022-01-01 00:37:36,30.346366281965476,-97.73138285597729,55.565838723523676
cfcd208495d565ef66e7dff9f98764da,iOS,2022-01-01 00:47:04,30.346372413166982,-97.7316461323459,36.26985794887487
cfcd208495d565ef66e7dff9f98764da,iOS,2022-01-01 00:56:20,30.346384675569983,-97.73156056752609,144.42897717113405
cfcd208495d565ef66e7dff9f98764da,iOS,2022-01-01 01:05:23,30.34640511290831,-97.73140918361412,43.246590356759704
cfcd208495d565ef66e7dff9f98764da,iOS,2022-01-01 01:14:25,30.34634584462714,-97.73150791225235,46.35681270455376
cfcd208495d565ef66e7dff9f98764da,iOS,2022-01-01 01:23:57,30.346441900117316,-97.73162638661823,157.83382148809042

2、获取公开的位置数据

在美国选择一个位置并创建一个 x 米的 bbox（边界框）。 接下来，获取公共数据集：

• 居住地点 — bbox 内的建筑物示例 (ArcGIS Rest API)
• POI（兴趣点）— bbox 内的企业 （Kaggle数据集）
• Roads — bbox 内的所有道路，以图表形式表示（OSM Overpass API，使用 osmnx）

创建包围框：

import osmnx as ox
lat, lng, radius = 45.496920, -122.665803, 10000
bbox = ox.utils_geo.bbox_from_point((lat, lng), radius)
bbox

使用 arcgis_rest_url 获取 bbox 内建筑物的多边形，最多2000 个多边形样本。

import geopandas as gpdARCGIS_REST_URL = 'https://services.arcgis.com/P3ePLMYs2RVChkJx/ArcGIS/rest/services/MSBFP2/FeatureServer/0/query?f=json&returnGeometry=true&spatialRel=esriSpatialRelIntersects&geometry={"ymax":%s,"ymin":%s,"xmax":%s,"xmin":%s,"spatialReference":{"wkid":4326}}&geometryType=esriGeometryEnvelope&outSR=4326'residence_df = gpd.read_file(ARCGIS_REST_URL % bbox)
residence_df['lat'] = residence_df['geometry'].centroid.y
residence_df['lng'] = residence_df['geometry'].centroid.x
residence_df.head()

结果如下：

	BlockgroupID	geometry	lat	lng
0	410510005012	POLYGON ((-122.59120 45.48139, -122.59120 45.4...	45.481425	-122.591104
1	410510059001	POLYGON ((-122.67869 45.49746, -122.67869 45.4...	45.497409	-122.678765
2	410510074002	POLYGON ((-122.60844 45.56566, -122.60844 45.5...	45.565693	-122.608380
3	410050216022	POLYGON ((-122.59890 45.43605, -122.59890 45.4...	45.435959	-122.598965
4	410510017012	POLYGON ((-122.58212 45.52941, -122.58212 45.5...	45.529375	-122.582223

使用 Kaggle API 下载 POI 数据集。 然后解析它，将其加载到 geopandas 并将数据集过滤为仅 bbox 内的点。

def isin_box(lat, lng, bounds):ymax, ymin, xmax, xmin = boundswithin = Falseif xmin < lng < xmax:if ymin < lat < ymax:within = Truereturn withinimport os
import json
import pandas as pd
from shapely.geometry import Pointexport_path = '<your_export_path_here>'API_TOKEN = {"username": 'your_kaggle_username_here', "key": 'your_kaggle_key_here'}os.system('mkdir ~/.kaggle')
os.system('touch ~/.kaggle/kaggle.json')with open('/root/.kaggle/kaggle.json', 'w') as file:json.dump(API_TOKEN, file)os.system('chmod 600 ~/.kaggle / kaggle.json')
os.system(f'mkdir /{export_path}/pois_data/')os.system(f'kaggle datasets download -d timmofeyy/-pizza-hut-restaurants-locations-in-us -p /{export_path}/pois_data/pizza-hut --unzip')
os.system(f'kaggle datasets download -d timmofeyy/-walmart-stores-location -p /{export_path}/pois_data/walmart --unzip')
os.system(f'kaggle datasets download -d appleturnovers/dunkin-locations -p /{export_path}/pois_data/dunkin --unzip')
os.system(f'kaggle datasets download -d timmofeyy/-subway-locations-in-us -p /{export_path}/pois_data/subway --unzip')pois_dfs = []for chain_name in os.listdir('/content/pois_data'):chain_file = [i for i in os.listdir(f'/{export_path}/pois_data/{chain_name}') if i.endswith('.csv')][0]chain_df = pd.read_csv(os.path.join(f'/{export_path}/pois_data', chain_name, chain_file))chain_df['poi_name'] = chain_namechain_df = chain_df.rename(columns={'latitude': 'lat','longitude': 'lng','lon': 'lng','loc_lat': 'lat','loc_long': 'lng','Latitude': 'lat','Longitude': 'lng'})[['poi_name', 'lat', 'lng']]chain_df['lat'], chain_df['lng'] = chain_df['lat'].astype(float), chain_df['lng'].astype(float)pois_dfs.append(chain_df)pois_df = pd.concat(pois_dfs, ignore_index=True)
pois_df = pois_df.dropna(subset=['lat', 'lng'], how='any')pois_df = pois_df[pois_df.apply(lambda x: isin_box(x['lat'], x['lng'], bbox), axis=1)]
pois_df['id'] = pois_df.index
pois_df['geometry'] = pois_df.apply(lambda x: Point(x['lng'], x['lat']), axis=1)
pois_df = gpd.GeoDataFrame(pois_df)pois_df.head()    

结果如下：

	poi_name	lat	lng	id	geometry
2996	Pizanos Pizza	45.516327	-122.559344	2996	POINT (-122.55934 45.51633)
3626	Washington Street Palace	45.518799	-122.654747	3626	POINT (-122.65475 45.51880)
6131	Fred Meyer	45.586069	-122.738125	6131	POINT (-122.73813 45.58607)
7392	Moberi	45.510439	-122.682163	7392	POINT (-122.68216 45.51044)
7920	Ether	45.526996	-122.698438	7920	POINT (-122.69844 45.52700)

从 Overpass API 获取 OSM 道路：

G = ox.graph_from_bbox(*bbox, network_type="drive", retain_all=True, truncate_by_edge=True)

结果如下：

现在我们已经拥有了创建电话时间线所需的一切 — 居住位置（将用于在家中、家人和朋友家中的停留）、POI 位置（将用于商店访问）和道路（将用于 住宿期间开车）。 在生成实际的移动数据之前，我们将生成一个综合时间线，其中包含手机停留时间及其时间范围。

3、生成合成时间线

合成时间线逻辑将迭代开始日期和结束日期之间的所有日期，并随机化在工作场所、居住地点和 POI 的停留时间。 为了保证正常的人类行为，该逻辑将仅在工作日产生工作停留，并确保用户晚上回家。

在运行逻辑之前，请确保：

• 设置随机的家庭和工作地点
• 设置时间范围（开始日期和结束日期）
• 设置指定日期要访问的最大 POI 和最大居住地点

home_info = residence_df.sample(1).to_dict(orient='records')[0] # pick a random home location
work_info = pois_df.sample(1).to_dict(orient='records')[0] # pick a random work location
start_date = '2022-10-01' # timeline start time
end_date = '2022-11-01' # timeline end time
max_residence = 3 # max residence stay in a day
max_pois = 3 # max poi stay in a dayfrom datetime import timedelta
import numpy as npphone_timeline = []stays_counter = 0
for day in pd.date_range(start_date, end_date):n_residence = np.random.choice(range(1, max_residence))n_pois = np.random.choice(range(1, max_pois))levaing_hour = np.random.choice([6,7,8])start_time = dayend_time = day + timedelta(hours=int(levaing_hour))phone_timeline.append({'stay_id': stays_counter,"start_time": start_time,"end_time": end_time,'entity_id': 'home','entity_name': 'home','entity_type': 'home','lat': home_info['lat'],'lng': home_info['lng']})stays_counter+=1if day.weekday() not in [5,6]: # if working dayworking_hours = np.random.choice([7,8,9])start_time = end_timeend_time += timedelta(hours=int(working_hours))phone_timeline.append({'stay_id': stays_counter,"start_time": start_time,"end_time": end_time,"entity_id": 'work','entity_name': 'work','entity_type': 'work','lat': work_info['lat'],'lng': work_info['lng']})stays_counter+=1time_left = 24 - end_time.hourif n_residence !=0 or n_pois !=0:time_left -= 1 # leave 1 hour to stay at home at the end of the timelinetime_per_stay = np.floor(time_left/(n_residence + n_pois))if n_residence !=0:for row in residence_df.sample(n_residence).itertuples():start_time = end_timeend_time += timedelta(hours=int(time_per_stay))time_left -= time_per_stayphone_timeline.append({'stay_id': stays_counter,"start_time": start_time,"end_time": end_time,'entity_name': row.Index,'entity_type': 'building','entity_id': row.Index,'lat': row.lat,'lng': row.lng})stays_counter+=1if n_pois !=0:for row in pois_df.sample(n_pois).itertuples():start_time = end_timeend_time += timedelta(hours=int(time_per_stay))time_left -= time_per_stayphone_timeline.append({'stay_id': stays_counter,"start_time": start_time,"end_time": end_time,'entity_name': row.poi_name,'entity_type': 'venue','entity_id': row.Index,'lat': row.lat,'lng': row.lng})stays_counter+=1start_time = end_timeend_time += timedelta(hours=int(time_left))phone_timeline.append({'stay_id': stays_counter,"start_time": start_time,"end_time": end_time + timedelta(hours=abs(end_time.hour-24)),'entity_id': 'home','entity_name': 'home','entity_type': 'home','lat': home_info['lat'],'lng': home_info['lng']})stays_counter+=1
phone_timeline = pd.DataFrame(phone_timeline)
phone_timeline.head()

结果如下：

	stay_id	start_time	end_time	entity_id	entity_name	entity_type	lat	lng
0	0	2022-10-01 00:00:00	2022-10-01 08:00:00	home	home	home	45.410639	-122.552044
1	1	2022-10-01 08:00:00	2022-10-01 10:00:00	1290	1290	building	45.449074	-122.758561
2	2	2022-10-01 10:00:00	2022-10-01 12:00:00	249212	Al's Garden Center	poi	45.458765	-122.710870
3	3	2022-10-01 12:00:00	2022-10-01 14:00:00	502973	Lloyd Center	poi	45.532621	-122.653473
4	4	2022-10-01 14:00:00	2022-10-01 16:00:00	225659	B&B Barbershop	poi	45.483497	-122.581378

下面的 gif 显示了我们合成时间线中的第一天：

4、生成合成移动数据（信号）

我们的合成时间线已准备就绪，需要新的逻辑将其转换为合成信号。 我们时间线中的第一个事件是寄宿家庭 (00:00 -> 08:00)，所以让我们从为这次住宿生成信号开始。

4.1 静态模式信号

以下脚本将生成停留开始和停留结束之间的信号数据帧。 采样率（相邻信号之间的时间间隔）是一个可配置的参数。 我将其设置为 600 秒（5 分钟）。 每个信号的 lat,lng 将使用随机“噪声系数”进行噪声处理

def generate_static_signals(lat, lng, start_time, end_time, sampling_rate=600):signals_list = []signal_time = start_time + timedelta(seconds=int(np.random.choice([5,10,15]))) # start smaple couple of second after visits startnoise_list = np.linspace(0.9999997, 1.000003) # noise factorwhile end_time > signal_time:signals_list.append({'lat': lat*np.random.choice(noise_list), 'lng': lng*np.random.choice(noise_list),'timestamp': signal_time})signal_time += timedelta(seconds=sampling_rate - int(np.random.choice(range(0,59))))return pd.DataFrame(signals_list).sort_values('timestamp')

在第一次停留时应用逻辑将产生以下输出：

4.2 驾驶模式信号

我们时间线上的下一个事件是在“Residence 1290”停留，但在为本次停留生成信号之前，我们需要为将手机从出发地（家）带到目的地（“Residence 1290”）的驱动器生成信号 。

为此，我们将使用道路图并寻找从起点到目的地的最短路径。 然后，我们将以 60 秒的采样率在有序路段上随机生成信号。

from shapely.ops import linemerge, nearest_points, snap
from shapely.geometry import LineString, Point, MultiLineString
import taxicab as tcdef generate_route_signals(route_geo, start_time, sampling_rate=60):end_time = start_time + timedelta(minutes=int(np.random.choice(range(15,50))))N_POINTS = round((end_time - start_time).total_seconds() / sampling_rate)POINTS_PER_SEGMENT = 10noise_list = np.linspace(0.9999999, 1.000001) # add noise to pointssignal_time = start_time - timedelta(seconds=int(np.random.choice([5,10,15])))timestamps_list = []for i in  range(0,N_POINTS):timestamps_list.append(signal_time)signal_time += timedelta(seconds=sampling_rate + int(np.random.choice(range(0,9))))line_to_points = np.array([{'lat': y*np.random.choice(noise_list), 'lng': x*np.random.choice(noise_list)}for p1, p2 in zip(route_geo.coords, route_geo.coords[1:]) # iterate through line segmentsfor x, y in zip(np.linspace(p1[0], p2[0], POINTS_PER_SEGMENT),np.linspace(p1[1], p2[1], POINTS_PER_SEGMENT),)])indexes = np.sort(np.random.choice(range(0, len(line_to_points)),replace=False,size=min(N_POINTS, len(line_to_points)))) random_points = line_to_points[indexes]signals_df = pd.DataFrame([i for i in random_points])signals_df['timestamp'] = timestamps_list[:len(signals_df)]return signals_df, signal_timedef get_route_geometry(route, orig, dest):x, y  = [], []for u, v in zip(route[tc.constants.BODY][:-1], route[tc.constants.BODY][1:]):# if there are parallel edges, select the shortest in lengthdata = min(G.get_edge_data(u, v).values(), key=lambda d: d["length"])if "geometry" in data:# if geometry attribute exists, add all its coords to listxs, ys = data["geometry"].xyx.extend(xs)y.extend(ys)else:# otherwise, the edge is a straight line from node to nodex.extend((G.nodes[u]["x"], G.nodes[v]["x"]))y.extend((G.nodes[u]["y"], G.nodes[v]["y"]))final_route = []if route[2]:final_route.append(LineString([Point(orig[1], orig[0]),nearest_points(Point(orig[1], orig[0]),route[2])[1]]))final_route.append(route[2])else:final_route.append(LineString([Point(orig[1], orig[0]),Point(x[0],y[0])]))final_route.append(LineString([Point(lng,lat) for lng, lat in list(zip(x, y))]))if route[3]:final_route.append(route[3])final_route.append(LineString([nearest_points(Point(dest[1],dest[0]), route[3])[1],Point(dest[1],dest[0])]))else:final_route.append(LineString([Point(x[-1],y[-1]), Point(dest[1],dest[0])]))final_route_geo = linemerge(final_route)if isinstance(final_route_geo, MultiLineString):final_route_geo = linemerge([snap(i, j,0.000001) for i,j in zip(final_route, final_route[1:])])if isinstance(final_route_geo, MultiLineString):final_route_geo = list(final_route_geo)[np.argmax(np.array([len(i.coords) for i in list(final_route_geo)]))]return final_route_geo

这就是合成驾驶信号在地图上的样子：

4.3 完整的合成移动数据生成

在最后一步中，我们将迭代所有的合成时间线。 对于每次停留，我们将生成静态模式信号，并且在每两次停留之间，我们将生成驾驶模式信号。

def generate_signals_df(device_timeline, G):signals_dfs = []drive_end = Nonefor row in device_timeline\.join(device_timeline[['lat','lng']].shift(-1), lsuffix='_orig', rsuffix='_dest')\.dropna(subset=['lat_orig','lng_orig','lat_dest','lng_dest'],how='any').itertuples():if not drive_end:drive_end = row.start_timestay_signals = generate_static_signals(row.lat_orig, row.lng_orig, drive_end, row.end_time)signals_dfs.append(stay_signals)route = tc.distance.shortest_path(G, (row.lat_orig, row.lng_orig), (row.lat_dest,row.lng_dest))route_geo = get_route_geometry(route, (row.lat_orig, row.lng_orig), (row.lat_dest,row.lng_dest))drive_signals, drive_end = generate_route_signals(route_geo, row.end_time, 45)signals_dfs.append(drive_signals)device_signals = pd.concat(signals_dfs, ignore_index=True).sort_values('timestamp')return device_signals

好了！ 我们现在拥有由开源包和免费数据生成的完整的移动合成数据。

5、结束语

生成合成移动数据是可行的，不需要特殊资源。 我们需要的所有数据都在那里并且可以免费使用。 话虽如此，仍有改进的空间：

• 逻辑的幼稚：逻辑非常简单，因此输出并不能真正代表完整的人类行为（例如出国旅行、移民等）
• 逻辑的效率：逻辑需要一些时间来运行。 主要瓶颈是驾驶生成部分（最短路径计算）
• 逻辑的覆盖范围：仅支持美国，这是因为我使用了 USBuildingFootprints，它非常准确，但不幸的是仅覆盖美国。

---

原文链接：交通流动合成数据生成 - BimAnt