114 条人工同刻事件的 Hawkes 评分对比:pre-batch 事件项 183.475922 nats,逐行先更新得 396.770155 nats,多出 213.294233 nats;补偿积分均为 48.237046 nats。
图:同刻事件的两种评分约定;所有数字来自可独立运行的人工反例,并非 AAPL 模型拟合。

整理一份 NASDAQ ITCH 订单流样本时,我碰到一个不适合直接交给普通逐事件 Hawkes 循环的数字:114。

在 2019 年 7 月 30 日 AAPL 的常规交易时段数据里,有一个纳秒时间戳对应 114 条 execution 事件记录。还有各出现一次的 131 条、187 条批次。更重要的是,当天 49,172 条执行类消息中,有 30,080 条(61.17%)与其他记录共享纳秒时间戳。并列现象很常见,但上百条同刻记录属于极端尾部。

这里的“执行事件”是按研究定义归并的 ITCH 消息行,不是 114 个不同订单,更不是证明 114 笔成交在物理时间上同时完成。

问题是:如果程序每读一行就更新 Hawkes 激励状态,即使这一批记录的时间差全为零,第 114 条仍会被前 113 条“提前发生”的记录激励。若同刻批次混有不同事件类型,非对称激励矩阵还可能让评分依赖消息排列。说白了,逐行更新会让文件里的行序决定同一纳秒内谁先激励谁。真正需要决定的是:相同的记录时间戳,究竟允许哪种批内先后假设?

我们在拟合前就把研究的评分规则固定下来:同股票、同时间戳先成批;用批次之前的历史给整批评分;整批评分结束后,再一次性用实际事件数量更新历史。 这套 pre-batch 约定只规定当前研究怎样计算 Hawkes 事件项。交易所的消息发送顺序有价值,但不能未经建模就把它当作纳秒以下的物理时间。

实际检查代码,我先看三个地方:分组键是否同时包含股票和整数纳秒时间戳;状态更新是否等整批完成评分后才按真实事件数量进行;补偿积分是否覆盖整个有效会话。后面的公式和测试,就是对这三处做核对。

并列时间戳很常见,但 114 条是尾部大批次

Nasdaq TotalView-ITCH 5.0 规范将消息时间戳定义为“午夜以来的纳秒数”;E 为 Order Executed,C 为 Order Executed With Price。规范也说明,同一订单引用号可能有多条部分执行消息。

这次研究按股票、事件类型和整数纳秒时间戳聚合,并自行规定三类建模事件:add=A/F、cancel_delete=X/D、execution=E/C。根据 ITCH 5.0 规范,C 带有 Printable=Y/N 字段;这里的 execution=E/C 包含 C:Y 和 C:N 两种执行消息,统计的是订单执行记录,而非仅限可打印成交。P(非展示订单的 Trade 消息)和 Q(Cross Trade 消息)属于另外的消息类型,没有纳入这三类事件。这是一项研究事件定义,不是 Nasdaq 对所有成交的统一归并规则。

以下是 2019-07-30、AAPL、09:30–16:00 美东交易所时间 的 execution 聚合结果:

指标观测值
execution 消息记录数49,172
不同纳秒时间戳数29,379
含多条记录的时间戳数10,287
落在并列时间戳上的记录数30,080(61.17%)
并列批次大小的中位数/众数2 条/2 条
大小为 2–3 条的并列批次8,439(占并列批次 82.04%)
恰好 114 条记录的时间戳批次1
更大的 131 条、187 条批次各 1

这些数字来自当日完整常规交易时段的脱敏聚合;此前还做过 09:30–11:00 的跨交易日共同窗口对照,两种时间范围必须分开。现有汇总可以确认批次大小,却不能还原 114 条记录各自的订单身份或微观发生顺序。

把单记录时间戳也算进去,平均每个时间戳有 1.67 条执行记录;只看 10,287 个并列批次,平均是 2.92 条,中位数和众数均为 2 条。114 条约为并列批次均值的 39 倍;达到 114 条及以上的批次,当天一共只有 3 个(114、131、187)。因此,大批次适合用来揭示计算陷阱,却不能用它代表通常的并列规模。

并列记录多,还意味着这套 pre-batch 评分不利用批次内部的亚纳秒间隔来识别激励:同类型自激及跨类型激励主要靠不同批次之间的时间变化与类型组合获得信息。不能只凭 61.17% 的比例断言参数已经不可识别;能否识别还取决于真实间隔分布和激励矩阵。

同一个时间戳,谁有资格先激励谁?

三类型指数核 Hawkes 的条件强度可以写成:

lambda_k(t-) = mu_k(t)
             + sum_j [alpha[k,j] / tau
               * sum_{u<t} m_j(u) * exp(-(t-u)/tau)]

k 是目标事件类型,j 是来源类型,m_j(u) 是来源类型在时间戳 u 的实际记录数量。激励矩阵方向是 目标 × 来源;mu_k 单位为事件/秒,tau 单位为秒,alpha[k,j] 表示一个来源事件带来的积分激励量(无量纲),因此核在刚发生后的强度系数才是 alpha[k,j]/tau,单位为 1/秒。补偿积分中直接出现 alpha 而不是 alpha/tau,正是对指数核做了时间积分。

这里最重要的是求和中的 u<t:同一时间戳的当前批次,不应在本研究的 pre-batch 约定中提前进入自己的历史。

于是同时间戳、同股票的一批事件,对数似然的事件项为:

event_term(t) = sum_k m_k(t) * log(lambda_k(t-))

先计算全部类型的 lambda_k(t-),再乘各类型的 multiplicity。最后更新一次历史,历史里保存的是 m_add、m_cancel_delete、m_execution 三个真实计数,而不是把 114 条压成“一次事件”。

ITCH 的消息流具有发送顺序。按发送顺序逐条评分也可以成为另一种研究方案,但必须先交代发送顺序如何映射为事件时间与强度历史:它是另一套观测假设,不能伪装成相同模型下的代码优化。这里的逐行反例只检验时间差为零、却允许批内前序记录自激会带来的数值变化。

具体到强度计算,同批次所有事件类型必须读取同一份 t- 历史,批次评分结束后再按真实 multiplicity 更新。对应的合成回归测试还应覆盖:同批次换序不改变评分结果、不同股票的历史彼此隔离,以及跨交易日重置。

用 114 条人工记录,把两种算法的差算出来

下面不是 2019 年 AAPL 的参数估计,而是专门放大顺序问题的 Python 标准库教学反例:假设会话起点有 114 条相同类型事件,起点之前没有激励历史,基准强度 mu=5 次/秒,自激 alpha=0.6,tau=1 秒;会话观察到起点之后一秒。

import math

m, mu, alpha, tau, horizon = 114, 5.0, 0.6, 1.0, 1.0

# 同一批全部使用批次前的强度
pre_batch = m * math.log(mu)

# 错误对照:假装同一纳秒的前 i 条已激励第 i+1 条
row_ordered = math.fsum(
    math.log(mu + alpha * i / tau) for i in range(m)
)

# 完整似然还需扣除整个观测窗口的强度积分
baseline_area = mu * horizon
excitation_area = alpha * m * (-math.expm1(-horizon / tau))
compensator = baseline_area + excitation_area

print(f"prebatch_event_log_term: {pre_batch:.6f}")
print(f"row_ordered_event_log_term: {row_ordered:.6f}")
print(f"artificial_excess: {row_ordered - pre_batch:.6f}")
print(f"one_second_compensator: {compensator:.6f}")
print(f"one_second_full_LL: {pre_batch - compensator:.6f}")

在 Python 3.13.5 的隔离 Linux 环境中实际执行,退出码为 0,输出:

prebatch_event_log_term: 183.475922
row_ordered_event_log_term: 396.770155
artificial_excess: 213.294233
one_second_compensator: 48.237046
one_second_full_LL: 135.238876

仅改变“同刻记录是否逐行更新历史”这一计算约定,事件对数项就多出了 213.294233 nats。更重要的是:在这个示例中,所有记录的时间戳都为会话起点,两个评分约定在整批结束后保留完全相同的 114 条激励历史。对零时长批次而言,后续一秒的积分不依赖批次内处理顺序,因此两者补偿积分同为 48.237046。逐行评分的完整 LL 是 348.533109,pre-batch 是 135.238876,两者完整 LL 的差仍恰好为 213.294233。程序没有多看到一个事件,却让同一纳秒里前面的记录替后面的记录“刷”出了两百多个 nats。这不应被解释为某模型真实预测效果的增益,也不应把两种工作评分视为同一个连续时间观测模型的两个精确似然。

反过来,禁止批次内自激也不等于丢掉 114 条的影响。在一秒后的时刻,按这个教学参数,它们对自身类型强度的激励贡献为 114 × 0.6 × exp(-1) ≈ 25.163 次/秒。把这个批次错误地当成只有一条事件,下一时刻的强度同样会被低估。

事件项算对了,还不能叫完整 likelihood

另一个常见漏项比较隐蔽:程序把每个事件的 log(lambda) 累加完,就返回“对数似然”。实际上,点过程的完整条件对数似然还要扣掉整个有效会话的补偿积分:

LL = sum_k sum_batches m_k(t) * log(lambda_k(t-))
     - sum_k integral_session lambda_k(t) dt

我们在拟合前将日内基准固定为 78 个五分钟区间,覆盖 09:30–16:00 的 390 分钟。基准积分按每段与观测窗口的半开重叠时间精确求和;指数核历史项则按每批实际 multiplicity 积分:

baseline_area[k] = sum_r mu[k,r] * overlap_seconds(r, session)

excitation_area[k] = sum_j sum_{u in session}
    alpha[k,j] * m_j(u) * (1 - exp(-(end-u)/tau))

实现里使用 -expm1(-(end-u)/tau),避免极短时间差下 1-exp(x) 的浮点消减误差。

前面的一秒人工反例中,基准积分为 5,激励积分约为 43.237046,所以完整补偿项为 48.237046,完整 LL 才是 135.238876。如果只给出 183.475922,连这段会话的总暴露都没有算进去——只记了事件发生时的强度得分,没扣掉模型在整段观察时间里维持强度的账单。

真实批次直方图能告诉我们多大的情景差异?

114 条人工记录说明了单个大批次的计算机制。为了看看不同规模的并列批次会怎样累积评分差异,再把前面那份 2019-07-30 完整常规交易时段 AAPL execution 批次大小直方图用作输入,进行人工参数情景计算。这里使用的是已经保存的聚合结果,没有重新读取原始 ITCH 逐事件文件。

假设只考虑执行类型自身的正向激励,每个批次处理前的执行强度 lambda^- 不小于人为指定的 mu > 0,该批次包含 m 条记录。如果错误地在同一时间戳里按行增加自激,其事件项增量满足:

0 <= Delta(m) = sum_{i=0}^{m-1} log(1 + alpha*i/(tau*lambda^-))
   <= sum_{i=0}^{m-1} log(1 + alpha*i/(tau*mu))

右端是有条件的同类型情景上界,不是完整三类型 Hawkes 的误差上界。我们不知道真实批次前的强度,也没有把同时间戳的其他事件类型加入这项上界,因此它不能替代实际 AAPL likelihood 计算。

取 mu = 5 次/秒,扫一组人工 alpha 和 tau,将各个批次的上界乘以其真实出现次数后求和:

人工 alphatau=0.1 秒tau=1 秒tau=10 秒
0.18,392.651,436.32178.46
0.317,088.263,472.79500.41
0.625,186.195,837.61928.45

表中单位为 nats;每个数字是带真实批次大小分布的假设参数计算。例如中间情景 mu=5, alpha=0.3, tau=1 的累计条件上界为 3,472.79 nats,除以 49,172 条执行记录,得到 0.0706 nats/记录的条件上界。同一组参数下,一个典型的 2 条并列批次,上界仅约 0.0583 nats;3 条批次约 0.1716 nats。把人为强度下界从 5 改到 1 或 20,全体批次的累计上界则分别变为 11,044.42 和 1,125.70 nats。

如果只看表中的九组人工参数、固定 mu=5,计算出的条件情景上界落在 178.46–25,186.19 nats,除以相同的 49,172 条记录,对应 0.00363–0.51221 nats/记录。这只是九组人工参数下的条件上界范围,无法由此反推真实 AAPL 的似然偏差或已拟合的市场参数。实际批次前强度越高,同刻逐行自激带来的评分差异就越小。

对中间情景做“删去最大批次”的贡献核算:

排除的最大批次上界减少(nats)占原情景上界
187 条(top 1)321.549.26%
187、131、114 条(top 3)665.6219.17%
最大 10 个批次1,123.5732.35%

三个最大批次合计只有 432 条记录,却贡献了这个特定假设下上界的近五分之一,值得优先排查。现有直方图只存了批次大小,没留这些大批次发生在盘中哪个时刻。因此,现在还看不出它们是否扎堆在开盘或收盘,也判断不了是否与撮合或时间戳记录方式有关。

这次分析使用 Python 标准库完成,直方图计数与原聚合重新对账;新增 13/13 项测试通过,包含 114 条手算对照、参数单调性、top-k 删除恒等式与仓库实际目录布局下的复现。这 13 项与后文的 99 项合成研究测试分别运行,各有自己的执行记录。

实验到底验收了什么

真实 AAPL 汇总回答了“相同纳秒的执行记录是否存在、规模如何”;已知参数的合成实验则检验“批次评分、补偿积分与历史更新是否按约定实现”。两组证据解决不同的问题。

合成测试范围实际验证
批次评分(10 项)同纳秒统一评分、保留 multiplicity,输入行序不影响结果
完整似然(33 项)基准与指数核补偿积分、解析和数值积分交叉核验
M0/M1/M2 合同(34 项)齐次 Poisson、78 段季节性 NHPP、三类型 Hawkes 在已知参数下的评分与边界
合成事件流(22 项)固定种子的 Ogata thinning 模拟、独立历史求和对照、双进程确定性
合计99/99 项研究切片测试通过(Python 3.13.5 / 隔离 Linux x86_64)

其中一个专门制作的并列时间戳压力夹具共有 209 条合成事件、4 个时间戳批次:会话起点 90 条 Add,1 秒后 114 条 Execution,2 秒后 3 条 Cancel/Delete,3 秒后再出现同批次的 1 条 Add 和 1 条 Execution。四批分别为 90、114、3、2 条,合计 209 条。它与研究内另一条直接历史求和/积分计算路径的最大绝对差约为 1.42e-14;输入行顺序反转,结果不变。

另外两组用 Ogata thinning(薄化模拟法)生成的连续时间合成事件流分别有 285、265 条事件,独立似然对照的最大绝对差约为 7.11e-15。这些连续时间模拟并没有“自然产生”114 条精确并列记录。大批次是专门构造的压力输入,两者的来源清清楚楚。

这些结果支持已知参数下的计算一致性:我们使用另一条独立编写的直接历史求和/积分计算路径交叉核验,而非只让同一套算法重复计算。这 99 项测试在隔离 Linux 环境完成并保留了运行回执,研究提交尚无对应的 GitHub Actions 验证记录。至于真实交易日的参数拟合、模型选择与留出日评分,目前还没做;本文不提供“Hawkes 在真实市场上优于季节性 NHPP”的实证结论。

工作似然的边界

连续时间简单点过程通常不会以正概率生成完全相同的到达时刻。ITCH 的纳秒记录是观察/编码后的时间标记;消息有发送顺序,不等于我们知道每条记录的亚纳秒物理时间。因此,上文使用的是 pre-batch 约定下的工作似然(可称伪似然评分),并不是已推导出的时间戳量化观测数据的精确联合密度。也不能未经构造就随口称为 composite likelihood:复合似然有特定的组成密度与依赖处理要求。

如果只关心固定时间桶里的事件数量,可以研究离散时间计数模型;如果希望推断时间桶内无法观察的具体到达时刻,就需要考虑区间删失模型。小批次的潜在顺序比较则更适合作为敏感性分析。它们各自有观察假设,不能当作直接替换 pre-batch 的一行代码。

关于 log(m!),有一个可以直接说清楚的代数事实:同一数据、同一分组及同一观测模型定义下,如果两个候选模型的评分都遗漏了相同的、只依赖数据而不依赖模型参数的常数 C(D),那么 (LL_A + C(D)) - (LL_B + C(D)) = LL_A - LL_B。例如某些离散计数观测模型中的 -sum log(m!) 就可能是这样的项。它能在评分差值中抵消,不代表绝对 LL 不变,也不意味着只给 pre-batch 工作评分添个阶乘项就得到了严格的观测数据似然。换分组、换观察窗口或换时间戳约定时,不能假设还是同一个常数。

同样,即使两模型使用同一评分约定,也不能自动套用经典卡方似然比检验——非负激励参数的边界,以及零假设下可能不可识别的衰减参数,都可能破坏常规渐近条件。

分组时间残差该怎样检查

同时间戳的问题还会影响时间重标度诊断。如果把同类型、同纳秒的 114 条记录当成 114 个连续到达时间,相邻间隔会大量变为零。为此,我们在拟合前补充了诊断规则:按每种事件类型的唯一时间戳形成分组时间残差,同时单列每批 multiplicity。QQ 图、残差自相关和 KS 距离用于描述失配,不把普通连续无并列假设下的 KS p 值当作正式推断。如果将来使用 bootstrap 或置换检验,还需先定义零假设、生成机制和保留批次依赖的重采样单位。

这些检查比先挑一组漂亮的 Hawkes 参数更早决定结果是否可信。让文件中的记录顺序不再偷偷变成模型里的因果时间。

公开复现材料:Python 脚本、直方图、测试与 SHA-256 校验.