Re: show residuals of circle fit
- To: mathgroup at smc.vnet.net
- Subject: [mg79951] Re: [mg79900] show residuals of circle fit
- From: DrMajorBob <drmajorbob at bigfoot.com>
- Date: Thu, 9 Aug 2007 05:26:15 -0400 (EDT)
- References: <4787716.1186570751288.JavaMail.root@m35>
- Reply-to: drmajorbob at bigfoot.com
Data points:
angleRange = {-\[Pi]/4, \[Pi]};
points = Table[theta = RandomReal@angleRange;
{2*Cos[theta] + 1 + .05*Random[],
2*Sin[theta] - 3 + .05*Random[]}, {30}];
Nonlinear fit and plot:
Clear[func]
func[x_] = Fit[points, {1, x, x^2, x^3}, x]
-1.26986 + 0.743417 x - 0.366836 x^2 - 0.0257307 x^3
Show[Graphics[
Plot[func[x], Evaluate@Flatten@{x, Sort[points][[{1, -1}, 1]]}]],
ListPlot[points, PlotStyle -> {Red, PointSize[Small]}], Axes -> True,
AspectRatio -> Automatic]
(no need to copy/paste the limits!)
Residuals:
ListPlot[points /. {x_, y_} :> {x, y - func[x]}]
Circle fit and plot:
hypersphereFit[data_List] := Module[{rhs, qq, rr, x, y, r},
{qq, rr} = QRDecomposition[Append[#, 1] & /@ data];
rhs = #.# & /@ data;
{x, y, r} = {1/2, 1/2, 1} LinearSolve[rr, qq.rhs];
{{x, y}, Sqrt[r + x^2 + y^2]}]
{center, r} = hypersphereFit@points
Show[Graphics[Circle[center, r, angleRange]], ListPlot@points,
AspectRatio -> Automatic]
{{1.02827, -2.97834}, 2.00092}
Circular residuals:
ListPlot[points /. {x_, y_} :> {ArcTan[x, y],
Norm[{x, y} - center] - r}]
Directly minimizing the sum of squared residuals gives a nearly identica=
l =
circle:
Clear[error, r]
error[xc_, yc_, r_] = #.# &[
points /. {x_, y_} :> Norm[{x, y} - {xc, yc}] - r];
Minimize[{error[x, y, r], r > 0}, {x, y, r}]
{0.00679177, {r -> 2.00117, x -> 1.02822, y -> -2.97891}}
Bobby
On Wed, 08 Aug 2007 03:53:39 -0500, will <willpowers69 at hotmail.com> wrot=
e:
> Dear Mathforum,
>
> I am having a problem with displaying residuals after fitting a circle=
to
> a set of points. I am using some code that I found on the mathgroup
> site =
> (see:http://forums.wolfram.com/mathgroup/archive/2003/Jan/msg00019.htm=
l)
> to generate the points and fit the circle:
>
> points generation and display:
>
> points = Table[theta = Random[];
> {2*Cos[theta] + 1 + .05*Random[],
> 2*Sin[theta] - 3 + .05*Random[]}, {30}]
>
> plot = ListPlot[points, AspectRatio -> Automatic]
>
> circle fitting and display with points:
>
> fit circle:
>
> hypersphereFit[data_List] :=
> Module[{mat, rhs, qq, rr, params, cen},
> mat = Map[Append[#, 1] &, data];
> rhs = Map[#.# &, data];
> {qq, rr} = QRDecomposition[mat];
> params = LinearSolve[rr, qq.rhs];
> cen = Drop[params, -1]/2;
> {Sqrt[Last[params] + cen.cen], cen}]
>
> display circle (using Graphics[Circle[]] inputing the results from the=
=
> above
> function) and points:
>
> Show[Graphics[
> Circle[Part[hypersphereFit[points], 2],
> Part[hypersphereFit[points], 1], {0, 1}]], plot]
>
> now i really want to be able to display the residuals for the points, =
=
> but can't
> figure out how to.
>
> If it is any help i use the following code when displaying residuals f=
rom
> a quadratic curve fit:
>
> quadratic fit:
>
> func = Fit[points, {1, x, x^2}, x]
>
> define range of plot:
>
> Insert[Part[Sort[points], {1, -1}, 1], x, 1]
>
> out: {x, 2.1305, 3.04856}
>
> copy and paste this output into the range for the Plot function:
>
> Show[Graphics[
> Plot[func, {x, 2.1305, 3.04856}]],
> ListPlot[points, PlotStyle -> {Red, PointSize[Small]}], Axes -> True,=
> AspectRatio -> Automatic]
>
> this shows the quadratic fitted line with the points
>
> then calculate and display the residuals for the quadratic fit:
>
> res = Transpose[{First /@ points,
> Last /@ points - (func /. x -> First /@ points)}];
> ListPlot[res]
>
> I would really, really like to be able to do this residuals calculatio=
n
> and display for the circle fit.
>
> Any help would be much appreciated,
>
> William Parr
>
>
-- =
DrMajorBob at bigfoot.com