# fec_stan: Modelling of faecal egg count data (one-sample case) In eggCounts: Hierarchical Modelling of Faecal Egg Counts

## Description

Models faecal egg counts data in a one-sample case with (zero-inflated) Poisson-gamma model formulation using Stan modelling language. It is computationally several-fold faster compared to conventional MCMC techniques. For the installation instruction of Stan, please read https://github.com/stan-dev/rstan/wiki/RStan-Getting-Started.

## Usage

 ```1 2 3 4``` ```fec_stan(fec, rawCounts = FALSE, CF = 50, zeroInflation = TRUE, muPrior, kappaPrior, phiPrior, nsamples = 2000, nburnin = 1000, thinning = 1, nchain = 2, ncore = 1, adaptDelta = 0.95, saveAll = FALSE, verbose = FALSE) ```

## Arguments

 `fec` vector of faecal egg counts `rawCounts` logical. If true, `preFEC` and `postFEC` correspond to raw counts (as counted on equipment). Otherwise they correspond to calculated epgs (raw counts times correction factor). Defaults to `FALSE`. `CF` a positive integer or a vector of positive integers. Correction factor(s) `zeroInflation` logical. If true, uses the model with zero-inflation. Otherwise uses the model without zero-inflation `muPrior` a list with hyper-prior information for the group mean epg parameter μ. The default prior is `list(priorDist = "gamma",hyperpars=c(1,0.001))`, i.e. a gamma distribution with shape 1 and rate 0.001, its 90% probability mass lies between 51 and 2996 `kappaPrior` a list with hyper-prior information for the group dispersion parameter κ. The default prior is `list(priorDist = "gamma",hyperpars=c(1,0.7))`, i.e. a gamma distribution with shape 1 and rate 0.7, its 90% probability mass lies between 0.1 and 4.3 with a median of 1 `phiPrior` a list with hyper-prior information for zero-inflation parameter. The default prior is `list(priorDist = "beta",hyperpars=c(1,1))`, i.e. a uniform prior between 0 and 1 `nsamples` a positive integer specifying the number of samples for each chain (including burn-in samples) `nburnin` a positive integer specifying the number of burn-in samples `thinning` a positive integer specifying the thinning parameter, the period for saving samples `nchain` a positive integer specifying the number of chains `ncore` a positive integer specifying the number of cores to use when executing the chains in parallel `adaptDelta` the target acceptance rate, a numeric value between 0 and 1 `saveAll` logical. If TRUE, posterior samples for all parameters are saved in the `stanfit` object. If FALSE, only samples for μ, κ and φ are saved. Default to FALSE. `verbose` logical. If true, prints progress and debugging information

## Details

The first time each non-default model is applied, it can take up to 20 seconds for Stan to compile the model. Currently the function only support prior distributions with two parameters. For a complete list of supported priors and their parameterization, please consult the list of distributions in Stan http://mc-stan.org/documentation/.

The default number of samples per chain is 2000, with 1000 burn-in samples. Normally this is sufficient in Stan. If the chains do not converge, one should tune the MCMC parameters until convergence is reached to ensure reliable results.

## Value

Prints out summary of `meanEPG` as the posterior mean epg. The posterior summary contains the mean, standard deviation (sd), 2.5%, 50% and 97.5% percentiles, the 95% highest posterior density interval (HPDLow95 and HPDHigh95) and the posterior mode. NOTE: we recommend to use the 95% HPD interval and the mode for further statistical analysis.

The returned value is a list that consists of:

 `stan.samples` An object of S4 class `stanfit` representing the fitted results. For more information, please see the `stanfit-class` in `rstan` reference manual. `posterior.summary` A data frame that is the same as the printed posterior summary.

## Author(s)

Craig Wang

`simData1s` for simulating faecal egg count data with one sample
 ```1 2 3 4 5 6 7 8 9``` ```## Not run: ## load the sample data data(epgs) ## apply zero-infation model model <- fec_stan(epgs\$before, rawCounts=FALSE, CF=50) samples <- stan2mcmc(model\$stan.samples) ## End(Not run) ```