Searching for multivariate polynomials: restrict the maximum degree?

[ErikSenn]

I have an application where I know that the data is generated by a multivariate polynomial (e.g. $a \cdot x_0 + b \cdot x_0 \cdot x_0 + c \cdot x_0 \cdot x_1$) up to order $degree_{max}$ (2 in example).
I want to restrict the pysr search space to not include polynomials with higher order than $degree_{max}$ (e.g. for $degree_{max} = 2$, an expression including the term $x_0 \cdot x_0 \cdot x_1$ should not be searched).

Relevant pysr-settings:

binary_operators=["+", "*", "-"]
unary_operators=[]

Now, I do not understand how to set correct constraints / nested_constraints / complexities to restrict the search space.
My idea was to restrict the usage of the "*" operator (note I am not using "^") to only allow for at most max_degree factors (but a constant should be allowed).

Thank you for your help!

[FelixSchloms]

you could restrict the "*" operator to just use a constant or a variable as arguments.
It is documented at: https://astroautomata.com/PySR/options/#constraining-use-of-operators
with
constraints={'mult': (1, 1)}
you restrict the multiplication operator to use a max. complextiy of 1 at each side, which is by default just a variable or a constant.

[zzccchen]

Hi FelixSchloms, if I want to restrict the order of the same feature, how can I do it? For example, x1*x2*x3 is allowed, but x1*x1*x1 is not allowed.

[ErikSenn]

Thank you for the idea.

Setup
So I tested the solution proposed by @FelixSchloms.
Here, I expect that only polynomials up to $degree_{max} = 2$ are allowed.

model = pysr.PySRRegressor(
    niterations=15,  
    binary_operators=['+', '*', '-'],
    complexity_of_operators={'+': 1, '*': 0, '-':1},
    complexity_of_variables=1,
    complexity_of_constants=0,
    constraints={'mult': (1, 1)}, 
    maxsize = 60,
    populations = 3*num_cores, 
    turbo = True # speed up 
)

Note that I set complexity_of_constant = 0 because including the constant should not reduce the possible order of a monomial.
(This is not great in general, but should simplify the problem. Alternatively, I used a low complexity value for the constant and added it to the constraining complexities.)
I also set the complexity_of_operators = {'*' : 0} since the number of multiplications within one monomial should be limited by the constraint.

Issue
I somehow still find solutions on the pareto front that include higher order polynomials; e.g. including the degree 3 monomial
$x_0 \cdot [x_0 \cdot (x_0)] = x_0 \cdot x_0 \cdot x_0$
(I skipped the multiplicative constants for this post as complexity of '*' and constant are 0. The notation is including parenthesis as displayed inmodel.equations_['equation']).

As Miles pointed out, the total complexity in the example with my defintions from above is the sum:
C[x_0] + C[*] + C[x_0] + C[*] + C[x_0] = 1 + 0 + 1 + 0 + 1 = 3

I am now unsure how the constraint works: I would expect that complexities would need to be computed from 'bottom-up' of the computation tree (= resolving inner parts of parenthesis first). The operation $x_0 \cdot [x_0 \cdot x_0]$ would multiply complexity C[x_0] = 1 with complexity ( C[x_0] + C[*] + C[x_0] ) = 1 + 0 + 1 = 2.
From my understanding, this should not be allowed according to the constraint.

I might still make a mistake in the constraint settings, or must be misunderstanding how constraints are enforced.

PS: As I am new to these discussions, let me know if you need more details & if my response should be a comment instead of a reply.

[Edited after comment of Miles Cramner for clarification.]

[MilesCranmer]

Just a quick note – complexity is calculated as a sum. So with variables having complexity 1, and * having complexity 0, x0 * x0 * x0 would be complexity 3 (C[x0] + C[*] + C[x0] + C[*] + C[x0] = 1 + 0 + 1 + 0 + 1 = 3)