# Is there an "Algebra of basis-functions"

**URL:** https://discourse.julialang.org/t/is-there-an-algebra-of-basis-functions/131750
**Category:** General Usage
**Tags:** basis
**Created:** [August 21, 2025, 10:44am UTC](https://discourse.julialang.org/t/is-there-an-algebra-of-basis-functions/131750 "2025-08-21T10:44:51Z")
**Posts on this page:** 9
**Page:** 1

<div class="post-metadata">

### Author: ![JADekker](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jadekker/32/210281_2.png) [@JADekker](https://discourse.julialang.org/u/JADekker)
#### Post date: [August 21, 2025, 10:44am UTC](https://discourse.julialang.org/t/is-there-an-algebra-of-basis-functions/131750/1 "2025-08-21T10:44:51Z")

</div>

In my research, I often work with the following kind of problem that essentially reduces to approximating a conditional expectation, say \mathbb{E}[g(X)|X=x], by a linear combination of basis-functions, \mathbb{E}[g(X)|X=x] \approx \sum\_{i} a\_i \psi\_i(x), where the coefficients a\_i are e.g. estimated using ordinary least squares. To apply least-squares, we want to compute for a sample x^1,...,x^n the matrix \mathbf{X} = (\psi\_i(x^n))\_{i,n}, while for pseudo regression, we may wish to compute the matrix \Psi = (\mathbb{E}[\psi\_i(X)\psi\_j(X)])\_{i,j}, etc.

So that raises the question: is there a package that allow you to generate some `Basis` object based on simple rules (such as: take polynomials in X\_1 of degrees 1:d, take Laguerre polynomials in X\_2 of degrees 1:d, and take all products such that the degree of the resulting polynomials are at most of degree d) and then constructs and computes the functions \psi\_i (and possibly the matrix \Psi, if possible)? (This can then be seen as an `algebra of basis-functions` in the sense that we join together sets of basis-functions, take their products, etc.)

For my own research I’ve built quite a bit of code to do just this and although I have spent some time optimizing for speed, I’m not using any GPU things and the like, while I think some kernel abstractions are possible here. I’ve also written API for various uses, e.g.

- to automatically generate a `Basis` that consists of the orthonormal/orthogonal basis functions for a (product)distribution that you specify
- to check that a certain subset of the `Basis` satisfies a certain condition on their moments
- to split the basis to compute only the basis-functions corresponding to certain variables

If there already is a package that does this (ideally with high performance), I’d like to check if I can adapt my API for that and use it in my projects (and perhaps contribute my code there). If there is not yet a package for this, I’d like to consider making my code into a stand-alone package once my PhD-schedule quiets down a bit.

---

<div class="post-metadata">

### Author: ![Tarny\_GG\_Channie](https://avatars.discourse-cdn.com/v4/letter/t/3bc359/32.png) [@Tarny\_GG\_Channie](https://discourse.julialang.org/u/Tarny_GG_Channie)
#### Post date: [August 21, 2025, 11:59am UTC](https://discourse.julialang.org/t/is-there-an-algebra-of-basis-functions/131750/2 "2025-08-21T11:59:23Z")

</div>

The problem is with the sampling. Normal least-square systems were decided with a finite data point. Maybe the function properties offer a nice solution. Should you sample points or do some fancier calculus tricks to get a more accurate least-square with less cost? That depends on a lot of factors like the function chosen. Numerical integration on its own, the art of sampling points and perhaps doing other fancier tricks to get a descriptive picture of a fixed function, let alone fitting the functions with a least-square. I think that there would probably be some nicely constructed functions where you could perform a nice integration and then perform a minimization on the integral itself instead of having to sample. I’m not that good at calculus though.

It depends on whether or not you want a quick “good enough” solution or you want the best and so on.

---

<div class="post-metadata">

### Author: ![JADekker](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jadekker/32/210281_2.png) [@JADekker](https://discourse.julialang.org/u/JADekker)
#### Post date: [August 21, 2025, 12:06pm UTC](https://discourse.julialang.org/t/is-there-an-algebra-of-basis-functions/131750/3 "2025-08-21T12:06:37Z")

</div>

I see your point. My question is actually concerned less with the regression part (which serves more as the motivation for wanting this) and more with the `algebra of basis-functions` part of my post: just like we can make fancy models by composing parts in MTK and nice plots using AlgebraOfGraphics, I was wondering if there is a way to make complicated sets of basis functions from easier ones (and then, if you have a unified API for this, you can think about fitting OLS (simply sampling points and regressing), pseudo-regression (computing the necessary expectations, if possible), or other things).

---

<div class="post-metadata">

### Author: ![langestefan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/langestefan/32/207923_2.png) [@langestefan](https://discourse.julialang.org/u/langestefan)
#### Post date: [August 21, 2025, 2:21pm UTC](https://discourse.julialang.org/t/is-there-an-algebra-of-basis-functions/131750/4 "2025-08-21T14:21:45Z")

</div>

I think StatsModels.jl comes with this out of the box?

```julia
julia> using StatsModels, DataFrames

julia> using StableRNGs; rng = StableRNG(1);

julia> f = @formula(y ~ 1 + a + b + c + b&c)
FormulaTerm
Response:
  y(unknown)
Predictors:
  1
  a(unknown)
  b(unknown)
  c(unknown)
  b(unknown) & c(unknown)

julia> df = DataFrame(y = rand(rng, 9), a = 1:9, b = rand(rng, 9), c = repeat(["d","e","f"], 3))
9×4 DataFrame
 Row │ y a b c
     │ Float64 Int64 Float64 String
─────┼────────────────────────────────────
   1 │ 0.585195 1 0.236782 d
   2 │ 0.0773379 2 0.943741 e
   3 │ 0.716628 3 0.445671 f
   4 │ 0.320357 4 0.763679 d
   5 │ 0.653093 5 0.145071 e
   6 │ 0.236639 6 0.021124 f
   7 │ 0.709684 7 0.152545 d
   8 │ 0.557787 8 0.617492 e
   9 │ 0.05079 9 0.481531 f

julia> f = apply_schema(f, schema(f, df))
FormulaTerm
Response:
  y(continuous)
Predictors:
  1
  a(continuous)
  b(continuous)
  c(DummyCoding:3→2)
  b(continuous) & c(DummyCoding:3→2)

julia> resp, pred = modelcols(f, df);

julia> pred
9×7 Matrix{Float64}:
 1.0 1.0 0.236782 0.0 0.0 0.0 0.0
 1.0 2.0 0.943741 1.0 0.0 0.943741 0.0
 1.0 3.0 0.445671 0.0 1.0 0.0 0.445671
 1.0 4.0 0.763679 0.0 0.0 0.0 0.0
 1.0 5.0 0.145071 1.0 0.0 0.145071 0.0
 1.0 6.0 0.021124 0.0 1.0 0.0 0.021124
 1.0 7.0 0.152545 0.0 0.0 0.0 0.0
 1.0 8.0 0.617492 1.0 0.0 0.617492 0.0
 1.0 9.0 0.481531 0.0 1.0 0.0 0.481531

julia> coefnames(f)
("y", ["(Intercept)", "a", "b", "c: e", "c: f", "b & c: e", "b & c: f"])

```

---

<div class="post-metadata">

### Author: ![JADekker](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jadekker/32/210281_2.png) [@JADekker](https://discourse.julialang.org/u/JADekker)
#### Post date: [August 21, 2025, 3:41pm UTC](https://discourse.julialang.org/t/is-there-an-algebra-of-basis-functions/131750/5 "2025-08-21T15:41:01Z")

</div>

Thanks, I didn’t think about that yet; so you could do all the logic for generating a set of functions (adding sets of functions together, adding cross-products, etc.) on the `FormulaTerm`s, and then parse things along the `FormulaTerm` to do fancy things (e.g. compute expectations when the data follows a product-distribution, etc.).

I’ll keep this in mind and revisit this once I have time to check if this would support all the things that I would want to do and how I can write API for those things!

---

<div class="post-metadata">

### Author: ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)
#### Post date: [August 21, 2025, 5:30pm UTC](https://discourse.julialang.org/t/is-there-an-algebra-of-basis-functions/131750/6 "2025-08-21T17:30:52Z")

</div>

> [@JADekker](#):
>
> So that raises the question: is there a package that allow you to generate some `Basis` object based on simple rules (such as: take polynomials in X\_1 of degrees 1:d, take Laguerre polynomials in X\_2 of degrees 1:d, and take all products such that the degree of the resulting polynomials are at most of degree d)

I use Chebyshev polynomials (with transformations) most of the time. I have coded them with derivatives, Smolyak combinations, etc in

> **[GitHub - tpapp/SpectralKit.jl: Building blocks of spectral methods for Julia.](https://github.com/tpapp/SpectralKit.jl)**
>
> Building blocks of spectral methods for Julia.

---

<div class="post-metadata">

### Author: ![danielwe](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/danielwe/32/35657_2.png) [@danielwe](https://discourse.julialang.org/u/danielwe)
#### Post date: [August 21, 2025, 9:20pm UTC](https://discourse.julialang.org/t/is-there-an-algebra-of-basis-functions/131750/7 "2025-08-21T21:20:16Z")

</div>

You might want to take a look at ApproxFun.jl if you haven’t already. Not sure if it supports everything you need, but it’s the first-stop shop for working with function spaces via orthogonal bases (similar to [chebfun.org](http://chebfun.org) for Matlab, whose thorough and pedagogical documentation you might find helpful even if you’re using Julia and ApproxFun).

---

<div class="post-metadata">

### Author: ![JADekker](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jadekker/32/210281_2.png) [@JADekker](https://discourse.julialang.org/u/JADekker)
#### Post date: [August 22, 2025, 8:59am UTC](https://discourse.julialang.org/t/is-there-an-algebra-of-basis-functions/131750/8 "2025-08-22T08:59:36Z")

</div>

Thank you, that is a good pointer. The way I understand it is that you should provide it with a function already, which is kind of the problem in my case: I’m interested in estimating a conditional expectation by sampled points, so I want to be able to do some form of regression, such as least-squares, and to be able to go up to higher dimensions. For my typical use-case I’m looking at roughly 6 dimensions and in the order of 100-1000 basis-functions. (I know this is not ideal, but the algorithms NEED basis functions with certain properties for this, so fancier approximation methods are off the table for me, unfortunately)

---

<div class="post-metadata">

### Author: ![danielwe](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/danielwe/32/35657_2.png) [@danielwe](https://discourse.julialang.org/u/danielwe)
#### Post date: [August 22, 2025, 9:23am UTC](https://discourse.julialang.org/t/is-there-an-algebra-of-basis-functions/131750/9 "2025-08-22T09:23:50Z")

</div>

If you’re doing least-squares, ApproxFun.jl can help build the Vandermonde matrix for your orthogonal function basis. See [Frequently Asked Questions · ApproxFun.jl](https://juliaapproximation.github.io/ApproxFun.jl/latest/faq/) (the example starts with interpolation, but progresses to regression by the third code block and multivariate regression by the fourth).
