同じ1ナノ秒に114件の執行記録：Hawkes尤度はどこで計算を誤るか
NASDAQ ITCH の注文イベントを整理していて、普通の「1 行読むたびに Hawkes の状態を更新する」実装には、そのまま渡せない数字に出会った。114 だ。
2019 年 7 月 30 日の AAPL 通常取引時間帯では、ある 1 ナノ秒のタイムスタンプに execution メッセージが 114 行も記録されていた。131 行と 187 行のバッチもそれぞれ 1 回ある。さらに重要なのは、49,172 行の執行系メッセージのうち 30,080 行（61.17%）がほかの行と同じタイムスタンプを持つこと。同時刻の記録そのものは珍しくないが、100 行を超えるバッチは極端な尾部だ。
ここでいう「執行イベント」は本研究で E/C と定義したメッセージ行であり、114 件の別々の注文や、物理的に同時だった 114 約定を意味しない。
問題は行単位の状態更新だ。時刻差がすべてゼロでも、114 行目は前の 113 行から自己励起を受けてしまう。励起行列が非対称なら、同時刻のイベント種別を並べる順番でもスコアが変わり得る。同じ記録時刻の内部に、どんな前後関係を仮定してよいのか？
そこで推定前にルールを固定した。銘柄と整数ナノ秒でまとめ、バッチ直前の履歴で全件を採点してから、実件数をまとめて履歴に追加する。送信されたメッセージの順序を、モデル化せずにサブナノ秒の物理時間へ読み替えないためだ。
同時刻の記録は多いが、114 行は尾部の巨大バッチ
Nasdaq TotalView-ITCH 5.0 仕様はタイムスタンプを「深夜からのナノ秒数」と定める。E は Order Executed、C は Order Executed With Price。同じ注文参照番号に複数の部分執行メッセージがある場合もある。
本研究のイベント定義は add=A/F、cancel_delete=X/D、execution=E/C の 3 種類。C の Printable=Y と Printable=N は両方含める。表示可能な約定だけの集計ではない。P（非表示注文の Trade）と Q（Cross Trade）は含めない。Nasdaq 全取引の普遍的な分類ではなく、研究上の選択だ。
以下は 2019-07-30、AAPL、米国東部の取引所時刻 09:30–16:00 に対応する通常取引時間帯全体の内部集計である。
| 指標 | 観測値 |
|---|---|
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 の日付間比較とは時間範囲が異なる。保存したヒストグラムで分かるのはバッチサイズまでで、個別注文の識別子やサブナノ秒の実際の順序は復元できない。
全 29,379 タイムスタンプの平均は 1.67 行。10,287 件の重複バッチに限定すると平均 2.92 行、中央値と最頻値はともに 2 行。114 行は重複バッチ平均の約 39 倍で、114 行以上は 114／131／187 の計 3 バッチしかない。巨大バッチは不具合をあぶり出す例として有効だが、通常の規模ではない。
また、61.17% という重複率だけでパラメータが識別不能だとは言えない。pre-batch 採点ではバッチ内の未観測間隔を使わず、異なるバッチ間の時間差やイベント型の組み合わせが重要になる。
同じ時刻に、誰が誰を励起できるのか
3 種類の指数核 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 にある型 j の実行数で、行列の向きは対象 × 発生源。mu_k はイベント／秒、tau は秒、alpha[k,j] は次元を持たない積分励起量。発生直後の核係数は alpha[k,j]/tau（1／秒）となる。補償積分に alpha が直接現れるのは、指数核を時間積分した結果だ。
重要なのは u<t。現在のバッチを、現在のバッチ自身を採点する履歴に入れない。同一銘柄・同一時刻のイベント項は次の形になる。
event_term(t) = sum_k m_k(t) * log(lambda_k(t-))まず全種別の lambda_k(t-) を共通の過去から計算し、実件数の m_add、m_cancel_delete、m_execution を掛けて採点。最後に一度だけ履歴を更新する。114 行を 1 イベントに潰してはいけない。合成テストでは同一バッチ内の行順変更、銘柄間の隔離、取引日をまたぐ状態リセットを検査する。
送信順に逐次採点するのも別の研究仮説にはなり得る。ただし観測時間との対応を宣言しない限り、同じ尤度の単なる高速化とは呼べない。
114 行の人工例でスコア差を計算する
ここからは AAPL の推定結果ではなく、Python 標準ライブラリによる教材用反例。同一型のイベント 114 件がセッション開始時刻にあり、それ以前の励起はゼロ。mu=5 件／秒、alpha=0.6、tau=1 秒で、開始後の 1 秒を観測する。
import math
m, mu, alpha, tau, horizon = 114, 5.0, 0.6, 1.0, 1.0
# バッチの全イベントを開始前の強度で評価する
pre_batch = m * math.log(mu)
# 反例：同一時刻の前の行が後続の行を励起したと仮定
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}")隔離 Linux 上の Python 3.13.5 で実行した標準出力（終了コード 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 イベントが残る。ゼロ長バッチの後に続く 1 秒の補償積分もともに 48.237046 nats。完全スコアは逐次更新が 348.533109、pre-batch が 135.238876で、差は同じ 213.294233 nats だ。
これは実市場で予測が良くなったことを示さず、同じ連続時間観測モデルに対する二つの厳密な尤度を比較したものでもない。pre-batch で後続の影響が消えるわけでもない。1 秒後の自己励起は、この人工設定なら 114 × 0.6 × exp(-1) ≈ 25.163 件／秒となる。
イベント項だけでは完全尤度にならない
イベントが生じた時点の log(lambda) を足し上げるだけでは足りない。観測した有効セッション全体の強度積分（補償項）を差し引く必要がある。
LL = sum_k sum_batches m_k(t) * log(lambda_k(t-))
- sum_k integral_session lambda_k(t) dt本研究では推定前に 09:30–16:00 の 390 分間を 5 分 × 78 区間の季節性ベースラインに固定した。各区間と観測窓の半開区間の重なりを正確に積分し、過去のバッチは各型の実際の多重度を使って積分する。
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 秒の人工例ではベースライン積分が 5、励起積分が約 43.237046、合計の補償項は 48.237046。完全 pre-batch 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))右辺はあくまで同型の条件付きシナリオ上界だ。未知の実際の強度や同時刻の他型イベントを無視しているため、3 型 Hawkes の完全尤度誤差の上界にはならない。
mu=5 件／秒の人工設定で、保存されたバッチサイズの出現回数に上界を掛けて合計すると：
| 仮定した alpha | tau=0.1 秒 | tau=1 秒 | tau=10 秒 |
|---|---|---|---|
| 0.1 | 8,392.65 | 1,436.32 | 178.46 |
| 0.3 | 17,088.26 | 3,472.79 | 500.41 |
| 0.6 | 25,186.19 | 5,837.61 | 928.45 |
単位はすべて nats、alpha と tau は実市場で推定した値ではない。mu=5, alpha=0.3, tau=1 の中央シナリオは 3,472.79 nats、49,172 行で割ると 0.0706 nats／行の条件付き上界。サイズ 2／3 の典型的なバッチ単独では 0.0583／0.1716 nats。人工の下限 mu を 1／20 に変えると合計は 11,044.42／1,125.70 nats になる。
表の 9 条件で得られる上界範囲は 178.46–25,186.19 nats（0.00363–0.51221 nats／行）。実際の AAPL 尤度誤差や市場パラメータを逆算することはできない。
中央シナリオから最大のバッチを取り除いてみる。
| 除外した最大バッチ | 上界の減少量（nats） | 元の上界に占める割合 |
|---|---|---|
| 187 行（top 1） | 321.54 | 9.26% |
| 187／131／114 行（top 3） | 665.62 | 19.17% |
| 上位 10 バッチ | 1,123.57 | 32.35% |
上位 3 バッチはわずか 432 行だが、この仮定上界のほぼ 5 分の 1 を占める。ただし保存されたヒストグラムにはそのバッチの日中の発生時刻がない。寄り付き・引け間際への集中やマッチング機構との因果関係は評価できない。
この部分の人工計算・集計照合は別に 13/13 のテストを通している。後述の 99 件とは異なる実行だ。
テストで確かめた範囲
AAPL の集計は同時刻記録の存在と規模を確認する。一方、既知パラメータでの合成テストはバッチ採点・補償積分・履歴更新が仕様どおりかを検査する。証拠の種類が異なる。
| 合成テスト群 | 実際に検証した内容 |
|---|---|
| バッチ採点（10 件） | 同時刻の共通強度・多重度の保持・行順を変えても不変 |
| 完全尤度（33 件） | 基準強度と指数核の補償積分、解析解と数値積分の照合 |
| M0/M1/M2 契約（34 件） | 同次 Poisson、78 区間の季節性 NHPP、三種 Hawkes の既知パラメータ評価 |
| 合成イベント列（22 件） | 固定シードの Ogata thinning、独立履歴和、二プロセスでの再現性 |
| 合計 | 99/99 件 PASS（Python 3.13.5、隔離 Linux x86_64） |
同時刻バッチ専用の合成ストレス入力は 209 イベント、4 タイムスタンプ群。開始時に Add が 90 件、+1 秒で Execution が 114 件、+2 秒で Cancel/Delete が 3 件、+3 秒で Add と Execution が各 1 件。サイズは 90／114／3／2。独立実装の直接履歴和／積分経路との差は最大約 1.42e-14で、入力行を逆順にしても結果は同じだった。
Ogata thinning で生成した別の連続時間合成列には 285 件と 265 件があり、独立尤度の差は最大約 7.11e-15。これは自然に 114 件の厳密な同時到着を生んだわけではない。114 件のバッチは意図的なストレス入力である。
99 件のテストは Python 3.13.5／隔離 Linux x86_64 で行った過去の実行記録であり、現行の編集 HEAD に対する GitHub Actions CI 成功ではない。実市場でのパラメータ推定、アウト・オブ・サンプル評価、季節性 NHPP に対する優位性はまだ検証していない。
pre-batch は厳密密度ではなく作業尤度である
単純な連続時間点過程では、正の確率で完全な同時到着は通常起こらない。ITCH のナノ秒値は記録・符号化されたタイムスタンプであり、送信順が分かっても、それより細かな物理的到着時刻を観測したことにはならない。ここでの **pre-batch 作業尤度（擬似尤度スコア）**は観測上の取り決めで、丸められた時刻データの厳密な同時確率密度を導出したわけではない。composite likelihood とも、構成密度と依存構造を定義せずには呼べない。
時間桶ごとの件数が主目的なら離散時間カウントモデル、桶内の潜在到着順を推測するなら区間打ち切りの観測モデルを検討できる。小さな同時刻バッチの並べ替えは感度分析であり、万能な代替実装ではない。
同じデータ・同じグループ化・同じ観測モデルの二つのスコアが、共通のパラメータ非依存定数 C(D) を省略しているなら、(LL_A+C(D))-(LL_B+C(D))=LL_A-LL_B と差分から消える。離散カウントモデルでの -sum log(m!) がその例になる場合もある。しかし絶対 LL は変わり、階乗項を足すだけで pre-batch スコアが厳密な尤度になるわけではない。グループ化や時刻規則を変えた場合には共通定数の前提が崩れる。
また、励起パラメータの非負制約境界や、帰無仮説下での減衰パラメータの識別不能性により、古典的なカイ二乗尤度比検定をそのまま適用できるとは限らない。
同じタイムスタンプの残差をどう診断するか
114 行を 114 回の逐次到着とみなすと到着間隔にゼロが大量発生し、時間再スケーリング診断を歪める。本研究ではイベント型ごとの一意なタイムスタンプで残差をまとめ、バッチ多重度を別に表示する。QQ プロット、残差自己相関、KS 距離は不一致の記述に使うが、同時刻のない連続到着を仮定した KS p 値を正式な推論にはしない。bootstrap／置換検定を行うなら、帰無仮説・生成過程・依存性を保つ再標本化単位を先に決める必要がある。
まず、ファイル内の行順がひそかに因果的な時刻へ変わらないことを確かめる。