Блог

Я пытаюсь получить сводную статистику ( summary() ) для линейной модели (ниже), которая использует 1000 перестановок исходного набора данных для создания 1000 случайных наборов данных (большая матрица).

 random_model <- rep(NA,1000)

for (i in c(1:1000)) {
  random_data <- final_data
  random_data$weighted_degree <- rowSums(node.perm_1000[i,,],na.rm=T)
  random_model[i] <- coef(lm(weighted_degree ~ age   sex   age*sex, data=random_data))
}
  

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

Попробуйте использовать tidy() из broom пакета. Он возвращает ожидаемые значения, подобные этому (пример):

 # A tibble: 2 x 5
  term             estimate std.error statistic  p.value
  <chr>               <dbl>     <dbl>     <dbl>    <dbl>
1 (Intercept)         6.53      0.479     13.6  6.47e-28
2 iris$Sepal.Width   -0.223     0.155     -1.44 1.52e- 1
  

В вашем случае предыдущий вывод будет сохранен для каждого элемента в списке цикла в соответствии с вашим определением:

 library(broom)
#Data
random_model <- rep(NA,1000)
#Loop
for (i in c(1:1000)) {
  random_data <- final_data
  random_data$weighted_degree <- rowSums(node.perm_1000[i,,],na.rm=T)
  random_model[i] <- broom::tidy(lm(weighted_degree ~ age   sex   age*sex, data=random_data))
}
  

Вы должны сохранить интересующие результаты (оцененные коэффициенты и t-значения) в списке.

Вот воспроизводимый пример, использующий 10 репликаций в mtcars наборе данных, который отбирается со скоростью 50% для каждой репликации.

Интересующие результаты извлекаются с использованием $coefficients атрибута summary() вывода lm объекта.

 # The data
data(mtcars)

# Define sample size of each replication
N <- nrow(mtcars)
sample_size <- floor(N/2)

# Number of replications (model fits) and initialization of the list to store the results 
set.seed(1717)
replications <- 10
random_model <- vector( "list", length=replications )
for (i in seq_along(random_model)) {
  shuffle = sample(N, sample_size)
  mtcars_shuffle = mtcars[shuffle, ]
  random_model[[i]] <- summary(lm(mpg ~ cyl   disp   cyl*disp, data=mtcars_shuffle))$coefficients
}
  

Например, модель, подходящая для повторений 1 и 10, является:

 > random_model[[1]]
               Estimate  Std. Error   t value     Pr(>|t|)
(Intercept) 48.26285335 8.219065181  5.872061 7.573836e-05
cyl         -3.33999161 1.366231326 -2.444675 3.089262e-02
disp        -0.12941685 0.063269362 -2.045490 6.337414e-02
cyl:disp     0.01394436 0.007877833  1.770076 1.020931e-01

> random_model[[10]]
               Estimate  Std. Error   t value     Pr(>|t|)
(Intercept) 54.27312267 7.662593317  7.082866 1.277746e-05
cyl         -4.40545653 1.586392001 -2.777029 1.674235e-02
disp        -0.15330770 0.047932153 -3.198431 7.654790e-03
cyl:disp     0.01792561 0.006707396  2.672514 2.031615e-02