{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "#"
      ],
      "id": "28b3f807-ab68-49a9-a4ca-7623651ecce6"
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "metadata": {},
      "outputs": [],
      "source": [
        " # Get set up, but we'll start on the board today\n",
        "import polars as pl\n",
        "import scipy as sp\n",
        "import numpy as np"
      ],
      "id": "b69b2efe-61be-4d29-a60a-f8bd06455017"
    },
    {
      "cell_type": "code",
      "execution_count": 6,
      "metadata": {},
      "outputs": [],
      "source": [
        "# change penguins to carnivora\n",
        "url = \"https://roualdes.sfo3.digitaloceanspaces.com/data/carnivora.csv\"\n",
        "df = pl.read_csv(url, null_values = [\"NA\"]).select(\"SuperFamily\", \"BW\").drop_nulls()"
      ],
      "id": "1689e6bb-b9ac-453e-9aeb-c27a7733cf78"
    },
    {
      "cell_type": "code",
      "execution_count": 13,
      "metadata": {},
      "outputs": [],
      "source": [
        "x1 = df.select(pl.col(\"SuperFamily\") == \"Feliformia\").cast(pl.Float64).to_numpy()\n",
        "X = np.c_[np.ones_like(df[:, 1]), x1]"
      ],
      "id": "b2e89148-f2e8-4086-b7a7-22500d2c7e91"
    },
    {
      "cell_type": "code",
      "execution_count": 18,
      "metadata": {},
      "outputs": [],
      "source": [
        "def ll_lm(theta, data):\n",
        "    X = data[\"X\"]\n",
        "    y = data[\"y\"]\n",
        "    N = np.size(y)\n",
        "    lm = np.zeros_like(y)\n",
        "    for n in range(N):\n",
        "        lm[n] = np.sum(X[n] * theta)\n",
        "    d = y - lm\n",
        "    return np.sum(d * d)"
      ],
      "id": "7caf3fd6-daa0-4345-8a45-3c538d7a1950"
    },
    {
      "cell_type": "code",
      "execution_count": 22,
      "metadata": {},
      "outputs": [
        {
          "output_type": "display_data",
          "metadata": {},
          "data": {
            "text/plain": [
              "np.float64(351.7162524773699)"
            ]
          }
        }
      ],
      "source": [
        "init = sp.stats.norm().rvs(size = 2)\n",
        "d = {\"X\": X, \"y\": df[\"BW\"].to_numpy()}\n",
        "o = sp.optimize.minimize(ll_lm, init, args = (d,), method = \"L-BFGS-B\")\n",
        "np.sum(o[\"x\"])"
      ],
      "id": "67f2912f-5088-41bc-9fb4-85b94a1394b3"
    },
    {
      "cell_type": "code",
      "execution_count": 24,
      "metadata": {},
      "outputs": [
        {
          "output_type": "display_data",
          "metadata": {},
          "data": {
            "text/html": [
              "</div>"
            ]
          }
        }
      ],
      "source": [
        "(df\n",
        " .group_by(\"SuperFamily\")\n",
        " .agg(pl.col(\"BW\").mean().alias(\"mean BW\"))\n",
        ")"
      ],
      "id": "ac601ce7-c55b-4a4d-9142-20c4073d6cbc"
    }
  ],
  "nbformat": 4,
  "nbformat_minor": 5,
  "metadata": {
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3 (ipykernel)",
      "language": "python"
    },
    "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.12.16"
    }
  }
}