format long x=[0.6 0.7 0.8 0.9]; e=[cos(0.6) cos(0.7) cos(0.8) cos(0.9)]; y=x-e; n=size(x)(2)-1; Q=zeros(n+1); Q(:,1)=x'; z=0; for i=1:n for j=1:i Q(i+1,j+1)=(((z-y(i+1-j))*Q(i+1,j))-((z-y(i+1))*Q(i,j)))/(y(i+1)-y(i+1-j)); end end Q(n+1,n+1)