add binomial coverage example

example/15-binomial-coverage.ipynb

@@ -0,0 +1,130 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "id": "49fd1559-64b2-41b8-a985-7af99528571f",
+   "metadata": {},
+   "source": [
+    "Compute frequentist coverages for the binomial distribution using Jeffreys and Laplace prior with various values of p and n."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "id": "94098bc9-95f3-4ad8-9975-5c43bd77f09e",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from scipy.stats import beta\n",
+    "import math"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "id": "a533ddc1-235c-4088-96b5-886326d70c19",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "def coverage_test(n, theta_true, prior_term):\n",
+    "    alpha = 0.95\n",
+    "    low = (1.0 - alpha) / 2.0\n",
+    "    high = 1.0 - low\n",
+    "    res = 0.0\n",
+    "    for y in range(0, n+1):\n",
+    "        dist = beta(y + prior_term, n - y + prior_term)\n",
+    "        t = dist.cdf(theta_true)\n",
+    "        if t > low and t < high:\n",
+    "            res += math.comb(n, y) * theta_true ** y * (1 - theta_true) ** (n - y)\n",
+    "    return res"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "id": "6db34fdb-7a3f-4360-84a3-d79374a39595",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "thetas = [1.0e-4, 1.0e-3, 1.0e-2, 1.0e-1, 0.25, 0.5]\n",
+    "nx = [5, 10, 20, 100]"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 4,
+   "id": "028bbb5b-d214-40aa-bf46-b3671b21622d",
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "n theta cov_laplace cov_jeffreys\n",
+      "5 0.0001 0.0 0.9995000999900006\n",
+      "5 0.001 0.0 0.995009990004999\n",
+      "5 0.01 0.9509900498999999 0.9509900498999999\n",
+      "5 0.1 0.9185400000000001 0.9914400000000001\n",
+      "5 0.25 0.984375 0.984375\n",
+      "5 0.5 0.9375 0.9375\n",
+      "10 0.0001 0.0 0.9990004498800211\n",
+      "10 0.001 0.0 0.9900448802097482\n",
+      "10 0.01 0.9043820750088044 0.9043820750088044\n",
+      "10 0.1 0.9298091736000003 0.9872048016000002\n",
+      "10 0.25 0.9802722930908203 0.9239587783813477\n",
+      "10 0.5 0.978515625 0.978515625\n",
+      "20 0.0001 0.0 0.9980018988604845\n",
+      "20 0.001 0.0 0.9801888648295347\n",
+      "20 0.01 0.8179069375972308 0.9831406623643482\n",
+      "20 0.1 0.9568255047155371 0.9568255047155371\n",
+      "20 0.25 0.9347622074283208 0.9347622074283208\n",
+      "20 0.5 0.9586105346679688 0.9586105346679688\n",
+      "100 0.0001 0.0 0.9900493386913719\n",
+      "100 0.001 0.9047921471137089 0.9047921471137089\n",
+      "100 0.01 0.920626797747819 0.9816259635553496\n",
+      "100 0.1 0.9363983902254425 0.9556901071912257\n",
+      "100 0.25 0.9512948142448159 0.9512948142448159\n",
+      "100 0.5 0.9431120663590193 0.9431120663590193\n"
+     ]
+    }
+   ],
+   "source": [
+    "print('n', 'theta', 'cov_laplace', 'cov_jeffreys')\n",
+    "for n in nx:\n",
+    "    for theta in thetas:\n",
+    "        cov_laplace = coverage_test(n, theta, 1.0)\n",
+    "        cov_jeffreys = coverage_test(n, theta, 0.5)\n",
+    "        print(n, theta, cov_laplace, cov_jeffreys)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "27150810-3dc6-4636-932e-613b886da79d",
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3 (ipykernel)",
+   "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.11.4"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}