#
# Epage 
#
# Written by:
# -- 
# John L. Weatherwax                2009-04-21
# 
# email: wax@alum.mit.edu
# 
# Please send comments and especially bug reports to the
# above email address.
# 
#-----

param_bootstrap <- function(x,y,B=999){

    # compute a the parametric fit predictions of y ~ x: 
    #
    m_param         <- lm( y ~ x )
    y_hat_param     <- m_param$fitted.values
    residuals_param <- m_param$residuals
    var_hat         <- summary(m_param)$sigma^2 

    # compute the loess fit predictiosn of y ~ x :
    #
    m_loess     <- loess( y ~ x )
    y_hat_loess <- lm$fitted 

    G = sum( ( y_hat_param - y_hat_loess )^2 )/ var_hat

    n <- length( residuals_param )
    all_boots <- NULL
    for (i in 1:B) {

	e_hat_star = sample( residuals_param, n, replace=TRUE )    # sample the parametric residuals
	y_hat_star = y_hat_param + e_hat_star                      # compute bootstrapped responses

	# do a bootstrapped parametric (linear) fit
	b_m_param <- lm( y_hat_star ~ x ) 
	# do a loess parametric fit
	b_m_loess <- loess( y_hat_star ~ x )

	#b_G <- sum( ( b_m_param$fitted.values - b_m_loess$fitted )^2 )/ (summary(b_m_param)$sigma^2) 
	b_G <- sum( ( b_m_param$fitted.values - b_m_loess$fitted )^2 )/ var_hat 

	all_boots <- rbind(all_boots,b_G)

    }

    # using these bootstrapped samples of G to compute the significance level of the G statistic: 
    p <- mean( all_boots > G ) 
    #print(p)

    p
}