集成到 R

Question

我正在尝试在 R 中计算 。我有以下值。

nb <- 100
tb <- 25
ns <- 90
ts <- 15
A0 <- 1
S_norm <- 0.4

R <-tb/ts
y_meas <- (ns-nb/R)/A0
sigma_meas = sqrt(ns+(nb+1)/R^2)/A0

我很困惑如何积分 L(psi)，比如从 -10 到 10。因为我正在积分 A。

Answer 1

您可以用logA 代替固定值psi，您可以按如下方式进行积分：

psi <- 5
integrate(function(x) exp(-0.5*(((x/A0)/S_norm)^2 + ((psi-y_meas*A0/exp(x))/sigma_meas)^2)),
                                                                                   -10, 10)
# 0.1775989 with absolute error < 6.6e-05

Answer 2

除了@SandipanDey 的出色回答之外，如果您可以将限制扩展到 -Infinity...+Infinity，还有一种更好的方法可以将函数与 e^{-x2 集成} 内核：Gauss-Hermite quadrature，并且有一个 R 包。

简单示例：

library(gaussquad)

n.quad <- 128 # integration order

# get the particular (weights,abscissas) as data frame
# with 2 observables and n.quad observations
rule <- ghermite.h.quadrature.rules(n.quad, mu = 0.0)[[n.quad]]

# test function - integrate 1 over exp(-x^2) from -Inf to Inf
# should get sqrt(pi) as an answer
f <- function(x) {
    1.0
}

q <- ghermite.h.quadrature(f, rule)
print(q - sqrt(pi))