这不是一个优雅的解决方案，但在编程方面很简单：

set.seed(581745)
dat <- MASS::mvrnorm(100, rep(0, 7), diag(7))
colnames(dat) <- c("y", paste0("m", 1:3), paste0("x", 1:3))
head(dat)
#>                y          m1         m2         m3          x1          x2
#> [1,]  0.19461669 -0.64146102 -0.7205590  0.7927056 -0.28878053  2.43030745
#> [2,] -1.04997220  1.48515283 -1.2392574  0.3022706 -0.12727050 -0.69668095
#> [3,]  0.31000671 -0.11000903  0.3791038 -1.1362097 -0.39751624 -0.32798805
#> [4,]  0.03023118  0.02059596  0.6318862 -0.7585338 -0.04914222  1.52437138
#> [5,]  0.20720111 -0.28528315 -0.5581268  0.3280841  0.68032932  0.48519911
#> [6,] -1.81659430 -0.11340674  0.3485680  0.3389938  1.03419571  0.09771754
#>               x3
#> [1,] -0.83107742
#> [2,] -1.64153599
#> [3,] -1.17252108
#> [4,]  1.11641882
#> [5,]  0.33163842
#> [6,]  0.08189787

x <- c("x1", "x2", "x3")
m <- c("m1", "m2", "m3")

mods <- paste0(
"# Path c (direct effect)
y ~ c*", x, "\n",
"# Path a", "\n",
m, " ~ a*", x, "\n",
"# Path b
y ~ b*", m, "\n",
"#Indirect effect (a*b)
ab := a*b")

cat(mods[1])
#> # Path c (direct effect)
#> y ~ c*x1
#> # Path a
#> m1 ~ a*x1
#> # Path b
#> y ~ b*m1
#> #Indirect effect (a*b)
#> ab := a*b
cat(mods[2])
#> # Path c (direct effect)
#> y ~ c*x2
#> # Path a
#> m2 ~ a*x2
#> # Path b
#> y ~ b*m2
#> #Indirect effect (a*b)
#> ab := a*b
cat(mods[3])
#> # Path c (direct effect)
#> y ~ c*x3
#> # Path a
#> m3 ~ a*x3
#> # Path b
#> y ~ b*m3
#> #Indirect effect (a*b)
#> ab := a*b

library(lavaan)
#> This is lavaan 0.6-15
#> lavaan is FREE software! Please report any bugs.
fits <- lapply(mods, sem, data = dat, se = "boot", iseed = 1234, bootstrap = 50)
lapply(fits, parameterEstimates)
#> [[1]]
#>   lhs op rhs label    est    se      z pvalue ci.lower ci.upper
#> 1   y  ~  x1     c  0.170 0.081  2.111  0.035   -0.012    0.336
#> 2  m1  ~  x1     a  0.072 0.086  0.839  0.402   -0.133    0.226
#> 3   y  ~  m1     b -0.121 0.087 -1.387  0.165   -0.310    0.067
#> 4   y ~~   y        0.974 0.119  8.187  0.000    0.703    1.185
#> 5  m1 ~~  m1        0.862 0.114  7.555  0.000    0.609    1.099
#> 6  x1 ~~  x1        0.980 0.000     NA     NA    0.980    0.980
#> 7  ab := a*b    ab -0.009 0.017 -0.524  0.600   -0.056    0.027
#> 
#> [[2]]
#>   lhs op rhs label   est    se     z pvalue ci.lower ci.upper
#> 1   y  ~  x2     c 0.029 0.090 0.319  0.750   -0.131    0.296
#> 2  m2  ~  x2     a 0.101 0.108 0.932  0.352   -0.084    0.385
#> 3   y  ~  m2     b 0.033 0.084 0.389  0.697   -0.117    0.244
#> 4   y ~~   y       1.010 0.117 8.656  0.000    0.728    1.228
#> 5  m2 ~~  m2       1.198 0.214 5.592  0.000    0.840    1.732
#> 6  x2 ~~  x2       1.103 0.000    NA     NA    1.103    1.103
#> 7  ab := a*b    ab 0.003 0.019 0.172  0.864   -0.019    0.084
#> 
#> [[3]]
#>   lhs op rhs label    est    se      z pvalue ci.lower ci.upper
#> 1   y  ~  x3     c -0.065 0.097 -0.675  0.499   -0.240    0.119
#> 2  m3  ~  x3     a -0.054 0.097 -0.558  0.577   -0.279    0.182
#> 3   y  ~  m3     b -0.177 0.099 -1.779  0.075   -0.378    0.062
#> 4   y ~~   y        0.986 0.121  8.150  0.000    0.720    1.245
#> 5  m3 ~~  m3        0.718 0.112  6.420  0.000    0.492    0.929
#> 6  x3 ~~  x3        1.118 0.000     NA     NA    1.118    1.118
#> 7  ab := a*b    ab  0.010 0.018  0.534  0.594   -0.041    0.052

首先，形成x变量和m变量的向量。其次，使用paste()创建模型语法字符串的向量。对cat()的调用是可选的，用于验证生成的模型语法字符串。第三，简单地使用lapply()来拟合每个生成的模型。如果需要bootstrap CI，则需要设置iseed，以确保在所有模型中使用相同的bootstrap样本集。
lavaan有semList()用于将模型拟合到数据集列表。但是，我不知道有类似的函数可以将一系列模型拟合到同一个数据集。如果是，则可以使用该函数代替lapply()。