knitr::opts_chunk$set( message = FALSE, warning = FALSE, error = FALSE, tidy = FALSE, cache = FALSE )

Version 2.0 integrates two well-known cubature libraries in one place:

- The cubature C library of Steven G. Johnson.
- The Cuba C library of Thomas Hahn.

It also provides a single function `cubintegrate`

that allows one to
call all methods in a uniform fashion, as I explain below.

**N.B.** One has to be aware that there are cases where one library
will integrate a function while the other won't, and in some cases,
provide somewhat different answers. That still makes sense and depends
on the underlying methodology used.

Following a suggestion by Simen Guare, we now have a function
`cubintegrate`

that can be used to try out various integration methods
easily. Some examples.

library(cubature) m <- 3 sigma <- diag(3) sigma[2,1] <- sigma[1, 2] <- 3/5 ; sigma[3,1] <- sigma[1, 3] <- 1/3 sigma[3,2] <- sigma[2, 3] <- 11/15 logdet <- sum(log(eigen(sigma, symmetric = TRUE, only.values = TRUE)$values)) my_dmvnorm <- function (x, mean, sigma, logdet) { x <- matrix(x, ncol = length(x)) distval <- stats::mahalanobis(x, center = mean, cov = sigma) exp(-(3 * log(2 * pi) + logdet + distval)/2) }

First we try the scalar invocation with `hcubature`

.

cubintegrate(f = my_dmvnorm, lower = rep(-0.5, 3), upper = c(1, 4, 2), method = "pcubature", mean = rep(0, m), sigma = sigma, logdet = logdet)

We can compare that with Cuba's `cuhre`

.

cubintegrate(f = my_dmvnorm, lower = rep(-0.5, 3), upper = c(1, 4, 2), method = "cuhre", mean = rep(0, m), sigma = sigma, logdet = logdet)

The Cuba routine can take various further arguments; see for example,
the help on `cuhre`

. Such arguments can be directly passed to
`cubintegrate`

.

cubintegrate(f = my_dmvnorm, lower = rep(-0.5, 3), upper = c(1, 4, 2), method = "cuhre", mean = rep(0, m), sigma = sigma, logdet = logdet, flags = list(verbose = 2))

As there are many such method-specific arguments, you may find the
function `default_args()`

useful.

str(default_args())

`cubintegrate`

provides vector intefaces too: the parameter `nVec`

is
by default 1, indicating a scalar interface. Any value > 1 results in
a vectorized call. So `f`

has to be constructed appropriately, thus:

my_dmvnorm_v <- function (x, mean, sigma, logdet) { distval <- stats::mahalanobis(t(x), center = mean, cov = sigma) exp(matrix(-(3 * log(2 * pi) + logdet + distval)/2, ncol = ncol(x))) }

Here, the two underlying C libraries differ. The cubature library
manages the number of points used in vectorization dynamically and
this number can even vary from call to call. So any value of `nVec`

greater than 1 is merely a flag to use vectorization. The Cuba C
library on the other hand, will use the actual value of `nVec`

.

cubintegrate(f = my_dmvnorm_v, lower = rep(-0.5, 3), upper = c(1, 4, 2), method = "pcubature", mean = rep(0, m), sigma = sigma, logdet = logdet, nVec = 128)

cubintegrate(f = my_dmvnorm_v, lower = rep(-0.5, 3), upper = c(1, 4, 2), method = "cuhre", mean = rep(0, m), sigma = sigma, logdet = logdet, nVec = 128)

```
sessionInfo()
```

**Any scripts or data that you put into this service are public.**

Embedding an R snippet on your website

Add the following code to your website.

For more information on customizing the embed code, read Embedding Snippets.