Beginner trouble with bootstrapping in R

Question

I am doing a "Statistical learning in R" online course which requires me to do exercises in R and right now, i'm stuck at a part that requires me to use the bootstrap.

I'm not sure whether it is allowed to ask about homework problems on here, but i have been looking around different forums and tutorials for a while and i can't really put my head around this.

So the original question is

Now, use the (standard) bootstrap to estimate s.e.($ \hat{\beta_1} $) . To within 10%, what do you get?

A dataset is given and in previous questions i was required to estimate the standard error of $\beta_1$ from a linear model on the data. The set Xy consists of three variables with 1000 entries each (X1, X2 and Y).

In my attempt of a solution, i have tried to stay close to the original code used in the Video lesson.

library(ISLR) 
library(boot)
alpha.fn=function(data,index){
X=data$X[index]
Y=data$Y[index]
return((var(Y)-cov(X,Y))/(var(X)+var(Y)-2*cov(X,Y)))}
#"renaming my variables to fit the alpha function"
Xy$X=Xy$X1;
Xy$Y=Xy$y;
set.seed(1)
alpha.fn(Xy,sample(100,100,replace=T))
boot(Xy,alpha.fn,R=1000)

Now i don't know whether this is the "standard" bootstrap i intend to use or how to pass the function that i want to estimate the std. error to within 10%. does it mean i just want it to 0.1 accuracy?

Answer 1

First of all, Tung Nguyen is correct; you have to estimate the full model and extract the relevant estimate (for $\beta_1$), not just run a subset of the model using only $X_1$ on the r.h.s.

Following up on Michael M's suggestion, which is also how I would approach it, we get the following code, with some of R's printout removed for brevity:

# Run basic regression, calculate coefs and s.e.s
summary(lm(y~x1+x2, data=df))

* blah blah blah *

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  1.03242    0.03213   32.14   <2e-16 ***
x1           0.97044    0.03214   30.19   <2e-16 ***
x2           1.01010    0.03291   30.69   <2e-16 ***
---

# Bootstrap (1000 replications)
beta_hat <- rep(0,1000)
for (i in 1:length(beta_hat)) {
    local_df <- df[sample(1000, replace=TRUE),]
    beta_hat[i] <- coef(lm(y~x1+x2, data=local_df))["x1"]
}
sd(beta_hat)
[1] 0.03258774

Our bootstrap estimate of the standard error is, as one would expect given the large bootstrap replication count and Normally-distributed errors, very close to the calculation made by lm.

You could shorten the R code a little by rewriting the two statements inside the loop as one"

beta_hat[i] <- coef(lm(y~x1+x2, data=df[sample(1000, replace=TRUE),]))["x1"]

but this is likely less clear for expository purposes!