Find best fit curve from multiple curves of different size

Hi,

I have multiple curves with different lengths. When plotted they look like the image below. I want to find one best fit line that best describes all the curves. The best fit line should look something like the second image below. The data is attached. From the data Z is on the x axis, R is on the y axis.

Answer 1

[Torsten]

I parametrized the curves by normalized arc length. To compute the average curve, I took the mean of R- and Z-coordinates of the 9 curves that correspond to the same normalized arc length.

Maybe you find a better solution.

  load("R.mat")
  load("Z.mat")
  n = 100;
  I = linspace(0,1,n);
  for j = 1:9
    Zj = Z{j};
    Rj = R{j};
    Lj = zeros(numel(Zj),1);
    for i = 2:numel(Zj)
      Lj(i) = Lj(i-1) + sqrt((Zj(i)-Zj(i-1))^2+(Rj(i)-Rj(i-1))^2);
    end
    Ljnorm = Lj/Lj(end);
    Zjinter(j,:) = interp1(Ljnorm,Zj,I);
    Rjinter(j,:) = interp1(Ljnorm,Rj,I);
  end
  Zaverage = mean(Zjinter,1);
  Raverage = mean(Rjinter,1);
  hold on
  for j = 1:9
    plot(Z{j},R{j},'b')
  end
  plot(Zaverage,Raverage,'r')
  hold off

●

[Konvictus177]

This is actually a quite simple logic to understand and answered my question. I ended up doing something very similar to your solution above but your solution is more compact.

Lets take this a step further. Imagine now that I have this "average" curve and I want to create many more curves where each of them goes through the same point for example at Z=31.5 and R=15. Additionally the starting point of each curve is at the bottom at R=-15 and Z=user defined.

So basically now the common point is fixed at Z=31.5, R=15. The user enters multiple different starting points at R=-15, Z=user input. These lines then all run in the same shape as the "average" line but all of them cross at Z=31.5, R=15. I know that in theory there are unlimited options but how would I create these lines so that they are similar to the average line?

Something like the image below:

[Torsten]

As long as you don't quantify what "similar to the average line" means, I agree: you have unlimited options.

[Sam Chak]

Hi @Konvictus177

Somehow, my solver intuition tells me that this resembles a boundary value problem. What is shown in the image is probably a phase portrait. Each trajectory is generated from the same dynamical system but at different starting points. You want to find the initial at the point such that it converges to a point at and . Please clarify if my assessment is correct.

[Konvictus177]

What is shown in the image are multiple rotationary axes R and one linear Z axis of a robot. I want the rotational axes to start at different Z locations but then I want them all to stop at the same R and Z position (that the crossing point). The path to get there or shape of the curve should somehow be similar to the average line.

Answer 2

[John D'Errico]

On the right end, you have an infinite slope. On the left end, not so much of a problem. But it will be better to do a couple of things. First, these look almost like arcs of an ellipse. So lets pretend they are that, or at least something close.

What you need to understand is an ellipse is just a circle that has been stretched. So we can undo that.

First, I'll rescale the axes.

load Z.mat
load R.mat
whos

  Name      Size             Bytes  Class    Attributes

  R         1x9             433320  cell               
  Z         1x9             433320  cell               
  ans       1x29                58  char               

I'll set the center of the coordinate system at (Z,R) = (30,-15), and then rescale R and Z, so they will now be roughly circular. Then I could transform to polar coordinates. At least, if I needed to do so. But do I?

R0 = -15;
Z0 = 30;
Rtrans = @(r) (r - R0)/(max(r) - R0);
Ztrans = @(z) (z - Z0)/(max(z) - Z0);
for i = 1:9
  plot(Ztrans(Z{i}),Rtrans(R{i}))
  hold on
end
axis equal
xlim([0 1])
ylim([0,1])
hold off
grid on

Not too bad. Do you see that by simply rescaling the axes by the endpoints of the curves, we now have all the curves virtually on top of each other?

My gut tells me that this may be sufficient for your purposes. Essentially, all you need to know are those endpoints, and you can reconstruct any curve, from this almost circular common arc.

Now, how, you might ask, can you use that curve, as a common arc, and to then reconstruct the original curves from it? Simpler (maybe) than you think.

Rmax = cellfun(@max,R)

Rmax = 1×9

13.8180 16.2501 15.5263 15.4948 15.5853 15.5561 15.5919 15.7650 15.6824

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

Zmax = cellfun(@max,Z)

Zmax = 1×9

54.9765 62.7786 60.7609 59.5278 58.4727 57.2374 56.0746 55.0867 53.5833

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

Next, once all of the curves are in this globally consistent form, convert to polar and model rhe result.

Rcommon = [];
Zcommon = [];
for i = 1:9
  Rcommon = [Rcommon;Rtrans(R{i})];
  Zcommon = [Zcommon;Ztrans(Z{i})];
end

[thet,rad] = cart2pol(Zcommon,Rcommon);
plot(thet,rad,'.')

There we see that one of your curves is a little unlike the others. But the difference is not too immense, if you look at the scale on the y-axis in that last plot.

mdl = fit(thet,rad,'poly4')

mdl = 
     Linear model Poly4:
     mdl(x) = p1*x^4 + p2*x^3 + p3*x^2 + p4*x + p5
     Coefficients (with 95% confidence bounds):
       p1 =    -0.05061  (-0.05154, -0.04967)
       p2 =      0.2918  (0.2889, 0.2948)
       p3 =     -0.3361  (-0.3392, -0.3331)
       p4 =    0.004989  (0.003827, 0.006151)
       p5 =       1.001  (1, 1.001)

plot(mdl,thet,rad)

If you look carefully, the blue curve falls neatly in the middle of that band. Now, we can reconstruct each of the curves, as well as plot a common curve on top, with an average set of parameters.

The model is a function of polar angle theta, varying from 0 to pi/2 radians. So the global curve will be:

thetglob = linspace(0,pi/2,20)';
radglob = mdl(thetglob);
[Zglob,Rglob] = pol2cart(thetglob,radglob);
Zmaxbar = mean(Zmax);
Rmaxbar = mean(Rmax);
plot(Zglob*(Zmaxbar - Z0)+ Z0,Rglob*(Rmaxbar - R0) + R0)

I've used average values to undo the transformations.