Weighted Least-Squares for a Nearly Perfect Min-Max Fit

Show more

1. Introduction

Finding the min-max, or best ${L}_{\infty}$ , polynomial approximation to a function, in some standard interval, is of the greatest interest in numerical analysis [1] [2] . For a polynomial function the least error distribution is a Chebyshev polynomial [3] [4] [5] .

The usual procedure [6] [7] to find the best ${L}_{\infty}$ approximation to a general function is to start with a good approximation, say in the ${L}_{2}$ sense, easily obtained by the minimization of a quadratic functional for the coefficients, then iteratively improving this initial approximation by a Remez-like correction procedure [8] [9] that strives to produce an error distribution that oscillates with a constant amplitude in the interval of interest.

In this note, we bring ample and varied computational evidence in support of the novel, worthy of notice, empirical numerical observation that taking the error distribution of a least squares, ${L}_{2}$ , best polynomial fit to a function, squared, as weight in a second, weighted, least squares approximation, results in an error distribution that is remarkably close to the best ${L}_{\infty}$ , or uniform, approximation.

2. Fixing Ideas; The Best Quadratic in [−1, 1]

The monic Chebyshev polynomial

${T}_{2}\left(x\right)={x}^{2}-\frac{1}{2},\text{\hspace{0.17em}}-1\le x\le 1$ (1)

is the solution of the min-max problem

$\underset{a}{\mathrm{min}}\underset{x}{\mathrm{max}}e\left(x\right),\text{\hspace{0.17em}}e\left(x\right)={x}^{2}-a,\text{\hspace{0.17em}}-1\le x\le 1.$ (2)

This min-max solution, the least function in the ${L}_{\infty}$ sense, is a polynomial that has two distinct roots, and oscillates with a constant amplitude in $-1\le x\le 1,$ $e\left(-1\right)=-e\left(0\right)=e\left(1\right).$ Indeed, say ${e}_{1}={x}^{2}+{a}_{0}+{a}_{1}x$ is such a polynomial, and ${e}_{2}={x}^{2}+{p}_{0}+{p}_{1}x$ is another quadratic polynomial, then ${e}_{1}\le {e}_{2}$ in the interval, for otherwise ${e}_{1}$ and ${e}_{2}$ would intersect at two points, which is absurd; ${x}^{2}+{a}_{0}+{a}_{1}x={x}^{2}+{p}_{0}+{p}_{1}x$ is either an identity, or has but the one solution $x=-\left({p}_{0}-{a}_{0}\right)/\left({p}_{1}-{a}_{1}\right)$ .

Thus, the monic Chebyshev polynomial of degree n is the least, uniform, or pointwise, error distribution in approximating ${x}^{n}$ by a polynomial of degree $n-1$ .

To obtain a least squares, a best ${L}_{2}$ , approximation to ${T}_{2}\left(x\right)$ we first minimize $I\left(a\right)$

$I\left(a\right)={\displaystyle {\int}_{-1}^{1}}{\left({x}^{2}-a\right)}^{2}\text{d}x,\text{\hspace{0.17em}}{I}^{\prime}\left(a\right)={\displaystyle {\int}_{-1}^{1}}\left({x}^{2}-a\right)\text{d}x=0$ (3)

to have the value $a=1/3=0.3333$ .

Minimizing next $I\left(p\right)$ , under the weight ${\left({x}^{2}-a\right)}^{2},a=1/3$

#Math_28# (4)

now with respect to p, we obtain $p=11/21=0.5238$ , which is surprisingly much closer to the optimal value of one half.

We may replace the difficult
${L}_{\infty}$ measure by the computationally easier
${L}_{m}$ measure with an even
$m\gg 1$ . Let a_{0} be a good approximation, and
${a}_{1}={a}_{0}+\delta $ be an improved one. Minimization cum linearization produces the equation

${\int}_{-1}^{1}}{\left({x}^{2}-{a}_{0}\right)}^{n}\text{d}x-n\delta {\displaystyle {\int}_{-1}^{1}}{\left({x}^{2}-{a}_{0}\right)}^{n-1}\text{d}x=0$ (5)

where $n\gg 1$ is odd.

Starting with ${a}_{0}=11/21=0.5238$ , we obtain from the above equation, for $n=17$ , the value ${a}_{1}=0.495$ , as compared with the optimal $a=0.5$ .

3. Optimal Cubic in [−1, 1]

Seeking to reproduce the optimal monic Chebyshev polynomial of degree three

${T}_{3}\left(x\right)={x}^{3}-\frac{3}{4}x,\text{\hspace{0.17em}}-1\le x\le 1$ (6)

we start by minimizing $I\left({a}_{1}\right)$

$I\left({a}_{1}\right)={\displaystyle {\int}_{-1}^{1}}{\left({x}^{3}-{a}_{1}x\right)}^{2}\text{\hspace{0.05em}}\text{d}x,\text{\hspace{0.17em}}\text{\hspace{0.17em}}{I}^{\prime}\left({a}_{1}\right)={\displaystyle {\int}_{-1}^{1}}\text{\hspace{0.05em}}x\left({x}^{3}-{a}_{1}x\right)\text{d}x=0$ (7)

and have ${a}_{1}=3/5=0.6$ .

Then we return to minimize the weighted $I\left({p}_{1}\right)$ with respect to ${p}_{1}$

$\begin{array}{l}I\left({p}_{1}\right)={\displaystyle {\int}_{-1}^{1}}{\left({x}^{3}-{a}_{1}x\right)}^{2}{\left({x}^{3}-{p}_{1}x\right)}^{2}\text{d}x,\\ {I}^{\prime}\left({p}_{1}\right)={\displaystyle {\int}_{-1}^{1}}\text{\hspace{0.05em}}\text{\hspace{0.05em}}x{\left({x}^{3}-{a}_{1}x\right)}^{2}\left({x}^{3}-{p}_{1}x\right)\text{d}x=0\end{array}$ (8)

and obtain ${p}_{1}=195/253=0.770751$ , which is considerably closer to the optimal value of 0.75. See Figure 1.

We are ready now for a Remez-like correction to bring the error function closer to optimal. The minimum of $e\left(x\right)={x}^{3}-0.770751x$ occurs at m = 0.50687. We write a new tentative $e\left(x\right)={x}^{3}-{a}_{1}x$ and request that $-e\left(m\right)=e\left(1\right)$ , by which we have

${a}_{1}=\frac{1+{m}^{3}}{1+m}=0.750047$ (9)

as compared with the Chebyshev optimal value of ${a}_{1}=3/4=0.75$ .

4. Optimal Quartic in [0, 1]

Starting with

$e\left(x\right)={x}^{4}+{a}_{3}{x}^{3}+{a}_{2}{x}^{2}+{a}_{1}x+{a}_{0}$ (10)

we minimize

$I\left({a}_{0}\mathrm{,}{a}_{1}\mathrm{,}{a}_{2}\mathrm{,}{a}_{3}\right)={\displaystyle {\int}_{0}^{1}}\text{\hspace{0.05em}}\text{\hspace{0.05em}}e{\left(x\right)}^{2}\text{d}x$ (11)

and obtain the best, in the ${L}_{2}$ sense, $e\left(x\right)$ shown in Figure 2.

Then we return to minimize

$I\left({p}_{0},{p}_{1},{p}_{2},{p}_{3}\right)={\displaystyle {\int}_{0}^{1}}\text{\hspace{0.05em}}\text{\hspace{0.05em}}e{\left(x\right)}^{2}{\left({x}^{4}+{p}_{3}{x}^{3}+{p}_{2}{x}^{2}+{p}_{1}x+{p}_{0}\right)}^{2}\text{d}x$ (12)

weighted by the previous $e\left(x\right)$ squared, and obtain the new, nearly perfectly uniform $e\left(x\right)$ of Figure 3.

By comparison, the amplitude of the monic Chebyshev polynomial of degree four in [0,1] is 1/128 = 0.0078125.

Figure 1. (a) Least squares cubic. (b) Weighted least squares cubic.

Figure 2. Least squares quartic.

Figure 3. Weighted least squares quartic.

5. Best Cubic Approximation of e^{x} in [0, 1]

To facilitate the integrations we use the approximation

${\text{e}}^{x}=1+x+\frac{1}{2!}{x}^{2}+\frac{1}{3!}{x}^{3}+\frac{1}{4!}{x}^{4}+\frac{1}{5!}{x}^{5}+\frac{1}{6!}{x}^{6}+\frac{1}{7!}{x}^{7}$ (13)

and minimize

$I\left({a}_{0},{a}_{1},{a}_{2},{a}_{3}\right)={\displaystyle {\int}_{0}^{1}}\text{\hspace{0.05em}}\text{\hspace{0.05em}}e{\left(x\right)}^{2}\text{d}x,\text{\hspace{0.17em}}e\left(x\right)={\text{e}}^{x}+{a}_{0}+{a}_{1}x+{a}_{2}{x}^{2}+{a}_{3}{x}^{3}$ (14)

with respect to ${a}_{0}\mathrm{,}{a}_{1}\mathrm{,}{a}_{2}\mathrm{,}{a}_{3}$ . The best $e\left(x\right)$ obtained from this minimization is shown in Figure 4.

Then we use the square of the minimal $e\left(x\right)$ just obtained, as weight in the next minimization of

$I\left({p}_{0},{p}_{1},{p}_{2},{p}_{3}\right)={\displaystyle {\int}_{0}^{1}}\text{\hspace{0.05em}}\text{\hspace{0.05em}}e{\left(x\right)}^{2}{\left({\text{e}}^{x}+{p}_{0}+{p}_{1}x+{p}_{2}{x}^{2}+{p}_{3}{x}^{3}\right)}^{2}\text{d}x$ (15)

with respect to ${p}_{0}\mathrm{,}{p}_{1}\mathrm{,}{p}_{2}\mathrm{,}{p}_{3}$ .

The nearly perfect result of this last minimization is shown in Figure 5.

6. Best Cubic Approximation of sinx in [0, 1]

To facilitate the integrations we take

$\mathrm{sin}x=x-\frac{1}{3!}{x}^{3}+\frac{1}{5!}{x}^{5}-\frac{1}{7!}{x}^{7}+\frac{1}{9!}{x}^{9}$ (16)

and obtain the least squares error distribution as in Figure 6.

The subsequent nearly perfect weighted least squares error distribution is shown in Figure 7.

Figure 4. Least squares cubic fit to e^{x}.

Figure 5. Weighted least squares cubic fit to e^{x}.

Figure 6. Least squares cubic fit to sinx.

Figure 7. Weighted least squares cubic fit to sinx.

7. Best Quadratic Fit to $\sqrt{x}$ in [0, 1]

We start with

$e\left(x\right)=\sqrt{x}-\left({a}_{0}+{a}_{1}x+{a}_{2}{x}^{2}\right),\text{\hspace{0.17em}}0\le x\le 1$ (17)

under the condition

$e\left(0\right)=-e\left(1\right),\text{\hspace{0.17em}}{a}_{0}=\frac{1}{2}\left(1-{a}_{1}-{a}_{2}\right)$ (18)

and minimize

$I\left({a}_{1},{a}_{2}\right)={\displaystyle {\int}_{0}^{1}}{\left(\sqrt{x}-\frac{1}{2}-{a}_{1}\left(x-\frac{1}{2}\right)-{a}_{2}\left({x}^{2}-\frac{1}{2}\right)\right)}^{2}\text{d}x$ (19)

with respect to ${a}_{1}$ and ${a}_{2}$ , to have

$e\left(x\right)=\sqrt{x}-\left(\frac{1}{10}+\frac{121}{70}x-\frac{13}{14}{x}^{2}\right),\text{\hspace{0.17em}}0\le x\le 1$ (20)

shown as curve a in Figure 8.

Next we minimize

Figure 8. (a) Least squares quadratic fit to $\sqrt{x}$ . (b) Weighted least squares quadratic fit to $\sqrt{x}$ .

$\begin{array}{l}I\left({p}_{1},{p}_{2}\right)={\displaystyle {\int}_{0}^{1}}{\left(\sqrt{x}-\frac{1}{2}-{p}_{1}\left(x-\frac{1}{2}\right)-{p}_{2}\left({x}^{2}-\frac{1}{2}\right)\right)}^{2}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\cdot {\left(\sqrt{x}-\left(\frac{1}{10}+\frac{121}{70}x-\frac{13}{14}{x}^{2}\right)\right)}^{2}\text{d}x\end{array}$ (21)

and obtain

$e\left(x\right)=\sqrt{x}-\left(0.064+1.949x-1.077{x}^{2}\right),\text{\hspace{0.17em}}0\le x\le 1$ (22)

shown as graph b in Figure 8, as compared with the optimal, in the ${L}_{\infty}$ sense

$e\left(x\right)=\sqrt{x}-\left(0.0674385+1.93059x-1.06547{x}^{2}\right),\text{\hspace{0.17em}}0\le x\le 1.$ (23)

8. Best Cubic Fit to x^{1/4} in [0, 1]

We start with

$e\left(x\right)={x}^{1/4}+{a}_{0}+{a}_{1}x+{a}_{2}{x}^{2}+{a}_{3}{x}^{3},\text{\hspace{0.17em}}0\le x\le 1$ (24)

under the restriction $e\left(0\right)=e\left(1\right)$ , or ${a}_{3}=-1-{a}_{1}-{a}_{2}$ , and minimize

$I\left({a}_{0},{a}_{1},{a}_{2}\right)={\displaystyle {\int}_{0}^{1}}{\left({x}^{1/4}-{x}^{3}+{a}_{0}+{a}_{1}\left(x-{x}^{3}\right)+{a}_{2}\left({x}^{2}-{x}^{3}\right)\right)}^{2}\text{d}x$ (25)

with respect to ${a}_{0}\mathrm{,}{a}_{1}\mathrm{,}{a}_{2}$ to have the minimal $e\left(x\right)$ shown in Figure 9.

Then we minimize

$I\left({p}_{0},{p}_{1},{p}_{2}\right)={\displaystyle {\int}_{0}^{1}}\text{\hspace{0.05em}}\text{\hspace{0.05em}}e{\left(x\right)}^{2}{\left({x}^{1/4}-{x}^{3}+{p}_{0}+{p}_{1}\left(x-{x}^{3}\right)+{p}_{2}\left({x}^{2}-{x}^{3}\right)\right)}^{2}\text{d}x$ (26)

and obtain the nearly optimal error distribution as in Figure 10.

9. Another Difficult Function

We now look at the error distribution

$e\left(x\right)=\mathrm{ln}\left(1.001+x\right)-\left({a}_{3}{x}^{3}+{a}_{2}{x}^{2}+{a}_{1}x+{a}_{0}\right),\text{\hspace{0.17em}}-1\le x\le 1$ (27)

under the condition that $e\left(1\right)=e\left(-1\right)$ , or ${a}_{3}=3.8007012-{a}_{1}.$

Least squares minimization of $e\left(x\right)$ yields the error distribution in Figure 11.

Next we minimize

#Math_100# (28)

Figure 9. Least squares cubic fit to ${x}^{1/4}$ .

Figure 10. Weighted least squares cubic fit to ${x}^{1/4}$ .

Figure 11. Least squares cubic fit to $ln\left(1.001+x\right)$ .

Figure 12. Weighted least squares cubic fit to $ln\left(1.001+x\right)$ .

under the restriction that ${p}_{3}=3.8007012-{p}_{1}$ , and obtain the nearly perfect error distribution shown in Figure 12.

10. Conclusion

We experimentally demonstrate, on a variety of continuous, analytic and nonanalytic functions, the remarkable observation that if the least squares polynomial approximation is taken as weight in a repeated, now weighted, least squares approximation, then this new, second, approximation is nearly perfect in the sense of Chebyshev, barely needing any further correction procedure.

References

[1] Linz, P. and Wang, R. (2003) Exploring Numerical Methods: An introduction to Scientific Computing Using MATLAB. Jones and Bartlett, Boston.

[2] Meinardus, G. (1967) Approximation of Functions: Theory and Numerical Methods, Springer, Berlin.

https://doi.org/10.1007/978-3-642-85643-3

[3] Veidinger, L. (1960) On the Numerical Determination of the Best Approximations in the Chebychev Sense. Numerische Mathematik, 2, 99-105.

https://doi.org/10.1007/BF01386215

[4] Snyder, M.A. (1966) Chebyshev Methods in Numerical Approximation. Prentice-Hall, New Jersey.

[5] Cody, W.J. (1970) A Survey of Practical Rational and Polynomial Approximation of Functions. SIAM Review, 12, 400-423.

https://doi.org/10.1137/1012082

[6] Fraser, W. (1965) A Survey of Methods of Computing Minimax and Near-Minimax Polynomial Approximations for Functions of a Single Independent Variable. J. ACM, 12, 295-314.

https://doi.org/10.1145/321281.321282

[7] Psarakis, E.Z. and Moustakides G.V. (2003) A Robust Initialization Scheme for the Remez Exchange Algorithm. IEEE Signal Processing Letters, 10, 1-3.

https://doi.org/10.1109/LSP.2002.806701

[8] Huddleston, R.E. (1974) On the Conditional Equivalence of Two Starting Methods for the Second Algorithm of Remez. Mathematics of Computation, 28, 569-572.

https://doi.org/10.1090/S0025-5718-1974-0341804-6

[9] Pachon, R. and Trefethen, L.N. (2009) Barycentric-Remez Algorithms for Best Polynomial Approximation in the chebfun System. BIT Numerical Mathematics, 49, 721-741.

https://doi.org/10.1007/s10543-009-0240-1