I tried to program the algorithm for cdf for multidimensional t-distribution after Genz and Bretz, the reference package in R is mvtnorm.
When I tested my function, I found that my numbers did not match. In the following example, adjusted using the mvtnorm help, the multidimensional random variable t has independent components. So the integral should just be the product of three independent probabilities
> lower <- -1
> upper <- 3
> df <- 4
> corr <- diag(3)
> delta <- rep(0, 3)
> pmvt(lower=lower, upper=upper, delta=delta, df=df, corr=corr)
[1] 0.5300413
attr(,"error")
[1] 4.321136e-05
attr(,"msg")
[1] "Normal Completion"
Reported Error 4e-5, Error Compared to Independent Probability Product
> (pt(upper, df) - pt(lower, df))**3
[1] 0.4988254
there is
0.5300413 - 0.4988254 = 0.0312159
I get inconsistencies in my own code compared to R mvtnorm for different examples in about the same range.
I'm basically new to R. So what am I doing wrong or what's wrong?
( R-help, .)
: pchalasani, , , t.
, , .
% * 100 (10000 ) ( ).
([[26, 25, 24, 23],
[24, 23, 24, 25],
[24, 27, 24, 24],
[24, 23, 26, 25]])
t
([[29, 20, 22, 29],
[20, 31, 28, 21],
[20, 29, 29, 20],
[29, 18, 18, 29]])
.
(, R-, , R.)