<table>

        <tr><th>Library     <td>Scythe MCMC - A Scythe Markov Chain Monte Carlo C++ Framework

        <tr><th>Version     <td>0.1

    </table>

    Scythe MCMC is a C++ header library that eases the development of Markov Chain Monte Carlo (MCMC).  It is based on the

    <a href="http://scythe.wustl.edu/">Scythe Statistical Library</a>.

    allows development of models that cannot be sampled effectively using either Bugs or Jags, such

    nonparametric Bayesian models without truncation approximations.

    - Based on the <a href="http://scythe.wustl.edu/">Scythe Statistical Library</a>.

    - Eliminates commandline and option parsing boilerplate code using

       - Slice sampler for Dirichlet Process and Indian Buffet Process.

    - MIT License.

    The license text below is the boilerplate "MIT License" used from:

    http://www.opensource.org/licenses/mit-license.php

mersenne myrng; 

 * These options are generally specified at the command line and are generic to all

 * MCMC samplers.  Other options are passed through an .ini config file.

 */

/// around everwhere.

MCMCOptions mcmc_options;

  Source: http://www.techbytes.ca/techbyte103.html

 */

bool FileExists(string strFilename) {

  struct stat stFileInfo;

  return(blnReturn);

}

   Diplays:

<pre>

   ScytheMCMC v. 0.1

       << "\t   --burnin=0 --thin=1 --random-seed=12345" << endl << endl;

}

  Parses and stores command line options and stores them in MCMCOptions instance

  mcmc_options.  Initializes random mersenne instance myrng.

  @see ShowUsage.

}

 * The model is defined as many instances of parameter subclasses.  Different Step types require different methods

 * to be implemented in each Parameter subclass.

 *

public:

  /// Default class constructor.

  Parameter() {}

   *

  */

  Parameter(bool track, string label) : track_(track), label_(label) {}

  /// Function for deterministic nodes.

  virtual double Function() { return 0.0; };

  /// Starting value.

  virtual double StartingValue() { return 0.0; };

  /// Log of the probability density (plus constant)

  virtual double LogDensity(double value) { return 0.0; }

  /// Return a random draw from the posterior.

  /// Random draw from posterior is called by GibbsStep.

  virtual double RandomPosterior() { return 0.0; }

  /// Value of parameter.

  virtual double Value() { return 0.0; };

  /// Save a new value of the parameter.

  virtual void Save(double new_value) {};

  /// Parameter is tracked / saved.

  bool Track() {

    return track_;

  }

  /// String label of parameter.

  string Label() {

    return label_;

  }

};

 * The step (Gibbs, Metropolis-Hastings, Slice, etc) describes how the sampler 

 * calculates the next value of each parameter.  The key function is Step::DoStep, which is

 * equivalent to the Execute function in a <a href="http://en.wikipedia.org/wiki/Command_pattern">Command Pattern</a>.

 */

class Step {

public:

   * Executes the main step function, such as taking a MH step or deterministic step.

   */ 

  virtual void DoStep() {}

   *  Calls the parameter's Parameter::Starting function and then saves the result using Parameter::Save.

   */

  virtual void Start() {}

  /// Return parameter value

  /// @see Parameter::Value

  virtual double ParameterValue() { return 0.0; }

  /// Return parameter label

  /// @see Parameter::Label

  virtual string ParameterLabel() { return ""; }

  /// Return the tracking status

  /// @see Parameter::Track

  virtual bool ParameterTrack() { return true; }

};

 * A univariate Gibbs step draws from the Parameter::RandomPosterior function of a Parameter and then saves the result

 * using Parameter:Save. 

 */

public:

  /// Default constructor, for copy assignments, etc.

  GibbsStep() {}

   */

  GibbsStep(ParameterType parameter) : parameter_(parameter) {

  }

   *  Draws from the parameter's Parameter::RandomPosterior function and then saves the result using Parameter::Save.

   */

  void DoStep() {

  ParameterType parameter_;  

};

  * Deterministic nodes can be helpful to track summary statistics or to reduce computations through intermediate use

  * of sufficient statistics across multiple parameters.  FunctionStep calls the parameters Parameter::Function method

  * and saves the result using Parameter::Save.

public:

  /// Default constructor, for copy assignments, etc.

  FunctionStep() {}

   */

  FunctionStep(ParameterType parameter) : parameter_(parameter) {

  }

   *  Calculates the value of Parameter::Function and then saves the result using Parameter::Save.

   */

  void DoStep() {

    parameter_.Save(parameter_.Function());

  }

   * parameters in an order that makes the function feasible.  StartingValue functions are called in the order

   * they are added.  Any implemention of a StartingValue method for a FunctionStep will not be used.

   */

};

 *  Normal proposal draws from N(starting_value, standard_deviation).

 */

class NormalProposal {

public:

  NormalProposal() {}

   * \param standard_deviation MH tuning parameter for normal model

   */

  NormalProposal(double standard_deviation) : standard_deviation_(standard_deviation) {}

   * \param starting_value Starting value.

   */

  double Draw(double starting_value) {

    return myrng.rnorm(starting_value, standard_deviation_);

  }

   * \param starting_value Starting value.

   * \param new_value Proposed new value.

   */

};

  * Basic univariate Metropolis-Hastings step.  Requires both a Parameter and Proposal instance.

  */

template <typename ParameterType, typename ProposalType>

public:

  /// Default constructor, for copy assignments, etc.

  MetropStep() {}

   */

  MetropStep(ParameterType parameter, ProposalType proposal) : parameter_(parameter), proposal_(proposal) {}

   *  Calculates the value of Parameter::Function and then saves the result using Parameter::Save.

   */

  void DoStep() {

  ProposalType proposal_;

};

  */

template <typename ParameterType>

class SliceStep: public Step {

public:

  /// Default constructor, for copy assignments, etc.

  SliceStep() {}

   */

  SliceStep(ParameterType parameter, double w = 1.0, double lower=-dInf, double upper=dInf) : 

    parameter_(parameter), w_(w), lower_(lower), upper_(upper) {}

   */

  void DoStep() {

    // Find the log density at the initial point.

        break;

      }

      L = L - w_;

    }

    while(true) {

        break;

      }

      R = R + w_;

    }

    // Shrink interval to lower and upper bounds.

      else {

        L = x1;

      }

    }

    // Save the point sampled

  double upper_;

};

 *  The sampler is the main MCMC object that holds all the MCMC steps for each parameter.  Running the sampler

 *  performs MCMC sampling for the model, saving results, and displaying progress.  In the language of the Command Pattern, the

 *  sampler is the Invoker or Command Manager.

 *

 */

class Sampler {

public:

  Sampler(MCMCOptions options) {

    cout << "Creating Sampler... " << endl;

    sample_size_ = options.sample_size;

    thin_ = options.thin;

    out_file_ = options.out_file;

  }

   *  All parameters should have an associated Step that is added to the sampler.  The sampler stack

   *  defines one sampler iteration.

  void AddStep(Step* step) {

    // Add step to step stack.

    steps_.push_back(step);

      tracks_.push_back(steps_.size() - 1);

    }

  }

   */

  void Iterate(int number_of_iterations, bool progress = false) {

    if(progress) {

      }

    }

  }

   * Run the MCMC sampling procedure, including burning in, sampling, thinning, displaying progress,

   * and saving results.

   */

    cout << "Total elapsed time: " << timer.elapsed() << " seconds" << endl;

  }

   */

  int NumberOfSteps() {

    return steps_.size();

  }

   */

  int NumberOfTrackedSteps() {

    return tracks_.size();

  vector<int> tracks_;

};

 *  Checks if two doubles are within the numeric_limits<double>::epsilon() precision of each other.

 *  Useful if === is too strict.

 */