3
PDL::Fit::Linfit - routines for fitting data with linear combinations of functions.
7
This module contains routines to perform general curve-fits to a set (linear combination)
8
of specified functions.
12
(y0, y1, y2, y3, y4, y5, ...ynoPoints-1)
14
The fit routine tries to model y as:
16
y' = beta0*x0 + beta1*x1 + ... beta_noCoefs*x_noCoefs
18
Where x0, x1, ... x_noCoefs, is a set of functions (curves) that
19
the are combined linearly using the beta coefs to yield an approximation
22
The Sum-Sq error is reduced to a minimum in this curve fit.
30
This is your data you are trying to fit. Size=n
34
2D array. size (n, noCoefs). Row 0 is the evaluation
35
of function x0 at all the points in y. Row 1 is the evaluation of
36
of function x1 at all the points in y, ... etc.
38
Example of $functions array Structure:
40
$data is a set of 10 points that we are trying to model using
41
the linear combination of 3 functions.
43
$functions = ( [ 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 ], # Constant Term
44
[ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9 ], # Linear Slope Term
45
[ 0, 2, 4, 9, 16, 25, 36, 49, 64, 81] # quadradic term
52
$yfit = linfit1d $data, $funcs
61
1D Fit linear combination of supplied functions to data using min chi^2 (least squares).
65
Usage: ($yfit, [$coeffs]) = linfit1d [$xdata], $data, $fitFuncs, [Options...]
69
Signature: (xdata(n); ydata(n); $fitFuncs(n,order); [o]yfit(n); [o]coeffs(order))
71
Uses a standard matrix inversion method to do a least
72
squares/min chi^2 fit to data.
74
Returns the fitted data and optionally the coefficients.
76
One can thread over extra dimensions to do multiple fits (except
77
the order can not be threaded over - i.e. it must be one fixed
78
set of fit functions C<fitFuncs>.
80
The data is normalised internally to avoid overflows (using the
81
mean of the abs value) which are common in large polynomial
82
series but the returned fit, coeffs are in
88
# Generate data from a set of functions
89
$xvalues = sequence(100);
90
$data = 3*$xvalues + 2*cos($xvalues) + 3*sin($xvalues*2);
92
# Make the fit Functions
93
$fitFuncs = cat $xvalues, cos($xvalues), sin($xvalues*2);
95
# Now fit the data, Coefs should be the coefs in the linear combination
97
($yfit, $coeffs) = linfit1d $data,$fitFuncs;
103
Weights Weights to use in fit, e.g. 1/$sigma**2 (default=1)
108
package PDL::Fit::Linfit;
110
@EXPORT_OK = qw( linfit1d );
111
%EXPORT_TAGS = (Func=>[@EXPORT_OK]);
116
@ISA = qw( PDL::Exporter );
117
use PDL::Options ':Func';
118
use PDL::Slatec; # For matinv()
121
my $opthash = ref($_[-1]) eq "HASH" ? pop(@_) : {} ;
122
my %opt = parse( { Weights=>ones(1) }, $opthash ) ;
123
barf "Usage: linfit1d incorrect args\n" if $#_<1 or $#_ > 3;
124
my ($x, $y, $fitfuncs) = @_;
126
($y, $fitfuncs) = @_;
130
my $wt = $opt{Weights};
132
# Internally normalise data
134
my $ymean = (abs($y)->sum)/($y->nelem);
135
$ymean = 1 if $ymean == 0;
136
my $y2 = $y / $ymean;
140
my $M = $fitfuncs->xchg(0,1);
141
my $C = $M->xchg(0,1) x ($M * $wt->dummy(0)) ;
142
my $Y = $M->xchg(0,1) x ($y2->dummy(0) * $wt->dummy(0));
144
# Fitted coefficients vector
146
$a = matinv($C) x $Y;
150
$yfit = ($M x $a)->clump(2); # Remove first dim=1
152
$yfit *= $ymean; # Un-normalise
154
my $coeff = $a->clump(2);
155
$coeff *= $ymean; # Un-normalise
156
return ($yfit, $coeff);
163
*linfit1d = \&PDL::linfit1d;
added
removed
Lib/Fit/Linfit.pm
