您可以使用gtsummary库来完成此操作。
首先，我们创建一个虚拟数据集进行演示。

library(tidyverse)
library(gtsummary)

# Create dummy dataset (50 cases).

set.seed(123)
df <- data.frame(
  Sex = sample(c("Female", "Male"), 50, replace = TRUE),
  Age = sample(18:80, 50, replace = TRUE),
  Race = sample(c("White", "Black"), 50, replace = TRUE),
  Income = sample(c("High", "Low"), 50, replace = TRUE),
  Treatment1 = sample(0:1, 50, replace = TRUE),
  Treatment2 = sample(0:1, 50, replace = TRUE)
)

然后我们可以使用gtsummary中的tbl_strata函数来创建一个包含所需级别的汇总表。

table <- df %>%
  tbl_strata(
  strata = c("Treatment1", "Treatment2"),
  ~ .x %>%
    tbl_summary(
      type = list(
        Age ~ "continuous",
        Sex ~ "categorical",
        Race ~ "categorical",
        Income ~ "categorical"
      ),
      statistic = list(
          Age ~ "{mean} ({sd})",
          all_categorical() ~ "{n}"
        )
    )
  )

我发现这些页面很有用：Epidemiologist R Handbook (Descriptive Tables)，Stratified gtsummary Tables .
另一种方法是使用table1 package。这个例子的代码是：

df2 <- df %>%
  mutate(Sex = factor(Sex, levels = c("Male", "Female")),
         Race = factor(Race, levels = c("White", "Black")),
         Income = factor(Income, levels = c("High", "Low")),
         Treatment1 = factor(Treatment1, levels = c(0, 1)),
         Treatment2 = factor(Treatment2, levels = c(0, 1)))

units(df2$Age) <- "years"

table1(~ factor(Sex) + Age + factor(Race) + factor(Income) | Treatment1, data = df2, overall = FALSE)

table1还支持使用x*y格式的两个变量的分层，例如table1(~ factor(Sex) + Age + factor(Race) + factor(Income) | Treatment1*Treatment2)。