Jump to content

Polynomial regression: Difference between revisions

C
m (<lang>)
(C)
Line 51:
1.000 2.000 3.000
</pre>
 
=={{header|C}}==
{{libheader|libgsl}}
'''Include''' file (to make the code reusable easily) named <tt>polifitgsl.h</tt>
<lang c>#ifndef _POLIFITGSL_H
#define _POLIFITGSL_H
#include <gsl/gsl_multifit.h>
#include <stdbool.h>
#include <math.h>
bool polynomialfit(int obs, int degree,
double *dx, double *dy, double *store); /* n, p */
#endif</lang>
'''Implementation''' (the examples [http://www.gnu.org/software/gsl/manual/html_node/Fitting-Examples.html here] helped alot to code this quickly):
<lang c>#include "polifitgsl.h"
 
bool polynomialfit(int obs, int degree,
double *dx, double *dy, double *store) /* n, p */
{
gsl_multifit_linear_workspace *ws;
gsl_matrix *cov, *X;
gsl_vector *y, *c;
double chisq;
 
int i, j;
 
X = gsl_matrix_alloc(obs, degree);
y = gsl_vector_alloc(obs);
c = gsl_vector_alloc(degree);
cov = gsl_matrix_alloc(degree, degree);
 
for(i=0; i < obs; i++) {
gsl_matrix_set(X, i, 0, 1.0);
for(j=1; j < degree; j++) {
gsl_matrix_set(X, i, j, pow(dx[i], j));
}
gsl_vector_set(y, i, dy[i]);
}
 
ws = gsl_multifit_linear_alloc(obs, degree);
gsl_multifit_linear(X, y, c, cov, &chisq, ws);
 
/* store result ... */
for(i=0; i < degree; i++)
{
store[i] = gsl_vector_get(c, i);
}
 
gsl_multifit_linear_free(ws);
gsl_matrix_free(X);
gsl_matrix_free(cov);
gsl_vector_free(y);
gsl_vector_free(c);
return true; /* we do not "analyse" the result (cov matrix mainly)
to know if the fit is "good" */
}</lang>
'''Testing''':
<lang c>#include <stdio.h>
 
#include "polifitgsl.h"
 
#define NP 11
double x[] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10};
double y[] = {1, 6, 17, 34, 57, 86, 121, 162, 209, 262, 321};
 
#define DEGREE 3
double coeff[DEGREE];
 
int main()
{
int i;
 
polynomialfit(NP, DEGREE, x, y, coeff);
for(i=0; i < DEGREE; i++) {
printf("%lf\n", coeff[i]);
}
return 0;
}</lang>
'''Output''' of the test:
<pre>1.000000
2.000000
3.000000</pre>
 
=={{header|Fortran}}==
Cookies help us deliver our services. By using our services, you agree to our use of cookies.