Simplex optimisation using fminsearch

Please how can i use the simplex optimisation syntax called 'fminsearch' to optimise a circle fitted to a set of random points? I've tried using the equation of a circle as the function for the syntax but am not making any head way. I presume am not getting the function right.

[Walter Roberson]

https://www.mathworks.com/matlabcentral/fileexchange/5557-circle-fit

https://www.mathworks.com/matlabcentral/fileexchange/22643-circle-fit-pratt-method

https://www.mathworks.com/matlabcentral/fileexchange/22678-circle-fit-taubin-method

https://www.mathworks.com/matlabcentral/fileexchange/15060-fitcircle-m

https://github.com/horchler/circfit

[Ebube Ezi]

Thanks for the resources. I've been through some of these files but i need to understand how to optimise the system using fminsearch.

[Stephan]

Do you need help with the usage of fminsearch or with the mathematical model to describe your task?

[Ebube Ezi]

Please I need help with the usage of fminsearch. I understand I need to state an initial guess and develop a function but how do I call it.

I've developed random points, fitted a circle through these points which was a little off from the points, how can I perform optimisation of the system using fminsearch?

[Ebube Ezi]

A mathematical model could also helpful. Thank you

[John D'Errico]

Why do you need to use fminsearch? This is a problem simply solved without recourse to fminsearch.

[Bruno Luong]

https://matlab.fandom.com/wiki/FAQ#How_can_I_fit_a_circle_to_a_set_of_XY_data?

There are more robust methods https://www.emis.de/journals/BBMS/Bulletin/sup962/gander.pdf

fminsearch is the last thing (even that is unlikely) I would choose.

Answer 1

[Walter Roberson]

   x = vector_of_x_coordinates;
   y = vector_of_y_coordinates;
   obj = @(xyr) sum( (x-xyr(1)).^2 + (y-xyr(2)).^2 - xyr(3).^2 );
   xyr0 = [guess_x, guess_y, guess_r];
   [XYR, residue] = fminsearch(obj, xyr0);

[Ebube Ezi]

I will give this a try. Thank you very much.

[Ebube Ezi]

Do I just run it as a script or create a function for 'obj' before running the script?

[Walter Roberson]

You can run it as a script. You would need to substitute values for vector_of_x_coordinates and vector_of_y_coordinates

[Ebube Ezi]

Ok. Thank you very much.

[Biraj Khanal]

@Walter Roberson, I tried your script like this.

x = C1_w_noise(1,:);
y =  C1_w_noise(2,:);

obj = @(xyr) sum( (x-xyr(1)).^2 + (y-xyr(2)).^2 - xyr(3).^2 );
xyr0 = [5,5,1]; % center and radius for C2 as initial guess
[XYR, residue] = fminsearch(obj , xyr0);

But this just runs to a maximum iteration number and stops. I am not sure if I understand the multivariable fminsearch here. Can you help?

EDIT:

the sum function can return in negative value and my assumption is that fminsearch, as shown in the graph allows it to go as low as possible which results in the max iteration . So I added an absolute to the sum. Now, it does try to fit but does not result in anything correct.C3 is the circle from the result of fminsearch.

[Walter Roberson]

Ah, I see now that increasing the radius would decrease the sum. I will need to think about a better formula.

[Biraj Khanal]

Shouldn't that be resolved when the absolute value of the sum is used?

[Biraj Khanal]

So, this is my adaptation of a tutorial I found. The idea is to minimize the error as distance between the coordinates of guessed circle and the actual data points.

obj = @(hkr) sumerror(hkr,x,y);

[sol ,residue] = fminsearch(obj,[Guess_h, Guess_k, Guess_r]);

function E = sumerror(hkr,x,y)

t=linspace(0,2*pi,length(x));
xx=hkr(1) + hkr(3)*cos(t); % create x and y for the guess iteration
yy=hkr(2) + hkr(3)*sin(t);

Ex=sum((xx-x).^2);   
Ey=sum((yy-y).^2);
E=(Ex+Ey)/2; %average of the eror in x and y coordinates

end

It seems to work for me. If anyone gets a simpler idea to do this, it would be great.

Answer 2

[John D'Errico]

You are still looking for help on this?

There is ABSOLUTELY no need to use a simplex optimizer for this problem. NONE AT ALL. Let me make up some data.

n = 50;
XY = normalize(randn(n,2),2,'norm') + 2 + randn(n,2)/5;
plot(XY(:,1),XY(:,2),'.')
axis equal

The data is a circle, with center at [2,2], and a radius of 1.

[C,R,rmse] = circlefit(XY)

C = 1×2

2.0039 2.0386

R = 1.0565

rmse = 0.2043

theta = linspace(0,2*pi)';
XYhat = C + R*[cos(theta),sin(theta)];
hold on
plot(XYhat(:,1),XYhat(:,2),'r-')
plot(C(1),C(2),'gx')
hold off

So a circle at center pretty near the point [2,2], with radius estimated as also very close to 1.

Even if you have only half an arc of a circle, it still works.

k = XY(:,1) > 2;
plot(XY(k,1),XY(k,2),'.')
axis equal
[C,R,rmse] = circlefit(XY(k,:))

C = 1×2

2.1750 2.0898

R = 0.9377

rmse = 0.2168

theta = linspace(-pi/2,pi/2)';
XYhat = C + R*[cos(theta),sin(theta)];
hold on
plot(XYhat(:,1),XYhat(:,2),'r-')
plot(C(1),C(2),'gx')
hold off

As expected, the estimates are less good there, but that data is seriously noisy. Honestly, if someone gave me this data and asked to fit a circle to it, I would have guessed the data is almost useless to fit a circle.

I've attached circlefit to this answer. (circlefit also works on higher dimensional data, where it can fit a sphere to data in 3-d, etc.)

Circlefit uses an algebraic trick to remove the quadratic terms in the circle equations. This reduces the problem to finding the best intersection of many straight lines. But that is a problem easily solved using linear algebra. It can also use robustfit from the stats toolbox, if you have it.

Don't write your own code to solve problems like this.