Mercurial > octave-nkf

annotate scripts/polynomial/polyfit.m @ 2261:1b6e1629fb91

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
[project @ 1996-05-23 00:52:07 by jwe] Initial revision
author jwe
date Thu, 23 May 1996 00:52:07 +0000
parents
children 5cffc4b8de57
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
2261
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
1 function p = polyfit (x, y, n)
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
2
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
3 # usage: polyfit (x, y, n)
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
4 #
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
5 # Returns the coefficients of a polynomial p(x) of degree n that
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
6 # minimizes sumsq (p(x(i)) - y(i)), i.e., that best fits the data
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
7 # in the least squares sense.
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
8
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
9 # Written by KH (Kurt.Hornik@ci.tuwien.ac.at) on Dec 13, 1994
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
10 # Copyright Dept of Statistics and Probability Theory TU Wien
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
11
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
12 if (nargin != 3)
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
13 usage ("polyfit (x, y, n)");
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
14 endif
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
15
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
16 if (! (is_vector (x) && is_vector (y) && size (x) == size (y)))
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
17 error ("polyfit: x and y must be vectors of the same size");
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
18 endif
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
19
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
20 if (! (is_scalar (n) && n >= 0 && ! isinf (n) && n == round (n)))
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
21 error ("polyfit: n must be a nonnegative integer");
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
22 endif
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
23
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
24 l = length (x);
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
25 x = reshape (x, l, 1);
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
26 y = reshape (y, l, 1);
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
27
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
28 X = ones (l, 1);
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
29
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
30 if (n > 0)
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
31 tmp = (x * ones (1, n)) .^ (ones (l, 1) * (1 : n));
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
32 X = [X, tmp];
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
33 endif
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
34
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
35 # Compute polynomial coeffients, making returned value compatible
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
36 # with Matlab.
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
37
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
38 [Q, R] = qr (X, 0);
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
39
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
40 p = flipud (R \ (Q' * y));
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
41
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
42 if (! prefer_column_vectors)
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
43 p = p';
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
44 endif
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
45
1b6e1629fb91 [project @ 1996-05-23 00:52:07 by jwe]
jwe
parents:
diff changeset
46 endfunction