marginal_posterior {aghq}

Marginal Posteriors

Description

Compute the marginal posterior for a given parameter using AGHQ. This function is mostly called within aghq().

Usage

marginal_posterior(...)

## S3 method for class 'aghq'
marginal_posterior(
  quad,
  j,
  qq = NULL,
  method = c("auto", "reuse", "correct"),
  ...
)

## S3 method for class 'list'
marginal_posterior(
  optresults,
  k,
  j,
  basegrid = NULL,
  ndConstruction = "product",
  ...
)

Arguments

...

Additional arguments to be passed to optresults$ff, see ?optimize_theta.

quad

Object returned by aghq::aghq.

j

Integer between 1 and the dimension of the parameter space. Which index of the parameter vector to compute the marginal posterior for.

qq

Numeric vector of length >=1 giving the points at which to evaluate the marginal posterior. The default, NULL, chooses these points in a 'clever' way, see Details.

method

Method for computing the quadrature points used to approximate moment. One of 'reuse' (default) or 'correct'. See details. The default SHOULD be 'correct'; it is currently set to 'reuse' to maintain compatibility of results with previous versions. This will be switched in a future major release.

optresults

The results of calling aghq::optimize_theta(): see return value of that function.

k

Integer, the number of quadrature points to use. I suggest at least 3. k = 1 corresponds to a Laplace approximation.

basegrid

Optional. Provide an object of class NIGrid from the mvQuad package, representing the base quadrature rule that will be adapted. This is only for users who want more complete control over the quadrature, and is not necessary if you are fine with the default option which basically corresponds to mvQuad::createNIGrid(length(theta),'GHe',k,'product').

ndConstruction

Create a multivariate grid using a product or sparse construction? Passed directly to mvQuad::createNIGrid(), see that function for further details. Note that the use of sparse grids within aghq is currently experimental and not supported by tests. In particular, calculation of marginal posteriors is known to fail currently.

Details

If qq=NULL, then it is set to the unique values in an adapted GHQ grid computed assuming that j=1 (there is nothing special about this procedure, it's just a way to provide an apparently sensible default).

If method='reuse', then the parameter vector is reordered so j=1, and the approximate marginal is computed by first computing the whole AGHQ grid, then summing over the other indices. This is an outdated method that does not have any theory pertaining to it, and is included for backwards compatibility. It does not use qq if supplied.

If method='correct' then the theoretically-justified approximation from Section 2.4 of the 'Stochastic Convergence Rates...' paper is returned.

method='auto' currently chooses 'reuse' for backwards compatibility, but this will be changed in a future release.

Value

a data.frame containing the normalized log marginal posterior for theta_j evaluated at the original quadrature points. Has columns "thetaj","logpost_normalized","weights", where j is the j you specified.

See Also

Other summaries: compute_pdf_and_cdf(), compute_quantiles(), interpolate_marginal_posterior()

Examples

## A 2d example ##
logfteta2d <- function(eta,y) {
  # eta is now (eta1,eta2)
  # y is now (y1,y2)
  n <- length(y)
  n1 <- ceiling(n/2)
  n2 <- floor(n/2)
  y1 <- y[1:n1]
  y2 <- y[(n1+1):(n1+n2)]
  eta1 <- eta[1]
  eta2 <- eta[2]
  sum(y1) * eta1 - (length(y1) + 1) * exp(eta1) - sum(lgamma(y1+1)) + eta1 +
    sum(y2) * eta2 - (length(y2) + 1) * exp(eta2) - sum(lgamma(y2+1)) + eta2
}
set.seed(84343124)
n1 <- 5
n2 <- 5
n <- n1+n2
y1 <- rpois(n1,5)
y2 <- rpois(n2,5)
objfunc2d <- function(x) logfteta2d(x,c(y1,y2))
funlist2d <- list(
  fn = objfunc2d,
  gr = function(x) numDeriv::grad(objfunc2d,x),
  he = function(x) numDeriv::hessian(objfunc2d,x)
)
opt_sparsetrust_2d <- optimize_theta(funlist2d,c(1.5,1.5))

# Now actually do the marginal posteriors
marginal_posterior(opt_sparsetrust_2d,3,1)
marginal_posterior(opt_sparsetrust_2d,3,2)
marginal_posterior(opt_sparsetrust_2d,7,2)