2
$\begingroup$

I have some experimental data and I need to do a derivative. I looked at various answers how to approach that at e.g. here. But I am not happy with them (the curve is not smooth as needed plus the process is somewhat a black box, especially the part when one filters the data). I found another way to do a derivative that for the type of data that I have works quite well. The process is as follows:

  1. Pick a point ($x_0$) on the curve of the experimental data and then select $\Delta n$ points before this point and $\Delta n$ points following this point.
  2. Fit the selected points with a straight line
  3. Take its slope as the derivative in the originally selected point $x_0$
  4. Do for all the points on the curve.

The width of the window $\Delta n$ controls the "level" of smoothing.

The implementation in Mathematica is as follows:

\[CapitalDelta]n = 100;
(* this is to provide the correct range on the x-axis *)
n = Range[\[CapitalDelta]n + 1, Length@data - (\[CapitalDelta]n + 1)];
(* list of ranges where we will fit the line *)
nrng = Range[n - \[CapitalDelta]n, n + \[CapitalDelta]n, 1];
(* get the slope *)
slope = Fit[Transpose@{data[[#, 1]], data[[#, 4]]}, {1, x}, x, 
      "BestFitParameters"][[2]] & /@ nrng;
(* construct pairs of (x, y'(x)) *)
der = Transpose@{data[[n, 1]], -slope};

Is there a way to make this faster? The dataset that I am using is on Wolfram cloud:

data = CloudGet["https://www.wolframcloud.com/obj/cd96d10f-65a9-4504-8abd-\
f7dd6a2e4668"];
$\endgroup$
3
  • $\begingroup$ "the process is somewhat a black box, especially the part when one filters the data" - depending on the field, some of those approaches are actually pretty standard (e.g. the Savitzky-Golay filter to smoothen spectroscopic data), whereas the one you mention does not sound as common to me. Perhaps you could familiarize yourself with those approaches first. $\endgroup$
    – MarcoB
    Oct 7, 2022 at 13:50
  • $\begingroup$ @MarcoB I read another answer that actually implements S-G filter. But as soon as I saw that it is a convolution, I realised I would really need to study this in detail to understand its effect on data. Anyway, I tried it and Mathematica crashes immeditelly with it due to insufficient memory. Simply put, it cannot be used for this dataset in the version that is implemented in the post (linked in the question). $\endgroup$
    – atapaka
    Oct 7, 2022 at 18:28
  • $\begingroup$ I recently found the ListD resource function that implements the fitting approach to calculating numerical derivatives. You may find it interesting as well. $\endgroup$
    – MarcoB
    Oct 12, 2022 at 11:47

3 Answers 3

4
$\begingroup$

Not your technique, but perhaps a better way to calculate the derivative of data:

der2 = Transpose@{data[[All, 1]],
     -Divide[
       DerivativeFilter[data[[All, 4]], {1}, 50], 
       DerivativeFilter[data[[All, 1]], {1}, 50]]}; // AbsoluteTiming

(*  {0.044465, Null}  *)

ListLinePlot[{der, der2}, 
 PlotStyle -> {AbsoluteThickness[3.5], AbsoluteThickness[1.5]}, 
 PlotRange -> All]

Mathematica graphics

ListLinePlot[Select[#, 8500 < First[#] < 9500 &] & /@ {der, der2}, 
 PlotStyle -> {AbsoluteThickness[3.5], AbsoluteThickness[1.5]}, 
 PlotRange -> All]

Mathematica graphics

$\endgroup$
3
$\begingroup$

Assuming that you use the first and 4th column of the data, your code results in a timing of 9.5 sec on my machine.

We may reduce the time by using "Partition" to get the pieces over which to do the fit. And we may assemble {x,f'[x]} at he same time we calculate the fit:

tmpdat = Partition[data[[All, {1, 4}]], 201, 1];
Fit[#, {1, x}, x][[2, 1]] & /@ tmpdat;

This takes 4.7 sec on my machine.

$\endgroup$
3
  • $\begingroup$ Partition is a good idea, but I realised that it the approach is not very good when I used a larger dataset. When I increase the window size, i keep getting out of memory. So I am now doing this: nrng = {# - \[CapitalDelta]n, # + \[CapitalDelta]n} & /@ n; slp = Fit[data[[#[[1]] ;; #[[2]], {1, 4}]], {1, x}, x][[2, 1]] & /@ nrng; which is based on just providing the upper and lower bounds of indices where we fit. It takes 25 secs on my machine, your approach takes 21 seconds but I cannot use it for larger datasets with larger windows. But I hope it can still be improved? $\endgroup$
    – atapaka
    Oct 6, 2022 at 23:51
  • $\begingroup$ Just a comment, that this method is also implemented in one software where it does the operation in about 2 sec and it is running on the same machine as Mathematica under virtual windows. This is for dataset with 306000 points. That one in Mathetica with the above code takes on my machine 105 seconds, which is ridiculously slow. $\endgroup$
    – atapaka
    Oct 7, 2022 at 0:03
  • $\begingroup$ @atapaka I can pretty much guarantee that second software will be using linear algebra to calculate the results of linear fitting directly, avoiding the overhead in FindFit. FindFit is a much more general function that will carry out some check, format data, include weighting, and produce a "nice" result. A purpose-built implementation avoids most of that: see my answer below for a similar approach in MMA that takes roughly the same time. $\endgroup$
    – MarcoB
    Oct 7, 2022 at 14:17
0
$\begingroup$

Here is a direct implementation of the linear fitting method that skips FindFit and instead implements the explicit solution for the linear regression problem:

ClearAll[linfit]
linfit[data_] :=
 Module[{x, y, m, mt},
  {x, y} = Transpose[data];
  m = {1, #} & /@ x;
  mt = Transpose[m];
  Quiet@LinearSolve[mt . m, mt . y][[2]]
]

Let's import your data and try this out, using a 5-point smoothing window. Like others, below I assumed that you are interested in columns 1 and 4 of your dataset.

data = CloudGet["https://www.wolframcloud.com/obj/cd96d10f-65a9-4504-8abd-f7dd6a2e4668"]

parted = Partition[data[[All, {1, 4}]], UpTo[5], 1];
RepeatedTiming[results = linfit /@ parted;]

(* Out: {2.26477, Null} *)

As you can see, that takes roughly 2 seconds on my machine, probably comparable to the other program you are using. This is, however, not better than using the built in interpolation in my opinion:

ListLinePlot[
 {
   Style[int'["ValuesOnGrid"], Thickness[0.02], Black],
   Style[results, Red]
 },
 DataRange -> data[[{1, -1}, 1]],
 PlotLegends -> {"from Interpolation", "from fitting"},
 InterpolationOrder -> 2, PlotRange -> {Automatic, {-0.022, 0.005}}
]

Comparison of results from interpolation and from fitting

$\endgroup$
4
  • $\begingroup$ I noticed that Fit uses all cores on my machine (4 cores). When I use this dataset CloudObject["https://www.wolframcloud.com/obj/f487156f-b8ae-4d07-a728-\ 56f528ce1bd6"] the Fit take 2-3 minutes with 1201 point window, linfit does not end before my patience does. This exact dataset with this smoothing takes less than 3 sec in the software which runs in a virtual environment (only 2 cpus). Unless there is some obvious way to speed it up by 2 om (i need this type of derivative anyway to keep consistent with all analyses that were done in that software), i am giving up on optimization. $\endgroup$
    – atapaka
    Oct 7, 2022 at 19:03
  • $\begingroup$ @atapaka I can't get the link to the cloud object in the comment above to work. In any case, if you already have a solution through the other software, then I guess your best bet would be to use that. $\endgroup$
    – MarcoB
    Oct 7, 2022 at 19:07
  • $\begingroup$ I used it for most of the analyses but the issue is that I need to do some analysis with the derivative that is not available in the software. The link is wolframcloud.com/obj/f487156f-b8ae-4d07-a728-56f528ce1bd6 I dont know why the slash appeared there. It is an interesting problem though, I just wrote the code in matlab, hoping it would be faster than MA, but it is almost the same performance... The S-G method you mentioned is really more and more intriguing but how do you deal with the memory issue (when i run S-G implemented from the link in the post, i run out of mem.? $\endgroup$
    – atapaka
    Oct 7, 2022 at 19:58
  • $\begingroup$ @atapaka Thank you for the link. On this second data set the linfit approach took 7 minutes and did not run out of memory on my relatively limited laptop (screenshot). $\endgroup$
    – MarcoB
    Oct 7, 2022 at 20:51

Your Answer

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

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