Hi, Ive been trying to write a program for bootstrapping residuals in R but without much success.

A lecturer wants to predict the performance of students in an end-of-year physics exam, y. The lecturer has the students results from a mid-term physics exam, x, and a mid-term biology exam, z.

The lecturer proposes the following linear model for the end-of-year exam result

yi = α + βxi + γzi + qi, where q is the error.

Y is a matrix of the data and we have y=first column of the data and X=second 2 columns(the x & z data)

Now I need to write a program for obtaining bootstrap estimates, i have:

x=scan(data)

Y=matrix(x,ncol=3,byrow=T)

y=Y[,1]

X=Y[,2:3]

ls=lsfit(X,y)

beta=ls$coef

yest=beta[1]+beta[2]*X[,1]+beta[3]*X[,2]

res=y-yest

boot=function(X,res,beta,b)

{

n=24

output=matrix(0,ncol=2,nrow=b)

for(i in 1:b)

{

error=sample(res,n,replace=T)

ystar=beta[1]+beta[2]*X[,1]+beta[3]*X[,2]+error

ls=lsfit(X,ystar)

output[i,]=ls$coef

}

output

}

I think the first 8 lines are right but my function might be wrong?

Any help?