# How to guess the correct fitting function to some data?

Imagine we are plotting some date points $\left(x_i,f(x_i)\right)$ that we obtained experimentally, and that we want to know what $f(x)$ is. The way to do this is to use some software and try to fit the data to some guessed function. For example, if the behavior of the data points looks like exponential decay we then choose an exponential decaying function ..etc.

My question is: sometimes the data points are perfectly fitted to the exponential decaying trial function only on a certain region, but then the rest of the points show deviation away from the trial function.

1. How to proceed in this case? is there a catalog (something like mathematical tables) for functions and their plots that one can use as a guide?
-
Isn't this pretty much what science is? There are good software programs that will fit functions to data; I expect that there are some that even try multiple functions. Ultimately, though, the software is only as good as the user. If you think an exponential decay would fit over some domain and a different function elsewhere, try a piecewise curve, look at the fitting statistics, and then try something else. If someone else has examined the same phenomenon, you might start with their approach. Also, there is a good database from NIST (dlmf.nist.gov) but it's probably too much. – AdamRedwine Oct 7 '11 at 12:42
""The way to do this is to use some software and try to fit the data to some guessed function."" That is the problem. In good old days one had a basic math knowledge and some understanding of the theory behind the experiment..... – Georg Oct 7 '11 at 14:39
@AdamRedwine: "Isn't this pretty much what science is?" - no, it is not. Science is about understanding, not guessing. Phenomenological trial and error is a really a minor tool in the research arsenal... – Slaviks Oct 8 '11 at 7:23
Slaviks, yes, you are correct. I should have said "Isn't that pretty much what engineering is?" :) Sorry, I've been in engineering the last three years and there's more of a difference than typically assumed. – AdamRedwine Oct 8 '11 at 12:05
One should know what function to fit the data to, because one should have a model of what is going on. Without any clue of what you're expecting, even a perfect fit is completely meaningless. – Danu Oct 25 '14 at 17:50

For example, if the behavior of the data points looks like exponential decay we then choose an exponential decaying function

No, if your theory predicts an exponential decay, then you use an exponential decay function. Or if your theory predicts a linear relationship, you use a linear fit. You really shouldn't have to guess which kind of function to use, because the "proper" way to analyze data is to test its consistency with some particular model, and the model tells you what kind of curve to expect.

This is something I feel isn't emphasized enough (or at all) in lab classes and such: if you don't have a model, the value of your data analysis is significantly diminished. In other words, just noticing that your data fit e.g. an exponential curve doesn't mean much by itself.

That being said, picking a functional form from the data isn't always totally worthless. It might hint at what sort of theory you should be looking at, for example. Or as Alexey suggested in a comment, it might allow you to find a simpler, approximate method for describing the data (which is kind of a special case of hinting at the kind of theory to look at). But any time the best you can say is that "the data looks like X, so we fit it with X", there is something deeply unsatisfying about that analysis.

-
I totally agree. But sometimes you have to rely on the data and the best fit to some guessed function in case you do not have concrete underlying theory or in case you are using the data to get a phenomenological model that you can use elsewhere. – Revo Oct 7 '11 at 21:53
Excellent point David. I believe I was once told that Max Plank fit the black body radiation spectrum with an equation first and explained it later. I don't know whether or not this is true, but it would be a rare example of figuring out the math before the physics. – AdamRedwine Oct 8 '11 at 11:46
You are certainly right of course, but you leave out of consideration an interesting mathematical problem, which can be stated intuitively as "find the simplest analytical formula describing the data within a given precision". It actually may be a rather powerful way of compressing data in some cases. – Alexey Bobrick May 15 at 16:15
@AlexeyBobrick yes, that's intentional, because that interesting mathematical problem is not really within the scope of physics. In other contexts it's certainly possible that fitting an arbitrary function to data is useful. – David Z May 16 at 8:55
@DavidZ: I had encountered the method used a few times rather popularly in astrophysics, see e.g. Eddington formula. Idea is as follows: If you have an approximate method to solve a certain problem, it might provide you with a complex computationally expensive analytical expression. One might prefer in this case to find another simpler function, which fits the original expression within the same error bounds. – Alexey Bobrick May 16 at 9:58

The possibilities are endless of course. You have to weigh the goodness of the fit against the prediction power of the hypothesis. You can fit any set of data "perfectly" by using as the "fit function" the data itself. This is of course complete nonsense since you would not gain any insight. The more degrees of freedom you admit, the better the fit can be (since, they're more parameters to adjust). With a common $\chi^2$ test, one oftentimes quotes the $\chi^2/\#dof.$, that is the sum of quadratic deviations divided by the number of degrees of freedom.

Ideally, when taking data, you already have an idea what kind of shape the data will have, since you would know the underlying physics process.

If this is not the case there are a number of ways to slice and dice the data to get more insights. In you example of exponential decay, you could try to have a sum of exponentials with different decay constants and different normalizations $\sum_i a_ie^{-x/a_i}$ and see if with a minimal extension your data fits better. (But again, you're losing prediction power).

Another possibility to dig up trends in data is Fourier analysis. This is useful when you suspect that the data has multiple overlapping periodical processes (modes). If this is the case you will see those nicely in a fourier transform of the data.

Still another common possibility is the convolution of two functions. This is for example the case when you're trying to measure a quantity but have a intrinsic uncertainty stemming from the apparatus. In that case the measurement would be a convolution of the true measurement values and the resolution function of the apparatus (often a good approximation is gaussian).

Having a table of functions would not really be useful, I think. In the end you would want to compare your data to some kind of theory, and theories are the better the more they explain (that is "fit data well") with the fewest possible number of assumptions (i.e. "parameters)

EDIT:

A neat example (including code in python) for Fourier Analysis is given here: http://linuxgazette.net/115/andreasen.html (Example 2). There, sunspot data is analyzed and it is found that there is a dominant frequency in the occurence of sunspots.

I'm no expert for a reference, but any introduction to signal processing or data analysis should do. This book by seems to have good rating on Amazon and seems to touch the main points: http://www.amazon.com/Introduction-Signal-Processing-Sophocles-Orfanidis/dp/0132091720/ref=ntt_at_ep_dpt_1

-

Physics, physics, physics.

You have to understand what is going on (this is what georg is hinting about in the comments): you shouldn't be "just guessing". If you're lucky there will be enough structure in the data to give you hint.

It is hard to be very specific without know the particulars of your situation.

The kinds of things to think about include

• Might there be a random process that fakes your signal? Probably you need a underling constant (i.e. fit to C + f(x)) if fitting in time. Non-random background may have a different structure. What other backgrounds exist? Can you measure them; or failing that can you model them?
• Is your function subtly wrong? Everyone reaches for a Gaussian first when fitting a peak, but if the process is resonant perhaps a Lorentzian would be better.
• Are there multiple channels through which this process can proceed. Maybe one (or more!) that you have neglected should be included in the fit.
• Are you sure your instrument is properly calibrated? Non-lineariteis especially can introduce artifacts in fitting.
• luksen mentions the convolution of your resolution with the data. This will have the biggest effect if the data changes rapidly over the scale of your instruments resolution. Fitting to convolved functions is harder, but can be done.

Note that when you get to the end you should not have "the fit works if I add these fudge factors" but "the fit is clear after we account for all the back-grounds that we have discussed". That is you aren't done until you know why all the fudge factor are needed.

Finally, there is a bit of craft to this and you learn it by doing. Worse, most of the experience you will get is somewhat specific to the particular sub-discipline or measurement technique. That's life for you.

-
I totally agree. Could you read my reply to luksen above. Any help will be greatly appreciated. – Revo Oct 7 '11 at 22:03
A large part of graduate education in experiment is an apprenticeship in this business; I'm not away of any static resource that pretend to teach it because it is a process. – dmckee Oct 7 '11 at 22:23

The theoretical decay model should be derived before you fit, and then the probability of the fit is given by the Baysian modification of your prior distribution for the parameters of the model. In this, I agree with the other answers.

But there are decay models which are sufficiently general that they can be used to fit large classes of experimental data, without fitting everything. If your data falls of exponentially for a while, then crosses over to falling like a power law with a slowly decreasing power, finally like a reciprocal logarithm, then like log-log, and then like log-log-log, you aren't going to get a good fit from any of these blind methods. But this is extraordinarily rare when studying physical systems--- such a complex decay only occurs in contrived mathematical situations, or when a system is actively moved from one phase to another according to a complicated plan.

### Fitting by complete sets of functions

The generic way in which you fit arbitrary data that you feel should be approximated by a smooth curve is to run a best-fit polynomial. The polynomials are dense in the continuous functions, so you can always approximate anything, but you must use a polynomial of lowest order which fits the structure you believe is there. Best-fit lines are most common.

The standard measure of quality of fit is the sum of the squares of the deviation from the fit. Minimizing this gives the least-squares line, and it is also easy to find least squares polynomials of arbitrary order, so long as the order is less than the number of data points.

But this type of fit is wrong for decaying data. In this case you have options:

• If the measurements are very accurate when the quantity is small, you can take the logarithm, and fit the log-linear plot with a polynomial. This will give you the exponential decay as a first approximation, with higher order terms correcting this. The problem with this approach is that you are going to infinite time, and a polynomial always blows up at infinite time very quickly, so that its exponential probably disappears asymptotically faster than a physical function.
• You can try to use a Laplace transform. The Laplace transform of a pure decaying exponential is a delta function, and if you are using a Laplace transform for numerical data, you are just superposing exponentials to get a fit. You can start with a two-exponential fit of the form $f(t) = A e^{-at} + B e^{-bt}$, where you fit the four-constants by least squares (you can do this by steepest descent). This is also a complete expansion, because it is related to the Fourier transform by analytic continuation.
• You can fit to a falling power, which you should do anyway, because this is a common non-exponential decay. For example, if you have a pendulum decaying by air friction, the profile falls off as $1\over A+ Bt$, and this is always badly approximated by decaying exponentials. This is not a complete expansion, but if you find that the result is asymptotically a power, you can write $f(t) = g(t)\over A+Bt$, and fit g by Pade approximations (polynomials over polynomials).

There is no general absolute method, because each asymptotic limit is different, but once you know even just a little bit, you can extract the leading behavior and fit the rest using a complete expansion.

-
Have anybody collected such best practices in a single book yet? – Revo Oct 8 '11 at 7:55

This is, generally, difficult question. But usually, programs for fits includes already a number of predefined functions. I use free Java program called SCaVis (http://jwork.org/scavis) for data fits and it has a number of predefined functions, such as polynomial, exponential, logs, etc. I just look at my data an make the decision about what function can fit the best. If I cannot find a predefined function, I create a script describing my functions using SCaVis python scripting

-

Finding good formulas to fit data is a simple task:
From Cornell Creative Machines Lab download 'Eureqa' (version II, codenamed 'formulize' for LINUX, W$, Mac) and use it. It will give you all the formulas and graphs you need. If you need more processing power you can pay for cloud computing. An example with SunSpots (monthly smoothed),$dS=S_{-1}-S_{-2}$: Sunspots delta_S time S dS t 58 4.6 1 62.6 4.6 2 70 4.6 3 55.7 7.4 4 85 -14.3 5 83.5 29.3 6 81.6 -1.5 7 ... ... 50.6 6.9 3153 78 6.7 3154 88 27.4 3155 96.7 10 3156 73 8.7 3157  set TARGET function to be like this S = f(dS, dS^2, dS^3, t) it gives the formula ("Correlation Coefficient 0.99097013") S = 0.4687512939*sma(dS, 34) + 0.321053344*sma((dS^2), 131) + 0.2694580519*sma((dS^2), 406) + -0.04853427812/sma((dS^2), 153) + sma(dS, 73) + sma(dS, 110) - 0.4618182571 where sma(variable,length) is the simple moving average We can identify several important periods : 34, 73, 110, 131, 153, 406 months and both$dS$and$(dS)^2$-- set TARGET function to be like: compute the Derivative of S in function of S D(S) = f(S)  it gives the formula ("Correlation Coefficient" 0.96109605) (D(S)) = 0.02477032993 + 2.767538634*S + 0.8820551019*sma(S, 3) - 3.650621347*sma(S, 2) we see that we only need the last three values of monthly sunspot count to achive a good formula With a Correlation Coefficient 0.99097013 (so close to 1 ;) I expect the formula is predictive to a certain degree. note: to obtain the sma(variable,length) in excell I used, as an example: column C = variable (dS) column F, first row F1 = 76 (length) next rows in column F will have sma(dS,76) F2 =AVERAGE(OFFSET($C2;-CELL("contents";F$1)+1;0;CELL("contents";F$1);1))
and propagate to all rows ...

-
See physics.stackexchange.com/questions/20548/sunspots-formula for why and how this is wrong headed. – dmckee Mar 12 '13 at 2:40