Запуск нескольких линейных регрессий по нескольким столбцам фрейма данных в R

Question

У меня есть структурированный набор данных: введите описание изображения здесь

Я хотел бы запустить модели линейной регрессии и ANOVA, используя V1, V2 ... и т. Д.как независимые переменные и данный столбец как зависимая переменная в каждом случае (т. е. lm (V1 ~ g), lm (V2 ~ g) и т. д.).Это было бы просто, за исключением того, что эти линейные регрессии должны быть сгруппированы по уровню в столбце пары, так что, например, мой вывод содержит lm (V1 ~ g) для всех строк с парой 1.1 и lm (V1 ~ g) для всехпары 1.201 и т. д.

Я пробовал несколько подходов, использующих for loop, lapply и пакет data.table, и ничто не дает мне именно тот результат, который мне нужен.Кто-нибудь может дать мне какое-либо представление о наилучшем способе решения этой проблемы?

Редактировать: Мой полный набор данных содержит 7056 различных пар в столбце пар и столбцах 100 В (V1 ... V100).Моя последняя попытка решить эту проблему:

df$pair <- as.factor(df$pair)
out <- list()
for (i in 3:ncol(df)){
    out[[i]] <- lapply(levels(df$pair), function(x) {
    data.frame(df=x, g = coef(summary(lm(df[,i]~ df$g, data=df[df$pair==x,])),row.names=NULL))})
    }

Answer 1

Давайте получим здесь некоторую tidyverse мощность, наряду с broom, и воздержимся от всех этих циклов ...

Сначала я сделаю фиктивную таблицу:

df <- data.frame(
  g = runif(50), 
  pair = sample(x = c("A", "B", "C"), size = 50, replace = TRUE), 
  V1 = runif(50), 
  V2 = runif(50), 
  V3 = runif(50), 
  V4 = runif(50), 
  V5 = runif(50),
  stringsAsFactors = FALSE
)

Примерно так выглядит ваша структура данных.Теперь о сути кода:

library(tidyverse)
library(broom)

df %>% 
  as_tibble %>% 
  gather(key = "column", value = "value", V1:V5) %>%       # first set the data in long format
  nest(g, value) %>%                                       # now nest the dependent and independent factors
  mutate(model = map(data, ~lm(g ~ value, data = .))) %>%  # fit the model using purrr
  mutate(tidy_model = map(model, tidy)) %>%                # clean the model output with broom
  select(-data, -model) %>%                                # remove the "untidy" parts
  unnest()                                                 # get it back in a recognizable data frame

Что дает нам следующее:

# A tibble: 30 x 7
   pair  column term        estimate std.error statistic  p.value
   <chr> <chr>  <chr>          <dbl>     <dbl>     <dbl>    <dbl>
 1 C     V1     (Intercept)  0.470       0.142    3.31   0.00561 
 2 C     V1     value        0.125       0.265    0.472  0.645   
 3 B     V1     (Intercept)  0.489       0.142    3.45   0.00359 
 4 B     V1     value       -0.0438      0.289   -0.151  0.882   
 5 A     V1     (Intercept)  0.515       0.111    4.63   0.000279
 6 A     V1     value       -0.00569     0.249   -0.0229 0.982   
 7 C     V2     (Intercept)  0.367       0.147    2.50   0.0265  
 8 C     V2     value        0.377       0.300    1.26   0.231   
 9 B     V2     (Intercept)  0.462       0.179    2.59   0.0206  
10 B     V2     value        0.0175      0.322    0.0545 0.957   
# … with 20 more rows

Да, это красавица!Обратите внимание, что я использовал lm(g ~ value) вместо lm(value ~ g), так как именно на это ссылается ваше текстовое описание.

Answer 2

Использование пакета tidyverse для фильтрации вашего фрейма данных:

library(tidyverse)

lm(V1~g, data=filter(yourData, pair==1.1))
lm(V2~g, data=filter(yourData, pair==1.201))

Это гарантирует, что вы удаляете строки, которые не содержат желаемого значения pair для каждой модели регрессии.Возможно, вы могли бы создать цикл, чтобы сделать это, но я думаю, что проще вручную продолжать фильтровать значения pair.Если вы действительно хотите использовать цикл, вот довольно простой способ сделать это:

for (i in levels(yourData$pair)) {
  if (i==1.1) {
    mod1 <- lm(V1~g, data=filter(yourData, pair==i))
  }

  if (i==1.201) {
    mod2 <- lm(V2~g, data=filter(yourData, pair==i))
  }
}

Но это все еще вручную с циклическим прохождением уровней pair.Мне бы пришлось увидеть весь ваш набор данных для автоматизации процесса зацикливания.

Кроме того, если в столбце g содержатся ваши зависимые значения, вызов должен быть lm(g~V1), lm(g~V2) и т. Д. Это не должен быть lm(V1~g).