{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "# Waiting for the Second Trade"
      ],
      "id": "180d7d55-a9b9-4028-ac74-7a6ce216c3d0"
    },
    {
      "cell_type": "raw",
      "metadata": {
        "raw_mimetype": "tex"
      },
      "source": [
        "\\renewcommand{\\perp}{\\mathrel\\bot}\n",
        "\\newcommand{\\ev}{\\operatorname{E}}\n",
        "\\newcommand{\\var}{\\operatorname{V}}\n",
        "\\newcommand{\\cov}{\\operatorname{Cov}}\n",
        "\\newcommand{\\Normal}{\\mathcal{N}}\n",
        "\\newcommand{\\qedf}{\\qed{\\parfillskip=0pt \\par}}\n",
        "\\renewcommand{\\vec}[1]{\\mathbf{#1}}\n",
        "\\newcommand{\\gvec}[1]{\\pmb{#1}}\n",
        "\\newcommand{\\inner}[1]{\\langle{#1}\\rangle}\n",
        "\\newcommand{\\norm}[1]{\\lVert{#1}\\rVert}\n",
        "\\newcommand{\\prob}{\\operatorname{P}}\n",
        "\\newcommand{\\1}[1]{1\\kern-0.25em\\text{l}_{\\{#1\\}}}"
      ],
      "id": "5bc0f999-053d-4568-9cda-2e5145a8c38b"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Introduction\n",
        "\n",
        "If trades in a bond arrive as a Poisson process with rate $\\lambda$, the\n",
        "times between trades are independent exponential random variables, and\n",
        "the waiting time $W$ until the second trade, the sum of two of them, has\n",
        "a $\\text{Gamma}(2, \\lambda)$ distribution. This notebook takes the gamma\n",
        "model to one day of corporate bond transactions from Enhanced TRACE. We\n",
        "keep $\\alpha$ and $\\lambda$ the same for every bond, but we estimate the\n",
        "shape $\\alpha$ instead of fixing it at 2, so the estimate tells us how\n",
        "far the data are from the Poisson benchmark.\n",
        "\n",
        "The density of a $\\text{Gamma}(\\alpha, \\lambda)$ random variable is\n",
        "<span id=\"eq-gamma-density\">$$\n",
        "  f(w) = \\frac{\\lambda^\\alpha}{\\Gamma(\\alpha)} w^{\\alpha - 1} e^{-\\lambda w}, \\qquad w > 0,\n",
        " \\qquad(1)$$</span> with mean $\\operatorname{E}(W) = \\alpha / \\lambda$\n",
        "and variance $\\alpha / \\lambda^2$."
      ],
      "id": "d045c4e5-532d-47a7-a211-5cdd8093dd24"
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "metadata": {},
      "outputs": [],
      "source": [
        "import io\n",
        "import os\n",
        "from contextlib import redirect_stdout\n",
        "from pathlib import Path\n",
        "\n",
        "import matplotlib as mpl\n",
        "import matplotlib.pyplot as plt\n",
        "import numpy as np\n",
        "import pandas as pd\n",
        "from scipy.optimize import minimize\n",
        "from scipy.special import gammainc, gammaincc, gammaln\n",
        "\n",
        "# Make SVG figures byte-identical across reruns, so the frozen outputs\n",
        "# change in git only when a figure does\n",
        "os.environ.setdefault(\"SOURCE_DATE_EPOCH\", \"0\")\n",
        "mpl.rcParams[\"svg.hashsalt\"] = \"waiting-for-the-second-trade\"\n",
        "\n",
        "DATA_DIR = Path(\"../_data\")\n",
        "\n",
        "\n",
        "def load_wrds(name, sql):\n",
        "    file = DATA_DIR / f\"{name}.parquet\"\n",
        "    if not file.exists():\n",
        "        import wrds\n",
        "        with redirect_stdout(io.StringIO()):\n",
        "            conn = wrds.Connection()\n",
        "            df = conn.raw_sql(sql)\n",
        "            conn.close()\n",
        "        DATA_DIR.mkdir(exist_ok=True)\n",
        "        df.to_parquet(file)\n",
        "    return pd.read_parquet(file)"
      ],
      "id": "7c305f97"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Data\n",
        "\n",
        "We use all trades executed on Wednesday, October 15, 2025. The sample\n",
        "consists of corporate bonds with at least one customer trade in the\n",
        "month before that day, which keeps out the thousands of bonds\n",
        "outstanding that almost never trade. Bonds in the sample that do not\n",
        "trade on October 15 stay in, because they tell us how long the wait can\n",
        "be."
      ],
      "id": "6363f7c9-efec-4bfe-8835-1c41a0c6f4ef"
    },
    {
      "cell_type": "code",
      "execution_count": 2,
      "metadata": {},
      "outputs": [],
      "source": [
        "day = load_wrds(\"trace_2025-10-15\", \"\"\"\n",
        "    SELECT cusip_id, trd_exctn_tm, trc_st, msg_seq_nb, orig_msg_seq_nb,\n",
        "           entrd_vol_qt, rptd_pr, rpt_side_cd, cntra_mp_id, sub_prdct, wis_fl\n",
        "    FROM trace_enhanced.trace_enhanced\n",
        "    WHERE trd_exctn_dt = '2025-10-15'\n",
        "\"\"\")\n",
        "\n",
        "universe = load_wrds(\"trace_universe_2025-10-15\", \"\"\"\n",
        "    SELECT DISTINCT cusip_id\n",
        "    FROM trace_enhanced.trace_enhanced\n",
        "    WHERE trd_exctn_dt BETWEEN '2025-09-15' AND '2025-10-14'\n",
        "      AND cntra_mp_id = 'C'\n",
        "      AND trc_st = 'T'\n",
        "      AND sub_prdct = 'CORP'\n",
        "\"\"\")"
      ],
      "id": "0c92e575"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "TRACE records errors as separate rows that point back to the original\n",
        "trade report. We follow the filter of Dick-Nielsen (2014) for data\n",
        "reported after February 2012. A cancellation (`X`) or correction (`C`)\n",
        "repeats the original report with the same message sequence number, and\n",
        "both the original and the cancellation or correction are deleted. A\n",
        "correction also files a replacement report (`R`) with the correct terms,\n",
        "which stays in the data. A reversal (`Y`) points to the original through\n",
        "`orig_msg_seq_nb`, and both are deleted. Besides the sequence number, a\n",
        "report is matched to the original on CUSIP, volume, price, execution\n",
        "time, side and counterparty.\n",
        "\n",
        "We then keep customer trades only (`cntra_mp_id = 'C'`). A trade between\n",
        "two dealers is reported by both, and a dealer who buys from a customer\n",
        "often sells to another dealer seconds later, so counting interdealer\n",
        "trades would produce many artificially short waits. Time is measured in\n",
        "hours from 8:00 a.m. Eastern, and the day ends at 5:00 p.m., nine hours\n",
        "later."
      ],
      "id": "01ef9dba-d541-4e85-927e-1a542c4300e0"
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "metadata": {},
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Bonds in the sample:          31,603\n",
            "Customer trades on the day:   69,085\n",
            "Bonds with no trade:          65.7%\n",
            "Bonds with two or more:       25.4%"
          ]
        }
      ],
      "source": [
        "match = [\"cusip_id\", \"entrd_vol_qt\", \"rptd_pr\", \"trd_exctn_tm\", \"rpt_side_cd\",\n",
        "         \"cntra_mp_id\", \"msg_seq_nb\"]\n",
        "cancelled = day[day.trc_st.isin([\"X\", \"C\"])]\n",
        "reversed_ = day[day.trc_st == \"Y\"].drop(columns=\"msg_seq_nb\").rename(\n",
        "    columns={\"orig_msg_seq_nb\": \"msg_seq_nb\"})\n",
        "dead = pd.concat([cancelled[match], reversed_[match]]).drop_duplicates()\n",
        "\n",
        "trades = day[day.trc_st.isin([\"T\", \"R\"])].merge(dead, on=match, how=\"left\", indicator=True)\n",
        "trades = trades[trades._merge == \"left_only\"].drop(columns=\"_merge\")\n",
        "\n",
        "trades = trades[\n",
        "    (trades.cntra_mp_id == \"C\")\n",
        "    & (trades.sub_prdct == \"CORP\")\n",
        "    & (trades.wis_fl != \"Y\")\n",
        "    & trades.cusip_id.isin(universe.cusip_id)\n",
        "]\n",
        "\n",
        "hours = pd.to_timedelta(trades.trd_exctn_tm.astype(str)).dt.total_seconds() / 3600 - 8\n",
        "trades = trades.assign(t=hours)[(hours > 0) & (hours <= 9)]\n",
        "\n",
        "counts = trades.groupby(\"cusip_id\").size().reindex(universe.cusip_id, fill_value=0)\n",
        "print(f\"Bonds in the sample:          {len(counts):,}\")\n",
        "print(f\"Customer trades on the day:   {counts.sum():,}\")\n",
        "print(f\"Bonds with no trade:          {(counts == 0).mean():.1%}\")\n",
        "print(f\"Bonds with two or more:       {(counts >= 2).mean():.1%}\")"
      ],
      "id": "9f7557d7"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "For each bond that trades at least twice, the waiting time $w_i$ is the\n",
        "time of its second customer trade."
      ],
      "id": "856f1ec5-4450-410a-a81a-61ca4cf406ac"
    },
    {
      "cell_type": "code",
      "execution_count": 4,
      "metadata": {},
      "outputs": [],
      "source": [
        "w = trades.sort_values(\"t\").groupby(\"cusip_id\").t.nth(1).to_numpy()\n",
        "n_bonds = len(counts)"
      ],
      "id": "91ebf015"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Maximum likelihood with censoring\n",
        "\n",
        "We observe each bond until some time $T$. A bond whose second trade\n",
        "happens at $w_i \\leq T$ contributes the density\n",
        "(<a href=\"#eq-gamma-density\" class=\"quarto-xref\">1</a>) to the\n",
        "likelihood. A bond that has not traded twice by $T$ contributes the\n",
        "probability of that event, $\\operatorname{P}(W > T) = 1 - F(T)$. Writing\n",
        "$n_c$ for the number of censored bonds, the log-likelihood is\n",
        "<span id=\"eq-loglik\">$$\n",
        "  \\ell(\\alpha, \\lambda) = \\sum_{w_i \\leq T} \\Big( \\alpha \\ln \\lambda - \\ln \\Gamma(\\alpha) + (\\alpha - 1) \\ln w_i - \\lambda w_i \\Big) + n_c \\ln Q(\\alpha, \\lambda T),\n",
        " \\qquad(2)$$</span> where\n",
        "$Q(\\alpha, x) = \\Gamma(\\alpha, x) / \\Gamma(\\alpha)$ is the regularized\n",
        "upper incomplete gamma function, so that\n",
        "$1 - F(T) = Q(\\alpha, \\lambda T)$. Dropping the censored bonds would\n",
        "keep only the bonds that trade often and bias $\\lambda$ upward. Because\n",
        "every bond is observed over the same window, the censoring term is the\n",
        "only adjustment needed.\n",
        "\n",
        "Later we let the rate change during the day. The code is written for\n",
        "that case from the start, using the expected number of trades by time\n",
        "$t$, $\\Lambda(t) = \\int_0^t \\lambda(s) \\, ds$, which equals $\\lambda t$\n",
        "when the rate is constant. We take $\\lambda(t)$ to be constant between\n",
        "break points, so a single break at $\\tau$ gives\n",
        "$\\Lambda(t) = \\lambda_1 \\min(t, \\tau) + \\lambda_2 \\max(t - \\tau, 0)$. We\n",
        "maximize over $\\ln \\alpha$ and the log rates to keep all parameters\n",
        "positive, and compute standard errors from a numerical Hessian of the\n",
        "log-likelihood."
      ],
      "id": "1f70a31f-fbfc-4ad6-a169-91c0cae82367"
    },
    {
      "cell_type": "code",
      "execution_count": 5,
      "metadata": {},
      "outputs": [],
      "source": [
        "def clock(t, rates, breaks):\n",
        "    edges = np.concatenate([[0.0], breaks, [np.inf]])\n",
        "    widths = np.diff(edges)\n",
        "    return sum(r * np.clip(t - e, 0, d) for r, e, d in zip(rates, edges, widths))\n",
        "\n",
        "\n",
        "def intensity(t, rates, breaks):\n",
        "    return np.asarray(rates)[np.searchsorted(breaks, t)]\n",
        "\n",
        "\n",
        "def neg_loglik(theta, T, breaks):\n",
        "    alpha, rates = np.exp(theta[0]), np.exp(theta[1:])\n",
        "    obs = w[w <= T]\n",
        "    L = clock(obs, rates, breaks)\n",
        "    observed = (alpha - 1) * np.log(L) - L - gammaln(alpha) + np.log(intensity(obs, rates, breaks))\n",
        "    censored = (n_bonds - len(obs)) * np.log(gammaincc(alpha, clock(T, rates, breaks)))\n",
        "    return -(observed.sum() + censored)\n",
        "\n",
        "\n",
        "def hessian(f, x, h=1e-4):\n",
        "    k = len(x)\n",
        "    H = np.empty((k, k))\n",
        "    E = h * np.eye(k)\n",
        "    for i in range(k):\n",
        "        for j in range(k):\n",
        "            H[i, j] = (f(x + E[i] + E[j]) - f(x + E[i] - E[j])\n",
        "                       - f(x - E[i] + E[j]) + f(x - E[i] - E[j])) / (4 * h**2)\n",
        "    return H\n",
        "\n",
        "\n",
        "def fit(T, breaks=()):\n",
        "    breaks = np.asarray(breaks, dtype=float)\n",
        "    f = lambda theta: neg_loglik(theta, T, breaks)\n",
        "    res = minimize(f, np.zeros(len(breaks) + 2), method=\"Nelder-Mead\",\n",
        "                   options={\"xatol\": 1e-8, \"fatol\": 1e-8, \"maxiter\": 20000})\n",
        "    se = np.sqrt(np.diag(np.linalg.inv(hessian(f, res.x))))\n",
        "    lo, hi = np.exp(res.x[0] - 1.96 * se[0]), np.exp(res.x[0] + 1.96 * se[0])\n",
        "    return {\"alpha\": np.exp(res.x[0]), \"alpha 95% CI\": f\"[{lo:.2f}, {hi:.2f}]\",\n",
        "            \"rates\": np.exp(res.x[1:]), \"breaks\": breaks, \"loglik\": -res.fun}"
      ],
      "id": "ef469845"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## The opening hours\n",
        "\n",
        "Corporate bonds trade over the counter, so there is no opening bell, but\n",
        "almost no customer trades happen before 8:00 a.m. Eastern, an hour and a\n",
        "half before the stock market opens. We first fit the model to that\n",
        "window alone, so $T = 1.5$ and a bond is censored if it has not traded\n",
        "twice by 9:30."
      ],
      "id": "f4194613-30e6-4e23-bd54-a550292832cb"
    },
    {
      "cell_type": "code",
      "execution_count": 6,
      "metadata": {},
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "alpha = 1.66, 95% CI [1.51, 1.82]\n",
            "lambda = 0.074 per hour"
          ]
        }
      ],
      "source": [
        "opening = fit(T=1.5)\n",
        "print(f\"alpha = {opening['alpha']:.2f}, 95% CI {opening['alpha 95% CI']}\")\n",
        "print(f\"lambda = {opening['rates'][0]:.3f} per hour\")"
      ],
      "id": "e24318b4"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "The estimate of $\\alpha$ is closer to 2, the value implied by Poisson\n",
        "arrivals, than any estimate below, although the confidence interval\n",
        "excludes 2. Before the stock market opens, bond trades arrive roughly\n",
        "the way the Poisson model describes.\n",
        "\n",
        "## The whole day\n",
        "\n",
        "Extending the window to the close, $T = 9$, and keeping one rate for the\n",
        "whole day gives a very different shape."
      ],
      "id": "007dfc12-8e45-4c7b-8e14-ef1b45a286bc"
    },
    {
      "cell_type": "code",
      "execution_count": 7,
      "metadata": {},
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "alpha = 1.30, 95% CI [1.26, 1.33]\n",
            "lambda = 0.054 per hour"
          ]
        }
      ],
      "source": [
        "constant = fit(T=9)\n",
        "print(f\"alpha = {constant['alpha']:.2f}, 95% CI {constant['alpha 95% CI']}\")\n",
        "print(f\"lambda = {constant['rates'][0]:.3f} per hour\")"
      ],
      "id": "bb531491"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "<a href=\"#fig-constant\" class=\"quarto-xref\">Figure 1</a> compares the\n",
        "fitted distribution function with the empirical one. Since every bond is\n",
        "observed for the same nine hours, the empirical $F(t)$ is simply the\n",
        "fraction of all bonds in the sample whose second trade came by time $t$."
      ],
      "id": "8c738bd4-b561-479f-8949-6d10fd0daba9"
    },
    {
      "cell_type": "code",
      "execution_count": 8,
      "metadata": {},
      "outputs": [],
      "source": [
        "grid = np.linspace(0, 9, 400)\n",
        "empirical = np.searchsorted(np.sort(w), grid, side=\"right\") / n_bonds\n",
        "\n",
        "\n",
        "def plot_fit(model, label):\n",
        "    fitted = gammainc(model[\"alpha\"], clock(grid, model[\"rates\"], model[\"breaks\"]))\n",
        "    fig, ax = plt.subplots(figsize=(7, 4.5))\n",
        "    ax.plot(grid, empirical, color=\"black\", label=\"Empirical\")\n",
        "    ax.plot(grid, fitted, linestyle=\"--\", label=label)\n",
        "    ax.set_xticks(range(10), [f\"{h}:00\" for h in range(8, 13)] + [f\"{h}:00\" for h in range(1, 6)])\n",
        "    ax.set_xlabel(\"Time of day (Eastern)\")\n",
        "    ax.set_ylabel(r\"$P(W \\leq t)$\")\n",
        "    ax.set_xlim(0, 9)\n",
        "    ax.set_ylim(0)\n",
        "    ax.legend(frameon=False)\n",
        "    plt.tight_layout()\n",
        "    plt.show()"
      ],
      "id": "f998dd40"
    },
    {
      "cell_type": "code",
      "execution_count": 9,
      "metadata": {},
      "outputs": [
        {
          "output_type": "display_data",
          "metadata": {},
          "data": {}
        }
      ],
      "source": [
        "plot_fit(constant, rf\"Gamma, constant rate, $\\hat\\alpha$ = {constant['alpha']:.2f}\")"
      ],
      "id": "cell-fig-constant"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "With $\\hat\\alpha = 1.30$ the fitted curve rises too early and too\n",
        "slowly. The shape parameter governs how the curve starts, since\n",
        "$F(t) \\approx (\\lambda t)^\\alpha / \\Gamma(\\alpha + 1)$ for small $t$. A\n",
        "single rate has to match both the quiet opening and the busy late\n",
        "morning, and lowering $\\alpha$ is the only way the model can make the\n",
        "curve rise faster later in the day, at the cost of rising too fast at\n",
        "the start.\n",
        "\n",
        "## A rate that changes at 9:30\n",
        "\n",
        "The empirical curve bends when the stock market opens. The simplest way\n",
        "to capture this is to let the rate jump at $\\tau = 1.5$, so that the\n",
        "bond market has one rate before 9:30 and another after. With the clock\n",
        "$\\Lambda(t)$ in place of $\\lambda t$, the waiting time is still gamma\n",
        "distributed in the sense that\n",
        "$\\operatorname{P}(W \\leq t) = P(\\alpha, \\Lambda(t))$, where $P = 1 - Q$\n",
        "is the regularized lower incomplete gamma function. The model adds one\n",
        "parameter."
      ],
      "id": "83fe6ab5-ba4d-4c84-93bf-5d577fd1b358"
    },
    {
      "cell_type": "code",
      "execution_count": 10,
      "metadata": {},
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "alpha = 0.79, 95% CI [0.76, 0.81]\n",
            "lambda before 9:30 = 0.0033 per hour\n",
            "lambda after 9:30  = 0.0230 per hour\n",
            "Log-likelihood gain over a constant rate: 440.8"
          ]
        }
      ],
      "source": [
        "jump = fit(T=9, breaks=[1.5])\n",
        "print(f\"alpha = {jump['alpha']:.2f}, 95% CI {jump['alpha 95% CI']}\")\n",
        "print(f\"lambda before 9:30 = {jump['rates'][0]:.4f} per hour\")\n",
        "print(f\"lambda after 9:30  = {jump['rates'][1]:.4f} per hour\")\n",
        "print(f\"Log-likelihood gain over a constant rate: {jump['loglik'] - constant['loglik']:.1f}\")"
      ],
      "id": "4d81b0c9"
    },
    {
      "cell_type": "code",
      "execution_count": 11,
      "metadata": {},
      "outputs": [
        {
          "output_type": "display_data",
          "metadata": {},
          "data": {}
        }
      ],
      "source": [
        "plot_fit(jump, rf\"Gamma, rate jumps at 9:30, $\\hat\\alpha$ = {jump['alpha']:.2f}\")"
      ],
      "id": "cell-fig-jump"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "<a href=\"#fig-jump\" class=\"quarto-xref\">Figure 2</a> shows that the jump\n",
        "at 9:30 is enough to reproduce the shape of the empirical curve, and the\n",
        "log-likelihood rises by more than 400 for one extra parameter. The rate\n",
        "after the stock market opens is about seven times the rate before it.\n",
        "The levels of the rates cannot be compared with the opening-window\n",
        "estimate, because a different $\\alpha$ needs a different clock to\n",
        "produce the same curve; only their ratio is informative. The estimate of\n",
        "$\\alpha$, however, falls below 1, further from 2 than before.\n",
        "\n",
        "## What the estimates say\n",
        "\n",
        "<a href=\"#tbl-summary\" class=\"quarto-xref\">Table 1</a> collects the\n",
        "three fits. The log-likelihoods of the two whole-day models can be\n",
        "compared because they use the same data; the opening-window fit uses a\n",
        "shorter window and has no comparable log-likelihood."
      ],
      "id": "f9c99f89-29d2-450d-b53b-e8eb9d8abe0f"
    },
    {
      "cell_type": "code",
      "execution_count": 12,
      "metadata": {},
      "outputs": [],
      "source": [
        "models = {\n",
        "    \"Opening hours, 8:00 to 9:30\": opening,\n",
        "    \"Whole day, constant rate\": constant,\n",
        "    \"Whole day, rate jumps at 9:30\": jump,\n",
        "}\n",
        "pd.DataFrame(\n",
        "    {\n",
        "        \"alpha\": [f\"{m['alpha']:.2f}\" for m in models.values()],\n",
        "        \"95% CI\": [m[\"alpha 95% CI\"] for m in models.values()],\n",
        "        \"lambda before 9:30\": [f\"{m['rates'][0]:.4f}\" for m in models.values()],\n",
        "        \"lambda after 9:30\": [f\"{m['rates'][-1]:.4f}\" if m is not opening else \"\" for m in models.values()],\n",
        "        \"Log-likelihood\": [f\"{m['loglik']:,.1f}\" if m is not opening else \"\" for m in models.values()],\n",
        "    },\n",
        "    index=models.keys(),\n",
        ")"
      ],
      "id": "1d0caed7-32ec-47c9-a444-f8d5de5631e8"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Once the intraday pattern is removed, the remaining departure from the\n",
        "Poisson model is a shape parameter below 1. That is what a mixture of\n",
        "slow and fast bonds produces: a mixture of gamma distributions with\n",
        "different rates has more mass at both short and long waits than any\n",
        "single member of the family. The counts point the same way. With a\n",
        "common rate, the number of trades a bond makes during the day would be\n",
        "Poisson, with variance equal to its mean."
      ],
      "id": "4036aa73-ac97-48a0-99e5-1be75406e6d1"
    },
    {
      "cell_type": "code",
      "execution_count": 13,
      "metadata": {},
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Mean trades per bond: 2.19\n",
            "Variance:             194.31"
          ]
        }
      ],
      "source": [
        "print(f\"Mean trades per bond: {counts.mean():.2f}\")\n",
        "print(f\"Variance:             {counts.var():.2f}\")"
      ],
      "id": "ef891cda"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "The variance is about ninety times the mean: a few bonds trade hundreds\n",
        "of times a day while most never trade. Some of the low $\\hat\\alpha$ may\n",
        "also come from trades clustering in time within a bond, which the\n",
        "customer-trade filter weakens but does not remove. With one day of data\n",
        "and one rate for all bonds, the two explanations cannot be told apart.\n",
        "\n",
        "The natural next step is to let the rate depend on bond characteristics,\n",
        "$\\lambda_i(t) = \\lambda(t) \\exp(\\mathbf{x}_i' \\pmb{\\beta})$, with amount\n",
        "outstanding, age, maturity and rating in $\\mathbf{x}_i$, and to check\n",
        "whether $\\hat\\alpha$ moves back toward 2 once that heterogeneity is\n",
        "accounted for. Time variation within a bond, beyond the time of day,\n",
        "leads to the autoregressive conditional duration model of Engle and\n",
        "Russell (1998).\n",
        "\n",
        "Dick-Nielsen, Jens. 2014. *How to Clean Enhanced TRACE Data*. Working\n",
        "Paper. Copenhagen Business School.\n",
        "\n",
        "Engle, Robert F., and Jeffrey R. Russell. 1998. “Autoregressive\n",
        "Conditional Duration: A New Model for Irregularly Spaced Transaction\n",
        "Data.” *Econometrica* 66 (5): 1127–62."
      ],
      "id": "d2474a6e-352b-43a5-a0c9-9139ed180859"
    }
  ],
  "nbformat": 4,
  "nbformat_minor": 5,
  "metadata": {
    "kernelspec": {
      "name": "wrds",
      "display_name": "Python (wrds)",
      "language": "python",
      "path": "/home/lnaranjo/.local/share/jupyter/kernels/wrds"
    },
    "language_info": {
      "name": "python",
      "codemirror_mode": {
        "name": "ipython",
        "version": "3"
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.13.16"
    }
  }
}