Fitting a Curve Using Polyfit and Polyval

Donations

£1.00

Generating Model Data

Let’s generate some model x and y data. First we will create an independent variable x that is a column vector ranging from 0 to 20 in steps of 2.

x=[0:2:20]'

We will then create the dependent variable y which will be a model quadratic function of x.

y=0.1*x.^2+1.2*x+5

Note that we use .^ to perform element by element exponentiation opposed to ^ which would perform array exponentiation. It then follows that y is a 2nd order polynomial with coefficients

p2(1)=0.1
p2(2)=1.2
p2(3)=5

For convenience we can save these as the first and second column in a sheet rawdata

rawdata(:,1)=x

rawdata(:,2)=y


Plotting the Raw Data

To suppress the output in the command window, we can end each of the lines above in a semicolon ;

We can instead go ahead and plot the second column of rawdata (y) with respect to the first column of rawdata (x) as figure(2). We want blue filled markers of size 8, no line. For the axes we want a fontsize of 20, gridlines and x and y labels. If unsure of any of this code see Plotting 2D Figures in Octave/MATLAB.

figure(2)
plot(rawdata(:,1),rawdata(:,2),'o','MarkerSize',8,'MarkerEdge','b','MarkerFace','b');
grid on;
grid minor;
set(gca,'FontSize',20)
xlabel('x')
ylabel('y')


Using polyfit to Find the Polynomial Coefficients

We can use polyfit to evaluate the polynomial coefficients. Polyfit has the form

[Pn]=polyfit(x,y,n)

Where the input arguments are the x data, y data and n, the order of polynomial. The output arguments is a row vector Pn of the polynomial coefficients. In most cases we would try fitting from a first order polynomial upwards however as we know our data is a second order polynomial, as it is created from a second order polynomial we can try fitting to a second order polynomial first.

[P2]=polyfit(rawdata(:,1),rawdata(:,2),2)

[P2]=polyfit(rawdata(:,1),rawdata(:,2),2)

\displaystyle \text{P2=}\left[ {\begin{array}{*{20}{c}} {0.100000} & {1.200000} & {5.000000} \end{array}} \right]

These are the three polynomial coefficients for the quadratic equation

\displaystyle 0.100000{{x}^{2}}+1.200000x+5.000000

P2(1) is the x.^2 coefficient P2(2) is the x coefficient and P2(3) is the constant.

Which as expected is a perfect match to our original equation

y=0.1*x.^2+1.2*x+5

p2(1)=0.1
p2(2)=1.2
p2(3)=5


Using polyval to Evaluate the Polynomial

We can use polyval to evaluate the polynomial at specified x values. Let’s first create a new variable called fitteddata and to it’s first column we can specify new x values. For illustration purposes, to make the new x values slightly different to the original x data we can use all the odd numbers opposed to the even numbers. Thus

fitteddata(:,1)=[1:2:19]'

If you prefer you can use the same original x data to evaluate the polynomial.

Now that we have the x data to evaluate fitteddata(:,1) the polynomial with coefficients in P2 we can use Polyval. Polyval has the form

[evaluatedydata]=polyval(Pn,evaluatedxdata)

In this example I will created a new column in fitteddata, column 3 (as I will later use column 2 to look at a linear fit)

fitteddata(:,3)=polyval(P2,fitteddata(:,1))



Plotting the Fitted Data with Respect to the Raw Data

We can once again suppress the output into the command window and instead look at graphical output. For convenience, the graphing will be moved to the bottom, at the bottom of the graph we can add the hold on and hold off commands and plot the fitted data in the same figure as the original data

figure(2)
plot(rawdata(:,1),rawdata(:,2),'o','MarkerSize',8,'MarkerEdge','b','MarkerFace','b');
grid on;
grid minor;
set(gca,'FontSize',20)
xlabel('x')
ylabel('y')
hold on;
plot(fitteddata(:,1),fitteddata(:,3),'o','MarkerSize',8,'MarkerEdge','g','MarkerFace','g');
hold off;

It is quite common to plot the raw data as circular markers and the fitted data as a red dotted line. We can amend the above to do this:

figure(2)
plot(rawdata(:,1),rawdata(:,2),'o','MarkerSize',8,'MarkerEdge','b','MarkerFace','b');
grid on;
grid minor;
set(gca,'FontSize',20)
xlabel('x')
ylabel('y')
hold on;
plot(fitteddata(:,1),fitteddata(:,3),'r:','LineWidth',5);
hold off;


Underfitting – Fitting a 1st Order Polynomial to 2nd Order Polynomial Data

To fit a second order polynomial we used lines 6-20.


x=[0:2:20]';
y=0.1*x.^2+1.2*x+5;
rawdata(:,1)=x;
rawdata(:,2)=y;

[P2]=polyfit(rawdata(:,1),rawdata(:,2),2);

fitteddata(:,1)=[1:2:19]';
fitteddata(:,3)=polyval(P2,fitteddata(:,1))

figure(2)
plot(rawdata(:,1),rawdata(:,2),'o','MarkerSize',8,'MarkerEdge','b','MarkerFace','b');
grid on;
grid minor;
set(gca,'FontSize',20);
xlabel('x');
ylabel('y');
hold on
plot(fitteddata(:,1),fitteddata(:,3),'r:','LineWidth',5);
hold off

These can be copied and pasted:


x=[0:2:20]';
y=0.1*x.^2+1.2*x+5;
rawdata(:,1)=x;
rawdata(:,2)=y;

[P2]=polyfit(rawdata(:,1),rawdata(:,2),2);

fitteddata(:,1)=[1:2:19]';
fitteddata(:,3)=polyval(P2,fitteddata(:,1))

figure(2)
plot(rawdata(:,1),rawdata(:,2),'o','MarkerSize',8,'MarkerEdge','b','MarkerFace','b');
grid on;
grid minor;
set(gca,'FontSize',20);
xlabel('x');
ylabel('y');
hold on
plot(fitteddata(:,1),fitteddata(:,3),'r:','LineWidth',5);
hold off

[P2]=polyfit(rawdata(:,1),rawdata(:,2),2);

fitteddata(:,1)=[1:2:19]';
fitteddata(:,3)=polyval(P2,fitteddata(:,1))

figure(2)
plot(rawdata(:,1),rawdata(:,2),'o','MarkerSize',8,'MarkerEdge','b','MarkerFace','b');
grid on;
grid minor;
set(gca,'FontSize',20);
xlabel('x');
ylabel('y');
hold on
plot(fitteddata(:,1),fitteddata(:,3),'r:','LineWidth',5);
hold off

Lines 22-36 are now a duplicate of the original code. We can now go ahead and modify lines 6-20 to fit the data to a first order polynomial and to plot in figure 1 opposed to figure 2.

Line 6 is changed from

[P2]=polyfit(rawdata(:,1),rawdata(:,2),2);

to

[P1]=polyfit(rawdata(:,1),rawdata(:,2),1);

Line 9 is changed from

fitteddata(:,3)=polyval(P2,fitteddata(:,1))

to

fitteddata(:,2)=polyval(P1,fitteddata(:,1))

Line 11 is changed from

figure(2)

to

figure(1)

Line 19 is changed from

plot(fitteddata(:,1),fitteddata(:,3),'r:','LineWidth',5);

to

plot(fitteddata(:,1),fitteddata(:,2),'r:','LineWidth',5);

Here, in contrast to the second order polynomial fit we can immediately see that the 1st order polynomial linear fit does not follow all the data points.

\displaystyle \text{P1=}\left[ {\begin{array}{*{20}{c}} {3.20000} & {-1.00000} \end{array}} \right]

P1(1) is the x coefficient P1(2) is the constant.


Overfitting – Fitting a 3rd Order Polynomial to 2nd Order Polynomial Data

Once again we can copy and paste lines 6-20. This time to lines 38-52. In order to fit a third order polynomial we will need to amend the following.

Line 38 is changed from

[P2]=polyfit(rawdata(:,1),rawdata(:,2),2);

to

[P3]=polyfit(rawdata(:,1),rawdata(:,2),3);

Line 41 is changed from

fitteddata(:,3)=polyval(P2,fitteddata(:,1))

to

fitteddata(:,4)=polyval(P3,fitteddata(:,1))

Line 43 is changed from

figure(2)

to

figure(3)

Line 51 is changed from

plot(fitteddata(:,1),fitteddata(:,3),'r:','LineWidth',5);

to

plot(fitteddata(:,1),fitteddata(:,4),'r:','LineWidth',5);

This time the data is overfitted and hence follows the data points very closely

\displaystyle \text{P3=}\left[ {\begin{array}{*{20}{c}} {7.5775\text{e}-\text{18}} & {1.0000\text{e}-\text{01}} & {\text{1}\text{.2000e}+\text{00}} & {\text{5}\text{.0000e}+\text{00}} \end{array}} \right]

P3(1) is the x.^3 coefficient, P3(2) is the x.^2 coefficient P3(3) is the x coefficient and P3(4) is the constant.

In this case we have overfitted a model quadratic 2nd order polynomial data to a cubic function and as a consequence the coefficient of the additional parameter P3(1)≈0 i.e. we can say this additional parameter is not required to fit this data.

Subplots

It is perhaps useful to plot the linear, quadratic and cubic fits together in one figure using subplots. It may also be useful to print the equation on the graphs. To do this we can rearrange the code, so all the graphing commands are at the bottom


x=[0:2:20]';
y=0.1*x.^2+1.2*x+5;
rawdata(:,1)=x;
rawdata(:,2)=y;

[P1]=polyfit(rawdata(:,1),rawdata(:,2),1);
fitteddata(:,1)=[1:2:19]';
fitteddata(:,2)=polyval(P1,fitteddata(:,1))

[P2]=polyfit(rawdata(:,1),rawdata(:,2),2);
fitteddata(:,1)=[1:2:19]';
fitteddata(:,3)=polyval(P2,fitteddata(:,1))

[P3]=polyfit(rawdata(:,1),rawdata(:,2),3);
fitteddata(:,1)=[1:2:19]';
fitteddata(:,4)=polyval(P3,fitteddata(:,1))

figure(1)

plot(rawdata(:,1),rawdata(:,2),'o','MarkerSize',8,'MarkerEdge','b','MarkerFace','b');
grid on;
grid minor;
set(gca,'FontSize',20);
xlabel('x');
ylabel('y');
hold on
plot(fitteddata(:,1),fitteddata(:,2),'r:','LineWidth',5);

hold off

figure(2)
plot(rawdata(:,1),rawdata(:,2),'o','MarkerSize',8,'MarkerEdge','b','MarkerFace','b');
grid on;
grid minor;
set(gca,'FontSize',20);
xlabel('x');
ylabel('y');
hold on
plot(fitteddata(:,1),fitteddata(:,3),'r:','LineWidth',5);

hold off

figure(3)
plot(rawdata(:,1),rawdata(:,2),'o','MarkerSize',8,'MarkerEdge','b','MarkerFace','b');
grid on;
grid minor;
set(gca,'FontSize',20);
xlabel('x');
ylabel('y');
hold on
plot(fitteddata(:,1),fitteddata(:,4),'r:','LineWidth',5);

hold off

To modify the code above to make a subplot I’m going to add Line 19:

subplot(3,1,1)

Replace Line 31 with

subplot(3,1,2)

Replace Line 43 with

subplot(3,1,3)

A few lines will need to be added to print the equation on each subplot. The linear equation looks like

\displaystyle y=\text{P1(1)}x+\text{P1(2)}

We will need to put together the equation as a string. For this we use [ ] to make an array and we type the equation as a string enclosed in ' '. Opposed to typing in the coefficients we want to read them from P1. We can use these as input arguments for the num2str function, this function will convert these inputs into strings.

FitEqnP1=[y= num2str(P1(1)) 'x +' num2str(P1(2))]

Now we can create a label in the form of a cell array using the string variable FitEqnP1. The { } are required as the strings of text in each row are different sizes.

labelP1={'Fitted Equation:'; FitEqnP1}

Now we can use the function text, to add the label to the graph

text(x,y,stringarray)

Here x and y are the co-ordinates on the figure, in this case I will use x=2 and y=60. I’ll use the labelP1 and add two additional input arguments to increase the Font Size to 20.

text(2,60,labelP1,'FontSize',20)

The quadratic equation looks like

\displaystyle y=\text{P2(1)}{{x}^{2}}+\text{P2(2)}x+\text{P2(3)}

We will once again need to put together the equation as a string. For this we use [ ] to make an array and we type the equation as a string enclosed in ' '. Opposed to typing in the coefficients we want to read them from P2. We can use these as input arguments for the num2str function, this function will convert these inputs into strings.

FitEqnP2=[y= num2str(P2(1)) 'x^2 +' num2str(P2(2)'x +' num2str(P2(3))]

Now we can create a label in the form of a cell array using the string variable FitEqnP2. The { } are required as the strings of text in each row are different sizes.

labelP2={'Fitted Equation:'; FitEqnP2}

Now we can use the function text, to add the label to the graph

text(x,y,stringarray)

Here x and y are the co-ordinates on the figure, in this case I will use x=2 and y=60. I’ll use the labelP2 and add two additional input arguments to increase the Font Size to 20.

text(2,60,labelP2,'FontSize',20)

The cubic equation looks like

\displaystyle y=\text{P3(1)}{{x}^{3}}\text{+P3(2)}{{x}^{2}}+\text{P3(3)}x+\text{P3(4)}

We will once again need to put together the equation as a string. For this we use [ ] to make an array and we type the equation as a string enclosed in ' '. Opposed to typing in the coefficients we want to read them from P3. We can use these as input arguments for the num2str function, this function will convert these inputs into strings.

FitEqnP3=[y= num2str(P3(1)) x^3 + num2str(P3(2)) 'x^2 +' num2str(P3(3)'x +' num2str(P3(4))]

Now we can create a label in the form of a cell array using the string variable FitEqnP2. The { } are required as the strings of text in each row are different sizes.

labelP3={'Fitted Equation:'; FitEqnP3}

Now we can use the function text, to add the label to the graph

text(x,y,stringarray)

Here x and y are the co-ordinates on the figure, in this case I will use x=2 and y=60. I’ll use the labelP2 and add two additional input arguments to increase the Font Size to 20.

text(2,60,labelP2,'FontSize',20)

When ran this gives the subplot with the equations


x=[0:2:20]';
y=0.1*x.^2+1.2*x+5;
rawdata(:,1)=x;
rawdata(:,2)=y;

[P1]=polyfit(rawdata(:,1),rawdata(:,2),1);
fitteddata(:,1)=[1:2:19]';
fitteddata(:,2)=polyval(P1,fitteddata(:,1))

[P2]=polyfit(rawdata(:,1),rawdata(:,2),2);
fitteddata(:,1)=[1:2:19]';
fitteddata(:,3)=polyval(P2,fitteddata(:,1))

[P3]=polyfit(rawdata(:,1),rawdata(:,2),3);
fitteddata(:,1)=[1:2:19]';
fitteddata(:,4)=polyval(P3,fitteddata(:,1))

figure(1)
subplot(3,1,1)
plot(rawdata(:,1),rawdata(:,2),'o','MarkerSize',8,'MarkerEdge','b','MarkerFace','b');
grid on;
grid minor;
set(gca,'FontSize',20);
xlabel('x');
ylabel('y');
hold on
plot(fitteddata(:,1),fitteddata(:,2),'r:','LineWidth',5);
FitEqnP1=['y=' num2str(P1(1)) 'x +' num2str(P1(2))]
labelP1={'Fitted Equation:'; FitEqnP1}
text(2,60,labelP1,'FontSize',20)
hold off

subplot(3,1,2)
plot(rawdata(:,1),rawdata(:,2),'o','MarkerSize',8,'MarkerEdge','b','MarkerFace','b');
grid on;
grid minor;
set(gca,'FontSize',20);
xlabel('x');
ylabel('y');
hold on
plot(fitteddata(:,1),fitteddata(:,3),'r:','LineWidth',5);
FitEqnP2=['y=' num2str(P2(1)) 'x^2 +' num2str(P2(2)) 'x +' num2str(P2(3))]
labelP2={'Fitted Equation:'; FitEqnP2}
text(2,60,labelP2,'FontSize',20)
hold off

subplot(3,1,3)
plot(rawdata(:,1),rawdata(:,2),'o','MarkerSize',8,'MarkerEdge','b','MarkerFace','b');
grid on;
grid minor;
set(gca,'FontSize',20);
xlabel('x');
ylabel('y');
hold on
plot(fitteddata(:,1),fitteddata(:,4),'r:','LineWidth',5);
FitEqnP3=['y=' num2str(P3(1)) 'x^3 +' num2str(P3(2)) 'x^2 +' num2str(P3(3)) 'x +' num2str(P3(4))]
labelP3={'Fitted Equation:'; FitEqnP3}
text(2,60,labelP2,'FontSize',20)
hold off