Least-squares line fitting: Difference between revisions

From JSXGraph Wiki
← Older editNewer edit →
No edit summary
No edit summary
Line 6: Line 6:
<jsxgraph width="600" height="600">
<jsxgraph width="600" height="600">
var brd = JXG.JSXGraph.initBoard('jxgbox',{boundingbox:[-5,5,5,-5], keepaspectratio:true, axis:true});
var brd = JXG.JSXGraph.initBoard('jxgbox',{boundingbox:[-5,5,5,-5], keepaspectratio:true, axis:true});
var i, j, p = [], angle, xr, yr, delta = 0.1;
brd.suspendUpdate();
if (true) {
    var i, p = [], angle, xr, yr, delta = 0.1;


// Random points are constructed which lie roughly on a line
    // Random points are constructed which lie roughly on a line
// defined by y = 0.3*x+1.
    // defined by y = 0.3*x+1.
// delta*0.5 is the maximal distance in y-direction of the random
    // delta*0.5 is the maximal distance in y-direction of the random
// points from the line.
    // points from the line.
brd.suspendUpdate();
    brd.suspendUpdate();
for (i=0;i<5;i++) {
    for (i=0;i<100;i++) {
  yr = 10*(Math.random()-0.5);
        yr = 10*(Math.random()-0.5);
  xr = 0.*yr+delta*(Math.random()-0.5);
        xr = 0.*yr+delta*(Math.random()-0.5);
  p.push(brd.create('point',[xr, yr], {withLabel:false}));
        p.push(brd.create('point',[xr, yr], {withLabel:false}));
    }
} else {
    var i, p = [], angle, co, si, delta = 0.2;
    // Random points are constructed which lie roughly on a circle
    // of radius 4 having the origin as center.
    // delta*0.5 is the maximal distance in x- and y- direction of the random
    // points from the circle line.
    for (i=0;i<100;i++) {
        angle = Math.random()*2*Math.PI;
        co = 4*Math.cos(angle)+delta*(Math.random()-0.5);
        si = 4*Math.sin(angle)+delta*(Math.random()-0.5);
        p.push(brd.create('point',[co+2, si-1], {withLabel:false}));
    }
}
}
brd.unsuspendUpdate();
brd.unsuspendUpdate();


var r = [], rbar = [], x = [], y = [], z = [], A = [[0,0,0],[0,0,0],[0,0,0]], n, d;
//
// Ab hier wird der beste Kreis, bzw. die beste Gerade ermittelt.
var j, r = [], rbar = [], x = [], y = [], z = [], A = [[0,0,0],[0,0,0],[0,0,0]], n, d,
    eigen, minIndex, minE, ev, c, xm, ym, zm, radius;
n = p.length;
n = p.length;
for (i=0;i<n;i++) {
for (i=0;i<n;i++) {
Line 51: Line 71:
     A[2][2] += r[i][2]*r[i][2];
     A[2][2] += r[i][2]*r[i][2];
}
}
/*
for (i=0;i<3;i++) {
for (i=0;i<3;i++) {
     for (j=0;j<3;j++) {
     for (j=0;j<3;j++) {
Line 56: Line 77:
     }
     }
}
}
/*
console.log(A[0]);
console.log(A[1]);
console.log(A[2]);
var res = JXG.Math.Numerics.Jacobi(A);
console.log(res);
*/
*/
/*
// Now, the general linear least-square fitting problem
//    min_z || M*z - y||_2^2
// is solved by solving the system of linear equations
//    (M^T*M) * z = (M^T*y)
// with Gauss elimination.
MT = JXG.Math.transpose(M);
B = JXG.Math.matMatMult(MT, M);
c = JXG.Math.matVecMult(MT, y);
z = JXG.Math.Numerics.Gauss(B, c);


// Finally, we can read from the solution vector z the coordinates [xm, ym] of the center
eigen = JXG.Math.Numerics.Jacobi(A);
// and the radius r and draw the circle.
minIndex = 0;
//var zm = z[0]*0.5;
minE = eigen[0][0][0];
var xm = z[0]*0.5;
for (j=1;j<3;j++) {
var ym = z[1]*0.5;
    if (eigen[0][j][j]<minE) {
var r = Math.sqrt(z[2]+xm*xm+ym*ym);
        minIndex = j;
        minE = eigen[0][j][j];
    }
}
ev = [eigen[1][0][minIndex],eigen[1][1][minIndex],eigen[1][2][minIndex]];
c = -(rbar[0]*ev[0]+rbar[1]*ev[1]+rbar[2]*ev[2]);


brd.create('circle',[ [xm,ym], r]);  
xm = -ev[1];
//alert([xm,ym,r].toString());
ym = -ev[2];
// Having constructed the points, we can fit a line described
zm = 2.0*(c+ev[0]);
// by homogeneous coordinates
console.log(c, c+ev[0]);
// through the point set, consisting of n points.
if (Math.abs(c)<0.01) {
// The (n times 2) matrix consists of
    brd.create('line',[zm,xm,ym], {strokeColor:'green'});
//  x_1, y_1
   
//  x_2, y_2
}  else {
//      ...
     var radius = Math.sqrt((xm*xm+ym*ym-2*c*zm)/(zm*zm));
//  x_n, y_n
     brd.create('circle',[[zm,xm,ym],radius]);
// where x_i, y_i is the position of point p_i
// y is equal to the all-minus-one vector.
var M = [], y= [], MT, B, c, z, n;
n = p.length;
for (i=0;i<n;i++) {
     M.push([p[i].X(), p[i].Y()]);
     y.push(-1.0);
}
}
// Now, the general linear least-square fitting problem
//    min_z || M*z - y||_2^2
// is solved by solving the system of linear equations
//    (M^T*M) * z = (M^T*y)
// with Gauss elimination.
MT = JXG.Math.transpose(M);
B = JXG.Math.matMatMult(MT, M);
c = JXG.Math.matVecMult(MT, y);
z = JXG.Math.Numerics.Gauss(B, c);
brd.create('line',[1, z[0], z[1]]);
*/
</jsxgraph>
</jsxgraph>


Revision as of 18:00, 9 November 2010

This little JXSGraph application finds the line - described by homogeneous coordinates [a,b,c] - that minimizes

[math]\displaystyle{ \sum_{i=1}^n (ax_i+by_i+cz_i)^2. }[/math]

Coming soon...

Retrieved from "http://jsxgraph.org/wiki/index.php?title=Least-squares_line_fitting&oldid=4731"