# HW07: Solving Linear Least Squares Problems
**Yuanxing Cheng, 810925466**

This homework studies the polynomial fitting problem. Given a set of $m$ different points, $(x_i , y_i )$, $i = 1, 2, \dots, m$, in the $xy\text{-coordinate}$ plane, determine a polynomial of degree $n − 1$.

$$p(x) = a_1 x^{n-1} + a_2 x^{n-2} + \cdots + a_{n-1} x + a_n$$

where $n\leq m$, such that $\displaystyle \sum_{i=1}^{m} \left|p(x_i) - y_i\right|^2$ is minimized.

This leads to a linear least squares problem, where the coefficient matrix is given by the $m \times n$ Vandermonde matrix

$$\left[\begin{array}{ccccc}
x_{1}^{n-1} & x_{1}^{n-2} & \cdots & x_{1} & 1 \\
x_{2}^{n-1} & x_{2}^{n-2} & \cdots & x_{2} & 1 \\
x_{3}^{n-1} & x_{3}^{n-2} & \cdots & x_{3} & 1 \\
\vdots & \vdots & \ddots & \vdots & \vdots \\
x_{m}^{n-1} & x_{m}^{n-2} & \cdots & x_{m} & 1 \\
\end{array}\right]$$

the right hand side is the $m\text{−dimensional}$ column vector 
$\left[\begin{array}{cccc}
y_{1} & y_{2} & \cdots & y_{m} \\
\end{array}\right]'$, and the unknown is the $n\text{−dimensional}$ column vector $\left[\begin{array}{cccc}
a_{1} & a_{2} & \cdots & a_{n} \\
\end{array}\right]'$ containing the coefficients of p(x). Note that several MATLAB functions, e.g., `vander`, `polyval`, may be useful in this homework.

## Q1 
Take the following $m$ equally spaced points on the interval $−1 \leq x \leq 1$,

$$\begin{align}
x_i=-1+ {{2(i-1)}\over{m-1}} && y_i = \frac{1} {1 + 25x_i^2}
\end{align}$$

Here $i = 1,2,\cdots, m$. Use the following three different approaches to solve the above least squares problem, respectively, to obtain the polynomial $p(x)$ of degree $n − 1$. Here take $m = 50$ and $n = 20$.

($\text{a}$) The Matlab backslash function `\`, which solves the least squares problem when the matrix is an $m \times n$ matrix ($m > n$).

($\text{b}$) The reduced $QR$ factorization. Here the $QR$ factorization can be computed via the MATLAB built-in function `qr`.

($\text{c}$) Solution of the normal equation.

Compare the coefficients of $p(x)$ obtained from the above three approaches. Which two approaches generate results more close to each other. You need to use format long to observe more digits in the coefficients for comparison.

$Answer$

The distance between using Matlab backslash function and reduced QR factorization is the lowest.

$Code$

```MATLAB
clear
clc

m = 50;
n = 20;

x = linspace(-1,1,m);
x = x';
y = 1./(1+25*x.^2);
A = vander(x);

A1 = A(:,m-n+1:m);
x1 = A1\y;

[Q,R] = qr(A1,0);
y2 = Q'*y;
x2 = backwardSubeps(R,y2);

A3 = A1'*A1;
L = chol(A3,'lower');
y3 = L^(-1)*A1'*y;
x3 = backwardSubeps(L',y3);

d12 = norm(x1-x2);
d23 = norm(x2-x3);
d31 = norm(x3-x1);

s1 = 'Matlab backslash function';
s2 = 'reduced QR factorization';
s3 = 'normal equation';

if d12 < d23
    if d12 < d31
        s = [s1,' and ',s2];
    else
        s = [s3,' and ',s1];
    end
else
    if d23 < d31
        s = [s2,' and ',s3];
    else
        s = [s3,' and ',s1];
    end
end

disp(['The distance between using ',s,' is the lowest.'])

```

## Q2
For the case $m = 50$ and $n = 20$, draw the fitting polynomial $y = p(x)$ obtained from one approach in Part 1 and compare it in the same figure with the plot of the original function $y = 1/(1 + 25x^2 )$. Do they match well?

$Answer$

Yes, all three methods fit very well.



$Code$

```MATLAB
clear
clc

m = 50;
n = 20;

x = linspace(-1,1,m);
x = x';
y = 1./(1+25*x.^2);
A = vander(x);

A1 = A(:,m-n+1:m);
x1 = A1\y;

[Q,R] = qr(A1,0);
y2 = Q'*y;
x2 = backwardSubeps(R,y2);

A3 = A1'*A1;
L = chol(A3,'lower');
y3 = L^(-1)*A1'*y;
x3 = backwardSubeps(L',y3);

s1 = 'Matlab backslash function';
s2 = 'reduced QR factorization';
s3 = 'normal equation';

t = linspace(-1,1,100);

subplot(3,1,1)
Y1 = polyval(x1,t);
plot(x,y,'bo')
hold on
plot(t,Y1,'r')
title(s1)
```

```MATLAB
subplot(3,1,2)
Y2 = polyval(x2,t);
plot(x,y,'bo')
hold on 
plot(t,Y2,'r')
title(s2)

subplot(3,1,3)
Y3 = polyval(x3,t);
plot(x,y,'bo')
hold on
plot(t,Y3,'r')
title(s3)
```

## Q3
Repeat Part 2, but taking $m = n = 50$, which is now a polynomial interpolation problem. Does the interpolation polynomial match the function $y = 1/(1 + 25x^2)$ well in this case?

$Answer$

For the result from using Matlab backslash function and reduced QR factorization, they don't fit as well as before. And for the last method, it even doesn't work here because the matrix we now are using is no more full rank.



$Code$

```MATLAB
clear
clc

m = 50;
n = 50;

x = linspace(-1,1,m);
x = x';
y = 1./(1+25*x.^2);
A = vander(x);
```

```MATLAB
A1 = A(:,m-n+1:m);
x1 = A1\y;

[Q,R] = qr(A1,0);
y2 = Q'*y;
x2 = backwardSubeps(R,y2);

s1 = 'Matlab backslash function';
s2 = 'reduced QR factorization';

t = linspace(-1,1,100);

subplot(2,1,1)
Y1 = polyval(x1,t);
plot(x,y,'bo')
hold on
plot(t,Y1,'r')
title(s1)

subplot(2,1,2)
Y2 = polyval(x2,t);
plot(x,y,'bo')
hold on 
plot(t,Y2,'r')
title(s2)
```

20c7b85c5fd106ee58937cc454af39a8fd9dcb41
ipynb
Jupyter Notebook
Computational/Intro to Numerical Computing/HW_Chap07.ipynb
XavierOwen/Notes
d262a9103b29ee043aa198b475654aabd7a2818d
MIT
2018-11-27T10:31:08.000Z
2019-01-20T03:11:58.000Z
Computational/Intro to Numerical Computing/HW_Chap07.ipynb
XavierOwen/Notes
d262a9103b29ee043aa198b475654aabd7a2818d
MIT
Computational/Intro to Numerical Computing/HW_Chap07.ipynb
XavierOwen/Notes
d262a9103b29ee043aa198b475654aabd7a2818d
MIT
2020-07-14T19:57:23.000Z
2020-07-14T19:57:23.000Z
Qwen/Qwen-72B
1. YES
2. YES