We create two column vectors for data about the population of China. The first has the years of census data and the other has the population, in millions of people.
year = [1982; 2000; 2010; 2015];
pop = [1008.18; 1262.64; 1337.82; 1374.62];
It’s convenient to measure time in years since 1980.
t = year - 1980;
y = pop;
Now we have four data points ( t 1 , y 1 ) , … , ( t 4 , y 4 ) (t_1,y_1),\dots,(t_4,y_4) ( t 1 , y 1 ) , … , ( t 4 , y 4 ) , so n = 4 n=4 n = 4 and we seek an interpolating cubic polynomial. We construct the associated Vandermonde matrix:
To solve for the vector of polynomial coefficients, we use a backslash to solve the linear system:
A backslash \ is used to solve a linear system of equations.
The algorithms used by the backslash operator are the main topic of this chapter. As a check on the solution, we can compute the residual .
Using floating-point arithmetic, it is not realistic to expect exact equality of quantities; a relative difference comparable to ϵ mach \macheps ϵ mach is all we can look for.
By our definitions, the elements of c are coefficients in descending-degree order for the interpolating polynomial. We can use the polynomial to estimate the population of China in 2005:
p = @(t) polyval(c, t - 1980); % include the 1980 time shift
p(2005)The official population value for 2005 was 1303.72, so our result is rather good.
We can visualize the interpolation process. First, we plot the data as points.
The scatter function creates a scatter plot of points; you can specify a line connecting the points as well.
scatter(year, y)
xlabel("years since 1980")
ylabel("population (millions)")
title("Population of China")We want to superimpose a plot of the polynomial. We do that by evaluating it at a vector of points in the interval.
The linspace function constructs evenly spaced values given the endpoints and the number of values.
tt = linspace(1980, 2015, 500); % 500 times in the interval [1980, 2015]
yy = p(tt); % evaluate p at all the vector elements
yy(1:4)Now we use plot! to add to the current plot, rather than replacing it.
Use hold on to add to an existing plot rather than replacing it.
The plot function plots lines connecting the given x x x and y y y values; you can also specify markers at the points.
hold on
plot(tt, yy)
legend("data", "interpolant", "location", "northwest");