0 引言

手机信令数据已成为震后应急响应中摸清人口时空分布的有效手段, 特别是在地震黑箱期, 掌握人口动态分布对预判地震灾情显得尤为重要(李东平等, 2017；张小咏等, 2018)。随着GPS和卫星导航系统定位服务在手机终端的应用, 以及三大通信运营商在无线电通讯、网络社交软件定位、地图导航定位等海量数据的采集快捷性(李林侗等, 2022), 手机信令数据也在数次地震应急响应和备灾中发挥了重要作用。李东平等(2017)利用24小时连续定位手机数据, 对2017年8月8日四川九寨沟7.0级地震灾区人口进行了量化分析, 并分析了基站退服对判断估算极灾区的重要性；张小咏等(2018)利用震前和震后的手机定位数据, 采用标准差椭圆模型, 进行基于缺失手机定位数据的地震灾害影响场方向模拟；李林侗等(2022)利用手机位置信息分析了2016年5月22日辽宁朝阳市朝阳县4.6级地震人口分布和流动趋势, 并绘制了人口流动图和应急疏散路线图；徐敬海等(2023)以四川九寨沟7.0级地震为例, 利用手机位置数据, 采用冷热点空间分析方法, 分析地震前后人员分布增减明显区域, 识别出地震重灾区, 确定了人员安置区；陈伟等(2023)基于手机信令数据, 对天津市小白楼中央商务区进行了地震应急避难疏散模拟, 给出避难需求建议；尹晶飞等(2023)基于手机位置研究了2022年9月5日泸定6.8级地震震前人口分布状态以及震后的人口热力变化情况, 认为人口热力数据可为快速获取灾区人口动态变化提供有效途径。由此可见, 手机信令数据作为新兴的信息源, 在灾害管理和应急服务中展现出重要价值, 可为地震事件发生后获取第一手地震灾情提供有效的辅助判断依据。

人口热力分析应用可主要归纳为人口分布与密度分析, 比如事前采用动态的人口密度代替静态的人口普查结果, 提高地震灾害损失快速评估系统中人员伤亡计算的准确性(肖爽等, 2023)；事中用于判断交通拥堵情况(于磊, 2024)、城市交通通行情况、识别重灾区影响场范围以及确定发震断裂带方向；事后利用人口流动性与城市功能区进行匹配, 寻找目前应急避难场所建设中存在的问题。但对基于手机信令数据的人口热力分析在应急响应和灾情情景构建上的应用性, 尚未进行充分的特征发掘和归纳总结, 特别是在生命线工程、地质灾害、次生灾害等方面的匹配特征。

2023年12月18日23时59分积石山发生6.2级地震(35.70°N, 102.79°E), 震源深度10km, 本次地震是逆冲型的浅源地震, 释放的能量较大, 造成了强烈的地表震动, 造成约151人死亡, 983人受伤, 地震引发当地部分老旧房屋倒塌, 与以往同等地震相比, 其表现出“同震重灾”的特征, 也引发了特殊流滑地质灾害(蒋伟等, 2024), 特别是地震发生于夜间, 这对应急救援提出了更高的要求。本文基于手机信令数据, 利用KDE(Kernel Density Estimation, 核密度估计)方法, 对震前和震后的人口动态变化进行了量化分析。通过提取震中区域内的手机位置数据, 比较了震前1h和震后1h的人口热力分布, 以识别出受灾严重的区域。区别于以往地震事件的人口热力分析, 发现本次地震人口流动性与道路交通主干线、地震地质灾害以及通信基站退服的之间有高度相关性, 其表现特征与震后的应急排查结果符合程度较高, 特别是在重灾区, 人口热力变化显著, 高变化率的优势方向与本次地震烈度影响场的长轴方向一致性较好。这表明, 手机信令数据对识别Ⅷ度以上的重灾区具有较好的适应性。这些表现可为第一时间掌握重灾区、生命线工程(如道路交通主干线和通信网络)及地震引发的地质灾害等灾情提供判断依据, 补充地震应急快速评估产品, 为应急处置、续报和对策建议的产出提供有力的帮助。

1 数据获取 1.1 研究区概况

本次地震震中距积石山保安族东乡族撒拉族自治县约8km, 距离临夏市约40km。震中位于地震基本烈度Ⅶ度设防区, 但本次地震造成的最大烈度为Ⅷ度, 宏观震中主要位于大河家镇、刘集乡、柳沟乡一带, 震源破裂过程显示主震呈现双侧破裂特征, 以NNW向破裂为主, 等震线长轴呈NNW走向, 长轴124km, 短轴85km, Ⅵ度区及以上面积8364km²。其中, 甘肃省5232km², 涉及甘肃省3个市(州)9个县(市、区)88个乡镇(街道)以及太子山天然林保护区和盖新坪林场(应急管理部, 2023)。Ⅷ度区的大河家镇(甘肃)、刘集乡(甘肃)、中川乡(青海)、石塬镇(甘肃)等灾情表现严重(应急管理部, 2023；李亦纲, 2024)。

1.2 数据源

本文研究所用数据主要来源于第三方推送平台通过SDK(Software Development Kit，软件开发工具包)获取的用户授权的设备定位信息(朱战宏, 2007；庞晓克等, 2019), 设备定位信息以Geohash编码形式存在, 通过解码可获取以150m×150m为格网中心点的经纬度、终端数量和采样时间(李东平等, 2017)。本数据主要是从第三方平台数据合作方(如移动位置应用开发商)收集的位置数据(如GPS数据), 每个数据点仅包含用户标识符、观察时间以及经纬度信息, 收集的数据点也为聚合数据, 以保护用户的隐私。

本次数据采集是以震中为中心, 截取宽度约70km, 长度约90km的长方形区域, 获取了该区域范围内12月18日23时0分(震前)和12月19日01时0分(震后)约460万条手机位置数据和43000个空间格网位置数据, 采样点分布如图 1所示。

2 KDE方法

KDE(Kernel Density Estimation)即核密度分析法, 又称为核密度估计, 是一种用于估计未知概率密度函数的非参数统计方法。该方法在处理空间数据时尤为有效, 其从地理学的角度考虑了空间位置的影响, 通过计算点要素附近的密度来构建平滑的密度场表面, 从而对地理上离散的点对象进行数据处理, 展示出静态点要素在空间分布上的密集程度(郑国强等, 2023), 常被用于数据的可视化、聚类分析和异常检测(张香成, 2023), 比如用于兴趣点(Point of Interest, POI)分布研究、人口分布研究、手机信令数据、热点探测等(张香成, 2023；张欣玥等, 2024；王宇宁等, 2024)。

对于给定的数据集x₁, x₂, …x_n, 核密度估计可以表示为

$\hat{f}_h(x)=\frac{1}{n h} {\sum}_{i=1}^n K\left(\frac{x-x_i}{h}\right)$

(1)

其中, n为样本数据中的数据个数, h为平滑带宽, x为样本数据中的一个数据点, K为核密度函数, 决定了权重的分配。核函数的形状及其值域直接影响在估计某一点x时所考虑的数据点数量和利用程度。

Silverman(1986)指出, 在大样本情况下, 非参数估计对核函数的选择相对不敏感, 而平滑带宽h的选取对估计结果的影响显著。如果h设置得过小, 密度估计可能会集中于观测数据附近, 从而产生多个错误的峰值；如果h设置过大, 则会导致密度估计过于平滑, 忽略数据中的重要波动特征。

在地理信息系统(GIS)的应用中, 通常将h定义为搜索半径。通过设置合理的h生成圆, 离圆心越近, 核密度越集聚。搜索半径的选择依赖于研究的具体尺度：较小的搜索半径能够更好地捕捉到局部特征, 而较大的搜索半径则更能反映宏观规律。基于不同的地理需求一般设置不同的搜索半径, 如尹晶飞等(2023)在研究泸定6.8级地震震前人口分布状态以及震后的人口热力变化情况时，采用2km的搜索半径；郑国强等(2023)在研究基于POI数据的城市功能区识别时对不同的POI分类设置了1~3km的搜索半径；张香成(2023)在研究上海市商业与居住空间关系式时，认为1.5km的搜索半径能够较好地识别该空间分布特征；张欣玥等(2024)在研究良渚聚落时空分布特征时，认为半径设置为30km可较为直观地反映良渚聚落的密度分布；王宇宁等(2024)在研究天津市零售业的空间演化特征时，选择5km作为双城零售业核密度估计的统一带宽。在人口热力空间分布上，为保持局部与整体的平衡, 通过反复测试, 认为在人口热力时空分析时设置1km的搜索半径足以反映出与地震相关的空间变化, 可以有效识别受灾严重的局部区域, 同时不会过于限制数据的空间分布。在人口热力对地震影响场的识别时, 则更关注整体的空间性, 此时设置h为5km。

综上所述, 核密度分析为本研究提供了强有力的工具, 能够直观展示地震前后人员流动的动态变化, 为后续的应急响应和资源配置提供重要支持。

3 人口热力分析 3.1 震前人口热力分析

从研究数据中提取震前1h(12月18日23时)的手机信令数据, 对其进行核密度分析, 发现震中2km范围内的人口热力密度较低, 介于200~400人/km², 距离震中约8km处的积石山保安族东乡族撒拉族自治县的人口热力密度高, 最高可达到4000人/km²以上, 其次为距离震中约40km的临夏回族自治州市。本次受灾较严重的大河家镇镇中心区域范围内, 官亭镇中心区域范围最高也可达2000人/km², 刘集乡、中川乡、石塬镇则介于1000~2000人/km²。从震前人口热力空间分布可知, 震中附近各乡镇行政所在地人口密度较大, 聚集程度高, 而相对较远的乡村的人口密度相对较小。人口热力显示, 距离震中20km(约以中川乡至震中的震中距为半径)范围内的高聚集人口区域较多(图 2), 这也是造成人员伤亡损失惨重的原因之一(陈虹, 2024)。从震后热力分布(图 3)发现, 震中区域范围较近的石塬镇、刘集乡、积石山县城、大河家镇、官亭镇、中川乡及较远的临夏回族自治州的人口热力高密度区的区域范围变大, 而部分热力密度区域也从高向低发展, 这些具体的时空特征需进行进一步人口热力时空分析。

3.2 人口热力时空分析

人口热力时空变化是预判地震灾情的重要手段, 因此将震后与震前人口热力作差值进行时空特征分析(图 4), 并特别关注其在道路交通、地震地质灾害以及通信基站退服的表现特征。

从图 4可知, 大区域范围内, 人口热力呈现出了明显的升高现象, 最明显的是距离震中较远的村落, 如新寺乡附近, 对比图 2, 出现了人口流动的增加, 这可能与深夜人群被地震惊醒, 查看手机, 移动端SDK服务被唤醒有关, 这种“升高”现象普遍穿插在各乡镇的村落区域范围内。对于人口集聚区较高的大河家镇、刘集乡、石塬镇等行政所在区域, 虽然出现明显的人员增加, 但同样人口减少现象也穿插其中, 特别是离震区较近的积石山县城、临夏回族自治州等(县)级政府驻地以及地级城市的中心区域呈现出明显的人口减少, 这可能与人群的避险行为和基站退服有关。为深入挖掘人口热力在生命线工程(如道路交通主干线和通信网络)及地震引发的地质灾害中的特征, 选取了7处区域进行了相关性分析, 见图 4和图 5。

3.2.1 道路交通主干线人口热力分析

本次地震地质灾害(崩塌、滑坡)对道路交通造成了较大的影响, 如分析区A(图 5(a))位于民和回族土族自治县满坪镇的满坪村和集场村、前河乡的张家寺村附近, 该区域为川大高速公路、前峡公路、S201(五十—官亭)省道的三角区域。分析区B(图 5(b))位于民和县中川乡下的虎狼城村附近, 分析区C(图 5(c))位于民和县官亭镇鲍家村附近, 二者均位于S201(五十—官亭)省道边缘, 均出现了“不规则斑点式”的人流减少现象, 初步推测可能是发生山体崩塌或滑坡导致车辆无法前行外加通信基站退服所致。12月19日, 中华人民共和国公安部发布的消息显示S201张家寺村段、虎狼城村路段发生落石、塌方^①, 中华人民共和国青海省交通厅发布的消息显示省道201五官线K89+300、K92+350(官亭鲍家桥头)处发生塌方^②，验证了这一现象。分析区F位于G310连共线大清段狐跳峡大桥附近, 虽然人口减少现象不如A、B、C三处明显, 但呈现出斑点式加外向的流动减少现象, 初步推测也是因为出现了崩塌或滑坡, 对桥梁处造成了影响, 12月20日中国交通广播微博公众号发布“G310连共线大清段K8+400至K8+500(循化, 狐跳峡大桥)路段路面有落石, 交通中断”^③, 证实确实如此。分析区域外地震地质灾害(崩塌、滑坡)也随处看见, 如G310连共线大清段K17+000(循化县积石峡)处发生塌, 其表现特征与分析区F类似, 靠近分析区C的大河家桥桥头(青海岸)发生山体坍塌, 其表现特征与分析区C相似, 其中图 5(i)是甘肃省交通厅发布的救援队伍开展机械作业清除道路塌方的照片^④。由此可见人口流动性在道路交通上具有高度相关性, 虽然目前大部分易发生地质灾害的区域实现了“智治”, 但在强烈的地震作用下, 地质灾害点会呈现明显递增趋势。据统计, 本次地震区域范围内, 已有历史地质灾害隐患点2758处, 新增滑坡、崩塌灾害点77处, 复活隐患点88处(郭富赟等, 2024), 这些新增和复活的隐患点, 特别是道路交通主干线两侧山体的塌方和滑坡, 造成了交通中断, 如何在地震灾害发生后第一时间将其识别出来并进行抢修, 以方便救援力量的到达和伤员的运送, 显得尤为重要, 此时人口热力在道路交通主干线上的高相关性表现特征, 可为判断这些受伤害路段提供帮助, 可弥补深夜地震中的遥感影像滞后、无人机无法拍摄及人为巡查周期长的缺陷, 同时分析出的关键节点，也可为其部署人力排查故障路段的轻重缓急提供支持, 可为交通系统的应急响应辅助决策提供新思路。

① https://www.gov.cn/lianbo/bumen/202312/content_6921239.htm

② https://jtyst.Qinghai.gov.cn/jtyst/2023-12/19/article_2023121915184179620.html

③ https://weibo.com/2789194754/NxZLI6tG2

④ http://jtys.gansu.gov.cn/jtys/c106395/202312/173821053.shtml

3.2.2 村庄地震地质灾害人口热力分析

通过分析可知, 分析区D主要为中川乡滑坡—泥流区域, 该滑坡在流动过程中经过了草滩村和金田村, 造成了19条道路的中断和输电塔倒塌(陈虹, 2024；郭富赟等, 2024；王辽等, 2024；李亚军等, 2024), 特别是流滑的物源区后缘造成了两处高压电塔倒塌, 如图 5(h)^⑤, 因此从该区域可明显发现, 输电塔的倒塌造成了通信中断或断电, 引起基站退服, 致使区域上方即盘各村、中川乡行政所在地的斑点式人口减少, 而推测流滑的金田村则是基站退服引起了斑点式的减少, 这也明显表明人口热力分布变化与地震地质次生灾害建立了明显的关联。

⑤ https://zhuanlan.zhihu.com/p/673571203

3.2.3 村庄区域人口热力分析

离交通线主干线较远处出现了“规则斑点式”人口流动减少现象, 初步推测为地震作用导致的通信基站退服现象, 如位于分析区E的大河家镇梅坡村、周刘家、甘河滩村区域和分析区C的鲍家村东南方向, 均呈现出典型的“斑点式”人口流动较少状态(李环宇等, 2024)。这些区域位于Ⅷ区内, 震区主要为土木、砖木和未设防砖混结构, 其抗震性能较差, 外加此次地震地震动场地放大效应显著(杨永强等, 2024), 导致灾区大规模房屋倒塌和损坏的同时电力和通信基础设施也遭受了不同程度的破坏, 造成基站退服现象较为严重。据统计, 本次地震退服了314座基站, 其中甘肃积石山县134座(陈虹, 2024)。

3.2.4 通信人口热力分析

通过分析可知, 分析区G整个区域呈现了连续的类线性的人口流动减少状态, 此区域主要分布S202省道平大线(平安—大力加山)和G310(乐都—大力加山高速公路), 沿线还分布了各种村庄, 因此通过分析推测, 连续的人口减少不单单是地震地质灾害对道路中断造成的影响, 更有可能是周边发生了比较严重的光缆信号中断, 对系列基站造成了退服, 从12月22日搜狐新闻给出的“中国广电青海公司第一时间抢通地震受损光缆信号”信息中也得到了相关验证^⑥, 信息指出“位于循化县的狐跳峡大桥因山体滑坡导致青海循化至民和方向21km、27km、55km等多处省干光缆线路在地震中受损, 初步统计受损光缆长度约为50km”。同时12月18日青海省交通运输厅也发布信息显示“168公里至174公里500m(大力加山)左侧路段为傍山路段, 山体风化严重, 易发生落石、泥石流, 存在安全隐患”^⑦, 12月19日光明时政也指出“S202(循化至大力加山)线, K30+800至K44+816(循化县大力加山)路段多处塌方, 部分道路中断”^⑧。综上可知, 这种双重的作用导致了该区域呈现出了连续的人口流动减少, 同样这种现象也反映在了整个清水乡至街子镇方向上。

⑥ https://www.sohu.com/a/746296509_121106869

⑦ https://jtyst.Qinghai.gov.cn/jtyst/2023-12/18/article_2023121808560647135.html

⑧ https://politics.gmw.cn/2023-12/19/content_37037746.htm

综合A、B、C、D、E、F和G区域的情况可知, 人口热力时空分析可为生命线工程(道路交通主干线、通信网络)、地震地质灾害判断评估起到辅助决策的作用, 在地震作用下不同形式的灾害具有的流动人口减少特征反响不同, 交通线主干线上的不规则斑点式人口减少多与地质灾害和基站退服引起的线路中断有关；交通主干线线上但离村庄较远区域呈现大区域范围的线性人口流动减少现象, 可能是地震地质灾害和通信中断的双重作用所致；村庄行政区范围内的斑点式人口流动减少, 则更多与通信基站退服有关。

4 人口热力对地震影响场的识别

有学者研究发现, 人口空间分布变化可识别重灾区的地震影响场, 如李东平等(2017)利用人口变化率方法分析了基站退服对判断估算极灾区的重要性；张小咏等(2018)采用标准差椭圆模型, 进行了基于缺失手机定位数据的地震灾害影响场方向模拟；徐敬海等(2023)采用冷热点空间分析方法识别出地震重灾区；尹晶飞等(2023)利用人口变化率方法推测出地震影响场与人口热力变化的关系。因此，本文利用人口变化率(震后与震前人口变化的比率)进行核密度分析, 以对地震影响场进行识别, 人口变化率核密度识别结果见图 6, 其中变化率核密度值是根据变化情况进行的分段, 仅代表变化的剧烈程度, 值越高代表人口热力的变化性越大, 即

$R=(N-K) / K$

(2)

式中, N为震后某一时刻单位格网内人口数量, K为震前某一时刻单位格网内人口数量, 本研究中单位格网大小为5km×5km, R为人口变化率, R的绝对值越高则人口变化率越高。

总体上, 清晰可见本次地震人口变化率核密度与中华人民共和国应急管理部发布的《甘肃积石山6.2级地震烈度图》(应急管理部，2023)中Ⅷ度区的方向(NNW向)匹配结果较好, 经计算相差约25°, Ⅷ度区(除烈度异常区)内, 人口热力变化最剧烈, 而Ⅶ度区内人口变化程度不及Ⅷ度区, 这一结果也与前人研究较为一致, 即认为通常地震烈度Ⅷ度及以上即可导致通信基站大量退服(李帆, 2022；程明超等, 2017；毛晨曦等, 2018), 震后的统计结果显示, 本次地震共造成314座基站退服, 2条通信光缆中断(陈虹, 2024)。大河家镇、刘集乡、石塬镇至积石山县城区域等的人口变化最明显, 其优势方向也与地震烈度圈长轴方向大致一致, 这种特征表现与地震灾害考察结果相符, 可为快速修正快速评估中产出的地震影响场提供一定参考。但本次的烈度异常区的人口变化率反映不明显, 可能与烈度异常区内人口密度较低有关, 从侧面反映出人口的密度集聚性对人口热力分析影响较大。

5 结论

通过对甘肃积石山6.2级地震前后人口热力时空分布变化的追踪分析, 探讨了该变化与道路交通主干线、地震引发的地质灾害以及通信基站退服等灾害识别之间的相关性。利用变化率对识别重灾区进行了深度挖掘, 得到以下结论：

(1) 人口热力在交通主干线上但周边村庄相对密集的区域, 若观察到大范围的“线性”流动减少特征, 可推测为严重的通信故障造成；在村庄行政区范围内出典型的“斑点式”流动减少特征, 可推测为通信基站退服引起。这些表现特征可为排查通信网络故障, 布设通信保障设备提供第一手灾情资料。

(2) 人口热力表现与地震导致的次生地质灾害(如滑坡和崩塌)具有较强的相关性, 特别是在道路交通主干线路上, 会呈现出“不规则斑点式”减少的特征, 在地震地质灾害易发区附近, 则呈现出大范围的斑点式减少的特征。这些表现特征可为判断这些受伤害路段提供帮助, 特别可弥补深夜地震中遥感影像滞后、无人机无法拍摄、人为巡查周期长的缺陷, 可为交通系统的应急响应辅助决策提供新思路。

(3) 人口热力变化率对于识别重灾区具有良好的适应性。受灾严重的大河家镇和刘集乡表现明显, 其变化方向与烈度影响场的长轴方向基本一致, 可为修正快速评估中产出的地震影响场提供参考, 同时需要指出，人口热力变化率在人口密度较低的地区识别影响场存在一定局限性。

通过积石山地震事件, 对人口热力在生命线工程(如道路交通主干线和通信网络)及地震引发的地质灾害的初步特征进行了挖掘, 展示了这一技术在优化地震应急服务能力和丰富灾后快速评估产品中的潜力。同时, 认识到地震灾害情景复杂, 受到地震强度、地形条件、承灾体能力、发生环境及人口分布等多种因素的影响, 这也导致人口热力在不同地震事件中的表现具有一定的复杂性。因此, 如何更充分地挖掘人口热力在应急响应中的服务能力尚需进一步的实践积累。建议在未来的工作中继续探索这一技术在自然灾害中的应用, 以提升整体的灾害应对效率和科学性。