Как получить p-значения коэффициентов из регрессии начальной загрузки?

Из Quick-R Роберта Кабакова у меня есть

# Bootstrap 95% CI for regression coefficients 
library(boot)
# function to obtain regression weights 
bs <- function(formula, data, indices) {
  d <- data[indices,] # allows boot to select sample 
  fit <- lm(formula, data=d)
  return(coef(fit)) 
} 
# bootstrapping with 1000 replications 
results <- boot(data=mtcars, statistic=bs, 
     R=1000, formula=mpg~wt+disp)

# view results
results
plot(results, index=1) # intercept 
plot(results, index=2) # wt 
plot(results, index=3) # disp 

# get 95% confidence intervals 
boot.ci(results, type="bca", index=1) # intercept 
boot.ci(results, type="bca", index=2) # wt 
boot.ci(results, type="bca", index=3) # disp

Как я могу получить p-значения коэффициентов регрессии начальной загрузки?$H_0:\, b_j=0$

Еще один вариант, который несколько упрощен, но я думаю, что доставить сообщение без явного использования библиотеки, bootкоторая может ввести некоторых людей в заблуждение синтаксисом, который она использует.

У нас есть линейная модель: ,$y = X \beta + \epsilon$$\quad \epsilon \sim N(0,\sigma^2)$

Ниже приведен параметрический начальный загрузчик для этой линейной модели, это означает, что мы не пересматриваем наши исходные данные, но на самом деле мы генерируем новые данные из нашей подобранной модели. Кроме того, мы предполагаем, что начальное распределение коэффициента регрессии является симметричным и является трансляционно-инвариантным. (Очень грубо говоря, что мы можем перемещать его ось, влияя на его свойства). Идея заключается в том, что флуктуации в обусловлены и поэтому при достаточном количестве выборок они должны обеспечивать хорошее приближение к истинному распределению. из -х. Как и прежде, мы снова тестируем и определяем наши p-значения как$\beta$$\beta$$\epsilon$$\beta$$H_0 : 0 = \beta_j$«вероятность, учитывая нулевую гипотезу для распределения вероятности данных, что результат будет таким же экстремальным, как или более экстремальным, чем наблюдаемый результат» (где наблюдаемые результаты в этом случае - это, которые мы получили для нашей оригинальной модели). Итак, здесь идет:$\beta$

# Sample Size
N           <- 2^12;
# Linear Model to Boostrap          
Model2Boot  <- lm( mpg ~ wt + disp, mtcars)
# Values of the model coefficients
Betas       <- coefficients(Model2Boot)
# Number of coefficents to test against
M           <- length(Betas)
# Matrix of M columns to hold Bootstraping results
BtStrpRes   <- matrix( rep(0,M*N), ncol=M)

for (i in 1:N) {
# Simulate data N times from the model we assume be true
# and save the resulting coefficient in the i-th row of BtStrpRes
BtStrpRes[i,] <-coefficients(lm(unlist(simulate(Model2Boot)) ~wt + disp, mtcars))
}

#Get the p-values for coefficient
P_val1 <-mean( abs(BtStrpRes[,1] - mean(BtStrpRes[,1]) )> abs( Betas[1]))
P_val2 <-mean( abs(BtStrpRes[,2] - mean(BtStrpRes[,2]) )> abs( Betas[2]))
P_val3 <-mean( abs(BtStrpRes[,3] - mean(BtStrpRes[,3]) )> abs( Betas[3]))

#and some parametric bootstrap confidence intervals (2.5%, 97.5%) 
ConfInt1 <- quantile(BtStrpRes[,1], c(.025, 0.975))
ConfInt2 <- quantile(BtStrpRes[,2], c(.025, 0.975))
ConfInt3 <- quantile(BtStrpRes[,3], c(.025, 0.975))

Как уже упоминалось, вся идея заключается в том, что у вас есть загрузочное распределение , приближенное к их истинному. (Очевидно, этот код оптимизирован для скорости, но для удобства чтения. :))$\beta$

Сообщество и @BrianDiggs могут исправить меня, если я ошибаюсь, но я считаю, что вы можете получить p-значение для вашей проблемы следующим образом. Значение p для двустороннего теста определяется как

Таким образом, если вы упорядочите загрузочные коэффициенты по размеру, а затем определите пропорции больше и меньше нуля, минимальная пропорция, умноженная на два, даст вам значение p.

Я обычно использую следующую функцию в такой ситуации:

twosidep<-function(data){
  p1<-sum(data>0)/length(data)
  p2<-sum(data<0)/length(data)
  p<-min(p1,p2)*2
  return(p)
}

Начальный загрузчик можно использовать для вычисления значений, но это потребует существенных изменений в вашем коде. Поскольку я не знаком с RI, я могу лишь дать вам ссылку, в которой вы можете посмотреть, что вам нужно сделать: глава 4 (Davison and Hinkley 1997).$p$

Дэвисон, AC и Хинкли, Д.В. 1997. Методы начальной загрузки и их применение. Кембридж: издательство Кембриджского университета.