help understanding polyval delta

I am looking to understand the procedure behind the determination of delta which is produced from the polyval function. I am not a stats buff so I am looking for some type of introduction or tutorial, I was unable to find a reference anywhere online (I have spent quite some time on google trying to find a explanation/tutorial/introduction). I can examine the calculation by typing 'edit polyval' but that only goes so far. It appears that delta is determined as follows:

E = V/R;% V is Vandermonde matrix for new X

e = sqrt(1+sum(E.*E,2));

delta = normr/sqrt(df)*e;

It is this procedure that I don't understand. I did find a vague reference to Numerical Recipes in C section 15.6 but no equations in that document follow the methodology in the polyval function.

Apart from the lack of theoretical understanding I have tried 1 piece of code which produce curious results. In the code snippet I create some data with large amount of noise and then extended the prediction into the future. I don't understand why delta increases as x increases. Delta at first decreases then increases in an apparent linear fashion. Why would delta increase linearly, or at all? I understand that delta represents prediction error estimates of the standard deviation at the future observation but the mechanism of how delta is computed is lost.

randn('seed',1);

%create initial data

x=[0,1,2,3,4,5,6]';

y = 2.5 * x + 80 + 50 * randn(length(x),1);

[p s] = polyfit(x, y, 1);

%create extended range

newX = 0:100;

[y_fit delta] = polyval(p, newX, s);

figure;

%plot data in blue

plot(x,y,'b*')

hold on

%add 2 sigma bounds

plot(newX, y_fit + 2 * delta, '--', newX, y_fit - 2 * delta, '--')

legend('Data','2 Sigma Bound')

%plot delta

figure;

plot(newX,delta)

[dpb]

It's just a technique for evaluation, the basic idea is illustrated for a linear model. The var or a predicted value is

var(y-yhat) --> var(y)+var(yhat) --> var(y)+var(mx+b)

Expanding the above for estimated slope and intercept eventually one gets that

var(prediction) --> s^2[1 + 1/m + (x-xbar)^2/Sum(x-xbar)^2]

NB the x-xbar term; that is the one that is responsible for the error estimate to be a minimum at the mean of x and to expand in both directions away from that point.

See Draper & Smith, Applied Regression Analysis or any other text on regression analysis for the theory.