MATLAB Answers

Finding the area under a polynomial equation in MATLAB

20 views (last 30 days)
mr mo
mr mo on 30 Apr 2018
Commented: mr mo on 30 Apr 2018
Hi. I have two vectors like these.
X=[0 0.0300 0.0333 0.0469 0.0575 0.0649 0.0960 0.1262 0.1610 0.1717 0.2017 0.2411 0.2711 0.3000 0.3035 0.3368 0.3711 0.4068 0.4411 0.4758 0.5253 0.5553 0.5853 0.6000]
Y=[0 126.2766 140.0068 190.5906 211.6557 219.3601 236.8709 253.2271 266.4663 269.3085 273.6913 279.3095 283.5479 286.1191 286.2725 282.7424 277.6845 272.1010 266.4556 260.1458 248.5416 241.1381 233.2857 230.3781]
I use the polyfit command to fit a polynomial curve of degree 6 to my data-set.
B=polyfit(X,Y,6)
Now my question is how can I calculate the area under this polynomial curve?
Thanks a lot.

  0 Comments

Sign in to comment.

Accepted Answer

John D'Errico
John D'Errico on 30 Apr 2018
Edited: John D'Errico on 30 Apr 2018
How do you doit? Basic calculus?
X=[0 0.0300 0.0333 0.0469 0.0575 0.0649 0.0960 0.1262 0.1610 0.1717 0.2017 0.2411 0.2711 0.3000 0.3035 0.3368 0.3711 0.4068 0.4411 0.4758 0.5253 0.5553 0.5853 0.6000];
Y=[0 126.2766 140.0068 190.5906 211.6557 219.3601 236.8709 253.2271 266.4663 269.3085 273.6913 279.3095 283.5479 286.1191 286.2725 282.7424 277.6845 272.1010 266.4556 260.1458 248.5416 241.1381 233.2857 230.3781];
xinterp = linspace(min(X),max(X),100);
P = polyfit(X,Y,6);
plot(X,Y,'ro',xinterp,polyval(P,xinterp),'b-')
As models go, I'd call it your basic piece of crapola. And DON'T use a higher order polynomial, thinking it will do better. Just use a spline instead, if you ABSOLUTELY insist on a model. Actually, pchip (shape preserving interpolant) would be my choice here.
If I had to guess, that looks like essentially a piecewise linear curve, with 5 linear segments. But why do you need to integrate a poorly fit polynomial at all?
The area under the curve is...
trapz(X,Y)
ans =
149.16
In fact, this is arguably a better estimator of the area than the integral of that crappy polynomial model.
What do I get for the pchip integral? My SLM toolbox has a nice little utility to integrate a spline. But it is easily enough done without that tool.
pp = pchip(X,Y);
integral(@(x) ppval(pp,x),min(X),max(X))
ans =
149.28
That is not too far off from what trapz gave. Bother are about as good as you can get.
Can we integrate the polynomial model? Again, of course. As I said, basic calc suffices. (Or, you could use polyint.)
P = polyfit(X,Y,6);
Pint = [P./(7:-1:1),0];
diff(polyval(Pint,[min(X),max(X)]))
ans =
149.78
Or if you are feeling lazy,
integral(@(x) polyval(P,x),min(X),max(X))
ans =
149.78
In any event, I'd use what pchip or what trapz gave here as better estimators.

More Answers (1)

James Tursa
James Tursa on 30 Apr 2018
Just integrate the polynomial and evaluate at the endpoints. E.g.,
IB = polyint(B);
area_under_poly = polyval(IB,X(end)) - polyval(IB,X(1));

  3 Comments

mr mo
mr mo on 30 Apr 2018
Thanks but if I have these vectors, the calculated area must be 102.5 but this code return 112.5
x=[0 10 15]
y=[0 10 11]
B=polyfit(x,y,2)
James Tursa
James Tursa on 30 Apr 2018
By hand (using the fact that x(1)=0):
>> (1/3)*B(1)*x(end)^3 + (1/2)*B(2)*x(end)^2 + B(3)*x(end)
ans =
112.5000
Why do you think the answer must be 102.5?
mr mo
mr mo on 30 Apr 2018
Because it is obvious and can be seen by using the plot command
plot(x,y,'-rs')
grid on
if the vectors are like these
x=[0 10 15]
y=[0 10 15]
then the area will be equal to 112.5

Sign in to comment.