{"nbformat":4,"nbformat_minor":0,"metadata":{"colab":{"provenance":[],"authorship_tag":"ABX9TyMHTnESJ4IMPXqcUXtQoWzk"},"kernelspec":{"name":"python3","display_name":"Python 3"},"language_info":{"name":"python"}},"cells":[{"cell_type":"code","execution_count":null,"metadata":{"colab":{"base_uri":"https://localhost:8080/","height":598,"output_embedded_package_id":"1bE38b3nFyAprSfzG4khRgCU6xx4YrT2B"},"id":"b5lWm1WxWPHB","executionInfo":{"status":"ok","timestamp":1785219164390,"user_tz":-330,"elapsed":75460,"user":{"displayName":"Raunak Bhattacharyya","userId":"05940311940370151459"}},"outputId":"245c7a44-997a-47c7-b58c-ab724461b81a"},"outputs":[{"output_type":"display_data","data":{"text/plain":"Output hidden; open in https://colab.research.google.com to view."},"metadata":{}}],"source":["%matplotlib inline\n","\n","import numpy as np\n","import matplotlib.pyplot as plt\n","from matplotlib.animation import FuncAnimation\n","from IPython.display import HTML\n","\n","# ============================================================\n","# Rejection Sampling Animation\n","# Target: p(x) = 2x\n","# Proposal: q(x) = Uniform(0,1)\n","# Envelope: Mq(x) = 2\n","# Acceptance probability = x\n","# ============================================================\n","\n","num_frames = 200\n","\n","accepted_x = []\n","accepted_y = []\n","\n","rejected_x = []\n","rejected_y = []\n","\n","# ------------------------------------------------------------\n","# Figure\n","# ------------------------------------------------------------\n","fig = plt.figure(figsize=(12,5))\n","\n","ax_left = fig.add_subplot(121)\n","ax_right = fig.add_subplot(122)\n","\n","# Curves\n","x = np.linspace(0,1,400)\n","target = 2*x\n","envelope = np.full_like(x,2)\n","\n","# ------------------------------------------------------------\n","# LEFT PANEL\n","# ------------------------------------------------------------\n","ax_left.plot(x, target, lw=3, color='royalblue', label='Target $p(x)=2x$')\n","ax_left.plot(x, envelope, 'r--', lw=2, label='Envelope $Mq(x)=2$')\n","\n","accepted_scatter = ax_left.scatter([], [], s=25, color='green', label='Accepted')\n","rejected_scatter = ax_left.scatter([], [], s=25, color='red', label='Rejected')\n","\n","current_point, = ax_left.plot([], [], 'ko', markersize=8)\n","\n","ax_left.set_xlim(0,1)\n","ax_left.set_ylim(0,2.1)\n","ax_left.set_xlabel(\"x\")\n","ax_left.set_ylabel(\"Density\")\n","ax_left.set_title(\"Rejection Sampling\")\n","ax_left.legend(loc=\"upper left\")\n","\n","# ------------------------------------------------------------\n","# RIGHT PANEL\n","# ------------------------------------------------------------\n","bins = np.linspace(0,1,21)\n","\n","ax_right.set_xlim(0,1)\n","ax_right.set_ylim(0,2.2)\n","ax_right.set_xlabel(\"x\")\n","ax_right.set_ylabel(\"Density\")\n","ax_right.set_title(\"Accepted Samples\")\n","\n","# ------------------------------------------------------------\n","# Animation update\n","# ------------------------------------------------------------\n","def update(frame):\n","\n","    # Generate proposal\n","    x_prop = np.random.rand()\n","    y_prop = np.random.rand()*2\n","\n","    current_point.set_data([x_prop],[y_prop])\n","\n","    # Accept / reject\n","    if y_prop <= 2*x_prop:\n","        accepted_x.append(x_prop)\n","        accepted_y.append(y_prop)\n","    else:\n","        rejected_x.append(x_prop)\n","        rejected_y.append(y_prop)\n","\n","    # Update scatter plots\n","    if accepted_x:\n","        accepted_scatter.set_offsets(np.column_stack([accepted_x, accepted_y]))\n","    else:\n","        accepted_scatter.set_offsets(np.empty((0,2)))\n","\n","    if rejected_x:\n","        rejected_scatter.set_offsets(np.column_stack([rejected_x, rejected_y]))\n","    else:\n","        rejected_scatter.set_offsets(np.empty((0,2)))\n","\n","    # Redraw histogram\n","    ax_right.cla()\n","\n","    ax_right.plot(\n","        x,\n","        target,\n","        lw=3,\n","        color='royalblue',\n","        label='Target $2x$'\n","    )\n","\n","    if accepted_x:\n","        ax_right.hist(\n","            accepted_x,\n","            bins=bins,\n","            density=True,\n","            color='lightgreen',\n","            edgecolor='black',\n","            alpha=0.8\n","        )\n","\n","    ax_right.set_xlim(0,1)\n","    ax_right.set_ylim(0,2.2)\n","\n","    ax_right.set_xlabel(\"x\")\n","    ax_right.set_ylabel(\"Density\")\n","    ax_right.set_title(\"Accepted Samples\")\n","\n","    total = len(accepted_x) + len(rejected_x)\n","    rate = len(accepted_x)/total if total else 0\n","\n","    ax_right.text(\n","        0.02,\n","        0.95,\n","        f\"Accepted: {len(accepted_x)}\\n\"\n","        f\"Rejected: {len(rejected_x)}\\n\"\n","        f\"Acceptance Rate: {rate:.2f}\",\n","        transform=ax_right.transAxes,\n","        va=\"top\",\n","        fontsize=11\n","    )\n","\n","    ax_right.legend(loc=\"upper left\")\n","\n","    return current_point, accepted_scatter, rejected_scatter\n","\n","# ------------------------------------------------------------\n","# Create animation\n","# ------------------------------------------------------------\n","ani = FuncAnimation(\n","    fig,\n","    update,\n","    frames=num_frames,\n","    interval=250,\n","    repeat=False,\n","    blit=False\n",")\n","\n","plt.close(fig)\n","\n","HTML(ani.to_jshtml())"]}]}