{
 "nbformat": 4,
 "nbformat_minor": 5,
 "metadata": {
  "kernelspec": {
   "name": "python3",
   "display_name": "Python 3",
   "language": "python"
  },
  "language_info": {
   "name": "python"
  },
  "colab": {
   "name": "matrix.ipynb"
  }
 },
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": "# matrix \u2014 Python demo\n\nNumerical companion to the entry [matrix](https://dictionaryofml.org/terms/matrix.html) of the [Dictionary of Applied Machine Learning](https://dictionaryofml.org/): it recomputes what the entry states and prints one line per check.\n\nOne block per paragraph of the entry (marked [P...]): each block verifies numerically what the corresponding statement asserts. Self-contained (numpy/matplotlib only), fixed seed.\n\nRequires NumPy and Matplotlib only, and uses fixed seeds, so the printed numbers reproduce exactly. Generated from [`pythondemos/matrix.py`](https://dictionaryofml.org/terms/matrix.py); CC BY 4.0."
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": "# Notebook shim: the script resolves output paths relative to __file__,\n# which a notebook kernel does not define; everything lands in the\n# working directory instead.\nimport os\n__file__ = os.path.join(os.getcwd(), \"matrix.py\")\nos.makedirs(\"pythondemos\", exist_ok=True)"
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": "\"\"\"\nmatrix.py \u2014 numerical companion to the glossary entry 'matrix'.\n\nOne block per paragraph of the entry (marked [P...]): each block verifies\nnumerically what the corresponding statement asserts. Self-contained\n(numpy/matplotlib only), fixed seed.\n\nBlocks\n------\n[P-featuremtx] Stacking the feature vectors of m data points row-wise\n               yields the m x d feature matrix: row r of X is x^(r) and\n               entry X[r, j] is feature j of data point r.\n[P-linsys]     A matrix represents a system of linear equations A w = y;\n               the normal equations X^T X w = X^T y of least-squares\n               linear regression are one instance \u2014 their solution\n               matches the least-squares fit.\n[P-linearmap]  A matrix defines a linear map: the image of basis vector\n               u^(j) is the linear combination sum_r A[r, j] v^(r) of the\n               target basis, and the map is additive and homogeneous.\n[P-array]      A matrix is the order-2 special case of an array: its\n               numpy representation has exactly two axes, while a scalar,\n               a vector, an image, and a stack of images have order\n               0, 1, 3, and 4.\n\nOutputs\n-------\nmatrix.png : preview figure (checking only).\n\nData generated by pythondemos/matrix.py.\n\"\"\"\n\nimport numpy as np\nimport matplotlib\n\nmatplotlib.use(\"Agg\")\nimport matplotlib.pyplot as plt\nfrom pathlib import Path\n\nOUT_DIR = Path(__file__).parent\n\nrng = np.random.default_rng(42)\nreport = []\n\n\ndef check(name, ok):\n    report.append((name, bool(ok)))\n    print(f\"  [{'ok' if ok else 'FAIL'}] {name}\")"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": "**[P-featuremtx]** Stacking the feature vectors of m data points row-wise yields the m x d feature matrix: row r of X is x^(r) and entry X[r, j] is feature j of data point r."
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": "print(\"[P-featuremtx] feature matrix stacks feature vectors row-wise\")\nm, d = 6, 3\nfeature_vecs = [rng.normal(size=d) for _ in range(m)]\nX = np.stack(feature_vecs, axis=0)\ncheck(\"shape is m x d\", X.shape == (m, d))\ncheck(\"row r equals x^(r)\",\n      all(np.array_equal(X[r], feature_vecs[r]) for r in range(m)))\ncheck(\"entry X[r, j] is feature j of data point r\",\n      X[2, 1] == feature_vecs[2][1])\n# the printed A_{i,j} convention: entry in row i, column j\nA_conv = np.array([[10 * i + j for j in range(1, 4)] for i in range(1, 3)])\ncheck(\"A[i, j] sits in row i and column j (A_{1,2} = 12, A_{2,3} = 23)\",\n      A_conv[0, 1] == 12 and A_conv[1, 2] == 23)"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": "**[P-linsys]** A matrix represents a system of linear equations A w = y; the normal equations X^T X w = X^T y of least-squares linear regression are one instance \u2014 their solution matches the least-squares fit."
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": "print(\"[P-linsys] matrices represent linear systems (normal equations)\")\ny = X @ np.array([1.0, -2.0, 0.5]) + 0.1 * rng.normal(size=m)\nw = np.linalg.solve(X.T @ X, X.T @ y)          # normal equations\nw_lstsq = np.linalg.lstsq(X, y, rcond=None)[0]\ncheck(\"normal-equation solution equals the least-squares fit\",\n      np.allclose(w, w_lstsq))\ncheck(\"residual satisfies X^T (y - X w) = 0\",\n      np.max(np.abs(X.T @ (y - X @ w))) < 1e-10)"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": "**[P-linearmap]** A matrix defines a linear map: the image of basis vector u^(j) is the linear combination sum_r A[r, j] v^(r) of the target basis, and the map is additive and homogeneous."
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": "print(\"[P-linearmap] a matrix defines a linear map\")\nA = rng.normal(size=(4, 3))\nU = np.eye(3)                                   # basis of the domain\nVb = np.eye(4)                                  # basis of the codomain\nfor j in range(3):\n    img = A @ U[:, j]\n    combo = sum(A[r, j] * Vb[:, r] for r in range(4))\n    check(f\"image of u^({j+1}) is sum_r A[r,{j+1}] v^(r)\",\n          np.allclose(img, combo))\nu1, u2, a = rng.normal(size=3), rng.normal(size=3), 2.7\ncheck(\"additivity: A(u + u') = A u + A u'\",\n      np.allclose(A @ (u1 + u2), A @ u1 + A @ u2))\ncheck(\"homogeneity: A(a u) = a A u\", np.allclose(A @ (a * u1), a * (A @ u1)))"
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": "**[P-array]** A matrix is the order-2 special case of an array: its numpy representation has exactly two axes, while a scalar, a vector, an image, and a stack of images have order 0, 1, 3, and 4."
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": "print(\"[P-array] a matrix is an order-2 array\")\nscalar = np.float64(3.0)\nvector = rng.normal(size=5)\nimage = rng.normal(size=(32, 32, 3))            # H x W x C\nbatch = rng.normal(size=(8, 32, 32, 3))         # a stack of 8 images\ncheck(\"matrix has exactly two axes\", X.ndim == 2)\ncheck(\"scalar / vector / image / image stack have order 0 / 1 / 3 / 4\",\n      (scalar.ndim, vector.ndim, image.ndim, batch.ndim) == (0, 1, 3, 4))\n\n# ------------------------------------------------------------ preview\nfig, ax = plt.subplots(figsize=(4.5, 3.2))\nim = ax.imshow(X, cmap=\"gray\", aspect=\"auto\")\nax.set_xlabel(\"feature j\"); ax.set_ylabel(\"data point r\")\nax.set_title(\"[P-featuremtx] feature matrix X (m x d)\")\nfig.colorbar(im, ax=ax, shrink=0.8)\nfig.tight_layout()\nfig.savefig(OUT_DIR / \"matrix.png\", dpi=110)\nprint(f\"\\n{sum(ok for _, ok in report)}/{len(report)} checks passed\")\nassert all(ok for _, ok in report)"
  }
 ]
}