{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "calm-0",
   "metadata": {},
   "source": [
    "# Continuous optimization (CALM)\n",
    "\n",
    "[![Open in Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/github/Abel-ai-lab/andrey/blob/main/examples/calm.ipynb)\n",
    "\n",
    "CALM learns a weighted DAG by gradient descent instead of searching over graphs. It first estimates\n",
    "the moral graph, which leaves only the pairs of variables that can share an edge. It then\n",
    "minimizes the Gaussian likelihood plus a penalty per edge, under a constraint that the graph is\n",
    "acyclic. The fit runs on torch, so it needs the `torch` extra."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7fb27b941602401d91542211134fc71a",
   "metadata": {},
   "source": [
    "## Start\n",
    "\n",
    "Choose a CPU runtime and run this cell first. In Colab it installs `andrey-core==0.9.0` from PyPI;\n",
    "local kernels use the installed environment. CALM needs torch, which Colab provides; the next cell\n",
    "installs it with uv in a local kernel that lacks it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "acae54e37e7d407bbb7b55eff062a284",
   "metadata": {
    "tags": [
     "colab-setup"
    ]
   },
   "outputs": [],
   "source": [
    "import os\n",
    "import sys\n",
    "\n",
    "IN_COLAB = \"google.colab\" in sys.modules\n",
    "\n",
    "if IN_COLAB:\n",
    "    os.environ[\"ANDREY_DEVICE\"] = \"cpu\"\n",
    "    !pip install --quiet \"andrey-core[viz]==0.9.0\" pandas"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "calm-install",
   "metadata": {},
   "outputs": [],
   "source": [
    "from importlib.util import find_spec\n",
    "\n",
    "if not IN_COLAB and find_spec(\"torch\") is None:\n",
    "    !uv pip install --python {sys.executable} --torch-backend=auto \"torch>=2.2\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "calm-1",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from IPython.display import display\n",
    "\n",
    "import andrey\n",
    "import andrey.metrics as metrics\n",
    "import andrey.viz\n",
    "from andrey.data import benchmarks"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "calm-2",
   "metadata": {},
   "source": [
    "## Data\n",
    "\n",
    "`linear_gauss_er` draws an Erdos-Renyi graph with linear mechanisms and Gaussian noise: here 10\n",
    "variables and 100 samples, ten per variable, the size the 0.7.0 benchmark measured."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "calm-3",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "100 observations, 10 variables, seed 7\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<svg xmlns=\"http://www.w3.org/2000/svg\" viewBox=\"0 0 464 332\" width=\"100%\" font-family=\"Inter, system-ui, sans-serif\"><rect width=\"464\" height=\"332\" fill=\"#fafbfc\"/><line x1=\"312.82\" y1=\"44.82\" x2=\"242.84\" y2=\"91.09\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"321.16,39.30 310.06,40.65 315.57,48.99\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"144.93\" y1=\"163.00\" x2=\"409.13\" y2=\"134.59\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"419.07,133.52 408.60,129.62 409.67,139.56\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"113.24\" y1=\"177.70\" x2=\"42.60\" y2=\"227.83\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"121.40,171.91 110.35,173.62 116.14,181.78\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"112.82\" y1=\"151.70\" x2=\"42.84\" y2=\"105.43\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"121.16,157.22 115.57,147.53 110.06,155.87\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"409.32\" y1=\"128.29\" x2=\"244.82\" y2=\"100.43\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"419.18,129.96 410.16,123.36 408.49,133.22\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"410.18\" y1=\"124.86\" x2=\"144.33\" y2=\"36.24\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"419.67,128.02 411.76,120.11 408.60,129.60\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"410.11\" y1=\"139.19\" x2=\"344.37\" y2=\"160.40\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"419.63,136.12 408.58,134.43 411.65,143.95\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"139.66\" y1=\"224.85\" x2=\"218.45\" y2=\"116.84\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"224.34,108.76 214.41,113.89 222.49,119.79\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"142.92\" y1=\"242.40\" x2=\"212.67\" y2=\"287.40\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"221.08,292.82 215.38,283.20 209.96,291.60\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"109.00\" y1=\"235.35\" x2=\"45.00\" y2=\"235.35\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"119.00,235.35 109.00,230.35 109.00,240.35\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"212.82\" y1=\"85.57\" x2=\"142.84\" y2=\"39.30\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"221.16,91.09 215.57,81.40 210.06,89.74\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"242.84\" y1=\"105.43\" x2=\"312.82\" y2=\"151.70\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"321.16,157.22 315.57,147.53 310.06,155.87\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"209.00\" y1=\"98.26\" x2=\"45.00\" y2=\"98.26\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"219.00,98.26 209.00,93.26 209.00,103.26\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"112.82\" y1=\"44.82\" x2=\"42.84\" y2=\"91.09\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"121.16,39.30 110.06,40.65 115.57,48.99\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"210.11\" y1=\"292.81\" x2=\"44.37\" y2=\"239.35\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"219.63,295.88 211.65,288.05 208.58,297.57\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><circle cx=\"332.00\" cy=\"32.13\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"332.00\" y=\"32.13\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">0</text><circle cx=\"132.00\" cy=\"164.39\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"164.39\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">1</text><circle cx=\"432.00\" cy=\"132.13\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"432.00\" y=\"132.13\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">2</text><circle cx=\"132.00\" cy=\"235.35\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"235.35\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">3</text><circle cx=\"232.00\" cy=\"98.26\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"232.00\" y=\"98.26\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">4</text><circle cx=\"132.00\" cy=\"32.13\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"32.13\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">5</text><circle cx=\"232.00\" cy=\"299.87\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"232.00\" y=\"299.87\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">6</text><circle cx=\"32.00\" cy=\"235.35\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"32.00\" y=\"235.35\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">7</text><circle cx=\"332.00\" cy=\"164.39\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"332.00\" y=\"164.39\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">8</text><circle cx=\"32.00\" cy=\"98.26\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"32.00\" y=\"98.26\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">9</text></svg>"
      ],
      "text/plain": [
       "'<svg xmlns=\"http://www.w3.org/2000/svg\" viewBox=\"0 0 464 332\" width=\"100%\" font-family=\"Inter, system-ui, sans-serif\"><rect width=\"464\" height=\"332\" fill=\"#fafbfc\"/><line x1=\"312.82\" y1=\"44.82\" x2=\"242.84\" y2=\"91.09\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"321.16,39.30 310.06,40.65 315.57,48.99\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"144.93\" y1=\"163.00\" x2=\"409.13\" y2=\"134.59\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"419.07,133.52 408.60,129.62 409.67,139.56\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"113.24\" y1=\"177.70\" x2=\"42.60\" y2=\"227.83\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"121.40,171.91 110.35,173.62 116.14,181.78\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"112.82\" y1=\"151.70\" x2=\"42.84\" y2=\"105.43\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"121.16,157.22 115.57,147.53 110.06,155.87\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"409.32\" y1=\"128.29\" x2=\"244.82\" y2=\"100.43\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"419.18,129.96 410.16,123.36 408.49,133.22\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"410.18\" y1=\"124.86\" x2=\"144.33\" y2=\"36.24\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"419.67,128.02 411.76,120.11 408.60,129.60\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"410.11\" y1=\"139.19\" x2=\"344.37\" y2=\"160.40\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"419.63,136.12 408.58,134.43 411.65,143.95\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"139.66\" y1=\"224.85\" x2=\"218.45\" y2=\"116.84\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"224.34,108.76 214.41,113.89 222.49,119.79\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"142.92\" y1=\"242.40\" x2=\"212.67\" y2=\"287.40\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"221.08,292.82 215.38,283.20 209.96,291.60\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"109.00\" y1=\"235.35\" x2=\"45.00\" y2=\"235.35\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"119.00,235.35 109.00,230.35 109.00,240.35\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"212.82\" y1=\"85.57\" x2=\"142.84\" y2=\"39.30\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"221.16,91.09 215.57,81.40 210.06,89.74\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"242.84\" y1=\"105.43\" x2=\"312.82\" y2=\"151.70\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"321.16,157.22 315.57,147.53 310.06,155.87\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"209.00\" y1=\"98.26\" x2=\"45.00\" y2=\"98.26\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"219.00,98.26 209.00,93.26 209.00,103.26\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"112.82\" y1=\"44.82\" x2=\"42.84\" y2=\"91.09\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"121.16,39.30 110.06,40.65 115.57,48.99\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"210.11\" y1=\"292.81\" x2=\"44.37\" y2=\"239.35\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"219.63,295.88 211.65,288.05 208.58,297.57\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><circle cx=\"332.00\" cy=\"32.13\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"332.00\" y=\"32.13\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">0</text><circle cx=\"132.00\" cy=\"164.39\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"164.39\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">1</text><circle cx=\"432.00\" cy=\"132.13\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"432.00\" y=\"132.13\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">2</text><circle cx=\"132.00\" cy=\"235.35\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"235.35\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">3</text><circle cx=\"232.00\" cy=\"98.26\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"232.00\" y=\"98.26\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">4</text><circle cx=\"132.00\" cy=\"32.13\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"32.13\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">5</text><circle cx=\"232.00\" cy=\"299.87\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"232.00\" y=\"299.87\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">6</text><circle cx=\"32.00\" cy=\"235.35\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"32.00\" y=\"235.35\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">7</text><circle cx=\"332.00\" cy=\"164.39\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"332.00\" y=\"164.39\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">8</text><circle cx=\"32.00\" cy=\"98.26\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"32.00\" y=\"98.26\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">9</text></svg>'"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "seed = 7\n",
    "dataset = benchmarks.scm(\"linear_gauss_er\", d=10).sample(n=100, seed=seed, scale=\"standardize\")\n",
    "X, truth = dataset.data, dataset.observed_graph\n",
    "print(f\"{X.shape[0]} observations, {X.shape[1]} variables, seed {seed}\")\n",
    "display(andrey.viz.draw(truth))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "calm-4",
   "metadata": {},
   "source": [
    "## Fit CALM\n",
    "\n",
    "`andrey.calm` returns a weighted DAG. Its optimizer draws random numbers and is seeded on every\n",
    "call, so a run repeats on one machine. The fit takes 10 to 30 seconds on a CPU."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "calm-5",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "dag, 10 edges, seed 0\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<svg xmlns=\"http://www.w3.org/2000/svg\" viewBox=\"0 0 364 325\" width=\"100%\" font-family=\"Inter, system-ui, sans-serif\"><rect width=\"364\" height=\"325\" fill=\"#fafbfc\"/><line x1=\"110.32\" y1=\"39.55\" x2=\"44.25\" y2=\"62.99\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"119.75,36.20 108.65,34.83 112.00,44.26\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"242.60\" y1=\"214.65\" x2=\"313.24\" y2=\"164.52\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"321.40,158.73 310.35,160.44 316.14,168.60\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"245.00\" y1=\"222.18\" x2=\"309.00\" y2=\"222.18\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"319.00,222.18 309.00,217.18 309.00,227.18\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"309.00\" y1=\"151.21\" x2=\"245.00\" y2=\"151.21\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"319.00,151.21 309.00,146.21 309.00,156.21\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"309.64\" y1=\"145.80\" x2=\"144.64\" y2=\"105.88\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"319.36,148.15 310.82,140.94 308.47,150.66\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"144.12\" y1=\"185.23\" x2=\"210.55\" y2=\"159.51\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"219.88,155.90 208.75,154.85 212.36,164.18\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"141.05\" y1=\"199.26\" x2=\"216.00\" y2=\"276.63\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"222.95,283.81 219.59,273.15 212.41,280.10\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"44.25\" y1=\"71.69\" x2=\"110.32\" y2=\"95.13\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"119.75,98.48 112.00,90.42 108.65,99.84\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"211.30\" y1=\"141.19\" x2=\"143.70\" y2=\"108.48\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"220.30,145.55 213.47,136.69 209.12,145.69\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"145.00\" y1=\"293.15\" x2=\"209.00\" y2=\"293.15\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"219.00,293.15 209.00,288.15 209.00,298.15\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><circle cx=\"132.00\" cy=\"31.85\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"31.85\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">0</text><circle cx=\"232.00\" cy=\"222.18\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"232.00\" y=\"222.18\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">1</text><circle cx=\"332.00\" cy=\"151.21\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"332.00\" y=\"151.21\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">2</text><circle cx=\"132.00\" cy=\"189.92\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"189.92\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">3</text><circle cx=\"32.00\" cy=\"67.34\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"32.00\" y=\"67.34\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">4</text><circle cx=\"232.00\" cy=\"151.21\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"232.00\" y=\"151.21\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">5</text><circle cx=\"132.00\" cy=\"293.15\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"293.15\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">6</text><circle cx=\"232.00\" cy=\"293.15\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"232.00\" y=\"293.15\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">7</text><circle cx=\"132.00\" cy=\"102.82\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"102.82\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">8</text><circle cx=\"332.00\" cy=\"222.18\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"332.00\" y=\"222.18\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">9</text></svg>"
      ],
      "text/plain": [
       "'<svg xmlns=\"http://www.w3.org/2000/svg\" viewBox=\"0 0 364 325\" width=\"100%\" font-family=\"Inter, system-ui, sans-serif\"><rect width=\"364\" height=\"325\" fill=\"#fafbfc\"/><line x1=\"110.32\" y1=\"39.55\" x2=\"44.25\" y2=\"62.99\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"119.75,36.20 108.65,34.83 112.00,44.26\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"242.60\" y1=\"214.65\" x2=\"313.24\" y2=\"164.52\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"321.40,158.73 310.35,160.44 316.14,168.60\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"245.00\" y1=\"222.18\" x2=\"309.00\" y2=\"222.18\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"319.00,222.18 309.00,217.18 309.00,227.18\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"309.00\" y1=\"151.21\" x2=\"245.00\" y2=\"151.21\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"319.00,151.21 309.00,146.21 309.00,156.21\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"309.64\" y1=\"145.80\" x2=\"144.64\" y2=\"105.88\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"319.36,148.15 310.82,140.94 308.47,150.66\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"144.12\" y1=\"185.23\" x2=\"210.55\" y2=\"159.51\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"219.88,155.90 208.75,154.85 212.36,164.18\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"141.05\" y1=\"199.26\" x2=\"216.00\" y2=\"276.63\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"222.95,283.81 219.59,273.15 212.41,280.10\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"44.25\" y1=\"71.69\" x2=\"110.32\" y2=\"95.13\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"119.75,98.48 112.00,90.42 108.65,99.84\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"211.30\" y1=\"141.19\" x2=\"143.70\" y2=\"108.48\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"220.30,145.55 213.47,136.69 209.12,145.69\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><line x1=\"145.00\" y1=\"293.15\" x2=\"209.00\" y2=\"293.15\" stroke=\"#3a4658\" stroke-width=\"2.0\" stroke-linecap=\"round\"/><polygon points=\"219.00,293.15 209.00,288.15 209.00,298.15\" fill=\"#3a4658\" stroke-linejoin=\"round\"/><circle cx=\"132.00\" cy=\"31.85\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"31.85\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">0</text><circle cx=\"232.00\" cy=\"222.18\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"232.00\" y=\"222.18\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">1</text><circle cx=\"332.00\" cy=\"151.21\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"332.00\" y=\"151.21\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">2</text><circle cx=\"132.00\" cy=\"189.92\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"189.92\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">3</text><circle cx=\"32.00\" cy=\"67.34\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"32.00\" y=\"67.34\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">4</text><circle cx=\"232.00\" cy=\"151.21\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"232.00\" y=\"151.21\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">5</text><circle cx=\"132.00\" cy=\"293.15\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"293.15\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">6</text><circle cx=\"232.00\" cy=\"293.15\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"232.00\" y=\"293.15\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">7</text><circle cx=\"132.00\" cy=\"102.82\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"132.00\" y=\"102.82\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">8</text><circle cx=\"332.00\" cy=\"222.18\" r=\"13.00\" fill=\"#ffffff\" stroke=\"#182233\" stroke-width=\"1.8\"/><text x=\"332.00\" y=\"222.18\" text-anchor=\"middle\" font-size=\"11\" font-weight=\"600\" fill=\"#182233\" dominant-baseline=\"central\">9</text></svg>'"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "out = andrey.calm(X)\n",
    "print(f\"{out.structure.kind}, {out.metadata['n_edges']} edges, seed {out.metadata['seed']}\")\n",
    "display(andrey.viz.draw(out.structure))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "calm-6",
   "metadata": {},
   "source": [
    "## Check accuracy\n",
    "\n",
    "SHD counts errors against the true graph; lower is better. Gaussian data identify the graph only up\n",
    "to its equivalence class, so `mec_shd` compares the two CPDAGs. The empty graph's SHD is the number\n",
    "of true edges, so an estimate at or above it has not improved on returning no edges.\n",
    "One dataset can land on either side: on the 0.7.0 benchmark's five datasets of this size,\n",
    "CALM's median was 6 errors against the empty graph's 10."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "calm-7",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CALM SHD: 14\n",
      "Empty graph SHD: 15\n"
     ]
    }
   ],
   "source": [
    "d = truth.n_nodes\n",
    "empty = andrey.GraphStructure.from_numpy(np.zeros((d, d), dtype=np.uint8), kind=\"dag\")\n",
    "print(\"CALM SHD:\", metrics.score(out, truth)[\"mec_shd\"])\n",
    "print(\"Empty graph SHD:\", metrics.score(empty, truth)[\"mec_shd\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "calm-8",
   "metadata": {},
   "source": [
    "## Cost and limits\n",
    "\n",
    "Every optimizer step computes a determinant and a matrix power of the full weight matrix, so a\n",
    "step's cost grows with the cube of the variable count. The median fit in the `andrey.calm`\n",
    "docstring took 18 s at 10 variables and 59 s at 50, on 4 cores. The fit is nonconvex and can stop\n",
    "at a local optimum, and one dataset can favor any method: change `seed`, `d`, or `n` and rerun."
   ]
  }
 ],
 "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.14.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
