# RunMh: Metropolis Hasting Algorithm Constrained on a Simplex In SALTSampler: Efficient Sampling on the Simplex

## Description

This function runs the Metropolis Hasting algorithm constrained on a simplex. The function can be used with any target distribution on the simplex defined by the user. Alternatively, two common target distributions are built into the function and can be specifed by the user. The function is designed to continue to perform well in difficult cases, such as those in high dimensions or with parameters that differ by orders of magnitude. Care is also taken to ensure accuracy even when some coordinates are numerically close to 0 or 1.

## Usage

 `1` ```RunMh(center, B, concentration = 1, h, type = 'user', dat = NULL, pars = NULL) ```

## Arguments

 `center` Vector of numeric values summing to 1 that define the center of the distributional parameters of the posterior. For type `'dirichlet'`, the parameter a is defined such that a_{i} is the ith element of `center` times `concentration`. For type `'multinom'`, the multinomial distribution parameter, p_{i}, is the ith value of `center` `B` Number of iterations to run the chain `concentration` This argument specifies the concentration parameter where a is defined such that a_{i} is the ith element of `center` times `concentration`. This is typically used with type `'dirichlet'`, but can also be used in a user-defined function. This arguments defaults to 1, so has no effect if it is not specified. `h` Vector of step sizes. Length of vector must match length of `center` `type` Specifies the target distribution. Select type `'user'` if a target distribution has already been defined (see details). Select type `'dirichlet'` for a Dirichlet distribution and type `'multinom'` for a multinomial distribution `dat` A matrix or vector passing data to the sampler. For type `'multinom'`, this is a matrix giving data from repeated multinomial draws where the data is formatted in the same way as data obtained via `GenData`. The number of the items in the ith bin on the jth multinomial trial should be in the ith column and the jth row of the matrix. For type `'user'`, any matrix or vector of data can be used to match the form specified in the user's target function. If unspecified, this argument defaults to `NULL` `pars` A list of additional parameters that can be passed to the user-specified target function for type `'user'` if desired. Argument defaults to `NULL`

## Details

Any target distribution on the simplex can be used with this function by defining a target distribution function in the environment prior to running `RunMh`. The function should be named `Target` and should take in parameters `ycand` and `ycurrent`, which are the current and proposed samples on the logit scale, and parameter `a`, which is `center` times `concentration`. Parameters `dat` and `pars` can be set to `NULL`. Alternatively, `dat` can be used to provide data to the target function and/or `pars` can be used to provide a list of additional parameters to the the target function. The target function should output the ratio of the log-likelihood of the posterior distribution for the proposal, θ = `ycand`, to the log-likelihood of the posterior for the current value, θ = `ycurrent`. For simple cases, there are built-in target distributions. For type `'dirichlet'`, `RunMh` uses a Dirichlet distribution as a posterior distribution. For type `'multinomial'`, `RunMh` samples the distributional parameters of a multinomial distribution that would have generated the data inputted for `dat`.

## Value

An object of class `mhOut`. `mhOut` has 12 attributes.

 `Y` Matrix of MCMC samples on logit scale `S` Matrix of MCMC samples on true scale `runTime` Summary of the MCMC runtime. The first entry gives the total user CPU time, the second entry gives the system CPU time, and the third entry gives the true elapsed time `moveCount` Number of steps where the proposal value was accepted `p` Length of `center` vector `center` Vector of numeric values summing to 1 that help to define distributional parameters. For type `'dirichlet'`, the parameter a is defined such that a_{i} is the ith element of `center` times `concentration`. For type `'multinom'`, the multinomial distribution parameter, p_{i}, is the ith value of `center` `B` Number of iterations to run the chain `concentration` For type `'dirichlet'`, this argument specifies the concentration parameter where a is defined such that a_{i} is the ith element of `center` times `concentration`. Otherwise, this argument takes on its default value of 1 and has no effect `h` Vector of step sizes. Length of vector must match length of `center` `type` Specifies the target distribution. Select type `'user'` if a target distribution has already been defined (see details). Select type `'dirichlet'` for a Dirichlet distribution and type `'multinom'` for a multinomial distribution `dat` A matrix or vector passing data to the sampler. For type `'multinom'`, a matrix giving data from repeated multinomial draws where the data is formatted in the same way as data obtained via `GenData`. The number of the items in the ith bin on the jth multinomial trial should be in the ith column and the jth row of the matrix. For type `'user'`, any matrix or vector of data can be used to match the form specified in the user's target function. If unspecified, this argument defaults to `NULL` `a` Dirichlet distribution parameters, a, where a_{i}, is the ith element of `center` times `concentration`

## Examples

 ``` 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39``` ```###Dirichlet sampling in 3-simplex dir <- RunMh(center = c(0.7, 0.2, 0.1), B = 2e3, concentration = 10, h = c(2, 2, 2), type = 'dirichlet', dat = NULL) ####Multinomial sampling ## Not run: sampData <- GenData(center = c(0.2, 0.3, 0.5), n = 100, size = 10) multinom <- RunMh(center = c(0.2, 0.3, 0.5), B = 1e4, h = c(2,2,2), type = 'multinom', dat = sampData) ## End(Not run) ####User-defined target distribution for a calibration problem ## Not run: #Known function which we want to calibrate CalibFn <- function(y, logit = FALSE) { if (logit == TRUE) { y <- exp(LogPq(y)\$logp) } out <- 1e3*y[1]^3*y[2]^3/sqrt(20 + y[3]) return(out) } #Generate data z <- rnorm(n = 1000, mean = CalibFn(c(1/3, 1/3, 1/3), 2)) #User defined target distribution Target <- function(ycand, ycurrent, a, dat, pars = NULL) { out <- sum(dnorm(dat, CalibFn(ycand, logit = TRUE), 2, log = TRUE)) - sum(dnorm(dat, CalibFn(ycurrent, logit = TRUE), 2, log = TRUE)) + sum((a - 1)*(LogPq(ycand)\$logp - LogPq(ycurrent)\$logp)) return(out) } #Run sampler inputDist <- RunMh(center = c(1/3, 1/3, 1/3), B = 3e4, concentration = 3, h = c(0.2, 0.2, 0.2), type = 'user', dat = z) ## End(Not run) ```

SALTSampler documentation built on Aug. 11, 2017, 1:02 a.m.