{ "cells": [ { "cell_type": "markdown", "id": "9656d9fa", "metadata": {}, "source": [ "# Classical confidence intervals in one call\n", "\n", "`conf_int` turns a bootstrap run into a confidence interval for a statistic of one\n", "series. You pass the series, the statistic, and a method spec; it runs the bootstrap\n", "and reads back the interval you asked for. This tour builds percentile, basic, and\n", "studentized intervals on a dependent AR(1) series, shows why BCa is refused for block\n", "methods (and why that refusal protects you), and finishes with per-series intervals\n", "over a ragged panel.\n", "\n", "The four interval families, shortest to say:\n", "\n", "- `percentile` and `basic` are first-order correct and need only the replicate\n", " distribution.\n", "- `studentized` is second-order correct for smooth statistics when its per-replicate\n", " standard error is dependence-aware, which `conf_int` arranges with a block jackknife.\n", "- `bca` corrects for bias and skewness, but its acceleration is defined under\n", " independent sampling, so it is offered for the `IID` spec only.\n", "\n", "Coverage here is approximate or asymptotic under temporal dependence, not finite-sample\n", "distribution-free, the same honest ceiling as the rest of the uq layer." ] }, { "cell_type": "code", "execution_count": 1, "id": "aabd4250", "metadata": { "execution": { "iopub.execute_input": "2026-07-02T19:33:04.024845Z", "iopub.status.busy": "2026-07-02T19:33:04.024536Z", "iopub.status.idle": "2026-07-02T19:33:05.321049Z", "shell.execute_reply": "2026-07-02T19:33:05.318641Z" } }, "outputs": [], "source": [ "# On Colab or Binder, install tsbootstrap first (skipped if already present):\n", "try:\n", " import tsbootstrap # noqa: F401\n", "except ImportError:\n", " %pip install -q \"tsbootstrap[examples]\"" ] }, { "cell_type": "markdown", "id": "05740f01", "metadata": {}, "source": [ "## A dependent series\n", "\n", "We use an AR(1) series with positive autocorrelation and a true mean of zero. The\n", "dependence is what makes the interval choice matter: a naive i.i.d. interval would\n", "understate the variance of the sample mean, so it would be too narrow." ] }, { "cell_type": "code", "execution_count": 2, "id": "153d512f", "metadata": { "execution": { "iopub.execute_input": "2026-07-02T19:33:05.325499Z", "iopub.status.busy": "2026-07-02T19:33:05.325174Z", "iopub.status.idle": "2026-07-02T19:33:05.333410Z", "shell.execute_reply": "2026-07-02T19:33:05.332030Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "series length 200, sample mean 0.030\n" ] } ], "source": [ "import numpy as np\n", "\n", "\n", "def ar1(phi, n, seed):\n", " \"\"\"AR(1) started at its first innovation; true mean 0 for zero drift.\"\"\"\n", " rng = np.random.default_rng(seed)\n", " e = rng.standard_normal(n)\n", " x = np.empty(n)\n", " x[0] = e[0]\n", " for t in range(1, n):\n", " x[t] = phi * x[t - 1] + e[t]\n", " return x\n", "\n", "\n", "x = ar1(0.5, 200, seed=0)\n", "print(f\"series length {x.size}, sample mean {x.mean():.3f}\")" ] }, { "cell_type": "markdown", "id": "75dd6a09", "metadata": {}, "source": [ "## Percentile, basic, and studentized, one call each\n", "\n", "Each interval is a single `conf_int` call. We use a `MovingBlock` method so the\n", "resampling respects the AR(1) dependence, a small `n_bootstraps` to keep the tutorial\n", "quick, and the built-in `\"mean\"` statistic. The studentized path estimates a\n", "dependence-aware standard error with a block jackknife on its own; you do not have to\n", "supply one." ] }, { "cell_type": "code", "execution_count": 3, "id": "38188cb5", "metadata": { "execution": { "iopub.execute_input": "2026-07-02T19:33:05.336519Z", "iopub.status.busy": "2026-07-02T19:33:05.336328Z", "iopub.status.idle": "2026-07-02T19:33:05.703798Z", "shell.execute_reply": "2026-07-02T19:33:05.702022Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "point estimate (sample mean): 0.030\n" ] } ], "source": [ "from tsbootstrap import MovingBlock, conf_int\n", "\n", "ALPHA = 0.1 # 90% intervals\n", "B = 399\n", "method = MovingBlock(block_length=\"auto\")\n", "\n", "intervals = {}\n", "point = None\n", "for kind in (\"percentile\", \"basic\", \"studentized\"):\n", " lower, upper, point = conf_int(\n", " x, \"mean\", method=method, kind=kind, alpha=ALPHA, n_bootstraps=B, random_state=0\n", " )\n", " intervals[kind] = (float(lower), float(upper))\n", "\n", "print(f\"point estimate (sample mean): {float(point):.3f}\")" ] }, { "cell_type": "code", "execution_count": 4, "id": "96be6c38", "metadata": { "execution": { "iopub.execute_input": "2026-07-02T19:33:05.706933Z", "iopub.status.busy": "2026-07-02T19:33:05.706536Z", "iopub.status.idle": "2026-07-02T19:33:05.917057Z", "shell.execute_reply": "2026-07-02T19:33:05.916266Z" } }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
lowerupperwidthcovers 0
percentile-0.1908650.3014420.492307yes
basic-0.2417490.2505580.492307yes
studentized-0.242680.2831170.525797yes
\n", "
" ], "text/plain": [ " lower upper width covers 0\n", "percentile -0.190865 0.301442 0.492307 yes\n", "basic -0.241749 0.250558 0.492307 yes\n", "studentized -0.24268 0.283117 0.525797 yes" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import pandas as pd\n", "\n", "summary = pd.DataFrame(\n", " {\n", " kind: {\n", " \"lower\": lo,\n", " \"upper\": hi,\n", " \"width\": hi - lo,\n", " \"covers 0\": \"yes\" if lo <= 0.0 <= hi else \"no\",\n", " }\n", " for kind, (lo, hi) in intervals.items()\n", " }\n", ").T\n", "summary.round(4)" ] }, { "cell_type": "markdown", "id": "b32aaaf1", "metadata": {}, "source": [ "All three intervals bracket the true mean of zero. The percentile and basic bands have\n", "the same width because basic is the percentile band reflected through the point\n", "estimate. The studentized band is a little wider: it accounts for the extra uncertainty\n", "in the per-replicate standard error, which is the cost of its second-order correctness\n", "under dependence." ] }, { "cell_type": "markdown", "id": "05a5a875", "metadata": {}, "source": [ "## A picture of the three intervals\n", "\n", "The same three intervals drawn as error bars around the point estimate, with the true\n", "mean marked." ] }, { "cell_type": "code", "execution_count": 5, "id": "a48313ab", "metadata": { "execution": { "iopub.execute_input": "2026-07-02T19:33:05.919364Z", "iopub.status.busy": "2026-07-02T19:33:05.918951Z", "iopub.status.idle": "2026-07-02T19:33:07.243113Z", "shell.execute_reply": "2026-07-02T19:33:07.242100Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAscAAAFLCAYAAAA6bCEQAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjExLjAsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvlcelbwAAAAlwSFlzAAAPYQAAD2EBqD+naQAAUU1JREFUeJzt3XlcVFX/B/DPsO8goAgCbiiIC4oK7opLmmuS+5amZmplVmrWz8esTMsy9w3NNZdEH01zKTUXVBTFFREVRQVExA0ERJg5vz8Y7uPIgDPAMDB+3q+Xr5oz597zvedeZ76eOfdcmRBCgIiIiIiIYKTvAIiIiIiIygomx0RERERESkyOiYiIiIiUmBwTERERESkxOSYiIiIiUmJyTERERESkxOSYiIiIiEiJyTERERERkRKTYyIDkJSUhGnTpmHgwIHo06cPAODbb7/FlClTNNp+xYoVeO+993QZYqmYPXs2Jk2apO8wSk1JH29CQgKmTZuGQYMGSdeRPkybNg3Tpk3TW/slbfXq1fjqq6+KvP1XX32F1atXl2BE5Vdx+/JNwWumeEz0HQDR6yQnJ2P16tW4fPkyzM3N0aZNGwwaNAgmJvkvX4VCgU2bNuHff/9FTk4OGjdujJEjR8LKykql3osXL7B48WJcunQJDRs2xLhx4/Lt7+rVq5g+fTqWLFkCJycnnR5jcWRkZKBJkyZo0KABhgwZAjMzMwDA0aNH8eTJE432ERkZiZ07d+owytIRFhaG+Pj4Im//7bffIj09HT/++GMJRqU7xT3el6WmpsLf3x+NGzfG0KFDYWpqWiL7LUhhff3vv//qtG11Vq5ciX379qFLly4YNWqU2joLFizA0aNHpddmZmZwd3dHjx490Lp1a7Xb3L9/HxMmTMDy5cvzvZeZmYm9e/di165dSEtLw08//YQaNWrkq1evXj18+OGH6Nq1K1xcXIp4hOVfQX356nkxMjJC5cqV0ahRIwwaNAjm5ualHSpiYmKwfft2XLt2DdbW1qhbty6GDh0KGxubfHW1+Y7TtC6vmWISRGVYRESEcHR0FI0bNxZLly4Vc+bMEZUrVxatWrUS6enpKnVfvHghunTpIlxcXMTcuXPFihUrRN26dYWvr69ISUlRqduzZ0/RuHFjsWLFCuHn5yf69eun8r5CoRAtW7YUn332mc6PsbhCQ0MFABEdHa1SfvToUXHw4EGN9jFmzBhhb2+vg+hKV1hYmPjnn3+KvH2HDh1E48aNSzAi3erWrZvw8/MrkX1t3rxZ7XWkK4X1dcuWLUXLli1LJQ4hhMjOzhZubm7C2NhYuLm5iZycHLX1Bg8eLACIjRs3iq1bt4pVq1aJd999VwAQkydPVrvNhAkTRLVq1YRcLlcpDw0NFRUrVhTBwcGiV69eAoCIiIhQuw+5XC6qVq0qJkyYUKzjLO8K6stXz8vmzZvFDz/8IFxdXUX16tXFvXv3SjXOr7/+WjRp0kRMmzZNrF27Vvz000+iZs2awsnJSURGRqrU1eY7Tpu6vGaKh8kxlWl16tQRnp6eIjMzUyqLjo4WxsbGYurUqSp1582bJwCI8PBwqSw5OVnY2dmJkSNHSmWXL18WAMTly5eFEEKcPXtWABAxMTFSnYULF4oaNWrk+8Api3755RcBQDx48KDI+zCU5Li43uTkeM6cOcW+jrRRlpLjHTt2CADSZ8jOnTvV1stLwl7+PBJCiM6dOwuZTCZu3rypUp6eni7s7e3FtGnT8u3r7t27Ii0tTQghxKxZswpNjoUQYurUqcLBwaFcfCbpQmF9WdB5OXz4sACgdhtdio+Pz1eWlJQkLCwsROfOnVXKtfmO06auELxmioPTKqjMun79OqKjo/H555/DwsJCKvfx8UFgYCBCQkLw3XffwdjYGEDuz6J+fn4IDAyU6lasWBHBwcHYuHEjFixYACsrK1y9ehUmJibw9fUFAPj5+cHIyAjR0dGoXbs27ty5g6+++gr//e9/803HeJ2MjAxs3LgRp0+fRnZ2Nho1aoSRI0fC2tpaqpOUlITVq1fjypUrMDU1RbNmzTB06FBYWlpKdfLmW86YMQPr1q3D4cOHYW1tjSFDhqB58+ZSvffffx8RERHS/+dNqchTs2bNfD9bX7t2DatWrUJiYiLq1auHsWPHFng86enp2LBhA06fPo2srCzUq1cPo0aNgrOzs9axatNHmrSrzuzZs/Hw4UPMmTNH6/g+/vhjXLp0CdnZ2dJ8W0tLS6xfv75I/fHNN99g3bp1OHr0KOrVqweFQoGoqCisXr0aMplMJe4tW7YgNDQUy5cvh6OjI5YtW4YDBw4AyP2JuEKFCmjRogUGDhyY7xy/KjMzExs3bsTZs2eRmZkJb29vDBs2DG5ubgVuo+46euutt/DBBx8A0O6affW4P//883ztadLXQO40KU2uqaJeL3lWrlwJf39/TJgwAevWrcPKlSvRs2dPjbYFgA4dOmD//v2IjIxE9erVpfK9e/fi6dOnavfl7u6u8f4BoGfPnpg1axb27t2Ld999t9C6eedi+vTp+O2333D8+HG4urpiwoQJcHV1RXZ2NtauXYuwsDA4ODhg3LhxqF27dr79aNKvml6r2n5OvKqwviyIt7c3AODhw4cab1MSqlSpkq/MxcUFbm5uuH//vlSmzXectt+HgHbXDKlickxl1uPHjwFA7XxfJycnpKSk4Pr16/Dx8cGzZ88QFRWFESNG5KvbpEkTrFmzBhcuXEDz5s1haWkJuVyO7OxsmJmZISsrCwqFQkrOxowZg379+qFDhw5axRsTE4POnTsDAD744AO4ubnh/PnzaNKkCaKjowEAJ0+exNtvv43atWtj+PDhePbsGWbMmIF58+bh33//leaG5c25vH//Puzt7dG2bVts27YNrVu3xsGDB9G2bVsAQO/evfHixQtcvnwZvXv3hq2trRTPt99+i7i4OJUYDx8+jK5du6Jp06YYNGgQHjx4gF69eqFq1ar5juf69et46623YGpqijFjxsDOzg4bN27Er7/+isOHD6NOnTpaxappH2narjrq5uBqGl+3bt1w5MgRpKenY8CAAQCgMo9P2/64desWXFxcEBgYiJs3b6JBgwaYPHkyhg4dmu/a+v7772FjYwNHR0cAQNOmTaUEJDs7G9evX8fUqVOxevVqHDp0CEZG6u+lzsrKQps2bZCamooPP/wQzs7OiImJQVBQEHbu3AkfHx+126m7jmrWrAlA+2v21eNW53V9nefDDz987TVVnOsFyL0Jce/evdI81rFjx+LDDz9EQkKC2iRHnWfPngFAvvmkhw8fhrm5Ofz8/DTaT2H8/f1hbm6Of//997WJTt65uHPnDtzd3dG6dWusWLECmzdvRmRkJD788EPUrFkTrVu3RkhICNatW4fo6GiVuama9qum16o2nxPqFKUvt2zZAgBo167da+tu3rwZoaGhGu136dKlqFixosZxALnxx8XFYcaMGVKZNt9x2tTNo801Q6/Q99A1UUGSkpIEAPHee+/le6927doCgDS/NCYmRgBQ+9PS1q1bBQCxdetWIYQQ9+/fFxYWFiI0NFQIIcTvv/8urK2tRUpKili3bp1wdXUVjx8/FufPnxeffvqpGD9+vDh58mShsSoUCtGgQQNRtWpV8fTp03zHIYQQOTk5onr16sLX11c8f/5cev/WrVvCwsJCDBgwQCpr2bKlsLS0FJs2bZLK5HK5qFmzpujUqZPK/vN+kn315/BXf7bOyckRVatWFQEBASI7O1sqP3XqlLC0tMw3rcLf31/UqFFDPHv2TCWGgIAA0aJFC61j1aSPtGlXHXXTDLTpy8J+6temPywsLMSaNWuksszMTJGZmSkcHBzEoEGDVPZ76tQpAUCEhIQUemx504HyrmN1x3vkyBEBQISFhalsm5GRka/PX6XuOtL2mlV33AV53bQKTc9Zca4XIYT49ttvVX56Tk9PFw4ODuK7777LV1fdz/f3798Xnp6ewtnZWSWGvGP08vJ6bQyaTKsQQogaNWqIdu3avXZ/eefi5WslKSlJmJiYiCZNmqj0a3JysjA1Nc039aA4/aruWtXmnKpTWF/mnZfevXuLd999V7zzzjuiYcOGokKFCmL58uWv3bcQufOEAWj059atWxrtc+zYsSI4OFgEBgYKR0dHsWTJEqFQKKT3tfmO06buyzS9ZkgVR46pzHJxcUH37t2xZcsWfPLJJ/D39wcAhISE4Nq1awByR8qA3NUnAPWjT3l33T9//hwAUKlSJcyZMwcjRozA0qVLcfr0afz666+Qy+WYOHEiQkJCcOfOHbRo0QIff/wxHB0dERQUhAMHDqBly5ZqYz179iwuXryIuXPnws7OLt9xAMDp06dx69YtLFy4UOXu6WrVqiE4OBjbtm1DVlaW9J61tTX69+8v1TMyMkL79u2xY8cO7TpS6fTp07h9+zamT5+u0k8BAQHw8/OTRm4B4Ny5c4iMjMQvv/yiMt3ByMgII0aMwNixY3Hv3j24urpqHKsmfaRtu5oqbl9qG5e5uTkGDx4s1cv7GXTQoEH47bff8OTJEzg4OAAAVq1alS8+ADhw4AD27duHxMRE6fo2MjLCxYsXC1xmLa+dP//8E02aNJGuJUtLS5UpEJrS9pot6LiLQpNzVtzrRQiB3377DcOHD5emUFlZWeG9997Db7/9hq+//jrfFBgg9zwaGRkhNTUVp0+fRs2aNbFixQqVGAAgJSVF+jWgJDg5OeHBgwca1bW1tVUZLXRxcYG3tzdu3ryp0q8VK1aEr68vzp07J5Vp26+aXqvF+XuoSV/27dsXpqamEELg3r17WL16NX799Ve0bdtWmmJRkIEDB6Jhw4avjQPI/Q7RRI8ePZCamoq4uDisXLkSc+fORYsWLaTRb22+47Sp+zJtrhn6HybHVKatXr0aAwcORIsWLdC8eXM8e/YMiYmJGDRoEDZu3Cj9xJT3c2ZGRka+faSnpwOAypSDjz76CD179kR0dDTq1q0Ld3d3DBgwAEFBQejduzdGjBiBNm3aYPbs2QCA6Oho/Pjjj/jzzz/Vxnn9+nUAKPQn3LwpDrVq1cr3Xu3atZGVlYXExERpzmL16tXzfTFXrFgRDx48gEKhKPCn9YLcunULAKSfy1/m5eWlkhzHxMQAAHbs2IHTp09D5N68CyD3Z2gg9yfbvC9HTWLVpI+0bVdTxe1LbeOqWrWq2n+ojRo1CkuWLMHGjRsxbtw4ZGRkYPPmzejXr5/K9Tlq1Chs2LABo0ePRseOHWFrawuZTIYdO3YgNTW1wDgDAgLwySefYM6cOViyZAlatWqF9u3bY+DAgVrPcQW0v2YLOu6i0OScFfd6+eeffxAXF4czZ86oJHFJSUm4desWDhw4gE6dOuXbLi8JS0lJwaNHjxAfH682ic6btlVSsrKyNP4HR7Vq1fLF5OTkBDMzM7XlL8+F1aZftblWi/P3UJO+7N27t0r/jBw5ErVq1cLgwYNx5syZQretW7cu6tatW2gdbb399tvS/48dOxYNGjRAv379pP4FNP+O07ZuHm2uGfofJsdUpjk7O+Off/5BTEwMrl27BhsbGwQGBuL//u//YGpqivr16wMAPDw8YG5ujtu3b+fbR17Zq1/wnp6e8PT0BADs2rULf//9t5QgXr16VWUOXKNGjbBixYoC48wbOctLxNXJuzklMzMz33t52708Oqfuxqu8L5aiJMd5+1P3D4hXy/JG21u1aiWNULxs4sSJ8PLy0ipWTfpI23Y1Vdy+1DaulxPdlzVq1Aj+/v5YtWoVxo0bh9DQUKSmpqqsq3v16lWsWrUKc+bMwRdffCGVp6SkQC6XFxonAMyfPx9ffvkl/vnnH4SFheHnn3/G9OnTsWfPHo3mXr5M22u2oOMuCk3OWXGvl5CQELRr1w7jx4/P956FhQVCQkLUJscvJ2HDhg1DkyZN0KtXL1y5ckWlDypXrqwyIltcDx8+1Hh0s6D+K6hcoVBIrzXtV22v1eL8PSxKX1pbW6N58+bYvn07Hj58WOh69bqec2xnZ4cuXbpg+fLlSEpKQuXKlQFo/h2nbd082lwz9D9Mjqlc8Pb2ln4Wy8nJwX//+18EBwdLP/kZGxujXbt2CAsLg1wuV7lj99ChQ6hSpUqBI5apqakYO3Ys5s6dK/28b2lpqZIQZGZm5vvJ9GUBAQEwMjLCsWPHCrzxIe8D6tSpU3jnnXdU3jt16pR0N7OuNGrUCABw5swZdOnSRSpXKBSIjIxUqRsQEACZTAZTU9MSe1KaJn2ki3a1YWJiopIk6CKukSNHYvz48bhw4QJWrVoFHx8ftGjRQno/b2Tu1S+6Q4cOadyGq6srhg0bhmHDhuHHH3+Ep6cnli9frnVyrMtrtqC+1kZxzsuDBw/w559/YsmSJWq3zcjIwOjRo/HgwYNCEyErKyvMmzcPnTt3xqxZs/DDDz9I7zVt2hS7d+/G48ePUaFCBa3ie1VKSgoSEhKkFUR0SdN+PX/+PIDiXauaKmpfpqWlAYA08l2Qy5cvY9u2bRrt8+eff9Y6OQZyz6GJiYnaB4G87juuKHVL85oxNHx8NJVpYWFhuHHjhvT6xYsX+Oijj6SnSb1s8uTJuH//PubOnSuV/f333zh48CC+/PJLtT97AsCkSZPg6+uL4cOHS2WtW7fG/v37kZWVBblcjl27dhX4BCwgd1mmESNGYNmyZdi3b59UrlAopEd41qxZEz169MCyZctw4cIFqc7vv/+Oo0ePYuLEiZp1ShHVrFkTb7/9NhYuXCjNUQOAOXPmqPxjAsgdiR89ejR+/vnnfE8re/jwIZYsWaJ1+5r0kS7a1YaHhwfi4+PzjXqVZFyDBw+GpaUlvvrqKxw7dgwjR45Ueb9u3bowNTVVeWLhnTt3EBIS8tqn1p04cUJaVitPZmYm5HK5xsuavUyX12xBfa3tPop6XtauXYsXL16o/EPxZZ07d0Z2djbWrVv32jjeeusttG3bFvPnz1eZntCxY0cIIXDy5EkNj6hgx48fBwC1I9klTdN+Lc61qq2i9OXp06dx+PBhtGnT5rXX/8CBA7F161aN/hQ251ihUGDx4sX5poDs3LkTO3fuxIABA1SSY22+47SpC5TuNWNoOHJMZZqVlRX69OkDMzMzVKpUCWfOnIGnpyeOHj0qTYnI0759eyxbtgwTJ07Etm3bYGNjg2PHjuHzzz/HRx99pHb/hw8fxsaNG3Hp0iWV8okTJ2LPnj1o2LAhzMzMIJfL8fXXXxca6+LFi2FmZoaePXvCx8cHbm5uuHLlCt5++21pibm1a9di6NChCAwMRGBgINLT03Hx4kV88cUXmDRpUjF6SjNr167FO++8g0aNGiEwMBAPHz5Ely5d0LFjR2zevFml7qJFi2BnZ4euXbvCy8sLnp6eSExMxMOHDzFu3Lgita9JH+miXU2NHTsWGzduhJ+fH7y9vWFlZSWtvVtScdnb2+Pdd9/Fhg0bYGpqimHDhqm8X7lyZcyfPx8TJ07E8ePH4ezsjPj4ePz+++9o1apVofuuUKECZsyYgREjRsDHxwcymQzh4eFo27YtvvnmG637A9DdNVtYX2ujqOdl5cqVqFevXoHLtbm4uMDf3x8rV65Uu1bzq2bOnIlWrVrh+++/x8KFCwEAzZs3h6+vLzZv3oyuXbuq1L958yYmT54M4H9zfKdMmYIKFSrA2toaa9euVam/ZcsW+Pr6arQmcEnQpF+Lc61qq7C+zJN3o2TeDXkRERHo2rUrli5d+tr9l9ScY5lMhtjYWNSuXRtubm5wdHREbGws4uLiMHLkSMybN0+lvjbfcdrUBUr/mjEkMvG63xqI9Ewul+PChQu4d+8evLy8XnvXcWpqKk6dOoWcnBw0bNiw0Ju3wsPDYWRkhICAgHzvZWdn48yZM8jJyUHTpk01vqkhJSUFkZGRkMlkqF+/vjS37GU3b96UHqjQuHHjfKMahw8fBpB/fc6oqChER0fj3XfflUbCY2JicOnSJfTs2VNlTt+xY8eQnZ2N9u3bq+xDCIHIyEjcu3cPvr6+qFGjBs6dO4e7d++qXWA/LS0NkZGRSEtLg6enJ3x8fFTa0SZWbfrode2qc/z4cWRmZqJjx45Fju/x48c4c+YMUlNTYWxsnG86QVH742W3b99GREQEHBwcVGJ9WVJSEi5evAhzc3M0a9YM5ubm+O9//4saNWpId7urO14AuHfvHqKioiCTyVCrVi21X5yvKug6ylPUa7YwBfV1Ua4pba6XzMxM/PXXX6hatSqaNm1aYHxnz57FrVu30L17d1hYWCAiIgK3b99GcHCw2vmxe/bsgVwuR48ePaSylStX4tNPP8WdO3dUVlt4/PgxDh48qLbdvH9A5nn06BE8PDzw66+/avQTeUH9d+TIESgUCgQFBamUF/RZAWjWr5pcq0U5p68qqC/zzksemUwGZ2dneHl5abxWdUnLyMhAVFQUEhMT4ejoCD8/v3yr9OTR5jtO07raXjOkiskxERGRjsjlcvj5+aFnz54q85G1MXXqVOzatQsXLlzINwXqTVISffmm4DVTPJxWQUREpCPGxsbYtm2btJRiUbRp0wbDhw9/45OckujLNwWvmeLhyDERERERkRJXqyAiIiIiUmJyTERERESkxOSYiIiIiEiJyTERERERkRJXqyiDFAoFEhMTYWtr+9p1H4mIiIjo9YQQSEtLg5ubm9q1yvMwOS6DEhMT4eHhoe8wiIiIiAzO3bt34e7uXuD7TI7LIFtbWwC5J6+gJ+oQEb0sIzYCAGBVs+AnvhERvclSU1Ph4eEh5VkFYXJcBuVNpbCzs2NyTEQaSVjUFQBQeVWWniMhIirbXjdllTfkEREREREpMTkmIiIiIlIqk9MqVq1aBRMTE7z33ns6b2vx4sWoUKECBg0apLM2SvN4iIiI3iRCCOTk5CAnJ0ffoVAZYGZmBmNj42Lto8SS44ULF6JSpUro379/sff177//wsLCosSTSXUx/vPPP3B3d9dpcqyr4yEi0ie5XI5jx47h3r17cHV1RevWrYv9pUSkjaysLMTHxyM9PV3foVAZIZPJUK1atdfedFeYEkuO9+/fDy8vrxJJjnVFXYwfffQRrK2t9RgVEVH5s337dkyYMAHx8fFSmbu7O+bPn4/g4GA9RkZvCoVCgevXr8PY2Bienp4wMzPjswHecAqFAsnJyYiLi4Ovr2+R/7GuVXL87NkzrF+/HlevXoWbmxsGDhwIT09PrF69GhcvXkRcXBxGjRoFAPj++++xYsUKNGvWDG+99Za0jxUrVsDc3FxlFPXkyZPYunUrbG1t0aVLF7Vtx8fHY8OGDYiPj0eNGjUwdOhQVKxYUXp/7ty5qFatGhwdHXHgwAHI5XL07dsX/v7+AFBgjDExMahQoQKaN2+OGzduYPbs2fna9vb2xqRJkzSKQ9PjKSvkCoFHGS/0HQYRFZexGQDgwTPdr1axe+cOjBw6EEIIlfKEhAT06dMHq9ZvQvde7+g8DnqzZWc9R45cDs+q1TjIpQdGstev+qAPlSpVQlpaGl68eAFLS8si7UPj5FihUKBNmzawtrZG3759kZycjB49emD79u2oWbMmKlSogEqVKqFZs2YAAEtLS2zfvh1WVlYqyfGhQ4dgY2MjJcc7d+5Enz598MEHH8DFxQXDhw/H06dP0b17d2mbkydPokePHggODka9evUQERGBWbNm4dSpU6hRowYAYM+ePbh69Sq8vLzQq1cvREdHo1mzZggPD4e/v3+BMb48rcLW1lZ6DwBycnLw1VdfoV27dhrHocnxlBVbLyTio+2XkPyMyTFRuWe7Mfe/0//WbTsKObByPPBKYgxASpbfH/sRcMYcMOIUC9KdmvYmWN2lEh5n5uCZgksYljYjGWBvaQpL07L197ywJ99pSuPk+NatWzh37hySkpLg4uICAPjyyy+Rk5ODmjVrwsPDA15eXtKorCaEEPj888/x5Zdf4rvvvgMADBo0CF5eXir13n//fXz99deYOHGiVDZ06FBMnz4d69evl8ocHR1x8OBBaRj9xo0bWLt2Lfz9/dGmTZvXxuji4qLy3scffwxzc3MsXLhQozg0PZ5XZWVlISvrf3+xU1NTC61fUkb/cQFPn/MGBiLSQsIV4NnDwuukpeTW86hfOjERUalTCOBJZnaJJMfp6ekwNjaGhYXFa+tmZGQAAKysrIrdbkE0To4rV66MChUqYMaMGfj444/h4+MDGxubYjUeHx+P2NhY9O3bVyrz8PBAixYtpNexsbG4evUqTp48iZiYGAghIIRAbGws0tLSVPYXFBSkMr/E29tbZT6cNlatWoWQkBAcPnwYVapU0SgOTY5HnVmzZmHGjBlFipOIqFSlPyrZekT0RsvOzoa/vz+2bduGevXqvbb+6dOnMXnyZJw6dUpn0zo0To6tra1x5MgR/PTTT2jfvj0AYMiQIfj+++9hbm5epMaTk5MBAE5OTirlzs7O0v+npKQAAPz9/VXKAwIC8t2J+OrcEmNjY8jlcq3jOnHiBMaNG4fly5dL0yw0iUOT41Fn6tSp+Oyzz6TXeY831LWQfn6cVkFkIA49Gg4AaO+4RrcNWTuWbD0iA5eVlYWs589hZ2+v71BKVN60iuJas2YNateurTYxTktLg42NjUoS3K5dO8hkMuzYsQO9e/cudvvqaHVDXv369aVpDBEREejRowdcXFzwxRdfqM3ezc3N8eKFauL1+PFjacQ5LwG8e/cuqlSpItW5ffs2GjRoAABSeb169Yo9b1eTf2HEx8cjODgY48ePx/Dhw6VyTeLQ5HjUMTc3L/I/MIqjr58bguu78oY8IgPwYMIzAEDyjLdeU7N45PIOaHxiCe4lJua7IQ/I/Zx1dauCsyGfclk30qkXWc/xID4OTtZmsLQs/e9QTa3dtglLlizGqdMR+g6lRJXUDXnLli3DlClTVMp27tyJ8ePH48GDB7Czs8Ps2bMxcuRI6f0hQ4Zg2bJl+k+O4+LikJycjICAAABA06ZN4enpiYcPc+eeOTo64tEj1Z/RvLy8cOLECen1jRs3EBYWhoEDBwLIvaOwZcuWWLBgAQIDAyGTyXD8+HFERERIyaS7uzvatm2Lb7/9Fm3atIGdnR0A4OHDh7hw4YI0iq0JdTG+7Pnz5+jduzcaNGiAOXPmqLynSRyaHE9ZY2wkQ0WbsvuhQkSaSVH+tzT+Pi9csAB9+vSBTCZTSZDzvigXLpiPyva6mw9IBACZxgo8lMlgbJT7pyzKyclBRkY6cnJy8PhRbr7k4OAAuVyO58+fw97eHkIIpKWlwdTUFDk5OSq/iqelpcHExCTfL+PPnz9/7fzcrKwsqQ2FQoHs7GyVgbj09PRCV/lQ10Z2djaePn0KALC1tc03sPdqmzk5OTAzMyuwjfj4eERGRqJTp05S2c2bN9GvXz8sX74cw4YNw969e9G7d2/Uq1cPgYGBAICOHTvis88+w7Nnz4o9xVcdjW/pMzU1xYQJE9C0aVMMGTIETZs2RXJyMj788EMAQK9evRAaGooBAwZg1KhRSEpKwqRJkxAWFobWrVujb9++6NatG6pVq6ay3wULFmDfvn1o2rQpgoODMWjQINStW1elzoYNGyCXy+Ht7Y1+/fqhU6dOCAgI0PrGNXUxvmzz5s04c+YMbG1tMWbMGIwaNQqjRo2SEmVN4tDkeIiIyrPg4GCEhoaq/EIG5A4ihIaGcp1jIqWzZ89i6tSpiIqKgo+PD3x8fBAVFYX169ejTZs2GDFiBCpWrIhBgwZh/vz5GDp0qMr2Q4cOxS+//CK9Xr16NapVqwYnJydUrFgR8+fPL7Dt9evXo3Xr1ujfvz8qVqwIa2trjB49GidOnICPjw+cnZ3h7u6O8PBwle0Ka+PYsWPScTg4OCAgIACXL19WabNt27YYPXo0KlasCBsbG/Tt2xfZ2dlqYzx16hTc3d1RoUIFqez333+Ht7c3hg8fDiMjI3Tr1g0dOnTA6tWrpTo+Pj4wNTVFZGTka85A0Wg8clylShWcPHkSZ86cwbVr1zBixAi0adMGpqa5802Cg4Nx5swZXLhwAenp6bC0tETDhg0RExODsLAwWFtbo3nz5oiMjJS2AXLn8F6/fh0HDhyAra0tQkJCcOnSJZiY/C80d3d3nDlzBuHh4bhx4wbc3NwQGBio8q+Fzz//PN/c3sGDB6s8NUddjC8/BKR58+YICQnJd+yurq4ax6HJ8RARlXfBwcHo1asXn5BHVIjAwEAsXLgQixYtwpkzZ6TyiIgIXLx4ESNGjMDKlSthbGys9jkLL9u1axcmT56M3bt3IzAwENHR0Wjfvj28vb0LfKbCpUuXMHHiRGzevBlRUVFo3LgxwsPDsW/fPlSrVg1TpkzBp59+KiXIr2ujffv20j1YOTk5+OGHHzB06FCcO3dOavPChQsYP348li9fjvv376Np06bYvHlzvsQfAJKSkvI9K+LSpUto2LChSlnDhg1x7Ngx6bVMJoOTkxPu3btXaJ8VldYZW5MmTdCkSRO17/n6+sLX11elzMXFBe+++670OigoKN92Tk5OKk+ty1tX+GUymQzNmzdH8+bN1bb99ttv5ytTV/fVGDt27Cj9v7e3N7y9vdXuX9M4AM2Oh4iovDM2NubnG5U50SMLnlpU7T/hsKzaCABwd0Ewnl34S209557/h4q9pgEAHh/9DUlrx6q8X2dV8ddV9vDwwKeffqpx/eXLl+O9995DzZo1kZKSgooVK2LAgAEIDQ0tMDmuXbs2RowYASD3nilvb2+MGDFC+hX/3XffxeLFi7VuQ6FQIC0tDR988AG+/fZb3L9/X1rmt3bt2hg9ejSA3MHFoKAgREVFqY3PyMgo38IJGRkZ+RJmOzs7PHv2TKUsJydHZ/8Y53AmERERUSlzc3PTqv7du3dx5MgRrFmzRqVc3aBjnpenKwCAmZmZSpmZmRmeP3+ucRuxsbEYPXo0jh8/DgsLC5iamkIulyMxMVFKjl9t09zcXKWNl3l4eOD+/fsqZfb29nj8+LFK2aNHj+Dg4CC9zsnJwaNHj+Du7l7gsRcHk2MiIgPg8ekufYdAVCZoOqrr8cl2jepVaPM+KrR5v8jxGBsbq13Z5VXW1tYqU0GB3Eey5/1aX7NmTfTo0QPff/+9Sh1N9q2p17UxadIkVK9eHbt374aVlRWePXsGe3v7Ii2bCwDNmjVDcnKyNDULABo3biw9fC3PyZMnpZvxAODy5cswMjJCo0aNitTu6xT/GXtERKR3NvXfgk193S7jRkTac3d3R1xcHK5evYqUlBTk5Kh/Mm3Tpk0RFhaGffv24fbt25g1axbOnj0rvT9p0iQsWLAAq1evRlxcHM6dO4evv/4aS5cuLbFYX9eGEAJPnjzB/fv3ERUVheHDh0OhUBS5PWdnZ7Rt2xZ79+6VyoYOHYqnT59i6tSpiI2Nxdy5c3H27Fl88MEHUp09e/agR48eOlsGl8kxERERkY60bdsW/fr1Q9euXaXVKiwsLFSmCQC5o6gzZ87ExIkT0a1bN2RmZmLAgAHSY5KbN2+Offv2Ydu2bWjfvj1Gjx4NW1tblfV/X6auDQcHB5Xl2UxMTFQWM3hdG7/88gueP3+O1q1bY9CgQejSpQuqVq0qLbSgrk1bW9tCl1v76KOPVFaiqFixIv755x+cPn0arVu3RmhoKP766y/pnjAhBNavX4+PPvqowH0Wl0yU5Hg8lYjU1FTY29vj6dOn0nrKRESFubcmd1lN1+HL9BwJUenIzMzE9evXUatWrXzrAFP50r59e/z4449o2rTpa+vu27cPISEh2LZtm9r3C7suNM2vOOeYiMgAPDmWO/LC5JiIyptDhw5pXLdLly4Frs5RUjitgoiIiIhIickxEREREZESk2MiIiIiIiUmx0RERFRuFWcpMTI8JbHOBJNjIiIiKnfMzMwAIN+DM+jN9uLFCwC5y9QVFVerICIyABWCPtR3CESlytjYGI6OjkhKSgKQ+4Q5IyOO+b3JFAoF7t27B2traybHRERvuspD5us7BKJSV6VKFQCQEmQiIyMjuLu7QyaTFXkfTI6JiIioXJLJZHB3d4erq6v0czq92czNzYv9CwKTYyIiA5B2fjcAwLZhdz1HQlT6jI2N+ZQ8KjFMjomIDED8wncBAHVWZek5EiKi8o0z14mIiIiIlJgcExEREREpMTkmIiIiIlJickxEREREpMTkmIiIiIhIiatVEBEZAJMKVfQdAhGRQWByTERkAGr9fFPfIRARGQROqyAiIiIiUmJyTERERESkxOSYiMgARI80R/RIc32HQURU7jE5JiIiIiJSYnJMRERERKTE5JiIiIiISInJMRERERGREpNjIiIiIiIlJsdEREREREp8Qh4RkQGo9p9wfYdARGQQmBwTERkAy6qN9B0CEZFB4LQKIiIiIiIlJsdERAbg7oJg3F0QrO8wiIjKPU6rICIyAM8u/KXvEIiIDAJHjomIiIiIlJgcExEREREpMTkmIiIiIlJickxEREREpMTkmIiIiIhIiatVEBEZAOee/6fvEIiIDAKTYyIiA1Cx1zR9h0BEZBA4rYKIiIiISInJMRGRAXh89Dc8PvqbvsMgIir3OK2CiMgAJK0dCwCo0OZ9PUdCRFS+ceSYiIiIiEiJyTERERERkRKTYyIiIiIiJSbHRERERERKTI6JiIiIiJS4WgURkQEw92yo7xCIiAwCk2MiIgNQY/opfYdARGQQOK2CiIiIiEiJyTERkQHISX2AnNQH+g6DiKjc47QKIiIDcH2iOwCgzqosPUdCRFS+ceSYiIiIiEiJyTERERERkRKTYyIiIiIiJSbHRERERERKTI6JiIiIiJSYHBMRERERKXEpNyIiA+A1J1bfIRARGQQmx0REBsDU0V3fIRARGQROqyAiIiIiUmJyTERkAOJmtkHczDb6DoOIqNzjtAoiIgOQefOUvkMgIjIIHDkmIiIiIlJickxEREREpMTkmIiIiIhIickxEREREZESk2MiIiIiIiWuVkFEZABcBs3TdwhERAaByTERkQFw7DBW3yEQERkETqsgIiIiIlJickxEZABS/voJKX/9pO8wiIjKPU6rICIyAA+2TwMAOHebrOdIiIjKN44cExEREREpMTkmIiIiIlJickxEREREpMTkmIiIiIhIickxEREREZESV6sgIjIA1nU76TsEIiKDwOSYiMgAeH62W98hEBEZBE6rICIiIiJSYnJMRGQAshKuICvhir7DICIq9zitgojIANz8TyMAQJ1VWXqOhIiofOPIMRERERGREpNjIiIiIiIlJsdEREREREpMjomIiIiIlJgcExEREREpMTkmIiIiIlLiUm5ERAbAe1mqvkMgIjIITI6JiAyAkam5vkMgIjIInFZBRERERKTE5JiIyADcmOqLG1N99R0GEVG5x2kVREQGIDs5Vt8hEBEZBIMfOd65cydOnjyp0zb27NmDo0eP6rQNIiIiItI9gx85XrhwIZo0aYLmzZvrrI0VK1bA3d0dbdq00VkbRERljVwux7Fjx3Dv3j24urqidevWMDY21ndYRETFYvDJcWno1q0bHBwc9B0GEVGp2b59OyZMmID4+HipzN3dHfPnz0dwcLAeIyMiKp43Jjm+c+cOwsPDIZfL0bVrV9jb20vvHT58GJcvXwYAODk5oUmTJqhVq1a+fZw8eRJXr16Fm5sb2rZtCwsLCwBAlSpVYGNjo1I3JycHhw4dwt27d1G/fn0EBATo8OiKR64QeJTxQt9hEFExCOV/HzzL0nlbu3fuwMihAyGEUClPSEhAnz59sGr9JnTv9Y7O4yB6UzlamcHYSKbvMAyWTLz66WZgOnbsiIcPH+LRo0do3bo1oqKikJycjGPHjqFGjRoAgC1btuDYsWMAgPv372Pfvn34/vvvMWHCBGk//fr1w8mTJ9GpUyckJycjISEBe/bsgaurK9555x24u7tj0aJFAICbN2+ia9eukMvlaN68OWJiYtC4cWMsWbJEo5hTU1Nhb2+Pp0+fws7OroR7RNXWC4n4aPslJD9jckxUnp1LyR2tbeS8XbcNKeTAylHAs4cF17F1BkaGAEacYkGkC5VszLAouD76+rnpO5RyRdP86o0YOY6JicGVK1dQrVo15OTkoHPnzpgyZQq2bt0KAOjfvz/69+8v1Q8LC0OnTp0wZMgQODk54c6dO9i6dStiY2OlhPrGjRv5Rk3yDB06FJ6enti1axfMzXMX5s9LvtXJyspCVtb/RntSU0vvSVej/7iAp89zSq09ItKNybZflE5DCVcKT4wBIC0lt55H/dKJiegNk/zsBUb/cYHJsY68Eclxt27dUK1aNQCAiYkJxo4diyFDhkChUMDIKHfBjps3b+LMmTNISUmBQqFAdnY2oqOj0apVK9jY2MDMzAz79u3DqFGjYGZmBi8vL7Vt3b17FydOnMChQ4ekxBgAWrduXWB8s2bNwowZM0rugInojfOPeYvSaSj9UcnWIyIqYwx+KTcg9yaRl3l4eCArKwsPHjwAAMycORN+fn5Yv349Lly4gKtXr0Imk+Hhw9zREUdHR2zatAlLly6Fk5MTunTpgt9//11tW4mJiQCA6tWraxzf1KlT8fTpU+nP3bt3i3KYRRLSzw+VbMxKrT0iKuesHUu2HhFprZKNGUL6+ek7DIP1Rowcp6SkqLx+8OABTExM4OTkhIyMDEyfPh1//fUXOnfuDABIT0/H4sWLVaZNBAcHIzg4GPfu3cPu3bsxZswYZGdnY/jw4Sr7dnJyAgAkJSVJo9WvY25urjLKXJr6+rkhuL4rb8gjKueebZ8KALAJnqXTduTyDmh8YgnuJSaqnVomk8ng6lYFZ0M+5bJuRDrCG/J0641Ijvfs2YPU1FRp8vXvv/+O1q1bw8TEBA8fPoRcLkflypWl+uvWrVPZ/sGDB5DJZHB2doarqytGjx6NjRs3SitcvMzLywve3t5YtmwZmjVrJpXHxcVpnCyXNmMjGSra6Cc5J6KSkXJkMQCg+rC5Om9r4YIF6NOnD2QymUqCLJPJlO/PR2V7K53HQUSkC29EcmxtbY2WLVtiwIABOH/+PHbv3o0jR44AAFxcXBAUFISBAwfivffew82bN7Fz506YmPyvax4/foyuXbsiKCgI3t7euHr1Kk6fPo05c+aobW/16tXo0qULkpKS0LZtW8TExCA1NRXbt+v4LnIiolIQHByM0NBQtescz5s3j+scE1G5ZvDJ8TvvvANPT084ODjg8OHDqFOnDr777jv4+PhIdf766y+sWbMGN27cQK1atXDu3Dn8/PPP0soUtWvXRmRkJP744w9cu3YNderUQXR0NDw9PQHkfwhI3vJtmzZtQmJiIjp27IgBAwaU6nETEelScHAwevXqxSfkEZHBMfh1jsuj0lznmIgMQ/TI3KlRdVbp/iEgRETlkab51RuxWgURERERkSaYHBMRERERKRn8nGMiojeBXbOB+g6BiMggMDkmIjIAVUav0XcIREQGgdMqiIiIiIiUmBwTERmA9KtHkH71iL7DICIq9zitgojIANyZ8xYALuVGRFRcHDkmIiIiIlJickxEREREpMTkmIiIiIhIickxEREREZESk2MiIiIiIiWuVkFEZACMLGz1HQIRkUFgckxEZAC8F6foOwQiIoPAaRVEREREREpMjomIiIiIlJgcExEZgJjxzogZ76zvMIiIyj3OOSYiMgCK52n6DoGIyCBw5JiIiIiISInJMRERERGREpNjIiIiIiIlJsdEREREREpMjomIiIiIlLhaBRGRAfCc9Le+QyAiMghMjomIDIC1T1t9h0BEZBA4rYKIiIiISInJMRGRAUgIGY6EkOH6DoOIqNxjckxEZABSwzchNXyTvsMgIir3mBwTERERESkxOSYiIiIiUmJyTERERESkxOSYiIiIiEiJyTERERERkRIfAkJEZAAcO32i7xCIiAwCk2MiIgPgMmCOvkMgIjIInFZBRERERKTE5JiIyACkRmxDasQ2fYdBRFTucVoFEZEBSFg2CABg1zRLz5EQEZVvHDkmIiIiIlJickxEREREpMTkmIiIiIhIickxEREREZESk2MiIiIiIiWuVkFEZABMK9XUdwhERAaByTERkQHwmnVF3yEQERkETqsgIiIiIlJickxEZAAU2VlQZPMBIERExcVpFUREBiDmQzsAQJ1VTJCJiIqDI8dEREREREpMjomIiIiIlJgcExEREREpMTkmIiIiIlJickxEREREpMTkmIiIiIhIiUu5EREZgBrfntN3CEREBoHJMRGRATCv4qvvEIiIDAKnVRARERERKTE5JiIyAHfmdsedud31HQYRUbnHaRVERAYgPeoffYdARGQQOHJMRERERKTE5JiIiIiISInJMRERERGREpNjIiIiIiIlJsdEREREREpcrYKIyABUDP5O3yEQERkEJsdERAbAudtkfYdARGQQOK2CiIiIiEiJyTERkQF4dHApHh1cqu8wiIjKPU6rICIyAPc3fgoAcOwwVr+BEBGVcxw5JiIiIiJSYnJMRERERKTE5JiIiIiISInJMRERERGREpNjIiIiIiIlrlZBRGQALGsE6jsEIiKDwOSYiMgAVPv6qL5DICIyCJxWQURERESkxOSYiMgAZD+KR/ajeH2HQURU7nFaBRGRAbgxqSYAoM6qLD1HQkRUvnHkmIiIiIhIickxEREREZESk2MiIiIiIiUmx0RERERESkyOiYiIiIiUmBwTERERESlxKTciIgNQ61eucUxEVBKYHBMRGQATu4r6DoGIyCBwWgURERERkRKTYyIiA3BzRiBuzgjUdxhEROUep1UQERmArDvn9R0CEZFB4MgxEREREZESk2MiIiIiIiUmx0RERERESkyOiYiIiIiUmBwTERERESlxtQoiIgNQ+b2l+g6BiMggMDkmIjIAFdq8r+8QiIgMAqdVEBEREREpMTkmIjIAD3Z+hwc7v9N3GERE5R6TYy1du3YNcXFxBb4mItKHlD+/R8qf3+s7DCKico/JcSFiYmJw+/ZtlbLJkyfj559/LvA1EdGbQi6X4/Dhw9i0aRMOHz4MuVyu75CIiIqNyXEhPv/8c/z6668qZd7e3qhevbqeIiIiKhu2b9+OatWqISgoCIMGDUJQUBCqVauG7du36zs0IqJi0clqFffv30dSUhL8/PyQkpKChIQE1KlTB2ZmZvnqZmZm4tq1a3BycoK7u3uB+7l79y7u3buHxo0bw9jYGADw8OFDxMfHo1atWrCyssq37+TkZCQmJqJ69eqwt7dXeS86Oho2NjZwd3fH7du3IZfLUaNGDchkMgBAbGwsnjx5gqSkJISHhwMAGjZsiJEjR6o9Dm3aJt2RKwQeZbzQdxhEpU4o//vgWZbO29q9cwdGDh0IIYRKeUJCAvr06YNV6zehe693dB4HEZV/jlZmMDaS6TsMFTpJjrdu3Yoff/wRDRs2xPnz56FQKKBQKLB9+3Y0b95cqvfrr79ixowZqFKlCh4/fgw3Nzds3rwZXl5eKvvx9/fHuXPn4ObmhgMHDkAIgVGjRmHnzp2oUaMGkpOTMWPGDIwfPx4A8OzZM7z//vv4+++/4enpiVu3bqF///5YunQpTE1NAQAff/wxzM3NpfnCiYmJ8Pb2xt9//w07Ozts3LgRUVFRiI2Nleps2bIFkydPhru7OxYtWqT22DVpm3Rj64VEfLT9EpKfMTmmN8+5HAUAoNL0v3XbkEIOrBwPvJIYA5CS5ffHfgScMQeMjHUbCxGVe5VszLAouD76+rnpOxSJzqZV5I3o3r17F/Hx8ejVqxeGDRuGnJwcALk/yf344484deoUoqKiEB8fj8DAQAwZMiTffgICAnDnzh2Eh4fDxsYGH3zwAS5fvozr16/jypUriI+Pl0aTAWDcuHHIzMxEfHw8Ll68iLi4OISHh+OXX35R2ffRo0exefNmREVFIS4uDklJSViyZAkAYNq0aWjZsiX69++P8PBwhIeHo2rVqq89bk3bfllWVhZSU1NV/pD2Rv9xgYkxka4lXAGePSy8TlpKbj0iotdIfvYCo/+4oO8wVOgsOTY1NcWMGTMAADKZDDNnzkRsbCyOHj0KAFi0aBHeeustpKWlISIiAmfOnEHbtm1x6tQppKSkSPuxtrbGlClTpNcPHz7Eli1b8P3338PDwwMAYGFhgQ8//BAA8PjxY/z+++945513cPXqVURERCA2Nhbt2rXD7t27VWLs06cP6tevDwCwt7dH+/btcenSpSIfszZtv2zWrFmwt7eX/uQdFxGRpo6YNcERsya6byj9UcnWIyIqY3T2hDxXV1fY2tpKr52cnODs7IybN2+iffv2iImJQVxcHK5du6ayXWBgIJ48eQJnZ2cAgJubG0xM/hdmbGwshBDw8/NT2+6NGzegUCiwZMmSfNMY3NxUh+xdXFxUXltZWeHRo6J/oGvT9sumTp2Kzz77THqdmprKBLkIQvr5cVoFvbE+tfuqdBqydizZekT0RsubVlGW6Cw5Tk9PV1tmZ2cHIHe0d8CAAZg5c2ah+zEyUh3czrvxrqCpBxYWFgCAZcuWoWnTplrHXRxFbdvc3Bzm5ua6CuuN0dfPDcH1XXlDHpEOyeUd0PjEEtxLTMx3Qx6Q+0uhq1sVnA35VGW6GxGROm/MDXlA7vSH8+fPo2HDhgCA48ePIzMzE40aNQIAtGvXDqGhofjmm29URlnT09NhbW1d4H59fHxQqVIl7NixQ9o3ADx//hwWFhaoU6cOKlWqhI0bN+ZLUF+371dZWloiOztb4/ol2TYVjbGRDBVt+A8NevNk3j4HALCs2kjnbS1csAB9+vSBTCZTSZDzVvtZuGA+KtvnX0GIiKg80FlybGZmhv79+2PmzJlQKBSYPHkyhgwZglq1agEAvvnmGwQGBqJTp074+OOPYWVlhVOnTmH//v04efJkwQGbmOCXX37ByJEjkZ2djXbt2iEmJgb79+/Hrl27YGJiggULFmDo0KFQKBTo1q0bHj9+jD179qBKlSr44YcfND6GevXqYdOmTThw4ABsbGxUkvGCYiuptomItBH3bTMAQJ1Vul/KLTg4GKGhoZgwYQLi4+Olcnd3d8ybNw/BwcE6j4GISFd0lhxXr14dCxYswNq1a5GYmIjhw4fjyy+/lN738PDAuXPnsGDBAixduhQWFhZo1qwZ9u7dK9WpXLmy2oR0yJAh8PDwwIoVK/DDDz/Az88Pq1atkt7v378/qlevjuXLl2PmzJlwc3ND9+7dMXDgQKmOr68vPD0988X88ujuxIkTkZaWhlmzZiE9PR1btmyBt7c3KlWqJNV59bUmbRMRlXfBwcHo1asXjh07hnv37sHV1RWtW7fmVAoiKvdkQt2ksWJatGgRFi1ahKtXr5b0rt8IqampsLe3x9OnT6U52kREhYkemTudqDRGjomIyiNN8ys+PpqIiIiISEknyXFB0yGIiIiIiMoyncw57tOnD/r06aOLXRMRERER6YzObsijosubBs7HSBORpp694OcGEVFh8j4fX3e7HZPjMigtLQ0A+JQ8ItLeBnt9R0BEVKalpaXB3r7gz0qdrFZBxaNQKJCYmAhbW1tpUX3632O17969y1U8ygCej7KF56Ns4fkoe3hOyhZ9nA8hBNLS0uDm5pbvCcwv48hxGWRkZAR3d3d9h1Fm2dnZ8YOtDOH5KFt4PsoWno+yh+ekbCnt81HYiHEeLuVGRERERKTE5JiIiIiISInJMZUb5ubmmD59OszNzfUdCoHno6zh+ShbeD7KHp6TsqUsnw/ekEdEREREpMSRYyIiIiIiJSbHRERERERKTI6JiIiIiJSYHFO5cOPGDZw9exaZmZka1X/x4gUuXbqEW7duQaFQ6Di6N8/9+/cRERGBlJQUjbd59OgRwsLCkJycrMPIDFtaWhrOnDmDuLg4nW5DmsnJycGFCxcQFRX12sfRvuzcuXOIjIzUYWRvrmvXriEyMhJZWVka1c/KysLFixdx+/Ztrc4haSYxMRERERF4/PixxtvcuXMHFy5ckB71rBeCqAx7+PChaNWqlbCzsxO1atUS9vb2YuvWrQXWz8jIEF988YVwdHQU9evXF66ursLb21ucPHmyFKM2XAqFQowbN06Ym5sLX19fYW5uLqZMmVLoNteuXRPDhw8Xrq6uAoBYvXp16QRrYEJCQoSVlZXw9vYW1tbW4u233xbPnj0r8W1IMydOnBBVqlQRHh4eolKlSsLHx0fExMQUus2CBQtEnTp1hIODg/D29i6lSN8M9+7dEwEBAcLBwUHUqFFDODo6it27dxdYPzU1VXzyySeiQoUKokGDBsLFxUXUq1dPREZGlmLUhis7O1sMHTpUWFhYSN8VM2fOLHSbv//+W9SvX1/UrFlT1KtXT1haWoqJEycKhUJRSlH/D5NjKtMGDBgg/Pz8RFpamhBCiF9//VWYm5uL27dvq62fkJAg5syZIyUAOTk54v333xeVKlUSL168KLW4DdWyZcuEra2tuHz5shBCiFOnTgkzMzOxZcuWArfZuXOnWLVqlcjIyBDGxsZMjovg4sWLwsjISGzcuFEIIURycrKoVq2a+Pjjj0t0G9JMRkaGcHNzE+PGjRNC5H7OdO3aVfj7+xe63cSJE0VUVJT4+uuvmRyXsO7du4vAwECRkZEhhBDiu+++EzY2NuL+/ftq68fGxor58+dL9V+8eCEGDBggPD099ZKMGZrZs2cLZ2dncfPmTSGEEAcPHhRGRkZi//79BW7z22+/iRs3bkivIyIihJGRUaEDYrrC5JjKrKdPnwpTU1OxZs0aqSw7O1s4OjqKWbNmabyf8PBwAUBERUXpIsw3SkBAgBg+fLhKWffu3UXnzp012p7JcdF89tlnwsvLS6Vs9uzZwt7eXuTk5JTYNqSZ7du3C5lMJhITE6WysLAwAUCcO3futdszOS5ZSUlJQiaTidDQUKksIyNDWFtbi4ULF2q8nwMHDggABQ6+kOZq164tPv30U5WyVq1aif79+2u8D4VCIRwcHMS8efNKOrzX4pxjKrMuX76M7OxsNG7cWCozMTFBw4YNce7cOY33ExERAVNTU3h6euoizDeGEAIXLlxQOR8AEBAQoNX5IO2dO3dObb8/ffoUN2/eLLFtSDPnzp2Dm5sbXF1dpbKAgADpPSpdFy5cgBBC5Xq3tLRE3bp1tf6usLKyUjmvpL309HRcu3atSN8VaWlpCAsLw759+/Dee+/Bzc0NgwcP1mW4apmUeov0Rrt06RKePn1aaJ1WrVoByL2BCwCcnJxU3ndycpLee53r169j2rRpmDhxImxsbIoQsWG7c+cO7ty5U2id+vXrw97eHunp6cjKyirW+aCiefToEerUqaNSlnceCur7omxDmnn06FG+vwempqawtbVl3+pBSXxXXL58GTNnzsSXX34JU1PTEo/xTZJ3811Rzsfdu3fx5Zdf4tGjR4iPj8cPP/wAZ2dnncVaECbHVKpWrFjx2n85HjlyBMbGxtIH1PPnz1Xez8zMhJmZ2Wvbio+PR+fOndG2bVvMnDmz6EEbsAMHDuC3334rtM68efPQpEmTYp8PKjpTU1O1/Q6gwL4vyjakGXV9C+T+3WDflr6XP5tsbW2l8szMTJXXBYmNjUWXLl3Qo0cPfP311zqL801RnO8KX19fhIWFAcgdyQ8KCoKpqSnGjBmjm2ALwOSYStXChQs1rlu1alUAQEJCAqpXry6VJyQkIDAwsNBt4+Pj0a5dO9SrVw9//PEHTEx4qavz/vvv4/3339eorrm5OVxcXJCQkKBSnpCQwCkrOla1alW1/Q6gwL4vyjakmapVqyIpKQkKhQJGRrmzE5OTk5Gdnc2+1YOXvysqVqwolSckJKBJkyaFbnvz5k0EBQWhVatWWL9+vXQ+qegqVqwIKyurYn9XNG3aFEFBQdi7d2+pJ8e8CqjM8vb2hoeHB/7880+p7Pbt2zh//jw6deoklV2+fBlXrlyRXickJCAoKAh16tRBaGgoR3JKUKdOnbBr1y7ptVwux+7du1XOx927d3Hy5El9hGewOnXqhKNHj6pMSdq5cycaNWok/XT59OlThIWFIT09XeNtqGg6deqE1NRUHD58WCrbuXMnzMzM0KZNG6ksLCwM9+7d00OEb5aGDRvC2dlZ5bsiOjoa169fV/lsunjxImJiYqTXcXFxCAoKQvPmzbFhwwYYGxuXatyGysjICO3bt1c5Hy9evMDevXtVzkdcXBwiIiKk13mfXXkUCgVu3bqln8+rUr8FkEgLGzZsECYmJuKnn34S27ZtE/7+/qJZs2ZCLpdLdTp06CC6desmhBDi0aNHonbt2qJ27dri4MGD4tixY9KfJ0+e6OswDMbVq1eFra2tGDlypPjzzz/FgAEDhJOTk7hz545U57vvvhPW1tbS6ydPnkjnwNjYWHz11Vfi2LFj4tq1a/o4hHIpMzNT+Pr6ijZt2ogdO3aI//znP8LY2Fjs3btXqvPvv/+qrJagyTZUdEOHDhWenp7i999/FyEhIcLe3l785z//kd7Pzs4WAFRWS7hw4YI4duyYGDZsmPD09JT+XnCZyeJbvny5MDc3F/PmzRNbt24VdevWFe3bt1epExgYKK2WkJSUJKpVqybq1asnDh8+rPJdkZqaqo9DMChnz54VFhYW4pNPPhF//vmn6Nmzp3BzcxMPHjyQ6nz++eeiSpUq0ut69eqJn376Sezdu1eEhoaK7t27Czs7O2np0NIkE4KPhKGybffu3Vi9ejWePn2K5s2bY9KkSbCzs5PenzhxIszMzPDjjz8iJiYGI0eOVLufvLmzVDxRUVH45ZdfcPv2bdSsWROTJ0+Gl5eX9P66deuwYcMG/P333wBy797/+OOP8+2na9eu+Oqrr0ot7vIuJSUFP/74I86dOwcnJyd8+OGHCAoKkt7P6+c1a9ZI5+N121DRZWdnY9GiRdi/fz9MTEzw7rvvYvjw4ZDJZAByf1Vp27YtPvvsMwQHBwMAxowZg6ioqHz7+vPPP+Ho6Fiq8Rui7du3Y/369Xj27Blat26Nzz//HNbW1tL7Y8aMgaurK7755psCP5eA3HtjfH19Sytsg3X27FnMnz8fCQkJ8PHxwZQpU1SmVSxevBgHDx7E9u3bAeR+Xi1cuBARERGwsLBAgwYNMHbsWLi4uJR67EyOiYiIiIiUOOeYiIiIiEiJyTERERERkRKTYyIiIiIiJSbHRERERERKTI6JiIiIiJSYHBMRERERKTE5JiIiIiJSYnJMRPSGOXHiBK5du6bTNi5fvoyDBw/qtI2SUNw4z5w5g8uXL5dgRESkb0yOiYjeIAkJCejVqxcsLS112k5oaCimT5+u0zZKwqtxXrp0CYcOHdJ4+/T0dPTs2RPPnz/XRXhEpAcm+g6AiIhKz3fffYfg4GB4eHjotJ169erpdP8l5dU4t2zZgrCwMLRv316j7du2bYvKlSsjJCSkwMcRE1H5wuSYiKgIzp8/j6dPnyIgIADnz59HSkoKWrZsCUdHR2RkZOD48eMQQqBFixawsbHJt/3Vq1dx9epVuLq6wt/fH6amptJ7sbGxiIiIAADY2tqibt26qFatmtr2mzVrhvPnz+PRo0do2rQpnJ2dC4w5NTUV69evx4EDB1TKExMTcf78edjZ2aFx48b5RpULi/Xlfjhx4gRSU1PRu3dv+Pj4oEKFClodtyax5Dl16hSEEGjWrJlUFh0djbt37+Ktt97SuI9ejvPatWu4cuUKkpOTsXnzZgBAixYt4OnpWWhcQ4cOxfz585kcExkIJsdEREWwYcMG/Pnnn5DL5ahVqxYSEhIQHx+P+fPn49tvv0WtWrUQFxeHFy9eICIiAo6OjgCArKwsDBs2DEeOHEHTpk0RFxcHhUKBXbt2oUaNGgCAW7duYceOHQByE9pjx45h/PjxmD17tkr7f/31F0xMTODq6oonT57g2rVr+PvvvxEQEKA25kOHDkEmk6m8v2TJEkyZMgUtWrRAdnY2kpKSsGnTJvj5+WkU64YNG7B7924YGRnBzc0NVatWRe/evREaGooDBw6gQ4cOGh93YbG8aunSpcjJyVFJjnfu3InQ0FApOdakj16OMzY2FjExMUhJSZH6v2rVqti9e3ehcbVv3x7jxo3DjRs34OXlpcVVRERlEZNjIqIiunnzJsLDw9GkSRMoFAo0atQIo0ePRkREBBo0aICcnBzUqVMHK1euxOTJkwEA3377Le7evYubN2/CysoKAPDBBx9g3Lhx2LdvHwCgY8eO6Nixo9TOjRs34Ofnh759+6Jx48ZS+fXr13H8+HEEBgYCAPr27YuZM2di586dauONjIyEt7c3jI2NpbLp06djxYoVGDhwIADg7t27uH//vsaxAkBMTAwOHTqEoKCgAvtKk30VFktRadNHb7/9No4fP46wsDBp5BgAevbsWWhc3t7eMDMzQ0REBJNjIgPA5JiIqIgaN26MJk2aAACMjIwQEBAABwcHNGjQAABgYmKCJk2aqKwMsXr1avTs2RN79uyBEAJCCFSqVAnr16+HQqGAkVHufdIZGRmIjIxEUlIScnJy4OLigtOnT6skx/7+/lLSBwDt2rXDwoULC4w3JSUl31QHS0tLREVFIScnByYmJvDw8JDmI2saa/369QtNjDXdV2GxFJW2faSOJnE5ODggJSWlWLESUdnA5JiIqIheTTTNzc3VluWtZPDixQvcu3cP0dHRePLkiUq9Xr164fnz57CyssKBAwcwYMAAuLm5oUaNGrCwsEBaWhqSk5NVtsmbqqGuLXVsbGyQnp6uUrZmzRp89NFHWLx4Mdq0aYM+ffpg8ODByMnJ0ShWAHB1dS2wTW2Ou6BY8pLwotC2j9TRJK709HTY2toWOU4iKjuYHBMRlRJTU1NYWFigf//+GDduXIH1PvvsM4wfPx4zZsyQynx8fCCEKFb7tWvXxrp161TK2rdvjytXriAuLg779+/HpEmTEB0djZkzZ2oUKwDIZLJC39f0uAuK5YcffshX18jICAqFQqVMV8upvS6uBw8eID09Hd7e3jppn4hKF9c5JiIqJTKZDJ07d8bKlSshl8tV3ktISJD+PykpSSXRunjxIq5fv17s9jt06ID79+8jNjYWACCXy6W5s9WqVcOYMWMwePBghIeHaxyrJjTZV2GxqFOlShXcuHFDpezw4cNaxaWOjY2NSpKtSVzHjx+Hvb09mjZtWuz2iUj/OHJMRFSK5s6di7Zt26JFixYYMmQIhBA4evQozMzMsHHjRgDAO++8g6+++gpPnjxBeno6fv31V7XLwWmrevXq6NChAzZt2oT/+7//Q3Z2Nlq1aoXOnTujUaNGSElJwerVqzFz5kyNYy2p435dLK8aOHAgZs+ejXHjxqFx48b466+/cOXKFXh6eharj5o0aYJvvvkG8+fPh4uLC/z9/dGtW7dC49qyZQvee+89mJjwK5XIEPBvMhFRETRq1CjfmsKNGzdG9erVVcoCAwNVfv6vUaMGLl26hLVr1yIyMhIODg4YMmQIevXqJdVZtGgRVqxYgYiICNjZ2WHr1q34559/VB5Yoa79mjVrokePHoXG/Z///AdDhgzBpEmTYGFhgfPnz2Pt2rU4efIkbG1tsXXrVmn5NU1iVRcHkP/hGq/b1+tieZWvry9OnjyJDRs24Ny5cxg8eDCGDx+OU6dOadVHr8bZsWNHrFmzBv/++y/Cw8NRtWrVQuOKj4/H/v37cf78+UL7nYjKD5ko7iQ2IiIqV2bMmIGOHTuiZcuW+g6l3Pvjjz+QmpqKUaNG6TsUIiohTI6JiIiIiJR4Qx4RERERkRKTYyIiIiIiJSbHRERERERKTI6JiIiIiJSYHBMRERERKTE5JiIiIiJSYnJMRERERKTE5JiIiIiISInJMRERERGREpNjIiIiIiIlJsdEREREREr/DzaWLnq6G5tdAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "import matplotlib.pyplot as plt\n", "\n", "fig, ax = plt.subplots(figsize=(7, 3.2), layout=\"constrained\")\n", "names = list(intervals)\n", "for i, (lo, hi) in enumerate(intervals.values()):\n", " ax.plot([lo, hi], [i, i], color=\"#0072B2\", lw=3, solid_capstyle=\"round\")\n", " ax.plot(float(point), i, marker=\"o\", color=\"black\", zorder=5)\n", "ax.axvline(0.0, color=\"#D55E00\", ls=\"--\", lw=1.4, label=\"true mean (0)\")\n", "ax.set_yticks(range(len(names)))\n", "ax.set_yticklabels(names)\n", "ax.set_xlabel(\"mean (series units)\")\n", "ax.set_title(f\"90% confidence intervals for the AR(1) mean (B = {B})\")\n", "ax.legend(loc=\"upper right\", fontsize=9, framealpha=0.9)\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "d2eae5ca", "metadata": {}, "source": [ "## Why BCa is refused for block methods\n", "\n", "BCa corrects the percentile interval for median bias and skewness, which makes it a\n", "strong default on i.i.d. data. Its acceleration constant, though, is a delete-one\n", "jackknife estimate that is only defined under independent sampling (Efron 1987). There\n", "is no valid drop-in version of that constant for dependent data, so asking for\n", "`kind=\"bca\"` with a block method raises a typed error rather than quietly returning an\n", "interval built on an invalid constant. The refusal is doing its job: it protects the\n", "statistical validity of the result and points you at the studentized interval, which is\n", "the second-order-correct route for dependent data." ] }, { "cell_type": "code", "execution_count": 6, "id": "a0321213", "metadata": { "execution": { "iopub.execute_input": "2026-07-02T19:33:07.247635Z", "iopub.status.busy": "2026-07-02T19:33:07.247333Z", "iopub.status.idle": "2026-07-02T19:33:07.250843Z", "shell.execute_reply": "2026-07-02T19:33:07.250191Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "MethodConfigError raised, as intended:\n", "\n", "[TSB_UNSUPPORTED_MODEL_FEATURE] BCa intervals are only supported with the IID method spec; got MovingBlock. The acceleration constant is a delete-one jackknife estimate defined under independent observations (Efron 1987); for dependent data the second-order-correct route is the studentized interval (Gotze and Kunsch 1996). Hint: Use kind='studentized' or kind='percentile' with this method, or method=IID() when the observations are exchangeable.\n" ] } ], "source": [ "from tsbootstrap.errors import MethodConfigError\n", "\n", "try:\n", " conf_int(x, \"mean\", method=MovingBlock(block_length=\"auto\"), kind=\"bca\", alpha=ALPHA)\n", "except MethodConfigError as exc:\n", " print(\"MethodConfigError raised, as intended:\\n\")\n", " print(exc)" ] }, { "cell_type": "markdown", "id": "dae7133b", "metadata": {}, "source": [ "## BCa on i.i.d. data, where it belongs\n", "\n", "With the `IID` spec the acceleration is well defined, so BCa runs and returns a\n", "skew-corrected interval. On a right-skewed statistic this is where BCa earns its keep\n", "over the plain percentile band: both bands bracket the true mean of 1, but the BCa band\n", "shifts to correct the skew." ] }, { "cell_type": "code", "execution_count": 7, "id": "e8e9692d", "metadata": { "execution": { "iopub.execute_input": "2026-07-02T19:33:07.256622Z", "iopub.status.busy": "2026-07-02T19:33:07.256451Z", "iopub.status.idle": "2026-07-02T19:33:07.334413Z", "shell.execute_reply": "2026-07-02T19:33:07.332468Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "sample mean: 1.116\n", "percentile: [0.851, 1.432]\n", "BCa: [0.865, 1.473] (shifted right to correct the skew)\n" ] } ], "source": [ "from tsbootstrap import IID\n", "\n", "skewed = np.random.default_rng(1).exponential(1.0, size=60) # true mean 1, right-skewed\n", "lo_pct, hi_pct, _ = conf_int(\n", " skewed, \"mean\", method=IID(), kind=\"percentile\", alpha=ALPHA, n_bootstraps=999, random_state=1\n", ")\n", "lo_bca, hi_bca, _ = conf_int(\n", " skewed, \"mean\", method=IID(), kind=\"bca\", alpha=ALPHA, n_bootstraps=999, random_state=1\n", ")\n", "print(f\"sample mean: {skewed.mean():.3f}\")\n", "print(f\"percentile: [{lo_pct:.3f}, {hi_pct:.3f}]\")\n", "print(f\"BCa: [{lo_bca:.3f}, {hi_bca:.3f}] (shifted right to correct the skew)\")" ] }, { "cell_type": "markdown", "id": "f9f9cb94", "metadata": {}, "source": [ "## Intervals over a ragged panel\n", "\n", "`conf_int_panel` computes a per-series interval for a whole collection of series in one\n", "pass, even when they have different lengths. It builds on the panel bootstrap, so it\n", "takes the observation-resampling methods only. The studentized panel path needs an\n", "explicit block length: a replicate reducer sees one series at a time without its\n", "identity, so a per-series automatic block length cannot be resolved honestly. We pass\n", "`se_block_length` to supply it." ] }, { "cell_type": "code", "execution_count": 8, "id": "53aec782", "metadata": { "execution": { "iopub.execute_input": "2026-07-02T19:33:07.339447Z", "iopub.status.busy": "2026-07-02T19:33:07.339079Z", "iopub.status.idle": "2026-07-02T19:33:07.788114Z", "shell.execute_reply": "2026-07-02T19:33:07.786822Z" } }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
lengthpointlowerupper
series
01200.1328-0.02700.3817
1180-0.1221-0.26950.0518
290-0.0530-0.38990.2961
\n", "
" ], "text/plain": [ " length point lower upper\n", "series \n", "0 120 0.1328 -0.0270 0.3817\n", "1 180 -0.1221 -0.2695 0.0518\n", "2 90 -0.0530 -0.3899 0.2961" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from tsbootstrap import conf_int_panel\n", "\n", "panel = [ar1(0.4, n, seed=100 + i) for i, n in enumerate((120, 180, 90))]\n", "\n", "lower, upper, point = conf_int_panel(\n", " panel,\n", " \"mean\",\n", " method=MovingBlock(block_length=10),\n", " kind=\"studentized\",\n", " alpha=ALPHA,\n", " se_block_length=10,\n", " n_bootstraps=299,\n", " random_state=0,\n", ")\n", "\n", "panel_table = pd.DataFrame(\n", " {\"length\": [s.size for s in panel], \"point\": point, \"lower\": lower, \"upper\": upper}\n", ")\n", "panel_table.index.name = \"series\"\n", "panel_table.round(4)" ] }, { "cell_type": "markdown", "id": "f0207edd", "metadata": {}, "source": [ "## When to use which\n", "\n", "- Reach for `percentile` or `basic` as the distribution-free default; they assume the\n", " least and work with every method.\n", "- Use `studentized` when you want second-order accuracy under dependence and can afford\n", " the block-jackknife standard error. It is the recommended route for time series.\n", "- Use `bca` for i.i.d. data with a skewed statistic, where its bias and skew correction\n", " buys real accuracy. `conf_int` refuses it for dependent methods rather than return a\n", " constant that is not valid there.\n", "- Use `conf_int_panel` when you have many series and want an interval for each in one\n", " pass.\n", "\n", "Whatever the family, the coverage claim stays approximate or asymptotic under temporal\n", "dependence, so read these as calibrated guidance rather than exact guarantees." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.13.3" } }, "nbformat": 4, "nbformat_minor": 5 }