Linear regression prediction interval If the best linear approximation (using least squares) of my data points is the line $y=mx+b$, how can I calculate the approximation error? If I compute standard deviation of differences between observations and predictions $e_i=real(x_i)-(mx_i+b)$, can I later say that a real (but not observed) value $y_r=real(x_0)$ belongs to the interval $[y_p-\sigma, y_p+\sigma]$ ($y_p=mx_0+b$) with probability ~68%, assuming normal distribution?
To clarify:
I made observations regarding a function $f(x)$ by evaluating it a some points $x_i$. I fit these observations to a line $l(x)=mx+b$. For $x_0$ that I did not observe, I'd like to know how big can $f(x_0)-l(x_0)$  be. Using the method above, is it correct to say that $f(x_0) \in [l(x_0)-\sigma, l(x_0)+\sigma]$ with prob. ~68%?
 A: @whuber has pointed you to three good answers, but perhaps I can still write something of value.  Your explicit question, as I understand it, is:  
Given my fitted model, $\hat y_i=\hat mx_i + \hat b$ (notice I added 'hats'), and assuming my residuals are normally distributed, $\mathcal N(0, \hat\sigma^2_e)$, can I predict that an as yet unobserved response, $y_{new}$, with a known predictor value, $x_{new}$, will fall within the interval $(\hat y -\sigma_e, \hat y +\sigma_e)$, with probability 68%? 
Intuitively, the answer seems like it should be 'yes', but the true answer is maybe.  This will be the case when the parameters (i.e., $m, b,$ & $\sigma$) are known and without error.  Since you estimated these parameters, we need to take their uncertainty into account.  
Let's first think about the standard deviation of your residuals.  Because this is estimated from your data, there can be some error in the estimate.  As a result, the distribution you should use to form your prediction interval should be $t_\text{df error}$, not the normal.  However, since the $t$ converges rapidly to the normal, this is less likely to be a problem in practice.  
So, can we just use $\hat y_\text{new}\pm t_{(1-\alpha/2,\ \text{df error})}s$, instead of $\hat y_\text{new}\pm z_{(1-\alpha/2)}s$, and go about our merry way?  Unfortunately, no.  The bigger issue is that there is uncertainty about your estimate of the conditional mean of the response at that location due to the uncertainty in your estimates $\hat m$ & $\hat b$.  Thus, the standard deviation of your predictions needs to incorporate more than just $s_\text{error}$.  Because variances add, the estimated variance of the predictions will be:
$$
s^2_\text{predictions(new)}=s^2_\text{error}+\text{Var}(\hat mx_\text{new}+\hat b)
$$
Notice that the "$x$" is subscripted to represent the specific value for the new observation, and that the "$s^2$" is correspondingly subscripted.  That is, your prediction interval is contingent on the location of the new observation along the $x$ axis.  The standard deviation of your predictions can be more conveniently estimated with the following formula:
$$
s_\text{predictions(new)}=\sqrt{s^2_\text{error}\left(1+\frac{1}{N}+\frac{(x_\text{new}-\bar x)^2}{\sum(x_i-\bar x)^2}\right)}
$$
As an interesting side note, we can infer a few facts about prediction intervals from this equation.  First, prediction intervals will be narrower the more data we had when we built the prediction model (this is because there's less uncertainty in $\hat m$ & $\hat b$).  Second, predictions will be most precise if they are made at the mean of the $x$ values you used to develop your model, as the numerator for the third term will be $0$.  The reason is that under normal circumstances, there is no uncertainty about the estimated slope at the mean of $x$, only some uncertainty about the true vertical position of the regression line.  Thus, some lessons to be learned for building prediction models are: that more data is helpful, not with finding 'significance', but with improving the precision of future predictions; and that you should center your data collection efforts on the interval where you will need to be making predictions in the future (to minimize that numerator), but spread the observations as widely from that center as you can (to maximize that denominator).  
Having calculated the correct value in this manner, we can then use it with the appropriate $t$ distribution as noted above.  
