# linear regression example w/ GSL

Ok so this is kind of random, but I’ve been playing with some functionality that I thought people might like to know about, which is doing linear regression to get the best parameters for a mapping equation. In my case, it’s for mapping tracked coordinates “eye” to “screen” coordinates. The mapping is non linear, and I tried warping which didn’t work so well. this uses GSL (http://www.gnu.org/software/gsl/) gnu scientific library but there might be other simpler libraries (this one worked quite fast). The linear regression allows you to find the best fit for any equation, ie:

Ax + Bx^2 + Cy + Dy^2 + E

it finds the A,B,C,D,E that best fit the mapping, and you need to do it for X and Y.

I’ll post a mac version that works (with GSL) shortly, but I thought I’d put this up cause it’s pretty cool:

in .h:

``````

void calculateWeights(vector <ofPoint> eyePoints, vector <ofPoint>  screenPoints);

vector <ofPoint> follow;
vector <ofPoint> myMouse;

double	cxfit;
double	cyfit;
bool	bBeenFit;

``````

in .cpp

relevant includes:

``````

#include <gsl_math.h>
#include <gsl_test.h>
#include <gsl_fit.h>
#include <gsl_multifit.h>
#include <gsl_multifit_nlin.h>
#include <gsl_blas.h>
#include <gsl_ieee_utils.h>

``````

fitting function:

``````

void testApp::calculateWeights(vector <ofPoint> eyePoints, vector <ofPoint>  screenPoints){

int length = eyePoints.size();

int nTerms = 8;

gsl_matrix * x = gsl_matrix_alloc(length,nTerms);
gsl_vector * yx = gsl_vector_alloc(length);
gsl_vector * yy = gsl_vector_alloc(length);
gsl_vector * w = gsl_vector_alloc(nTerms);

double * ptr;
double * ptrScreenX;
double * ptrScreenY;

ptr = gsl_matrix_ptr(x,0,0);
ptrScreenX = gsl_vector_ptr(yx,0);
ptrScreenY = gsl_vector_ptr(yy,0);

for (int i = 0; i < length; i++){

float xPosEye = eyePoints[i].x;
float yPosEye = eyePoints[i].y;
*ptr++ = xPosEye;
*ptr++ = xPosEye*xPosEye;
*ptr++ = yPosEye;
*ptr++ = yPosEye*yPosEye;
*ptr++ = xPosEye*yPosEye;
*ptr++ = xPosEye*xPosEye*xPosEye;
*ptr++ = yPosEye*yPosEye*yPosEye;
*ptr++ = 1;

*ptrScreenX++ = screenPoints[i].x;
*ptrScreenY++ = screenPoints[i].y;

}

gsl_vector *cx = gsl_vector_calloc(nTerms);
gsl_vector *cy = gsl_vector_calloc(nTerms);

gsl_matrix *cov = gsl_matrix_calloc(nTerms, nTerms);
double chisq;

gsl_multifit_linear_workspace *work = gsl_multifit_linear_alloc(length, nTerms);

int res = gsl_multifit_linear (x,
yx,
cx,
cov,
&chisq,
work);

int res2 = gsl_multifit_linear (x,
yy,
cy,
cov,
&chisq,
work);

printf("-------------------------------------------- \n");

double * xptr = gsl_vector_ptr(cx,0);
double * yptr = gsl_vector_ptr(cy,0);

for (int i = 0; i < nTerms; i++){
printf("cx %i = %f \n", i, xptr[i]);
cxfit[i] =  xptr[i];
}

for (int i = 0; i < nTerms; i++){
printf("cy %i = %f \n", i, yptr[i]);
cyfit[i] =  yptr[i];
}

bBeenFit = true;

printf("-------------------------------------------- \n");

//return ;

}

``````

// how to use the remapping:

``````

float remappedx = (x) * cxfit + (x*x) * cxfit + (y) * cxfit + (y*y) * cxfit
+ (x*y) * cxfit + (x*x*x) * cxfit + (y*y*y) * cxfit + cxfit;

float remappedy = (x) * cyfit + (x*x) * cyfit + (y) * cyfit + (y*y) * cyfit
+ (x*y) * cyfit + (x*x*x) * cyfit + (y*y*y) * cyfit + cyfit;

``````

btw, it was the amazing Tim Hanson who sat with me for an hour to explain this and get it to work, who is part of this great project openprosthetics.org – please check it out, an effort to help support research and development of open source prosthetic limbs.

I misread the title and got super excited.

I thought it said “linear regression example w/ GLSL” Still, nice tip.

Also, OpenCV lets you do this too with CvStatModel.