2016-11-08 4 views
1

多変量正規desityに関連した二重積分を解く私は、既知の平均ベクトルと共分散行列を持つ多変量正規密度に関連した二重積分解決しようとしています

library(cubature) 

mu1 <- matrix(c(3,3), nrow=2) 
sigma1 <- rbind(c(4,-1), c(-1,6)) 

quadratic <- function(a,b) { 
    X <- matrix(c(a,b),nrow=2) 
    Q <- (-1/2)*t(X-mu1)%*%solve(sigma1)%*%(X-mu1) 
} 

NormalPDF <- function(x1,x2) { 
    f <- (1/(2*pi))*(1/sqrt(det(sigma1)))*exp(quadratic(x1,x2)) 
} 

# Solving for P(1 < X1 < 3, 1 < X2 < 3) 
P <- adaptIntegrate(NormalPDF(x1,x2), c(1,3), c(1,3)) 

しかし、それは続けて私にエラーを与えました:

Error in matrix(c(a, b), nrow = 2) : object 'x1' not found 

私のコードに明らかなエラーはありますか?

+1

は 'adaptIntegrate(NormalPDF、c(1,3)、c(1,3)) 'となります – HubertL

答えて

1

HubertLは、最初の引数は引数である関数呼び出しではなく、関数でなければならないことを指摘しています。この関数は "x"引数、長さ1の2つのベクトルを受け付けると想定されているので、NormalPDF関数をその引数とヘルパー関数の呼び出しで変更する必要があります。限界がどのように設定されているかという別のエラーがありました。

これを考慮してください

library(cubature) 

mu1 <- matrix(c(3,3), nrow=2) 
sigma1 <- rbind(c(4,-1), c(-1,6)) 

quadratic <- function(a,b) { 
    X <- matrix(c(a,b),nrow=2) 
    Q <- (-1/2)*t(X-mu1)%*%solve(sigma1)%*%(X-mu1) 
} 

NormalPDF <- function(x) { 
    f <- (1/(2*pi))*(1/sqrt(det(sigma1)))*exp(quadratic(x[1],x[2])) 
} 
# Solving for P(1 < X1 < 3, 1 < X2 < 3) 
P <- adaptIntegrate(NormalPDF, lowerLimit= c(1,1), upperLimit=c(3,3)) 
P 
#============== 
$integral 
[1] 0.09737084 

$error 
[1] 1.131395e-08 

$functionEvaluations 
[1] 17 

$returnCode 
[1] 0 

これは(3,3)で(1,1)と、「右上」に「左下」の正方形の上にその密度を積分します。ドメインは単一のポイントだったので、質問内の呼び出しは常に0を返していました。あなたが何か「数字」でそれをするつもりなら、P$integralでリストから抽出する必要があります。我々が(3,3)の最大値から1/4平面でのみ評価していたので、結果が0.25未満であることは妥当であると思われる。

+0

あなたの助けてくれてありがとう!それはそれを明らかにする。 – ywa136

+0

[積分]のSOタグは、単語の数学的使用とは明らかに関係していません。私は積分または積分のためのタグが表示されないので、それを挿入しませんでした。 –

関連する問題