0
$\begingroup$

I have this experimental data set which I know should follow the Complementary error function.

My goal is to find the σ of the associated gaussian distribution.

However, I cannot seem to do it in any way using NonlinearModelFit. Could you please give me some hints?

The data is in a list of lists, e.g.:

data = ToExpression@Import["http://pastebin.com/raw.php?i=BtcXHh9c"];

And when plotted it looks like this:

ListPlot[Transpose@data]

Improve this question
$\endgroup$
4
  • $\begingroup$ where is the data? $\endgroup$
    – Basheer Algohi
    Jul 7, 2014 at 8:15
  • $\begingroup$ Actually, I should have first of all put the data in ordered pairs: {{x1,y1},{x2,y2} ... }. $\endgroup$
    – Danno
    Jul 7, 2014 at 10:03
  • $\begingroup$ I think you mean something like this. $\endgroup$
    – Öskå
    Jul 7, 2014 at 10:04
  • $\begingroup$ To fix data download: data = ReplaceRepeated[ ToExpression[ ToString@Import["http://pastebin.com/raw.php?i=BtcXHh9c"]][[1, 1]], {a___, {b_, Null}, {c_}, d___} :> {a, {b, c}, d}]; ListPlot[data] $\endgroup$
    – Chris Degnen
    Jul 7, 2014 at 15:34

2 Answers 2

Reset to default
2
$\begingroup$

I don't think there is anything unusual here that can't be found in the documentation. Normalization parameters in the x and y dimensions (b, and a, respectively) can be added into the model, as can the x-axis shift (parameter c). Parameter estimates are taken from a look at the original data.

nlm = NonlinearModelFit[data, 
  a Erfc[b (x + c)], {{a, 8 10^-8}, {b, 1/10}, {c, -35000}}, x]
Plot[nlm[x], {x, 35000, 35400}, Epilog -> {Red, Point /@ data}]

Mathematica graphics

Improve this answer
$\endgroup$
6
  • $\begingroup$ agree with your comment. Originally the data posted was not such a beautiful fit...hence my hamfisted approach to force a function where NonlinearModelFit did not work...timing...but at least it led to provision of correct data...though this fit seems to good to be true $\endgroup$
    – ubpdqn
    Jul 7, 2014 at 20:35
  • $\begingroup$ @ubpdqn +1 for "spending" the time on the bad data! $\endgroup$
    – bobthechemist
    Jul 7, 2014 at 21:02
  • $\begingroup$ Thanks. I need to ask, are parameter estimates in NonlinearModelFit necessary? $\endgroup$
    – Danno
    Jul 8, 2014 at 15:50
  • $\begingroup$ @Danno I can't think of a time when knowing the uncertainty of your parameter model is not necessary. Get them with nlm["ParameterTable"] $\endgroup$
    – bobthechemist
    Jul 8, 2014 at 16:00
  • $\begingroup$ But, and it really may be a newbie question on many levels, isn't the fitting algorithm's job to find the parameter values? I am talking about the estimates you gave, e.g. the 8 10^-8besides $a$ etc. Notice I may not have fully understood your last comment. $\endgroup$
    – Danno
    Jul 8, 2014 at 17:15
2
$\begingroup$

I am a little unsure what the aim is here. Looking at the data only the first third of the plot could be fitted to Erfc like function.

I post this, perhaps, as motivation.

d = Transpose@data;
ListPlot[ds = SortBy[d, #[[1]] &], Joined -> True]

enter image description here

Extracting the part of plot I referred to:

fd = ds[[1 ;; 35]]
fdlp = ListPlot[fd]

enter image description here

Now as a rather ugly approach to trying to fit desired function. This will involve rescaling, removal of troublesome 0 or negative values, logarithmic transform and back transforming. This is NOT ideal and am sure expert users will have better approaches (assuming this is what is desired).

{min, max} = Through[{Min, Max}[fd[[All, 2]]]];
res = Rescale[fd[[All, 2]]];
trans = Transpose[{fd[[All, 1]], res}];
int = Interpolation[trans, InterpolationOrder -> 1];
{xmin, xmax} = {34800.`, 35137.`};
dif[u_] := D[int[x], x] /. x -> u
w = {#, -dif[#]} & @@@ trans;
wlp = ListPlot[w]
peak = Cases[w, {_, Max[w[[All, 2]]]}][[1, 1]]
shift = Cases[{#1 - peak, #2} & @@@ w, {_, _?Positive}];
lg = {#1, Log@#2} & @@@ shift;
lm = LinearModelFit[lg, {1, x^2}, x]
func[u_] := Exp[Normal@lm] /. x -> u
Show[Plot[func[t - peak], {t, xmin, xmax}], wlp]
ex = Rescale[v, {0, 1}, {min, max}]
h[u_] := ex /. v -> u;
soln[t_] := NIntegrate[func[x], {x, -Infinity, t}]
Show[fdlp, Plot[h[1 - soln[x - peak]], {x, xmin, xmax}]]

The derivative of interpolating function (and fit):

enter image description here

The final "fit":

enter image description here

The fit expression (with back transforming and rescaling):

enter image description here

Improve this answer
$\endgroup$
2
  • $\begingroup$ Thank you very much for your answer, however I didn't notice that the data I published was corrupted. I just edited the question with the actual data, that looks way different, and I am quite sure it fits a complementary error function (apart from some two constant pieces at the beginning and end of the data set). $\endgroup$
    – Danno
    Jul 7, 2014 at 14:03
  • $\begingroup$ @Danno...ok...bobthechemist has provided the very clear answer $\endgroup$
    – ubpdqn
    Jul 7, 2014 at 20:39

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.