Блог

Я хочу создать 1000 выборок из 200 двумерных нормально распределенных векторов

 set.seed(42)  # for sake of reproducibility
mu <- c(1, 1)
S <- matrix(c(0.56, 0.4,
              0.4, 1), nrow=2, ncol=2, byrow=TRUE)
bivn <- mvrnorm(200, mu=mu, Sigma=S)
 

чтобы я мог запускать регрессии OLS для каждой выборки и, следовательно, получать оценки 1000. Я попробовал это

 library(MASS)
bivn_1000 <- replicate(1000, mvrnorm(200, mu=mu, Sigma=S), simplify=FALSE)
 

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

Я был бы признателен за помощь, чтобы узнать, как выполнить эти 1000 регрессий, а затем извлечь коэффициенты.

Мы могли бы написать пользовательскую функцию регрессии.

 regFun1 <- function(x) summary(lm(x[, 1] ~ x[, 2]))
 

с помощью которого мы можем перебирать данные lapply :

 l1 <- lapply(bivn_1000, regFun1)
 

Коэффициенты сохраняются внутри списка и могут быть извлечены следующим образом:

 l1[[1]]$coefficients  # for the first regression
#              Estimate Std. Error   t value     Pr(>|t|)
# (Intercept) 0.5554601 0.06082924  9.131466 7.969277e-17
# x[, 2]      0.4797568 0.04255711 11.273246 4.322184e-23
 

Редактировать:

Если нам нужны только оценки без статистики, мы соответствующим образом корректируем вывод функции.

 regFun2 <- function(x) summary(lm(x[, 1] ~ x[, 2]))$coef[, 1]
 

Поскольку нам может понадобиться вывод в матричном виде, мы используем sapply next.

 m2 <- t(sapply(bivn_1000, regFun2))

head(m2)
#      (Intercept)    x[, 2]
# [1,]   0.6315558 0.4389721
# [2,]   0.5514555 0.4840933
# [3,]   0.6782464 0.3250800
# [4,]   0.6350999 0.3848747
# [5,]   0.5899311 0.3645237
# [6,]   0.6263678 0.3825725
 

где

 dim(m2)
# [1] 1000    2
 

гарантирует нам, что у нас есть наши 1000 оценок.