1
$\begingroup$

I tried to fit a fourier sin series onto this data with the two different following codes.

len = Subtract @@ data[[{-1, 1}, 1]];
func = Fit[data, Table[Sin[(\[Pi] n)/len x], {n, 1, 40}], x];
Show[ListPlot[data], 
Plot[func, {x, 0, len}, PlotStyle -> {{Thick, Red}}]]

and this one fits not necessarily a fourier sin series with interpolation

model1 = Sum[Subscript[a, i]*x^i, {i, 0, 5}];
sol = FindFit[data, {model1}, 
Flatten[{Table[{Subscript[a, i]}, {i, 0, 5}]}], x]
nOfControlPoints = 30;
controlPoints = Subdivide[#1, #2, nOfControlPoints - 1] & @@ 
MinMax[data[[;;, 1]]];
ListPlot[data, Epilog -> {Red, Point[{#, 0}] & /@ controlPoints}, 
Frame -> True, Axes -> False]
model[y : {__Real}] := Interpolation[Transpose[{controlPoints, y}]];
nlm = NonlinearModelFit[data, model[Array[y, nOfControlPoints]][x], 
Array[y, nOfControlPoints], x];
Plot[nlm["Function"][x], {x, 0, 3500}, Epilog -> {Point[data]}]

Whenever I try to take the fourier transform and plot it using below code (changing func in below code to nlm and viceversa I get a bunch of Nintegrate errors.

FT = FourierTransform[func[x], x, \[Omega]]
Plot[Evaluate@Abs@FT, {\[Omega], -10, 10}, PlotRange -> All]

How can I plot a Fourier Transform of these two fits that I have made? If it's not possible is there a method to fit similar curves like the ones above to the data so that I can take a Fourier transform plot of it afterwards? Should I be just using Fourier plot instead since I fit in a fourier sin series? Basically I'm hoping to get a fourier transform plot that shows some peaks.

Thanks.

$\endgroup$

1 Answer 1

3
$\begingroup$

I assign the data from your spreadsheet to dat. Then I delete elements 57, 85, and 152 which are duplicate data points. This may not be how you wish to handle the duplicates, but it's what I've done to create a quick example. Since the Fourier function is going to expect data to be evenly spaced, I use a linear interpolation to re-sample the data (whether this makes sense for your data is something that you would have to decide).

I sample the interpolation at intervals of 1 day to build up a new dataset, and take the Fourier transform of only the y-data. Then I take the square of the absolute value of the Fourier transform and plot the result. The peaks here indicate which frequencies contribute most to the creation of the dataset. To get the x- and y-axis ranges, I basically just messed around with them until I got something that looked nice.

dat = (* your data *);
dat2 = Delete[dat, {{57}, {85}, {152}}];
int = Interpolation[dat2, InterpolationOrder -> 1];
dat3 = Table[{x, int[x]}, {x, 0, 3528, 1}];
f = Fourier[dat3[[All, 2]]];
pwr = Abs[f]^2;
ListLinePlot[
  pwr,
  AxesLabel -> {"Fourier Component", "Power"},
  ImageSize -> 500,
  PlotRange -> {{0, 500}, {0, 2*^9}},
 ]

Power plot of Fourier transform of question asker's data.

One thing to note about the Fourier transform is that data points will be duplicated. If you plot the full x range, you will see a mirror image at high frequencies, so you rarely need to include more than half the range (in this case we could go from 1 to ~1764). Usually component #1 will have the highest power. This corresponds to a DC offset of the data from zero. If your data averages to zero, then the power of this component will be low.

Component 2 corresponds to a sine wave that oscillates exactly once across the entire data set. Component 3 corresponds to a sine wave that oscillates exactly twice across the entire data set (and so on). We can use FindPeaks to get the positions of these 2 peaks we see in the graph.

FindPeaks[pwr, Automatic, Automatic, 6*^7]

\begin{equation} \left( \begin{array}{cc} 1 & 7.52512\times 10^{10} \\ 117 & 1.3158\times 10^9 \\ 233 & 1.54078\times 10^8 \\ 3298 & 1.54078\times 10^8 \\ 3414 & 1.3158\times 10^9 \\ 3529 & 2.00314\times 10^9 \\ \end{array} \right) \end{equation}

We can see that components 117 and 233 have the highest power. These correspond to signals that oscillate 116 times and 232 times over the graph range (which is 3528 days). Since 232 is exactly twice 116, it's probably just the second harmonic of the 116 signal, meaning the signal with frequency 116/3528 days (0.0329/day) is the most interesting signal that I see. Period = 1/frequency, and 1/(116/3528 days) gets us 30.4138 days. This looks like basically a month, so perhaps that makes sense for whatever data you're using?

Keep in mind that this is the period of just the strongest frequency. We can also quickly check to see if we're in the right ballpark with this:

pks2 = dat2[[Floor@FindPeaks[dat2[[All, 2]]][[All, 1]]]];
Mean[Differences[pks2[[All, 1]]]]

That first line is a bit messy, but I'm using FindPeaks to find all of the peaks using default arguments. This returns an n x 2 list where the first column represents the position in the list dat2. I use these positions to extract the actual data points from dat2 and store them as pks2. Then I ask for the average difference between consecutive x values.

This gets me 35.9398 days. That's a bit longer than the Fourier transform. Why is that?

ListLinePlot[
  dat2,
  Epilog -> Point@pks2
]

OP's data with some peaks identified.

It looks like my hastily written peak finding code missed some of the peaks. We also have a fairly long stretch from day 2800 to day 3450 that is basically flat with no peaks. These things are probably driving up the average a bit. Based on that, I would say that the 30.4 day period the Fourier transform identified is probably an important and correct signal.

Depending on what kind of data this is, extracting the strongest frequency may be enough for you. If not, you could always find the peaks more accurately than I did and look at things like the standard deviation of the differences, and other statistical stuff.

One last thing: out of idle curiousity, I plotted the first 700 days with a sine wave with the frequency we found.

Show[
 ListLinePlot[
   dat2, 
   AxesLabel -> {"Days Between", "Amount"}, 
   ImageSize -> 500, 
   PlotRange -> {{0, 700}, Full}, 
   Epilog -> Point@pks2
  ],
 Plot[
   5000 Sin[2 \[Pi] 116/3528 t + 0.85 \[Pi]] + 5000, 
   {t, 0, 3528}, 
   PlotRange -> Full, 
   PlotStyle -> Red]
]

OP's data matched with a sine wave.

It matches pretty well! Of course, there are places later in the plot where it won't match up, and this isn't really of much use, but I thought it was neat.

$\endgroup$
3
  • $\begingroup$ This is a great detailed answer,thank you very much. BTW, these are invoice amounts of a customer however there are some irregular invoice dates (they are not duplicates but invoice of another product, but the frequency should be close to 30 as you found out and mentioned.) $\endgroup$
    – indiffer
    Feb 7, 2019 at 10:37
  • 1
    $\begingroup$ @indiffer Glad it helped a bit! If these are invoice dates, it might be worth taking the Differences of the peak days and making a histogram or something. You should see lots of invoices 30 to 31 days apart, and if there are any with an unusually long or short timing, they should show up to either side of that peak. $\endgroup$
    – MassDefect
    Feb 7, 2019 at 16:51
  • 1
    $\begingroup$ Doing this a bit differently I am getting a peak frequency .4127/(2Pi) suggesting a period of 15.2 or so. There is also a signal spike at half that frequency, but it does not go quite as high. I get this using irregular DFT and periodogram code from this prior MSE thread $\endgroup$ Feb 7, 2019 at 17:04

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.