Symbolic Regression

Question

7+ years ago I posted a question here regarding symbolic regression in Mathematica. At that time there did not seem to be much built-in support for this. I'm interested in multi-input models of the form y=f(x1,x2,x3,x4,x5) where I need to find the best fit "f" given a data set.

Would Mathematica accept this FindFormula[data,{x1,x2,x3,x4,x5}] ?

Here is the link to my original question: Does anyone have experience with performing symbolic regression using Mathematica?

Has any meaningful progress been made during the last 7 years in this regard in Mathematica ? I'm aware of TuringBot. Is that considered the state of the art in symbolic regression ?

Thanks for any knowledge on this subject.

Clarification and Concrete Example:

Here is a dataset consisting of 27 points. Each point has 3 independent variables, x1, x2 and x3 and 1 dependent variable, y.

The structure of this list is {{x1,x2,x3,y} . . .}

data = {{1000, 200, 0.00831, 34.5091}, {1000, 200, 0.01622, 21.5683}, {1000, 200, 0.02493, 12.8382}, {1000, 100, 0.00831, 5.1462}, {1000, 100, 0.01622, 4.3582}, {1000, 100, 0.02493, 3.4017}, {1000, 50, 0.00831, 0.6760}, {1000, 50, 0.01622, 0.6452}, {1000, 50, 0.02493, 0.5953}, {1100, 200, 0.00831, 37.9601}, {1100, 200, 0.01622, 23.7251}, {1100, 200, 0.02493, 14.1220}, {1100, 100, 0.00831, 5.6609}, {1100, 100, 0.01622, 4.7940}, {1100, 100, 0.02493, 3.7419}, {1100, 50, 0.00831, 0.7436}, {1100, 50, 0.01622, 0.7097}, {1100, 50, 0.02493, 0.6548}, {1200, 200, 0.00831, 41.4110}, {1200, 200, 0.01622, 25.8820}, {1200, 200, 0.02493, 15.4059}, {1200, 100, 0.00831, 6.1755}, {1200, 100, 0.01622, 5.2299}, {1200, 100, 0.02493, 4.0821}, {1200, 50, 0.00831, 0.8112}, {1200, 50, 0.01622, 0.7743}, {1200, 50, 0.02493, 0.7144}};

For datapoint #1 1000, 200 and 0.00831 combine algebraically in some unknown way to produce 34.5091 as an output.

y = f(x1,x2,x3)

or

34.5091 = f(1000,200,0.00831) etc.

The task is to find the simplest algebraic function f that explains the data. I'm not interested in an interpolating function, neural network or other non-algebraic result; an algebraic formula is desired.

Here is my attempt at this using Mathematica:

FindFormula[data, {x1, x2, x3}, 5, All]

And here is the Mathematica output:

FindFormula::wrgfmt: Argument {{1000.,200.,0.00831,34.5091},{1000.,200.,0.01622,21.5683},{1000.,200.,0.02493,12.8382},{1000.,100.,0.00831,5.1462},{1000.,100.,0.01622,4.3582},{1000.,100.,0.02493,3.4017},{1000.,50.,0.00831,0.676},{1000.,50.,0.01622,0.6452},<<12>>,{1200.,200.,0.02493,15.4059},{1200.,100.,0.00831,6.1755},{1200.,100.,0.01622,5.2299},{1200.,100.,0.02493,4.0821},{1200.,50.,0.00831,0.8112},{1200.,50.,0.01622,0.7743},{1200.,50.,0.02493,0.7144}} at position 1 does not have the right format. Data should be a numerical array of depth less or equal than 2.

From the help menus it appears that Mathematica can deal with just 1 independent variable.

Can this problem involving 3 independent variables be solved using Mathematica ?

If not, are there future plans for this capability ?

@ Mark Kotanchek Thank you for your efforts on my question but your answer is not turning out to be useful for us.

First of all, none of the 10 functions you supplied is equivalent to the theoretical function used to generate the 27 data points I supplied in my question.

Secondly, if I evaluate your first function at each of the 27 given points and compute the errors, I see that the maximum error your function produces is in excess of 172%, error = ((incorrect - correct)/correct)*100. This error magnitude (at a data point known to your software) far exceeds what would be considered acceptable and brings into question expected error magnitudes at unknown data points where such a function would need to be used and trusted. This error takes place at datapoint #7 {1000, 50, 0.00831} where the correct value is 0.67601 and your model prediction is 1.84302

Thirdly, any function found that truly describes the theoretical underpinnings of the data should be able to accurately extrapolate beyond the confines of the data set used to find the function. For example, the gravitational force model found in the referenced YouTube video should be able to return reasonable values when asked to extrapolate. When I use your first model to extrapolate a reasonable amount (no negative or zero valued inputs for example) errors exceeding 1000% can occur. I’d be happy to give you this 28th datapoint if you think it would be useful. But keep in mind, in my real world of engineering acquiring 27 accurate datapoints, equally spaced in a matrix is a luxury. Real-world scenarios (for us) are much messier, data starved, incomplete and irregular.

Could it be that your software requires more than 30 seconds of computation time to return a robust model ? Or do you think this is beyond the current state of the art and I should re-ask in another 7 years ?

Here is the error plot for your model #1. I see a kind of symmetry here, perhaps that holds clues for model improvement ?

Mark, here is the error plot for your 2nd model. Although it's greatly improved, the errors are still too large at the known points. I'm not aware of too many industries that would consider a 24% error acceptable. What your models are in competition with is the Mathematica function Interpolation that returns a data model with maximum error less than 10^-13, essentially zero. I realize that Interpolation does not return algebra and is not doing symbolic regression but error-prone algebra is of little value to those working on life-critical systems. Why don't you refrain from restricting yourself and utilize all available functions in your software and let it run not for 10 minutes but 10 hours (or more) so we can really see what it can do on this problem ? If you think 3 additional data points contained within the input space would help, tell me where in the volume you would like them. This would give your software a total of 30 to work with.

Answer 1

There is a commercial solution from www.evolved-analytics.com which also includes all of the surrounding functionality to actually extract the value from the data.

The next version will be released in the next few days. Free trials are available as are academic versions.

DataModeler is 100% Mathematica and one of the versions is an add-on package so all of the 500+ functions and symbols are fully documented and available for custom workflows. For those not interested in such, the GUI is also 100% Mathematica (with a little linux under the hood to allow easy export to word or powerpoint documentation).

CAVEAT: I am the CTO and original developer so I am biased.

Answer 2

If you want to see how Symbolic Regression can be used from within Mathematica (Wolfram Language), you can check out my video here: https://www.youtube.com/watch?v=eMcBMvy_tlk. It's an EXTREMELY powerful product.

Answer 3

I have been a user of DataModeler from Evolved Analytics for more than 10 years that is mentioned just above. It will be one of products that match your need. It has a function called SymbolicRegression that creates approximate functions with genetic programming for optimization. I have been using it for the analysis in the medical field.

I compared its ability of fitting to several datasets with other machine learning models, including Lasso, Ridge, Logistic, SVM, Random Forest, XGBoost and so on. Probably due to the high flexibility of symbolic regression, it achieved as much accuracy/AUC as or even better than these of other models.

Answer 4

Just for fun, I devoted 30 seconds to modeling the data set ...

Each dot represents a model and the ones in red define the ParetoFront exploring the trade-off of model complexity vs accuracy. Somewhat arbitrarily, I decided to focus on those models having a complexity less than 100 and an R2 better than 0.98. The models along the ParetoFront in this region are:

From the contenders, I chose a diverse model set from which to form an ensemble,

The attraction of a model ensemble is that it defines a trustable model since the models will agree where constrained by data but will diverge (since we chose them for their diversity) when asked to venture into unknown regions of parameter space or if the system has undergone some sort of fundamental change.

For real-world analysis, understanding the relative importance of variables, identification of metavariables as well as the creation of trustable models is very important.

Real-world data here is defined as up to tens of thousands of variables and millions of data records. Of course, things get easier with smaller data sets. For the big data sets, more than 30 seconds of model search will be required.

Per JimB's comment, here is a graphic using the standard Mathematica machine learning techniques relative to a SymbolicRegression ModelEnsemble. The blue line represents the model prediction and the envelope the trust for a simple ball trajectory problem. The left side shows the prediction if we stop the available data at the halfway point and the right if we only have the extremal values and are looking to interpolate. (The red points are known whereas the yellow are the unknown reality.) The gray lines in the case of SymbolicRegression are the trajectories of the constituent models which comprise the ensemble.