Hypothesis testing with non-parametric bootstrap on beta parameter of linear model [closed]

I'd like to test the value of a $$\beta$$ parameter estimated on a linear model. To be specific I want to test $$H_0: \beta_1 = -3$$ but I'm not sure if the procedure I used in R is correct.

This is the code I used:

attach(mtcars)
require(boot)
set.seed(5)

# Define the linear model
fit = lm(mpg ~ hp)
sum = summary(fit)
sum

# Define the null hypothesis for beta1
H0 = 1


Then, I calculate the test statistic with

$$T_{oss} = \frac{\hat{\beta_1} -\beta_1}{\sqrt{\hat{V}(\hat{\beta_1})}}$$

toss = (coefficients(fit)[2]-H0)/sum$coefficients[4]; toss pval = 2*(1-pt(toss, length(hp)-2)); pval # bootstrap procedure e = rstandard(fit) fitted = fit$$fitted.values toss.boot <- function(x, idx){ y = fitted + e[idx] fit = lm(y~hp) sum = summary(fit) toss1 = (coefficients(fit)[2]-H0)/sum$$coefficients[4] return (toss1) } out.boot = boot(mtcars, statistic = toss.boot, R=1000) # pvalue as proportion of values more extreme than the one observed pval = sum(out.boot$$t>toss)/length(out.boot$$t) pval  I'm pretty sure there is something wrong but I would like to get your feedback and hopefully suggestions on other/better approaches. • Can you give a better example for us to consider? The cyl variable in the mtcars variable is in fact a categorical variable with 3 levels, so it wouldn't make sense to test a null hpothesis like the one you stated for it. Why not choose a numerical variable like hp (horse power) and then re-do your example for it? Also, it is odd to state a null hypothesis the way you did - by convention, the null would be stated as Ho: beta1 = 3, etc. Commented May 8, 2018 at 18:19 • You're right. I apologize if it's not well written, hopefully it's better now, although far from perfect! Commented May 8, 2018 at 18:28 1 Answer Thanks for your edits. There are still some issues with your code, as the cyl variable continues to appear throughout the code, even though you replaced it with the hp variable. I suspect the reason your code doesn't produce the intended results is because you are not simulating data according to the null hypothesis for the computation of your bootstrapped p-value. If you look at the post on http://dwoll.de/rexrepos/posts/resamplingBootALM.html, you will find an example of such simulation in the context of performing an ANOVA F-test rather than a t-test. However, the principle is similar and should look more or less like this: # Test null hypothesis that coefficient of hp is equal to 1 # using the original data modelBase <- lm(I(mpg-1*hp) ~ hp, data=mtcars) tBase <- summary(modelBase)$coefficients["hp","t value"]
tBase  ## test statistic
pBase <- summary(modelBase)$coefficients["hp","Pr(>|t|)"] pBase ## p-value fit0 <- lm(I(mpg-1*hp) ~ 1, data=mtcars) ## fit null model E <- residuals(fit0) ## residuals from null model Er <- E / sqrt(1-hatvalues(fit0)) ## rescaled residuals from null model Yhat <- fitted(fit0) ## fitted values from null model ## function for performing model-based resampling gett <- function(dat, idx) { Ystar <- Yhat + Er[idx] model <- lm(Ystar ~ hp, data=dat) summary(model)$coefficients["hp","t value"]
}

library(boot)
nR <- 999
bootres <- boot(mtcars, statistic=gett, R=nR)
bootres

bootres$t hist(bootres$t)
abline(v=tBase, col="red", lty=2)

## Be careful with the computation of the bootstrap p-value, as it
## needs to reflect what alternative hypothesis you are testing:
## Ha: beta_hp > 1   OR
## Ha: beta_hp < 1   OR
## Ha: beta_hp != 1
##
## pvalue as proportion of values more extreme than the one observed
## pval <- sum(bootres$t > tBase)/length(bootres$t)
## pval


In the above, I used the fact that it is possible to test a null hypothesis like "coefficient of hp is equal to 1" by fitting the model:

lm(I(mpg-1*hp) ~ hp, data=mtcars)


and then looking at the t-test corresponding to the hp predictor variable. (See here for details: Test model coefficient (regression slope) against some value.)