f20xvals = -2:0.04:2;
f20fvals = sin(f20xvals) + 2*exp(-30*(f20xvals.^2));
fsd = std(f20fvals);
epssd = fsd/50;
numobs = size(f20xvals, 2);
eps = normrnd(0, epssd, 1, numobs);
f20yvals = f20fvals + eps;

knots = [-1.8, -0.4, -0.2, 0, 0.2, 0.4, 1.8];
designmat = zeros(numobs, size(knots, 2));
designmat(:, 1) = 1;
designmat(:, 2) = f20xvals;
designmat(:, 3) = f20xvals.^2;
designmat(:, 4) = f20xvals.^3;
for i = 1:size(knots, 2)
    designmat(:, 4+i) = (subplus(f20xvals-knots(i))).^3;
end

lmspline = regress(f20yvals', designmat);
yhat = designmat*lmspline;

plot(f20xvals, yhat, '-k', 'LineWidth', 2)
hold on;
plot(f20xvals, f20yvals, '.k', 'MarkerSize', 15)
set(gca, 'FontSize', 22, 'XLim', [-2.1, 2.1], 'YLim', [-1.1, 2.1], 'XTick', -2:1:2, 'YTick', -1:0.5:2, 'TickDir', 'out', 'Box', 'off')
xlabel({'';'x'}, 'FontSize', 22)
ylabel({'y';''}, 'FontSize', 22)
set(gcf, 'Position', [200, 100, 1000, 800])