[R] No convergence using ADAPT

Philip Turk Philip.Turk at NAU.EDU
Sat Jul 7 21:20:05 CEST 2007


I am trying calculate a probability using numerical integration. The first 
program I ran spit out an answer in a very short time. The program is below:

## START PROGRAM

trial <- function(input)

{
pmvnorm(lower = c(0,0), upper = c(2, 2), mean = input, sigma = matrix(c(.1, 0, 
0, .1), nrow = 2, ncol = 2, byrow = FALSE)) 
}

require(mvtnorm)
require(adapt)

bottomB <- -5*sqrt(.1)
topB <- 2 + 5*sqrt(.1)
areaB <- (topB - bottomB)^2

unscaled.Po.in.a <- adapt(2, lo = c(bottomB, bottomB), up = c(topB, topB), 
minpts = 1000, eps = 1e-4, functn = trial)

(1/areaB)*unscaled.Po.in.a$value

## FINISH PROGRAM

I tried to run the program again changing a.) sigma in the trial function, b.) 
upper in the trial function, and c.) the bounds of integration; that is, 
bottomB and topB.  The new program is below:

## START PROGRAM

trial <- function(input)

{
pmvnorm(lower = c(0,0), upper = c(10, 10), mean = input, sigma = matrix(c(.01, 
0, 0, .01), nrow = 2, ncol = 2, byrow = FALSE)) 
}

require(mvtnorm)
require(adapt)

bottomB <- -5*sqrt(.01)
topB <- 10 + 5*sqrt(.01)
areaB <- (topB - bottomB)^2

unscaled.Po.in.a <- adapt(2, lo = c(bottomB, bottomB), up = c(topB, topB), 
minpts = 1000, eps = 1e-4, functn = trial)

(1/areaB)*unscaled.Po.in.a$value

## FINISH PROGRAM

Now, the program just runs and runs (48 hours at last count!).  By playing 
around with the program, I have deduced the program is highly sensitive to 
changing the upper option in the trial function.  For example, using a vector 
like (4, 4) causes no problems and the program quickly yields an answer.  I 
have a couple of other programs where I can easily obtain a simulation-based 
answer, but I would ultimately like to know what's up with this program before 
I give up on it so I can learn a thing or two.  Does anyone have any clues or 
tricks to get around this problem?  My guess is that it will simply be very 
difficult (impossible?) to obtain this type of relative error (eps = 1e-4) and 
I will have no choice but to pursue the simulation approach.

Thanks for any responses (philip.turk at nau.edu)!

-- Phil

Philip Turk
Assistant Professor of Statistics
Northern Arizona University
Department of Mathematics and Statistics
PO Box 5717
Flagstaff, AZ 86011
Phone: 928-523-6884
Fax: 928-523-5847
E-mail: Philip.Turk at nau.edu
Web Site: http://jan.ucc.nau.edu/~stapjt-p/



More information about the R-help mailing list