Some time back, I implemented a function that duplicated the functionality of the old Splines` package, but returned an InterpolatingFunction[] instead (original here ):
splinefit[dat_?MatrixQ] :=
Module[{n = Length[dat], del = Differences[dat]},
del = LinearSolve[SparseArray[
{Band[{2, 1}] -> 1,
Band[{1, 1}] -> ArrayPad[ConstantArray[4, n - 2], 1, 2],
Band[{1, 2}] -> 1}],
ListCorrelate[{{3}, {3}}, del, {-{1, 1}, {1, 1}}, {0}]];
BSplineFunction[#1, SplineDegree -> 3, SplineKnots -> #2] & @@
{Join[{dat[[1]]},
Riffle[Most[dat] + Most[del]/3, Rest[dat] - Rest[del]/3],
{dat[[-1]]}],
ArrayPad[Riffle[#, #] &[Range[0, n - 1]], 2, "Fixed"]}]
For your data, then:
data = {{-2, 1}, {-1, 4}, {0, 5}, {1, 3}, {2, 2}};
fun1 = splinefit[data];
ParametricPlot[fun1[t], {t, 0, 4}, AspectRatio -> 1/GoldenRatio]
However, I do not recommend using this old functionality. Instead, you can use a "not-a-knot" spline interpolation like what Michael describes. Here's a compact implementation (original here):
notAKnotSpline[pts_?MatrixQ] := Module[{dy, h, p1, p2, sl, s1, s2, tr},
h = Differences[pts[[All, 1]]]; dy = Differences[pts[[All, 2]]]/h;
s1 = Total[Take[h, 2]]; s2 = Total[Take[h, -2]];
p1 = ({3, 2}.Take[h, 2] h[[2]] dy[[1]] + h[[1]]^2 dy[[2]])/s1;
p2 = (h[[-1]]^2 dy[[-2]] + {2, 3}.Take[h, -2] h[[-2]] dy[[-1]])/s2;
tr = SparseArray[{Band[{2, 1}] -> Append[Rest[h], s2],
Band[{1, 1}] -> Join[{h[[2]]}, ListCorrelate[{2, 2}, h], {h[[-2]]}],
Band[{1, 2}] -> Prepend[Most[h], s1]}];
sl = LinearSolve[tr, Join[{p1},
3 Total[Partition[dy, 2, 1]
Reverse[Partition[h, 2, 1], 2], {2}],
{p2}]];
Interpolation[MapThread[{{#1[[1]]}, #1[[2]], #2} &, {pts, sl}],
InterpolationOrder -> 3, Method -> "Hermite"]]
To try it out:
spl = notAKnotSpline[data];
Plot[spl[x], {x, -2, 2}, Epilog -> {Directive[AbsolutePointSize[6], ColorData[116, 4]], Point[data]}]