n = 4;
f[t_] := E^(-5*t);
x[0] = x[1] = x[2] = 0;
x[3] = x[4] = 1;
DividedDifferences[f_, from_, to_] := If[from == to, f[x[from]],

   If[x[from] ==
     x[to], (D[f[t], {t, to - from}] /. t -> x[from])/(to - from)!,
    													(DividedDifferences[f, from + 1, to] -
       DividedDifferences[f, from, to - 1])/(x[to] - x[from])]];
myProduct[t_, to_] := Product[t - x[k], {k, 0, to - 1}];
myPolynomial[f_] :=
  Sum[DividedDifferences[f, 0, k]*myProduct[x, k], {k, 0, n}];
p = N[Expand[myPolynomial[f]]]
Plot[Abs[f[x] - p], {x, 0, 1}, PlotRange -> All]

+++++++++++++++++++++++++++++++++++++++++++++++++++++++

f[t_] := E^(-5*t);
InterpolatingPolynomial[{{0, {1, -5, 25}}, {1, {E^(-5), -5 E^(-5)}}},
  x];
N[Expand[%]]