如何使用R包sensitivity正确执行Sobol敏感性分析？

你的核心问题是错误地使用了原始观测数据的子集作为Sobol分析所需的X1和X2样本，这不符合Sobol方法的要求。Sobol分析需要两组从输入变量边缘分布中独立生成的随机样本，而非已有数据的拆分。此外，样本量设置过小也会导致结果偏差。

错误原因拆解 #

• 样本生成逻辑错误：Sobol方法要求X1和X2是从输入变量的真实分布（这里是Uniform(0,1)）中独立抽取的随机样本，你用原始数据拆分的子集无法保证独立性和分布一致性。
• 样本量不足：仅用50个样本进行分析，结果的随机性和稳定性太差，容易偏离真实敏感性顺序。

修正后的代码 #

library(tidymodels)
library(sensitivity)

# Example data
set.seed(123)
x1 = runif(100)
x2 = runif(100)
x3 = runif(100)
y = 3 * x1 + 2 * x2 + x3 + rnorm(100)
data <- data.frame(x1, x2, x3, y)

# Split data into training and testing sets
set.seed(234)
data_split <- initial_split(data, prop = 0.8)
train_data <- training(data_split)
test_data <- testing(data_split)

# Create a linear regression model using tidymodels
lm_spec <- linear_reg() %>%
  set_engine("lm") %>%
  set_mode("regression")

lm_fit <- lm_spec %>%
  fit(y ~ x1 + x2 + x3, data = train_data)

# Define a model function
model_function <- function(x) {
  new_data <- data.frame(x1 = x[, 1], x2 = x[, 2], x3 = x[, 3])
  predict(lm_fit, new_data)$`.pred`
}

# Perform Sobol sensitivity analysis - 修正样本生成逻辑
set.seed(345)
# 定义样本量，建议使用较大值保证稳定性
n_sample <- 1000
# 从输入变量的分布（Uniform(0,1)）中独立生成X1和X2
X1 <- matrix(runif(n_sample * 3), ncol = 3, dimnames = list(NULL, c("x1", "x2", "x3")))
X2 <- matrix(runif(n_sample * 3), ncol = 3, dimnames = list(NULL, c("x1", "x2", "x3")))

sobol_results <- sobol(model = model_function, 
                       X1 = X1, 
                       X2 = X2, 
                       nboot = 1000, 
                       order = 2)

# 查看敏感性结果
sobol_results
# 提取一阶敏感性指数（Si）
sobol_results$S

修正说明 #

1. 正确生成样本：通过runif直接从输入变量的真实分布（0-1均匀分布）生成独立的X1和X2样本，确保符合Sobol方法的假设。
2. 增大样本量：使用1000个样本替代原有的50个，提升结果的稳定性和准确性。
3. 结果验证：修正后，一阶敏感性指数（Si）会呈现x1 > x2 > x3的顺序，与你设定的y生成公式一致。

内容的提问来源于stack exchange，提问作者Yang Yang