R中的三次積分（如何指定域）

我想計算R中三個變量函數的三重積分f(x,y,z)我使用包cubature和函數adaptIntegrate()。被積函數只在某個域（x<y<z，否則爲0）中等於1，我不知道如何指定。我嘗試了2種不同的函數實現，但它們都不起作用：R中的三次積分（如何指定域）

#First implementation 
fxyz <- function(w) { 
x <- w[1] 
y <- w[2] 
z <- w[3] 
x*y*z*(x < y)&(y < z) 
} 

#Second implementation 
fxyz <- function(w) { 
x <- w[1] 
y <- w[2] 
z <- w[3] 
if(x<y&y<z) 
    out<-1 
else 
    out<-0 
out 
} 

#Computation of integral 
library(cubature) 
lower <- rep(0,3) 
upper <- rep(1, 3) 
adaptIntegrate(f=fxyz, lowerLimit=lower, upperLimit=upper, fDim = 3)

有關如何正確指定域的任何想法？

來源

2015-04-04 Egodym

在你的第一個函數中，返回值是錯誤的。它應該是as.numeric(x<=y)*as.numeric(y<=z)。在你的第二個功能中，你還應該使用<=而不是<，否則`adapIntegrate將無法正常工作。您還需要指定最大數量的評估。試試這個

library(cubature) 
lower <- rep(0,3) 
upper <- rep(1,3) 

# First implementation (modified) 
fxyz <- function(w) { 
    x <- w[1] 
    y <- w[2] 
    z <- w[3] 
    as.numeric(x <= y)*as.numeric(y <= z) 
} 

adaptIntegrate(f=fxyz,lowerLimit=lower,upperLimit=upper,doChecking=TRUE, 
      maxEval=2000000,absError=10e-5,tol=1e-5) 
#$integral 
#[1] 0.1664146 
#$error 
#[1] 0.0001851699 
#$functionEvaluations 
#[1] 2000031 
#$returnCode 
#[1] 0

來源

2015-04-05 09:08:39 Bhas

太好了。使用'<='而不是'<'是否正確？我的域名是'x Egodym 2015-04-05 13:57:17

是的。但在這裏也是必要的。您不能在由<<構成的域中評估積分。想想吧！ – Bhas 2015-04-06 09:30:05

我不知道cubature包，但可以通過重複應用基R的integrate功能進行一維集成。

f.xyz <- function(x, y, z) ifelse(x < y & y < z, 1, 0) 
f.yz <- Vectorize(function(y, z) integrate(f.xyz, 0, 1, y=y, z=z)$value, 
        vectorize.args="y") 
f.z <- Vectorize(function(z) integrate(f.yz, 0, 1, z=z)$value, 
       vectorize.args="z") 

integrate(f.z, 0, 1) 
# 0.1666632 with absolute error < 9.7e-05

您可能想要使用控制參數來設置數值公差;內部整合中的小錯誤可以變成外部的大錯誤。

來源

2015-04-04 15:55:12

R中的三次積分（如何指定域）

回答

相關問題