AGL - Algebraic Geometric Learning: Grobner Bases, Ideals and Orthogonal Polynomials Inside a Neural Network

In standard neural networks, every hidden layer applies the same two operations:

1. An affine map $x \mapsto Wx + b$
2. A nonlinear activation (ReLU, tanh, sigmoid)

Both are chosen for optimization convenience, not for structure. The weights are opaque tensors — you train them, you use them, and you can’t say what the network means in any mathematical sense. The polynomial neural network asks a different question: what if the network itself were built out of the machinery of commutative algebra — so that its weights are elements of a coordinate ring, its hypothesis space is an algebraic variety, and its learned function comes out as a closed-form equation you can read?

That is what AGL sets out to build.

1. The idea in one sentence

Treat the network’s weights as polynomials over the coordinate ring $\mathbb{R}[x_1, \dots, x_n]$, then let the a system of algebraic constraints — Gröbner bases, ideals, varieties — do the heavy lifting during training.

Because polynomials are closed under composition, any multi-layer polynomial network collapses to one single high-degree polynomial map from inputs to outputs. And that function isn’t hidden inside the tensor stack anymore: AGL can expand it out into a symbolic equation at the end of training.

2. What AGL is

AGL is a PyTorch library that physically embeds algebraic geometry into polynomial neural networks. Instead of “fit a curve”, you ask the network to learn a function inside a constrained algebraic structure — and it hands the structure back to you as arithmetic.

You get:

• Exact symbolic extraction — the global map of the network printed as an equation (e.g. rediscover Kepler’s $T^2 \propto a^3$ as -0.667 + 1.5·log_a).
• Hard constraints, not soft penalties — Gröbner layers reduce weights modulo a constraint ideal to machine precision; the constraint becomes a property of the architecture, not of the optimizer.
• Stability at high degree — orthogonal bases (Chebyshev, Legendre, Hermite, Laguerre) keep the condition number orders of magnitude below the monomial basis, so degree-10+ network still train without divergence.
• Physics-grade interpretability — on the chaotic Lorenz system, an OrthoPolyNN both outperforms a deep ReLU MLP at 1000-step rollout and returns dz/dt ≈ 0.0099·xy − 0.0145·x + 0.0087·y + 0.9737·z − 0.0081.

flowchart LR
    DATA["real-world data<br/>(physics, engineering, finance)"] --> PX["polynomial layer<br/>x -> W(x)·x + b"]
    DATA --> OX["orthogonal layer<br/>Chebyshev / Legendre / Hermite"]
    OX --> GB["Groebner / variety layer<br/>exact ideal constraints"]
    GB --> EQ["symbolic extraction<br/>get_global_equation()"]
    EQ --> YOU["a readable closed-form law"]
    style GB fill:#4a148c,color:#fff

3. Who is this for?

4. The journey — documentation and tutorials

The docs develop the mathematics in dependency order, from definitions to working code, then apply the library to six real datasets:

5. Project links

6. Quick start — from zero to a closed-form equation

Requires: Python 3.10+, PyTorch, and the usual scientific stack. Works with a conda env (e.g. conda activate aisystem).

git clone https://github.com/ajeetkbhardwaj/agl && cd agl
pip install -r requirements.txt
pytest            # run the test suite

Step 1 — Symbolic regression / law discovery

python examples/symbolic_regression.py

Fits a monomial-basis PNN to a known cubic and prints the recovered polynomial — e.g. 3.269x³ − 1.722x² + 2.097x − 0.929 against the true 3.5x³ − 1.5x² + 2x − 1.

Step 2 — Exact coefficient extraction

python examples/equation_extraction.py

Recovers Chebyshev coefficients exactly to printed precision.

Step 3 — Chaotic dynamics (the fun one)

python examples/lorenz_attractor.py

Trains a degree-1 multilinear OrthoPolyNN on the Lorenz map and a deep ReLU MLP, rolls both forward 1000 steps, and prints the extracted equations alongside the mean trajectory error:

A symbolic Lorenz system, recovered from data, that beats the black-box MLP at long-horizon extrapolation.

Step 4 — More real-data examples

python examples/california_housing.py   # OrthoPolyNN R² 0.769 vs MLP 0.739
python examples/solve_ode.py            # solve an ODE with polynomial layers

7. Under the hood

flowchart TB
    subgraph AG["src.alggeom — the commutative algebra"]
        POL["polynomial / monomial orders"] --> GB2["GroebnerBasis (Buchberger)"]
        GB2 --> ID["Ideal (membership, normal form)"]
        ID --> V["AlgebraicVariety (dimension, Jacobian)"]
    end
    subgraph NN["src.nn + src.layers — the networks"]
        N["PolynomialNeuralNetwork<br/>global equation extraction"] --> O["OrthoPolyNetwork<br/>Chebyshev · Legendre · Hermite · Laguerre"]
        O --> R["RationalPolyLayer<br/>Pade-style poles"]
        O --> G["GroebnerLayer / VarietyLayer<br/>hard constraints"]
    end
    subgraph TR["src.training — the losses"]
        IDM["IdealMembershipLoss"] --> SA["Optimizers (sparsity, spectral reg)"]
    end
    AG --> NN
    NN --> TR

Everything is PyTorch-native and autograd-compatible; the symbolic side (SymPy) is called only at extraction time, so training stays fast and the exact equation comes out at the end.

8. What the experiments actually show

We were careful to report results as they are, not as a marketing story.

• Where AGL wins clearly — chaotic dynamics (Lorenz rollout, ~40% lower error than a deep ReLU MLP plus closed-form equations), exact constraint enforcement (Gröbner layer satisfies the constraint to ~1e-7, a soft penalty cannot), and stability at degree ≥ 10 where the monomial basis diverges.
• Where AGL is competitive — tabular regression on Airfoil (R² 0.953 vs 0.887) and Power Plant; exact recovery on symbolically structured problems.
• Where AGL loses, honestly — Concrete and California Housing (MLP still wins), and it trains slower than a plain MLP (the algebraic bookkeeping has a real cost). Scalability of the algebraic layers is open work.

The takeaway, stated plainly: use AGL where the problem has polynomial or constraint-bearing structure — it pays off there. For generic tabular regression it is not a drop-in replacement for an MLP.

9. Where to go next

1. Read the theory — the Algebraic Geometry and Neural Networks parts (§4) build everything from definitions to code.
2. Run the tutorials — each of the six real-data tutorials ships a runnable script under experiments/.
3. Try the chaos example — examples/lorenz_attractor.py is the best single demonstration of what algebraic structure buys you.
4. Contribute — the most valuable direction is making the Groebner/variety layers faster; the framework’s correctness is there, scalability is the gap.