voigt (version 2.0)

evoigt: Estimation of Voigt distribution parameters

Description

The function provides both point and interval estimates for the Voigt distribution parameters in a Bayesian way (see Details).

Usage

evoigt(data, hyper = NULL, init = NULL, S = 10000, burn = S/2,
thin = 10, fix.arg = NULL, chain = FALSE, ...)

Value

evoigt returns a list with the following components:

posterior mean

the vector of posterior means

posterior median

the vector of posterior medians

HPD interval

the highest probability density intervals for the parameters

chain

the thinned chain values after burn-in (only if argument chain is TRUE.)

Arguments

data

A numeric vector of length at least one containing only finite values.

hyper

A numeric vector of length six giving the hyperparameters for mu, sigma and gamma priors (in this order). If NULL all hyperparameters are set to 0.1.

init

The starting values for the chain parameters. If NULL (default) all chains start at 1. If scalar it is used as the common starting value. If character it can be set either to "rand" for random starting values or to "datadriven" for quantile based estimation of starting values.

S

The number of iterations in the Gibbs sampler. Default to \(10000\).

burn

The number of initial chain values to be discarded (burn-in period). Default to \(S/2\).

thin

This parameter allows the user to specify if and how the chain should be thinned after burn-in. By default thin = 5 is used, which corresponds to keeping 1/5 of the chain values.

fix.arg

An optional vector of length 3 giving the values of fixed parameters of the Voigt distribution, in this order: \(\mu,\sigma, \gamma\) with NA corresponding to an unknown parameter. By default all parameters are estimated. Parameters with fixed value are thus NOT estimated by the Gibbs sampler procedure.

chain

logical; if TRUE the output contains also the (thinned) chains values after burn-in. The first value of the chains is set according to the init argument.

...

Additional arguments to be passed to functions of the package coda. A prob argument can be used to specify the required probability of the credibility intervals via the HPDinterval function (default is 0.95; see the examples). Note that if S is too small the final probability can be different from the required one .

Author

tools:::Rd_package_author("voigt")

Details

The function runs a Gibbs sampler for estimating parameters using the following prior distributions: \(\mu\sim N(\mu_0, \nu_0^2)\), \(\sigma\sim \Gamma(a, b)\) and \(\gamma\sim \Gamma(c, d)\).

References

Cannas, M. and Piras, N. (2025) Estimation of Voigt Distribution Parameters: A Bayesian Approach. (conference paper) https://link.springer.com/chapter/10.1007/978-3-031-96303-2_53

Cannas, M. and Piras, N. Mixture representation and parameter estimation for the Voigt profile (submitted)

See Also

See HPDinterval for details on how the highest posterior density interval is built.

Examples

x = rvoigt(500, mu=0, sigma=1,gamma=1)
# point estimates and (default) 95% credibility intervals
evoigt(data=x)
# point estimates and 90% credibility intervals
evoigt(data=x, prob = 0.9)
# chain values
res = evoigt(data=x, prob = 0.9, chain=TRUE)
mu.chain = res$mu.chain
plot(0: (length(mu.chain)-1), type = "l", mu.chain, xlab="",ylab=c(expression(mu)) )
points(0,mu.chain[1],pch=1,col="red")
# if a parameter is known its value can be fixed using "fix.arg"
# e.g. set sigma =1:
res = evoigt(data=x, fix.arg=c(NA, 1, NA))
res["posterior mean"]
res["posterior median"]
res["HPD interval"]

# Inverse sampling from histogram
if (FALSE) {
x = seq(-15000,15000, by=0.1)
y = dvoigt(x, sigma=2, gamma=3, fast=FALSE) 
int<- 1:length(x)
breaks <-  x[- length(x)] + (x[-1] - x[-length(x)] ) /2
breaks = c( x[1]- (x[2]-x[1])/2 ,breaks, x[length(x)]+ (x[length(x)]-x[length(x)-1])/2)
bins <- sample(int, size=5000, prob=y, replace=TRUE)
samplehist <- vector()
for (j in 1:length(bins)){
samplehist[ j ] = runif(1, min= breaks[ bins[ j ] ], max= breaks[ bins[ j ] +1 ])
}
# estimation using sample histogram
evoigt(samplehist)}