R：如何计算数据框中列组合的多个线性模型？

Question

好吧，这个问题可能有点复杂，所以让我尽力解释我的场景。

我有一个数据框，类似于下面设置的示例，其中有一个预测列 (

xpred

) 和三个响应列 (

y1:y3

)。所有这些都有一个关联的对数转换列，因为我希望使用对数转换数据运行

lm()

的线性模型。

suppressWarnings(library(tidyverse))

# Create example data. 
set.seed(10)
df <- data.frame(
  subject = rep(paste("Subject", LETTERS[1:5]), each = 10),
  xpred = rep(1:10, 5),
  y1 = sort(runif(10, min = 130, max = 220), decreasing = TRUE),
  y2 = sort(runif(10, min = 10, max = 90), decreasing = TRUE),
  y3 = sort(runif(10, min = 2, max = 5), decreasing = TRUE)
  )

# Log transfrom pred_x:y3 columns.
df <- df %>%
  group_by(subject) %>%
  mutate(across(.cols = xpred:y3, .fns = ~log(.x), .names = "{col}_log")) %>%
  ungroup()

head(df)
#> # A tibble: 6 × 9
#>   subject   xpred    y1    y2    y3 xpred_log y1_log y2_log y3_log
#>   <chr>     <int> <dbl> <dbl> <dbl>     <dbl>  <dbl>  <dbl>  <dbl>
#> 1 Subject A     1  192.  76.9  4.59     0       5.26   4.34   1.52
#> 2 Subject A     2  185.  62.1  4.51     0.693   5.22   4.13   1.51
#> 3 Subject A     3  176.  57.7  4.33     1.10    5.17   4.05   1.46
#> 4 Subject A     4  169.  55.4  4.31     1.39    5.13   4.01   1.46
#> 5 Subject A     5  168.  44.3  4.12     1.61    5.13   3.79   1.42
#> 6 Subject A     6  158.  41.9  3.85     1.79    5.06   3.74   1.35

^{创建于 2024 年 10 月 25 日，使用 reprex v2.1.1}

现在，我可以轻松地单独运行模型（下面的示例），但是当我想在

subject

列中的每个独特主题的以下列组合上运行模型时，这会变得有点重复：

 y1_log ~ xpred_log

、

y2_log ~ xpred_log

和

y3_log ~ xpred_log

。

如果我单独做这件事，这就是我的代码对于一个主题的样子。

# Subset data for one subject.
df_subset <- df %>% filter(subject == "Subject A")

# Run model.
model <- lm(df_subset$y1_log ~ df_subset$xpred_log)

# Extract coefficients for prediction calculation.
c <- unname(exp(model$coefficients["(Intercept)"]))
n <- unname(model$coefficients["df_subset$xpred_log"])

正如您在我的数据中注意到的那样，

xpred

的范围从 1 到 10，因为

y1:y3

之前已针对这些值计算过。我想要为这个特定项目做的是完成

xpred

，使其范围从 1 到 30，并根据公式

y1:y3

预测

c*x^n

的值 11 到 30，其中

x

是

xpred

。例如，如果我想预测

y1

的值为 11，则它将是

c*11^n

。

希望你还在我身边！

所以，本质上，我想要一种干净/有效的方法：

1. 对每个主题的数据框中的上述对数转换列组合执行

   lm()

   。我假设这将依赖于某种形式的

   group_by()

   调用

   dplyr

   。

2. 使用系数（保存为变量

   c

   和

   n

   ）根据

   y1:y3

   值 11 到 30（已经知道 1 到 10）来预测

   xpred

   。

希望有一个解决方案，其中每个主题包含

xpred

值 1 到 30 的最终输出存储在整洁的数据框中以供以后分析。

感谢对此提供的任何帮助，如果有任何需要澄清的地方，请告诉我，因为这是我的第一篇文章。

Answer 1

这是一种基于嵌套 tibbles 和按行操作的方法，是对 dplyr 的按行操作小插图的建模示例的修改。

library(tidyverse)

# Pivot longer and then log transform;
# create nested tibble with data tibble for every subject-response combination (grouped rowwsie);
# fit models and extract coefficients from all combinations
df_coef <- 
  df |> 
  pivot_longer(y1:y3, names_to = "response", values_to = "response_value") |> 
  mutate(across(.cols = c(xpred, response_value), .fns = ~log(.x), .names = "{col}_log")) |> 
  nest_by(subject, response) |> 
  # creating unnamed frame in a mutate call adds multiple columns (c,n)
  mutate(
    lm(response_value_log ~ xpred_log, data = data) |> 
      coef() |> 
      bind_rows() |> 
      rename(`c` = "(Intercept)", n = "xpred_log")
    )
# long nested tibble:
df_coef
#> # A tibble: 15 × 5
#> # Rowwise:  subject, response
#>    subject   response               data     c      n
#>    <chr>     <chr>    <list<tibble[,4]>> <dbl>  <dbl>
#>  1 Subject A y1                 [10 × 4]  5.30 -0.132
#>  2 Subject A y2                 [10 × 4]  4.61 -0.621
#>  3 Subject A y3                 [10 × 4]  1.65 -0.226
#>  4 Subject B y1                 [10 × 4]  5.30 -0.132
#>  5 Subject B y2                 [10 × 4]  4.61 -0.621
#>  6 Subject B y3                 [10 × 4]  1.65 -0.226
#>  7 Subject C y1                 [10 × 4]  5.30 -0.132
#>  8 Subject C y2                 [10 × 4]  4.61 -0.621
#>  9 Subject C y3                 [10 × 4]  1.65 -0.226
#> 10 Subject D y1                 [10 × 4]  5.30 -0.132
#> 11 Subject D y2                 [10 × 4]  4.61 -0.621
#> 12 Subject D y3                 [10 × 4]  1.65 -0.226
#> 13 Subject E y1                 [10 × 4]  5.30 -0.132
#> 14 Subject E y2                 [10 × 4]  4.61 -0.621
#> 15 Subject E y3                 [10 × 4]  1.65 -0.226

# First data tibble (SubjectA, response y1)
df_coef$data[[1]]
#> # A tibble: 10 × 4
#>    xpred response_value xpred_log response_value_log
#>    <int>          <dbl>     <dbl>              <dbl>
#>  1     1           192.     0                   5.26
#>  2     2           185.     0.693               5.22
#>  3     3           176.     1.10                5.17
#>  4     4           169.     1.39                5.13
#>  5     5           168.     1.61                5.13
#>  6     6           158.     1.79                5.06
#>  7     7           155.     1.95                5.04
#>  8     8           155.     2.08                5.04
#>  9     9           150.     2.20                5.01
#> 10    10           138.     2.30                4.92

# Replace nested data tibbles with calculated values for xpred = 11:30;
# reverse nesting & pivoting
df_coef |> 
  mutate(data = tibble(xpred = 11:30, response_value = c*xpred^n) |> list(), .keep = "unused") |> 
  unnest(data) |> 
  pivot_wider(names_from = response, values_from = response_value)

结果：

#> # A tibble: 100 × 5
#> # Groups:   subject [5]
#>    subject   xpred    y1    y2    y3
#>    <chr>     <int> <dbl> <dbl> <dbl>
#>  1 Subject A    11  3.86 1.04  0.963
#>  2 Subject A    12  3.81 0.986 0.944
#>  3 Subject A    13  3.77 0.938 0.927
#>  4 Subject A    14  3.74 0.896 0.912
#>  5 Subject A    15  3.70 0.858 0.897
#>  6 Subject A    16  3.67 0.824 0.885
#>  7 Subject A    17  3.64 0.794 0.872
#>  8 Subject A    18  3.61 0.766 0.861
#>  9 Subject A    19  3.59 0.741 0.851
#> 10 Subject A    20  3.56 0.718 0.841
#> # ℹ 90 more rows

---

示例数据：

# Create example data. 
set.seed(10)
df <- data.frame(
  subject = rep(paste("Subject", LETTERS[1:5]), each = 10),
  xpred = rep(1:10, 5),
  y1 = sort(runif(10, min = 130, max = 220), decreasing = TRUE),
  y2 = sort(runif(10, min = 10, max = 90), decreasing = TRUE),
  y3 = sort(runif(10, min = 2, max = 5), decreasing = TRUE)
)