API

EquationSearch(X, y[; kws...])

Perform a distributed equation search for functions f_i which describe the mapping f_i(X[:, j]) ≈ y[i, j]. Options are configured using SymbolicRegression.Options(...), which should be passed as a keyword argument to options. One can turn off parallelism with numprocs=0, which is useful for debugging and profiling.

Arguments

• X::AbstractMatrix{T}: The input dataset to predict y from. The first dimension is features, the second dimension is rows.
• y::Union{AbstractMatrix{T}, AbstractVector{T}}: The values to predict. The first dimension is the output feature to predict with each equation, and the second dimension is rows.
• niterations::Int=10: The number of iterations to perform the search. More iterations will improve the results.
• weights::Union{AbstractMatrix{T}, AbstractVector{T}, Nothing}=nothing: Optionally weight the loss for each y by this value (same shape as y).
• varMap::Union{Array{String, 1}, Nothing}=nothing: The names of each feature in X, which will be used during printing of equations.
• options::Options=Options(): The options for the search, such as which operators to use, evolution hyperparameters, etc.
• numprocs::Union{Int, Nothing}=nothing: The number of processes to use, if you want EquationSearch to set this up automatically. By default this will be 4, but can be any number (you should pick a number <= the number of cores available).
• procs::Union{Array{Int, 1}, Nothing}=nothing: If you have set up a distributed run manually with procs = addprocs() and @everywhere, pass the procs to this keyword argument.
• multithreading::Bool=false: Whether to use multithreading. Otherwise, will use multiprocessing. Multithreading uses less memory, but multiprocessing can handle multi-node compute.
• runtests::Bool=true: Whether to run (quick) tests before starting the search, to see if there will be any problems during the equation search related to the host environment.
• saved_state::Union{StateType, Nothing}=nothing: If you have already run EquationSearch and want to resume it, pass the state here. To get this to work, you need to have stateReturn=true in the options, which will cause EquationSearch to return the state. Note that you cannot change the operators or dataset, but most other options should be changeable.
• addprocs_function::Union{Function, Nothing}=nothing: If using distributed mode (multithreading=false), you may pass a custom function to use instead of addprocs. This function should take a single positional argument, which is the number of processes to use, as well as the lazy keyword argument. For example, if set up on a slurm cluster, you could pass addprocs_function = addprocs_slurm, which will set up slurm processes.

Returns

• hallOfFame::HallOfFame: The best equations seen during the search. hallOfFame.members gives an array of PopMember objects, which have their tree (equation) stored in .tree. Their score (loss) is given in .score. The array of PopMember objects is enumerated by size from 1 to options.maxsize.

Options(;kws...)

Construct options for EquationSearch and other functions. The current arguments have been tuned using the median values from https://github.com/MilesCranmer/PySR/discussions/115.

Arguments

• binary_operators: Tuple of binary operators to use. Each operator should be defined for two input scalars, and one output scalar. All operators need to be defined over the entire real line (excluding infinity - these are stopped before they are input). Thus, log should be replaced with log_abs, etc. For speed, define it so it takes two reals of the same type as input, and outputs the same type. For the SymbolicUtils simplification backend, you will need to define a generic method of the operator so it takes arbitrary types.
• unary_operators: Same, but for unary operators (one input scalar, gives an output scalar).
• constraints: Array of pairs specifying size constraints for each operator. The constraints for a binary operator should be a 2-tuple (e.g., (-1, -1)) and the constraints for a unary operator should be an Int. A size constraint is a limit to the size of the subtree in each argument of an operator. e.g., [(^)=>(-1, 3)] means that the ^ operator can have arbitrary size (-1) in its left argument, but a maximum size of 3 in its right argument. Default is no constraints.
• batching: Whether to evolve based on small mini-batches of data, rather than the entire dataset.
• batchSize: What batch size to use if using batching.
• loss: What loss function to use. Can be one of the following losses, or any other loss of type SupervisedLoss. You can also pass a function that takes a scalar target (left argument), and scalar predicted (right argument), and returns a scalar. This will be averaged over the predicted data. If weights are supplied, your function should take a third argument for the weight scalar. Included losses: Regression: - LPDistLoss{P}(), - L1DistLoss(), - L2DistLoss() (mean square), - LogitDistLoss(), - HuberLoss(d), - L1EpsilonInsLoss(ϵ), - L2EpsilonInsLoss(ϵ), - PeriodicLoss(c), - QuantileLoss(τ), Classification: - ZeroOneLoss(), - PerceptronLoss(), - L1HingeLoss(), - SmoothedL1HingeLoss(γ), - ModifiedHuberLoss(), - L2MarginLoss(), - ExpLoss(), - SigmoidLoss(), - DWDMarginLoss(q).
• npopulations: How many populations of equations to use. By default this is set equal to the number of cores
• npop: How many equations in each population.
• ncyclesperiteration: How many generations to consider per iteration.
• ns: Number of equations in each subsample during regularized evolution.
• topn: Number of equations to return to the host process, and to consider for the hall of fame.
• complexity_of_operators: What complexity should be assigned to each operator, and the occurrence of a constant or variable. By default, this is 1 for all operators. Can be a real number as well, in which case the complexity of an expression will be rounded to the nearest integer. Input this in the form of, e.g., [(^) => 3, sin => 2].
• complexity_of_constants: What complexity should be assigned to use of a constant. By default, this is 1.
• complexity_of_variables: What complexity should be assigned to each variable. By default, this is 1.
• alpha: The probability of accepting an equation mutation during regularized evolution is given by exp(-delta_loss/(alpha * T)), where T goes from 1 to 0. Thus, alpha=infinite is the same as no annealing.
• maxsize: Maximum size of equations during the search.
• maxdepth: Maximum depth of equations during the search, by default this is set equal to the maxsize.
• parsimony: A multiplicative factor for how much complexity is punished.
• useFrequency: Whether to use a parsimony that adapts to the relative proportion of equations at each complexity; this will ensure that there are a balanced number of equations considered for every complexity.
• useFrequencyInTournament: Whether to use the adaptive parsimony described above inside the score, rather than just at the mutation accept/reject stage.
• fast_cycle: Whether to thread over subsamples of equations during regularized evolution. Slightly improves performance, but is a different algorithm.
• migration: Whether to migrate equations between processes.
• hofMigration: Whether to migrate equations from the hall of fame to processes.
• fractionReplaced: What fraction of each population to replace with migrated equations at the end of each cycle.
• fractionReplacedHof: What fraction to replace with hall of fame equations at the end of each cycle.
• shouldOptimizeConstants: Whether to use NelderMead optimization to periodically optimize constants in equations.
• optimizer_nrestarts: How many different random starting positions to consider when using NelderMead optimization.
• hofFile: What file to store equations to, as a backup.
• perturbationFactor: When mutating a constant, either multiply or divide by (1+perturbationFactor)^(rand()+1).
• probNegate: Probability of negating a constant in the equation when mutating it.
• mutationWeights: Relative probabilities of the mutations, in the order: MutateConstant, MutateOperator, AddNode, InsertNode, DeleteNode, Simplify, Randomize, DoNothing.
• annealing: Whether to use simulated annealing.
• warmupMaxsize: Whether to slowly increase the max size from 5 up to maxsize. If nonzero, specifies how many cycles (populations*iterations) before increasing by 1.
• verbosity: Whether to print debugging statements or not.
• bin_constraints: See constraints. This is the same, but specified for binary operators only (for example, if you have an operator that is both a binary and unary operator).
• una_constraints: Likewise, for unary operators.
• seed: What random seed to use. nothing uses no seed.
• progress: Whether to use a progress bar output (verbosity will have no effect).
• probPickFirst: Expressions in subsample are chosen based on, for p=probPickFirst: p, p(1-p), p(1-p)^2, and so on.
• earlyStopCondition: Float - whether to stop early if the mean loss gets below this value. Function - a function taking (loss, complexity) as arguments and returning true or false.
• timeout_in_seconds: Float64 - the time in seconds after which to exit (as an alternative to the number of iterations).
• max_evals: Int (or Nothing) - the maximum number of evaluations of expressions to perform.
• skip_mutation_failures: Whether to simply skip over mutations that fail or are rejected, rather than to replace the mutated expression with the original expression and proceed normally.
• enable_autodiff: Whether to enable automatic differentiation functionality. This is turned off by default. If turned on, this will be turned off if one of the operators does not have well-defined gradients.
• nested_constraints: Specifies how many times a combination of operators can be nested. For example, [sin => [cos => 0], cos => [cos => 2]] specifies that cos may never appear within a sin, but sin can be nested with itself an unlimited number of times. The second term specifies that cos can be nested up to 2 times within a cos, so that cos(cos(cos(x))) is allowed (as well as any combination of + or - within it), but cos(cos(cos(cos(x)))) is not allowed. When an operator is not specified, it is assumed that it can be nested an unlimited number of times. This requires that there is no operator which is used both in the unary operators and the binary operators (e.g., - could be both subtract, and negation). For binary operators, both arguments are treated the same way, and the max of each argument is constrained.
• deterministic: Use a global counter for the birth time, rather than calls to time(). This gives perfect resolution, and is therefore deterministic. However, it is not thread safe, and must be used in serial mode.

string_tree(tree::Node, options::Options; kws...)

Convert an equation to a string.

Arguments

• varMap::Union{Array{String, 1}, Nothing}=nothing: what variables to print for each feature.

eval_tree_array(tree::Node, cX::AbstractMatrix{T}, options::Options)

Evaluate a binary tree (equation) over a given input data matrix. The options contain all of the operators used. This function fuses doublets and triplets of operations for lower memory usage.

This function can be represented by the following pseudocode:

function eval(current_node)
    if current_node is leaf
        return current_node.value
    elif current_node is degree 1
        return current_node.operator(eval(current_node.left_child))
    else
        return current_node.operator(eval(current_node.left_child), eval(current_node.right_child))

The bulk of the code is for optimizations and pre-emptive NaN/Inf checks, which speed up evaluation significantly.

Returns

• (output, complete)::Tuple{AbstractVector{T}, Bool}: the result, which is a 1D array, as well as if the evaluation completed successfully (true/false). A false complete means an infinity or nan was encountered, and a large loss should be assigned to the equation.

node_to_symbolic(tree::Node, options::Options;
            varMap::Union{Array{String, 1}, Nothing}=nothing,
            index_functions::Bool=false)

The interface to SymbolicUtils.jl. Passing a tree to this function will generate a symbolic equation in SymbolicUtils.jl format.

Arguments

• tree::Node: The equation to convert.
• options::Options: Options, which contains the operators used in the equation.
• varMap::Union{Array{String, 1}, Nothing}=nothing: What variable names to use for each feature. Default is [x1, x2, x3, ...].
• index_functions::Bool=false: Whether to generate special names for the operators, which then allows one to convert back to a Node format using symbolic_to_node. (CURRENTLY UNAVAILABLE - See https://github.com/MilesCranmer/SymbolicRegression.jl/pull/84).

calculate_pareto_frontier(X::AbstractMatrix{T}, y::AbstractVector{T},
                        hallOfFame::HallOfFame, options::Options;
                        weights=nothing, varMap=nothing) where {T<:Real}

Compute the dominating Pareto frontier for a given hallOfFame. This is the list of equations where each equation has a better loss than all simpler equations.

calculate_pareto_frontier(dataset::Dataset{T}, hallOfFame::HallOfFame,
                        options::Options) where {T<:Real}