How to bootstrap prediction intervals for “customized” regression models in R

Question

Are there functions in R that could help me do the following?

We have a special type of regression which is called Geometric Mean Regression.

We have done some search and found the following:

https://stat.ethz.ch/pipermail/r-help/2011-July/285022.html

The question is: how to do the statistical inference on GMR results?

More specifically, we are looking for the prediction interval:

Lets say we regress y1, y2, ..., yn onto x1, x2, ..., xn:

we would like to know what's the prediction interval for a new data point:

x_new=x1+x2+x3

(i.e. the new data point is the sum of the existing first three data points)

In ordinary linear regression, we could derive prediction interval for an in-sample data point as well as a new data point...

For our x_new=x1+x2+x3, we can derive formulas for the prediction interval.

But for the above customized regression,

how do we obtain the prediction intervals?

---

Are there functions in R that can help us do this?

We are thinking of using bootstrapping, etc. Are there functions in R help us on this?

Thanks a lot!

---

I did some googling and researching... Reading the following article,

http://www.ecd.bnl.gov/pubs/BNL-79819-2008-JA.pdf

It seems that once we estimate the parameters of the bivariate normal distribution,

then we can plug into the formula of conditional distribution of Y|X=x1+x2+x3 ?

http://en.wikipedia.org/wiki/Multivariate_normal_distribution

My question is:

Is it a correct procedure to do the following:

Step 1: estimate the parameters of the bivariate normal distribution; Step 2: plug the estimated parameters into the Y|X=x1+x2+x3 formula and get the 95% quantile of it?

Do I need to repeat Step 2 many times following the bootstrapping procedure?

Or one shot of Step 2 is enough?

I got very much confused...

Any thoughts?

Thanks a lot!