Least-Squares Polynomials

• Fitting polynomials to data that do not plot linearly is common.
• It will turn out that the normal equations are linear for this situation (an added advantage).
• $n$ as the degree of the polynomial
• N as the number of data pairs.
• If $N = n + 1$, the polynomial passes exactly through each point and the methods discussed earlier apply,
• so we will always have $N > n+ 1$.
• We assume the functional relationship
  

  $$y=a_0+a_1x+a_2x^2+\ldots+a_nx^n$$ (1)

  
  

• With errors defined by
  
  $$e_i=Y_i - y_i=Y_i- a_0-a_1x_i-a_2x_i^2-\ldots-a_nx_i^n$$
  
  

• We again use $Y_i$ to represent the observed (experimental) value corresponding to $x_i$ (it is assumed that $x_i$ free of error for the sake of simplicity).
• We minimize the sum of squares;
  
  $$S=\sum_{i=1}^N e_i^2=\sum_{i=1}^N (Y_i - a_0 - a_1x_i - a_2x_i^2 - \ldots - a_nx_i^n)^2$$
  
  

• At the minimum, all the partial derivatives $\partial S/\partial a_0,\partial S/\partial a_n$ vanish.
• Writing the equations for these gives $n + 1$ equations:
  
  $$\begin{array}{l} \frac{\partial S}{\partial a_0}=0=\sum_{i=1... ..._i- a_0-a_1x_i-a_2x_i^2-\ldots-a_ix_i^n)(-x_i^n)\\ \end{array}$$
  
  

• Dividing each by $-2$ and rearranging gives the $n + 1$ normal equations to be solved simultaneously:
  

  $$\begin{array}{rl} a_0N+a_1\sum x_i+a_2 \sum x_i^2+\ldots+a_n... ...i^{n+2}+\ldots+a_n\sum x_i^{2n}= & \sum x_i^nY_i\\ \end{array}$$ (2)

  
  

• Putting these equations in matrix form shows the coefficient matrix (B).

$$\overbrace{\left[ \begin{array}{rrrrrl} N & \sum x_i & \sum x_i^2... ...Y_i\\ \sum x_i^2Y_i\\ \vdots \\ \sum x_i^n Y_i\\ \end{array} \right]$$ (3)

All the summations in Eqs. 2 and 3 run from 1 to $N$.

• Equation 3 represents a linear system.
• However, you need to know that if this system is ill-conditioned and round-off errors can distort the solution: the $a$'s of Eq. 1.
• Up to degree-3 or -4, the problem is not too great.
• Special methods that use orthogonal polynomials are a remedy.
• Degrees higher than 4 are used very infrequently.
• It is often better to fit a series of lower-degree polynomials to subsets of the data.
• Matrix $B$ of Eq. 3 is called the normal matrix for the least-squares problem.

• There is another matrix that corresponds to this, called the design matrix.
• It is of the form;
  
  $$A=\left[ \begin{array}{rrrrr} 1 & 1 & 1 & 1 & 1 \\ x_1 & ... ... x_1^n & x_2^n & x_3^n & \ldots & x_N^n\\ \end{array} \right]$$
  
  

• $AA^T$ is just the coefficient matrix of Eq. 3.
• It is easy to see that $Ay$, where $y$ is the column vector of $y$-values, gives the right-hand side of Eq. 3.
  

  $$\overbrace{\left[ \begin{array}{rrrrr} 1 & 1 & 1 & 1 & 1 \\ x... ...Y_i\\ \sum x_i^2Y_i\\ \vdots \\ \sum x_i^n Y_i\\ \end{array} \right]$$ (4)

  
  

 $$

We can rewrite Eq. 3 in matrix form, as

$$AA^Ta = Ba= Ay$$

1

$AA^T=B$. To find the solution (with MATLAB) $>> a = Ay\textbackslash A*transpose(A)$

$$\overbrace{\left[ \begin{array}{rrrrr} 1 & 1 & 1 & \ldots &... ...^n & x_2^n & x_3^n & \ldots & x_N^n\\ \end{array} \right]}^A*$$

$$\overbrace{\left[ \begin{array}{rrrrr} 1 & x_1 & x_1^2 & \l... ...& x_N & x_N^2 & \ldots & x_N^n \\ \end{array} \right]}^{A^T}=$$

$$\underbrace{\left[ \begin{array}{rrrrrl} N & \sum x_i & \sum... ... x_i^{n+3} & \ldots & \sum x_i^{2n}\\ \end{array} \right]}_B$$

2

$A^Ta=y$

$$\overbrace{\left[ \begin{array}{rrrrr} 1 & x_1 & x_1^2 & \l... ...& x_N & x_N^2 & \ldots & x_N^n \\ \end{array} \right]}^{A^T}*$$

$$\overbrace{\left[ \begin{array}{r} a_0\\ a_1\\ a_2\\ \ldots\\ a_n\\ \end{array} \right]}^{a}=$$

$$\overbrace{\left[ \begin{array}{r} y_1\\ y_2\\ y_3\\ \ldots\\ y_N\\ \end{array} \right]}^{y}$$

 $$

That is

$$\begin{array}{l} a_0+a_1x_1+a_2x_1^2+\ldots+a_nx_1^n=y_1\\ ... ...dots\\ a_0+a_1x_N+a_2x_N^2+\ldots+a_nx_N^n=y_N \\ \end{array}$$

 $$

Least-squares polynomials with all $x$-values (from given $xy$-pair data) are inserted.

• It is illustrated the use of Eqs. 2 to fit a quadratic to the data of Table 1.
  
  

  Table 1: Data to illustrate curve fitting.

  
  
  
  
• To set up the normal equations, we need the sums tabulated in Table 1. The equations to be solved are:
  
  $$\begin{array}{rl} 11a_0+6.01a_1+4.6545a_2& =5.905 \\ 6.01... ...839 \\ 4.6545a_0+4.1150a_1+3.9161a_2& =1.3357 \\ \end{array}$$
  
  

• The result is $a_0= 0.998$, $a_2=-1.018$,$a_3 = 0.225$, so the least- squares method gives
  
  $$y = 0.998 - 1.018x + 0.225x^2$$
  
  

• which we compare to $y=1- x +0.2x^2$.
• Errors in the data cause the equations to differ.

• Figure 1 shows a plot of the data.
• The data are actually a perturbation of the relation $y=1- x +0.2x^2$.

Figure 1: Figure for the data to illustrate curve fitting.

• Example: The following data:

  R/C: 0.73, 0.78, 0.81, 0.86, 0.875, 0.89, 0.95, 1.02, 1.03, 1.055, 1.135, 1.14, 1.245, 1.32, 1.385, 1.43, 1.445, 1.535, 1.57, 1.63, 1.755.

   
   
  

  $V_\theta /V_\infty$: 0.0788, 0.0788, 0.064, 0.0788, 0.0681, 0.0703, 0.0703, 0.0681, 0.0681, 0.079, 0.0575, 0.0681, 0.0575, 0.0511, 0.0575, 0.049, 0.0532, 0.0511, 0.049, 0.0532,0.0426.

   
   
  

• Let $x = R/C$ and $y=V_\theta/V_\infty$,
• We would like our curve to be of the form
  
  $$g(x)=\frac{A}{x}(1-e^{-\lambda x^2})$$
  
  

• and our least-squares equation becomes
  
  $$S = \sum_{i=1}^{21}(Y_i-\frac{A}{x_i}(1-e^{-\lambda x_i^2}))^2$$
  
  

• Setting $S_\lambda=S_A=0$ gives the following equations:
  
  $$\begin{array}{l} \sum_{i=1}^{21}(\frac{1}{x_i})(1-e^{-\lambd... ...^2})(Y_i-\frac{A}{x_i}(1-e^{-\lambda x_i^2}))=0 \\ \end{array}$$
  
  

• When this system of nonlinear equations is solved, we get
  
  $$g(x)=\frac{0.07618}{x}(1-e^{-2.30574x^2})$$
  
  

• For these values of $A$ and $\lambda, S = 0.00016$.
• The graph of this function is presented in Figure 2.

Figure 2: The graph of $V_\theta /V_\infty$ vs $R/C$.