How do I apply FindFit to data describing a spiral?

I have a list of numbers that, taken in order, describe a spiral in the complex plane. Here is a small subset of it, sufficient to illustrate the problem:

list = {0.01263537001717932 + 0.07067349025848388*I,
0.047250267751479655 + 0.05209549861304419*I,
0.06596386316698546 + 0.02009477948782494*I,
0.06598609040888194 - 0.014961756210224754*I,
0.049803550815569794 - 0.04397149975599973*I,
0.023258224576407013 - 0.06099441316314905*I,
-0.006595117201610087 - 0.06383941766600312*I,
-0.03328082763577616 - 0.053650337718323365*I,
-0.05206972138609968 - 0.03391359218040549*I,
-0.060484300582113376 - 0.009278027786894835*I,
-0.05825268030262476 + 0.015525647453509116*I,
-0.04688894972912606 + 0.03649893985640911*I,
-0.029084475019386297 + 0.05086222796961101*I,
-0.008065197460458237 + 0.05723002594987413*I,
0.012975410250310286 + 0.05553656470450358*I,
0.031305061954391124 + 0.046792399011505864*I,
0.04492053689497534 + 0.03275238482001683*I,
0.05266179241274281 + 0.015563156240077383*I,
0.054197917274563315 - 0.002559520264756074*I,
0.04991942530288897 - 0.019591820266882193*I,
0.04077290783116659 - 0.03388637260674089*I,
0.028070627915692578 - 0.044278439184983934*I,
0.01330140973390842 - 0.05012128281294888*I,
-0.002038137442440774 - 0.051263817289083816*I,
-0.0165802068618294 - 0.04798609535183286*I,
-0.02918463847368017 - 0.04090829084508081*I,
-0.03899948965681603 - 0.030887325316889896*I,
-0.0454864418481062 - 0.01891284287902329*I,
-0.04841647916489085 - 0.006011387028409941*I,
-0.04784247104545147 + 0.0068352246318077806*I,
-0.044055588207469146 + 0.018753932199973615*I,
-0.03753212093089012 + 0.029025665281149138*I,
-0.028876474184390682 + 0.03711227916123358*I,
-0.018765076191573443 + 0.04266692897277542*I,
-0.007894807565850093 + 0.04553040163855704*I,
0.0030615501673451687 + 0.04571633159705858*I,
0.013491358125476596 + 0.043388326822319356*I,
0.022869033045436377 + 0.038831902057961665*I,
0.030773906534906574 + 0.032423827007817124*I,
0.0368996370049279 + 0.02460111005525975*I,
0.041056149623619684 + 0.015831403761543628*I,
0.04316530554941611 + 0.006586176629337748*I,
0.043251606732296555 - 0.0026824244302409914*I,
0.04142925049656888 - 0.011560475046661134*I,
0.03788678231239426 - 0.01968671709393703*I,
0.03287047720230047 - 0.026762590899545816*I,
0.026667429059332358 - 0.0325581087282407*I,
0.019589158683998575 - 0.03691397289036705*I,
0.011956378358853229 - 0.03974044128341975*I,
0.004085383055940278 - 0.04101349823952391*I,
-0.0037236170298382623 - 0.04076890659052196*I};

ListLinePlot[ReIm[list], AspectRatio -> 1]


The full data set is thousands of entries long, and produces the following ListLinePlot:

I'd like to use FindFit to generate a model of the form a*x^b*E^(c*x) (assuming that is the most appropriate method and model). But the problem is that

Table[Arg[list[[j]]], {j, Length[list]}]


{1.393880369, 0.8341309359, 0.2957019528, -0.222970991, -0.723286234,
-1.206498742, -1.67373906, -2.126029227, -2.564296235, -2.989383456,
2.881124903, 2.480154188, 2.090243925, 1.710800332, 1.341276188,
0.9811660967, 0.6300023216, 0.2873511408, -0.0471903747,
-0.3739972084, -0.6934189241, -1.005781912, -1.311391395,
-1.610533215, -1.903475446, -2.190469834, -2.471753094, -2.747548088,
-3.018064883, 2.999683597, 2.739139467, 2.483310929, 2.232030327,
1.985138787, 1.742485618, 1.503927758, 1.269329268, 1.038560868,
0.8114995072, 0.5880279702, 0.3680345156, 0.1514125389,
-0.06193973712, -0.2721195498, -0.4792198625, -0.6833296112,
-0.8845339342, -1.082914385, -1.278549132, -1.47151314, -1.661878347}

produces arguments that do not respect the 'winding number' of the spiral. (Why it produces some negative arguments is also baffling to me, but secondary). Since the values of Arg are restricted to an interval of 2*Pi, and are not in regular increments, I can't apply FindFit to more than a single circuit of the origin.

I could of course use FindFit[list, a*x^b*Cos[c*x] - I*a*x^b*Sin[c*x], {a, b, c}, x] as the model, but this produces an unevaluated result - and even if it did produce a result, it would be harder to work with.

How do I go about fitting list to a meaningful model?

• Should we assume that the sample points are equally spaced in terms of x? If so, then modelling the real and imaginary components separately might be the first step to find the form of the model. Maybe first taking a look at dataRe = Transpose[{Range[Length[list]], Re[list]}]; dataIm = Transpose[{Range[Length[list]], Im[list]}]; ListLinePlot[{dataRe, dataIm}]. (And a link to the complete dataset might be helpful.)
– JimB
Commented Apr 23, 2021 at 16:30
• Can we assume that each turn of the spiral has at least one data point?
– A.G.
Commented Apr 23, 2021 at 17:51
• Do you have a fit criteria in mind? (least squares, max distance, etc.)
– A.G.
Commented Apr 23, 2021 at 17:52

It will be easier to fit the data if it does not wind around the origin. So, let's unwind the spiral. Applying the Arg function to the data list gives

ListLinePlot[Arg/@list]


We see that arguments are decreasing, except where they jump up by $$2\pi$$. We want to unwind the arguments by removing the jumps. Here is a function that will unwind the arguments

ClearAll[unwind]
unwind[arg_List] :=
n = First@
FirstPosition[xlist, SelectFirst[xlist, Last[#] > First[#] &]];
If[NumericQ[n],
Join[arg[[;; n]], unwind[-2 \[Pi] + arg[[n + 1 ;;]]]],
arg]
]
ListLinePlot[unwind[Arg /@ list], ImageSize -> Small]


Now we do the fit and plot the results

arg = unwind[Arg /@ list];
abs = Abs /@ list;
fit = FindFit[Thread[{arg, abs}], b x + c, {b, c}, x]

(* {b -> 0.00192386, c -> 0.0671975}  *)

ParametricPlot[ReIm[(b x + c /. fit) Exp[I x]], {x, -\[Pi], 5 \[Pi]},
Epilog -> {PointSize[.02], Red, Point[ReIm[list]]},
ImageSize -> Small]


This model uses a linear $$b x + c$$ relation to relate the radial distance to the angular coordinate. Other relationships may be possible.

Note that the unwind function only works with decreasing data that has $$2\pi$$ jumps. If we Reverse the data list we would have increasing data with $$-2\pi$$ jumps. unwind requires 2 modifications to handle increasing data (instead of decreasing data). First, change the > comparison to < and second, change $$-2\pi$$ to $$2\pi$$.

• Very nice. If one uses a quadratic rather than just a linear model, the fit looks perfect to the eye. But using $AIC_c$ as a guide, a 14-degree polynomial is even better and supported statistically.
– JimB
Commented Apr 24, 2021 at 4:47
• Thanks to all of you. @JimB, any chance you could share that model? Commented Apr 24, 2021 at 11:15

@LouisB supplied the answer and this is just an extended comment illustrating how other models using that answer might be selected.

First, rather than a linear model with a polynomial of order 1, consider a quadratic (2nd degree polynomial).

arg = unwind[Arg /@ list];
abs = Abs /@ list;
lm2 = LinearModelFit[Transpose[{arg, abs}], {x, x^2}, x];

ParametricPlot[ReIm[lm2[x] Exp[I x]], {x, -5 π, π},
Epilog -> {PointSize[.02], Red, Point[ReIm[list]]}, ImageSize -> Medium]


An absolute measure of the fit is the root mean square error for predicting abs from arg:

lm2["EstimatedVariance"]^0.5
(* 0.0000451409 *)


If by this measure the fit is good enough, one can stop. But it can pay to look at the residuals to see if there are any additional patterns that might be better fit with additional terms in the model.

ListPlot[Transpose[{lm2["PredictedResponse"], lm2["FitResiduals"]}]]


Clearly there is more that could be explained possibly by adding in a 3rd order term.

lm3 = LinearModelFit[Transpose[{arg, abs}], {x, x^2, x^3}, x];
lm3["EstimatedVariance"]^0.5
(* 3.69339*10^-6 *)


This seems a much better fit.

This process of adding polynomial terms can continue until the $$AIC_c$$ is minimized. For this particular dataset, a 14-degree polynomial has the miniminum $$AIC_c$$ value. This is unlikely to hold for the larger dataset, however. (I mostly deal with biological data and I have never ever seen anything more that a 6-th degree polynomial that might be appropriate.)

aicc[p_] := LinearModelFit[Transpose[{arg, abs}], Table[x^i, {i, p}], x]["AICc"]
ListPlot[Table[{p, aicc[p]}, {p, 20}], Frame -> True, FrameLabel -> {"Degree of polynomial", "AICc"}]