在R中正确使用integrate函数及解决混合高斯模型拟合与积分问题

1. 原始数据 #

x <- data.frame(
  V1 = c(15260.14, 15260.16, 15260.17, 15260.19, 15260.20, 15260.22, 15260.23,
         15260.25, 15260.26, 15260.28, 15260.29, 15260.31, 15260.32),
  V2 = c(0.04629, 0.22787, 0.68676, 0.89477, 0.50650, 0.13612, 0.07962,
         0.14235, 0.43131, 0.73034, 0.55780, 0.19124, 0.06062)
)

（数据对应双峰混合高斯模型，横轴为V1，纵轴V2为对应位置的密度值）

2. 混合高斯模型拟合参数异常 #

问题根源 #

你错误地将纵轴V2作为待拟合的样本数据，而实际应该对横轴V1建模，V2是每个V1点的权重/密度值。直接拟合V2导致参数完全偏离图示的峰位置。

修正方案 #

使用带权重的混合高斯拟合，将V2作为V1的权重：

library(mixtools)
# k=2指定双成分混合模型
x2_weighted <- normalmixEM(x$V1, w = x$V2, k = 2)
lambda_weighted <- x2_weighted$lambda
mu_weighted <- x2_weighted$mu
sigma_weighted <- x2_weighted$sigma

此时得到的mu会对应V1轴上的两个峰位置（约15260.18和15260.28），符合预期。

3. Gaussianmix函数输出全为0的问题 #

问题根源 #

1. 高斯密度公式错误：指数部分遗漏平方项，正确公式应为((x-mu)/sigma)^2
2. 参数不匹配：原拟合用V2的均值（0.12/0.62）匹配V1的数值（15260+），量级差过大导致指数项为极大负数，exp()后趋近于0

修正方案 #

用R内置的dnorm()简化公式，同时使用权重拟合后的参数：

# 定义单成分高斯密度函数
Gaussianmix <- function(lambda, mu, sigma, x) {
  lambda * dnorm(x, mean = mu, sd = sigma)
}
# 向量化处理多成分
Gaussianmix_vec <- Vectorize(Gaussianmix, vectorize.args = c("mu", "sigma", "lambda"))

# 测试：使用修正后的参数和V1区间
x_seq <- seq(15260.14, 15260.32, by = 0.005)
result <- Gaussianmix_vec(lambda_weighted, mu_weighted, sigma_weighted, x_seq)
head(result)

4. integrate函数使用错误 #

问题根源 #

integrate()要求传入的函数仅接受积分变量（x）作为参数，直接传入所有参数会导致函数提前执行，变成数值而非函数对象。

修正方案 #

用匿名函数固定lambda/mu/sigma，仅保留x作为变量：

# 定义混合高斯总密度函数（多成分求和）
mix_gaussian_total <- function(x, lambda, mu, sigma) {
  sum(lambda * dnorm(x, mean = mu, sd = sigma))
}

# 执行积分
integral_result <- integrate(
  f = function(x) mix_gaussian_total(x, lambda_weighted, mu_weighted, sigma_weighted),
  lower = 15260.14,
  upper = 15260.32
)
print(integral_result)

5. trapz结果过小的问题 #

trapz(x$V1, x$V2)的结果0.0698是正确的：V1的区间长度仅为0.18（15260.32-15260.14），V2峰值不到1，曲线下面积自然在0.07左右。如果认为结果不符合预期，可能是混淆了两个概念：

• 若要计算拟合的混合高斯分布在该区间的积分，使用上述integrate()方法
• 若要计算原始数据曲线下的面积，trapz的结果就是准确值

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