Define Multiple Outputs with TemplateExpressionSpec

[janzenchoi]

Hello,

I would like to define TemplateExpressionSpec for two different outputs with coupled parameters.

For instance, with

y_0 = A*x_0^n
y_1 = A*x_0-n

What would I add here to define templates for both outputs?

template = TemplateExpressionSpec(
    expressions=["f",],
    variable_names=["x_0"],
    parameters={"p": 2}
    combine="""
        A, n = p[1], p[2]
        y_0 = A * x_0 ^ n
        y_1 = A*x_0-n
        <???>
    """
)

Thanks.

[MilesCranmer]

You can set it up so that you pass y_0 and y_1 as features, and have the template returns the per-row loss:

template = TemplateExpressionSpec(
    expressions=["f",],
    variable_names=["x_0", "y_0", "y_1"],
    parameters={"p": 2}
    combine="""
        A, n = p[1], p[2]
        y_0_predicted = A * x_0 ^ n
        y_1_predicted = A * x_0 - n
        
        (y_0 - y_0_predicted)^2 + (y_1 - y_1_predicted)^2
    """
)

then, your PySRRegressor would have a custom loss that simply returns the predicted

model = PySRRegressor(
    expression_spec=template,
    elementwise_loss="my_custom_loss(predicted, target) = predicted"
)

model.fit(X, np.zeros(len(X)))

and pass y_0 and y_1 as the 2nd and 3rd columns of X basically. And then simply pass whatever for y.

[janzenchoi]

Thanks @MilesCranmer! That seems to work well.

Unfortunately though, I keep coming across this error.

juliacall.JuliaError: Element type of `x` is Float64 is different from element type of `y` which is Float32.

The following is a stripped down version of my code, which produces the same error message.

from pysr import PySRRegressor, TemplateExpressionSpec
import numpy as np

expression_spec = TemplateExpressionSpec(
    expressions    = ["f0", "f1"],
    variable_names = ["x0", "x1", "x2", "y0", "y1"],
    combine        = """
        A = 6.3757e-19; n = 7.643720763; M = 6.461e-20; phi = 17.01786425; chi = 5.9171705;
        z0 = A*x1^n * ((1-(phi+1)*M*x1^chi*x0)^((phi+1-n)/(phi+1))-1) / (M*x1^chi*(n-phi-1)) + f0(x0,x1,x2);
        z1 = 1/((phi+1)*M*x1^chi) + f1(x1,x2);
        (y0-z0)^2+(y1-z1)^2
    """
)
regressor = PySRRegressor(
    expression_spec  = expression_spec,
    populations      = 32,
    population_size  = 32,
    maxsize          = 32,
    niterations      = 64,
    binary_operators = ["+", "*", "/"],
    elementwise_loss = "take_first(p, t, w) = p"
)

x_data = np.array([
    [2.9520000e+03, 7.0000000e+01, 8.0000000e+02, 4.1000000e-05, 9.9026640e+06],
    [1.2022200e+06, 7.0000000e+01, 8.0000000e+02, 4.8268000e-02, 9.9026640e+06],
    [2.5709760e+06, 7.0000000e+01, 8.0000000e+02, 8.7537000e-02, 9.9026640e+06],
    [3.9732840e+06, 7.0000000e+01, 8.0000000e+02, 1.2346300e-01, 9.9026640e+06],
    [5.3164080e+06, 7.0000000e+01, 8.0000000e+02, 1.5799800e-01, 9.9026640e+06],
    [6.6475800e+06, 7.0000000e+01, 8.0000000e+02, 1.9575300e-01, 9.9026640e+06],
    [7.9175880e+06, 7.0000000e+01, 8.0000000e+02, 2.4064200e-01, 9.9026640e+06],
    [8.9542080e+06, 7.0000000e+01, 8.0000000e+02, 2.9397400e-01, 9.9026640e+06],
    [9.5383080e+06, 7.0000000e+01, 8.0000000e+02, 3.4689400e-01, 9.9026640e+06],
    [9.9026640e+06, 7.0000000e+01, 8.0000000e+02, 4.6810600e-01, 9.9026640e+06],
], dtype=np.float32)
y_data = np.zeros(x_data.shape[0], dtype=np.float32)
regressor.fit(x_data, y_data)

Do you know why this might be? I've tried forcing it's data type, but it doesn't seem to work.

[MilesCranmer]

You can just set precision=64 in the PySRRegressor. Though eventually I’d like to have this conversion happen automatically.

The issue is that by default, data is fit as Float32. And when you write the number 6.3757e-19, this is actually a Float64. The quick fix is to write it 6.3757f-19 which creates a Float32 instead. (Or just write Float32(6.3757e-19) and the compiler will do it)

But this isn’t obvious to anyone and especially Python native users so there is a patch (1 month+ until available via the normal PyPI installation, since we’re busy with v2.0) that will fix this

[janzenchoi]

That works! Thank you!
Incidentally, is there a way to define the bounds of the parameters during the optimisation?

[MilesCranmer]

For that I would just write a custom loss_function_expression to encourage values to be in bounds. There are some other discussions for that you could check out

[MilesCranmer]

I have since put up an example here: https://ai.damtp.cam.ac.uk/pysr/examples/#13-vector-valued-expressions