function Y=newton(x,y)
syms X;
n=length(x);
b=zeros(n,n);
b(:,1)=y(:);
%finite difference matrix
for j=2:n
for i=1:n-j+1
b(i,j)= (b(i+1,j-1)-b(i,j-1))/(x(i+j-1)-x(i));
end
end
%interpolation polynom
yint=b(1,1);
xt=1;
for j=1:n-1
xt=xt*(X-x(j));
yint=yint+b(1,j+1)*xt;
end
Y=yint
Sorry, it was supposed to come out like this.