NAG CL Interface
e02adc (dim1_​cheb_​arb)

1 Purpose

e02adc computes weighted least squares polynomial approximations to an arbitrary set of data points.

2 Specification

#include <nag.h>
void  e02adc (Integer m, Integer kplus1, Integer tda, const double x[], const double y[], const double w[], double a[], double s[], NagError *fail)
The function may be called by the names: e02adc, nag_fit_dim1_cheb_arb or nag_1d_cheb_fit.

3 Description

e02adc determines least squares polynomial approximations of degrees 0 , 1 , , k to the set of data points x r , y r with weights w r , for r=1,2,,m.
The approximation of degree i has the property that it minimizes σ i the sum of squares of the weighted residuals ε r , where
ε r = w r y r - f r  
and f r is the value of the polynomial of degree i at the r th data point.
Each polynomial is represented in Chebyshev series form with normalized argument x ¯ . This argument lies in the range -1 to +1 and is related to the original variable x by the linear transformation
x ¯ = 2 x - x max - x min x max - x min .  
Here x max and x min are respectively the largest and smallest values of x r . The polynomial approximation of degree i is represented as
1 2 a i + 1 , 1 T 0 x ¯ + a i + 1 , 2 T 1 x ¯ + a i + 1 , 3 T 2 x ¯ + + a i + 1 , i + 1 T i x ¯ ,  
where T j x ¯ is the Chebyshev polynomial of the first kind of degree j with argument x ¯ .
For i = 0 , 1 , , k , the function produces the values of a i + 1 , j + 1 , for j=0,1,,i, together with the value of the root mean square residual s i = σ i / m-i - 1 . In the case m = i + 1 the function sets the value of s i to zero.
The method employed is due to Forsythe (1957) and is based upon the generation of a set of polynomials orthogonal with respect to summation over the normalized dataset. The extensions due to Clenshaw (1960) to represent these polynomials as well as the approximating polynomials in their Chebyshev series forms are incorporated. The modifications suggested by Reinsch and Gentleman (Gentleman (1969)) to the method originally employed by Clenshaw for evaluating the orthogonal polynomials from their Chebyshev series representations are used to give greater numerical stability.
For further details of the algorithm and its use see Cox (1974), Cox and Hayes (1973).
Subsequent evaluation of the Chebyshev series representations of the polynomial approximations should be carried out using e02aec.

4 References

Clenshaw C W (1960) Curve fitting with a digital computer Comput. J. 2 170–173
Cox M G (1974) A data-fitting package for the non-specialist user Software for Numerical Mathematics (ed D J Evans) Academic Press
Cox M G and Hayes J G (1973) Curve fitting: a guide and suite of algorithms for the non-specialist user NPL Report NAC26 National Physical Laboratory
Forsythe G E (1957) Generation and use of orthogonal polynomials for data fitting with a digital computer J. Soc. Indust. Appl. Math. 5 74–88
Gentleman W M (1969) An error analysis of Goertzel's (Watt's) method for computing Fourier coefficients Comput. J. 12 160–165
Hayes J G (ed.) (1970) Numerical Approximation to Functions and Data Athlone Press, London

5 Arguments

1: m Integer Input
On entry: the number m of data points.
Constraint: m mdist 2 , where mdist is the number of distinct x values in the data.
2: kplus1 Integer Input
On entry: k+1 , where k is the maximum degree required.
Constraint: 0 < kplus1 mdist , where mdist is the number of distinct x values in the data.
3: tda Integer Input
On entry: the stride separating matrix column elements in the array a.
Constraint: tdakplus1 .
4: x[m] const double Input
On entry: the values x r of the independent variable, for r=1,2,,m.
Constraint: the values must be supplied in nondecreasing order with x[m-1] > x[0] .
5: y[m] const double Input
On entry: the values y r of the dependent variable, for r=1,2,,m.
6: w[m] const double Input
On entry: the set of weights, w r , for r=1,2,,m. For advice on the choice of weights, see the E02 Chapter Introduction.
Constraint: w[r] > 0.0 , for r=0,1,,m - 1.
7: a[kplus1×tda] double Output
On exit: the coefficients of T j x ¯ in the approximating polynomial of degree i . a[i×tda+j] contains the coefficient a i + 1 , j + 1 , for i=0,1,,k and j=0,1,,i.
8: s[kplus1] double Output
On exit: s[i] contains the root mean square residual s i , for i=0,1,,k, as described in Section 3. For the interpretation of the values of the s i and their use in selecting an appropriate degree, see the E02 Chapter Introduction.
9: fail NagError * Input/Output
The NAG error argument (see Section 7 in the Introduction to the NAG Library CL Interface).

6 Error Indicators and Warnings

NE_2_INT_ARG_GT
On entry, kplus1=value while the number of distinct x values, mdist = value. These arguments must satisfy kplus1 mdist .
NE_2_INT_ARG_LT
On entry, tda=value while kplus1=value .
The arguments must satisfy tdakplus1 .
NE_ALLOC_FAIL
Dynamic memory allocation failed.
NE_INT_ARG_LT
On entry, kplus1 must not be less than 1: kplus1=value .
NE_NO_NORMALISATION
On entry, all the x[r] in the sequence x[r] , r = 0 , 1 , , m - 1 are the same.
NE_NOT_NON_DECREASING
On entry, the sequence x[r] , r = 0 , 1 , , m - 1 is not in nondecreasing order.
NE_WEIGHTS_NOT_POSITIVE
On entry, the weights are not strictly positive: w[value] = value.

7 Accuracy

No error analysis for the method has been published. Practical experience with the method, however, is generally extremely satisfactory.

8 Parallelism and Performance

e02adc is not threaded in any implementation.

9 Further Comments

The time taken by e02adc is approximately proportional to m k+1 k+11 .
The approximating polynomials may exhibit undesirable oscillations (particularly near the ends of the range) if the maximum degree k exceeds a critical value which depends on the number of data points m and their relative positions. As a rough guide, for equally spaced data, this critical value is about 2 × m . For further details see page 60 of Hayes (1970).

10 Example

Determine weighted least squares polynomial approximations of degrees 0, 1, 2 and 3 to a set of 11 prescribed data points. For the approximation of degree 3, tabulate the data and the corresponding values of the approximating polynomial, together with the residual errors, and also the values of the approximating polynomial at points half-way between each pair of adjacent data points.
The example program supplied is written in a general form that will enable polynomial approximations of degrees 0 , 1 , , k to be obtained to m data points, with arbitrary positive weights, and the approximation of degree k to be tabulated. e02aec is used to evaluate the approximating polynomial. The program is self-starting in that any number of datasets can be supplied.

10.1 Program Text

Program Text (e02adce.c)

10.2 Program Data

Program Data (e02adce.d)

10.3 Program Results

Program Results (e02adce.r)