Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Better handling of adaptation and algorithm parameters #24

Merged
merged 2 commits into from Aug 11, 2017
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
37 changes: 36 additions & 1 deletion ChangeLog
@@ -1,3 +1,38 @@
2017-08-11 Leah South <leah.south@hdr.qut.edu.au>

* inst/include/adaptMethods.h: Changing the name (from algParams.h) and
removing the parameter object so it is not a member of the class.
Instead, it is passed by reference to the adaptation functions.
* inst/include/moveset.h: Passing a reference to the algorithm
parameters to all functions.
* inst/include/sampler.h: Adding the parameter object as a member of
this class and updating accordingly. Passing the algorithm parameters
by reference to the adaptation and initialisation/move/MCMC functions.
parameters to all functions.

* inst/include/LinReg_LA_adapt.h: Switching to a pointer for the
sampler object and passing the parameter objects by reference to
the adaptation and move functions.
* src/LinReg_LA_adapt.cpp: Idem. Also adding the adaptation object
type to the template arguments for moveset.

* inst/include/blockpfgaussianopt.h: Adding params to the moveset
function arguments with an empty class as the type.
* inst/include/LinReg.h: Idem.
* inst/include/LinReg_LA.h: Idem.
* inst/include/nonLinPMMH.h: Idem.
* inst/include/pflineart.h: Idem.
* inst/include/pfnonlinbs.h: Idem.

* src/blockpfgaussianopt.cpp: Adding params to the moveset
function arguments and template arguments for sampler and moveset.
Using an empty class as the type.
* src/LinReg.cpp: Idem.
* src/LinReg_LA.cpp: Idem.
* src/nonLinPMMH.cpp: Idem.
* src/pflineart.cpp: Idem.
* src/pfnonlinbs.cpp: Idem.

2017-08-09 Leah South <leah.south@hdr.qut.edu.au>

* R/LinReg.R: Accessing radiata data through RcppSMC:: namespace rather
Expand Down Expand Up @@ -62,7 +97,7 @@
2017-08-04 Leah South <leah.south@hdr.qut.edu.au>

* R/LinReg.R: Silencing R CMD check --as-cran by not loading
data(radiata) into the global workspace.
data(radiata) into the global workspace.
* R/LinRegLA.R: Idem.
* inst/NEWS: Adding a line for each of the recent major changes.

Expand Down
6 changes: 3 additions & 3 deletions inst/include/LinReg.h
Expand Up @@ -41,7 +41,7 @@ namespace LinReg {

double logWeight(long lTime, const rad_state & value);
double logPosterior(long lTime, const rad_state & value);
void fInitialise(rad_state & value, double & logweight);
void fMove(long lTime, rad_state & value, double & logweight);
bool fMCMC(long lTime, rad_state & value, double & logweight);
void fInitialise(rad_state & value, double & logweight, smc::nullParams & param);
void fMove(long lTime, rad_state & value, double & logweight, smc::nullParams & param);
bool fMCMC(long lTime, rad_state & value, double & logweight, smc::nullParams & param);
}
6 changes: 3 additions & 3 deletions inst/include/LinReg_LA.h
Expand Up @@ -46,9 +46,9 @@ namespace LinReg_LA {
double logLikelihood(const rad_state & value);
double logPrior(const rad_state & value);

void fInitialise(rad_state & value, double & logweight);
void fMove(long lTime, rad_state & value, double & logweight);
bool fMCMC(long lTime, rad_state & value, double & logweight);
void fInitialise(rad_state & value, double & logweight, smc::nullParams & param);
void fMove(long lTime, rad_state & value, double & logweight, smc::nullParams & param);
bool fMCMC(long lTime, rad_state & value, double & logweight, smc::nullParams & param);

double integrand_ps(long,const rad_state &, void *);
double width_ps(long, void *);
Expand Down
17 changes: 9 additions & 8 deletions inst/include/LinReg_LA_adapt.h
Expand Up @@ -48,19 +48,20 @@ namespace LinReg_LA_adapt {
double logLikelihood(const rad_state & value);
double logPrior(const rad_state & value);

void fInitialise(rad_state & value, double & logweight);
void fMove(long lTime, rad_state & value, double & logweight);
bool fMCMC(long lTime, rad_state & value, double & logweight);
void fInitialise(rad_state & value, double & logweight, smc::staticModelAdapt & param);
void fMove(long lTime, rad_state & value, double & logweight, smc::staticModelAdapt & param);
bool fMCMC(long lTime, rad_state & value, double & logweight, smc::staticModelAdapt & param);

double integrand_ps(long,const rad_state &, void *);
double width_ps(long, void *);

//A derived class for adaptation which adapts parameters of type smc::staticModelAdapt
class rad_adapt:
public smc::algParam<rad_state,smc::staticModelAdapt>
public smc::adaptMethods<rad_state,smc::staticModelAdapt>
{
public:

void updateForMove(const smc::population<rad_state> & pop) {
void updateForMove(smc::staticModelAdapt & param, const smc::population<rad_state> & pop) {
unsigned long N = pop.GetNumber();
arma::vec loglike(N);
for (unsigned int i=0; i<N; i++){
Expand All @@ -69,7 +70,7 @@ namespace LinReg_LA_adapt {
param.ChooseTemp(pop.GetLogWeight(),loglike,rho*N);
}

void updateForMCMC(const smc::population<rad_state> & pop, double acceptProb, int nResampled, int & nRepeats) {
void updateForMCMC(smc::staticModelAdapt & param, const smc::population<rad_state> & pop, double acceptProb, int nResampled, int & nRepeats) {
if (nResampled == 0) {
nRepeats = 0;
} else {
Expand All @@ -87,6 +88,6 @@ namespace LinReg_LA_adapt {

};

smc::algParam<rad_state,smc::staticModelAdapt> * myParams;

smc::sampler<rad_state,smc::staticModelAdapt> * Sampler;
smc::adaptMethods<rad_state,smc::staticModelAdapt> * myAdapt;
}
27 changes: 10 additions & 17 deletions inst/include/algParam.h → inst/include/adaptMethods.h
@@ -1,6 +1,6 @@
// -*- mode: C++; c-indent-level: 4; c-basic-offset: 4; indent-tabs-mode: nil; -*-
//
// algParam: A class that holds all of the algorithm parameters that can be adapted
// adaptMethods: A class that adapts the algorithm parameters
//
// Copyright (C) 2017 Dirk Eddelbuettel, Adam Johansen and Leah South
//
Expand All @@ -21,41 +21,34 @@


//! \file
//! \brief A base class containing algorithm parameters and virtual
//! functions to adapt them.
//! \brief A base class with virtual functions to adapt parameters
//!

#ifndef __SMC_ALGPARAM_H
#define __SMC_ALGPARAM_H 1.0
#ifndef __SMC_ADAPTMETHODS_H
#define __SMC_ADAPTMETHODS_H 1.0

#include <population.h>


namespace smc {

/// A base class which contains the algorithm parameters and virtual functions to adapt them.
template <class Space, class Params> class algParam {

protected:
Params param;

template <class Space, class Params> class adaptMethods {

public:

/// Free the workspace allocated for the algorithm parameters.
virtual ~algParam() {
virtual ~adaptMethods() {
}

/// Holder function for updates to be done before the move step.
virtual void updateForMove(const population<Space> & pop) {}
virtual void updateForMove(Params &, const population<Space> & pop) {}

/// Holder function for updates to be done before the MCMC step.
virtual void updateForMCMC(const population<Space> & pop, double acceptProb, int nResampled, int & nRepeats) {}
virtual void updateForMCMC(Params &, const population<Space> & pop, double acceptProb, int nResampled, int & nRepeats) {}

/// Holder function for updates to be done at the end of each iteration.
virtual void updateEnd(const population<Space> & pop) {}

/// Returns the parameters.
const Params & GetParams(void) const {return param;}
virtual void updateEnd(Params &, const population<Space> & pop) {}
};
}

Expand Down
6 changes: 3 additions & 3 deletions inst/include/blockpfgaussianopt.h
Expand Up @@ -22,7 +22,7 @@

#include "smctc.h"

namespace BSPFG {
void fInitialise(arma::vec & value, double & logweight);
void fMove(long lTime, arma::vec & value, double & logweight);
namespace BSPFG {
void fInitialise(arma::vec & value, double & logweight, smc::nullParams & param);
void fMove(long lTime, arma::vec & value, double & logweight, smc::nullParams & param);
}
92 changes: 46 additions & 46 deletions inst/include/moveset.h
Expand Up @@ -35,68 +35,68 @@
namespace smc {

/// A template class for a set of moves for use in an SMC samplers framework.
template <class Space> class moveset
template <class Space, class Params> class moveset
{
private:
///The number of moves which are present in the set
long number;
///The function which initialises a single particle.
void (*pfInitialise)(Space &, double &);
void (*pfInitialise)(Space &, double &, Params &);
///The function which selects a move for a given particle at a given time.
long (*pfMoveSelect)(long , const Space &, const double &);
long (*pfMoveSelect)(long , const Space &, const double &, Params &);
///The functions which perform actual moves on a single particle.
void (**pfMoves)(long, Space &, double &);
void (**pfMoves)(long, Space &, double &, Params &);
///One iteration of a Markov Chain Monte Carlo move for a single particle.
bool (*pfMCMC)(long, Space &,double &);
bool (*pfMCMC)(long, Space &,double &, Params &);

public:
///Create a completely unspecified moveset
moveset();
///Create a reduced moveset with a single move
moveset(void (*pfInit)(Space &, double &),
void (*pfNewMoves)(long, Space &,double &),
bool (*pfNewMCMC)(long,Space &,double &));
moveset(void (*pfInit)(Space &, double &, Params &),
void (*pfNewMoves)(long, Space &,double &, Params &),
bool (*pfNewMCMC)(long,Space &,double &, Params &));
///Create a fully specified moveset
moveset(void (*pfInit)(Space &, double &),long (*pfMoveSelector)(long , const Space &, const double &),
long nMoves, void (**pfNewMoves)(long, Space &,double &),
bool (*pfNewMCMC)(long,Space &, double &));
moveset(void (*pfInit)(Space &, double &, Params &),long (*pfMoveSelector)(long , const Space &, const double &, Params &),
long nMoves, void (**pfNewMoves)(long, Space &,double &, Params &),
bool (*pfNewMCMC)(long,Space &, double &, Params &));

///Initialise the population of particles
void DoInit(population<Space> & pFrom, long N);
void DoInit(population<Space> & pFrom, long N, Params &);
///Perform an MCMC move on the particles
bool DoMCMC(long lTime, population<Space> & pFrom, long N, int nRepeats, int & nAccepted);
bool DoMCMC(long lTime, population<Space> & pFrom, long N, int nRepeats, int & nAccepted, Params &);
///Select an appropriate move at time lTime and apply it to pFrom
void DoMove(long lTime, population<Space> & pFrom,long N);
void DoMove(long lTime, population<Space> & pFrom,long N, Params &);

///Free the memory used for the array of move pointers when deleting
~moveset();

/// \brief Set the initialisation function.
/// \param pfInit is a function which returns a particle generated according to the initial distribution
void SetInitialisor( void (*pfInit)(Space &, double &))
void SetInitialisor( void (*pfInit)(Space &, double &, Params &))
{pfInitialise = pfInit;}

/// \brief Set the MCMC function
/// \param pfNewMCMC The function which performs an MCMC move
void SetMCMCFunction(bool (*pfNewMCMC)(long,Space &,double &)) {pfMCMC = pfNewMCMC;}
void SetMCMCFunction(bool (*pfNewMCMC)(long,Space &,double &, Params &)) {pfMCMC = pfNewMCMC;}

/// \brief Set the move selection function
/// \param pfMoveSelectNew returns the index of move to perform at the specified time given the specified particle
void SetMoveSelectionFunction(long (*pfMoveSelectNew)(long , const Space &, const double &))
void SetMoveSelectionFunction(long (*pfMoveSelectNew)(long , const Space &, const double &, Params &))
{pfMoveSelect = pfMoveSelectNew;};

///Set the individual move functions to the supplied array of such functions
void SetMoveFunctions(long nMoves, void (**pfNewMoves)(long, Space &, double &));
void SetMoveFunctions(long nMoves, void (**pfNewMoves)(long, Space &, double &, Params &));

///Moveset assignment should allocate buffers and deep copy all members.
moveset<Space> & operator= (moveset<Space> & pFrom);
moveset<Space,Params> & operator= (moveset<Space,Params> & pFrom);
};


/// The argument free smc::moveset constructor simply sets the number of available moves to zero and sets
/// all of the associated function pointers to NULL.
template <class Space>
moveset<Space>::moveset()
template <class Space, class Params>
moveset<Space,Params>::moveset()
{
number = 0;

Expand All @@ -111,10 +111,10 @@ namespace smc {
/// \param pfInit The function which should be used to initialise a particle when the system is initialised
/// \param pfNewMoves An functions which moves particles at a specified time to a new location
/// \param pfNewMCMC The function which should be called to apply an MCMC move (if any)
template <class Space>
moveset<Space>::moveset(void (*pfInit)(Space &, double &),
void (*pfNewMoves)(long, Space &,double &),
bool (*pfNewMCMC)(long,Space &,double &))
template <class Space, class Params>
moveset<Space,Params>::moveset(void (*pfInit)(Space &, double &, Params &),
void (*pfNewMoves)(long, Space &,double &, Params &),
bool (*pfNewMCMC)(long,Space &,double &, Params &))
{
SetInitialisor(pfInit);
SetMoveSelectionFunction(NULL);
Expand All @@ -130,10 +130,10 @@ namespace smc {
/// \param nMoves The number of moves which are defined in general
/// \param pfNewMoves An array of functions which moves particles at a specified time to a new location
/// \param pfNewMCMC The function which should be called to apply an MCMC move (if any)
template <class Space>
moveset<Space>::moveset(void (*pfInit)(Space &, double &),long (*pfMoveSelector)(long ,const Space &, const double &),
long nMoves, void (**pfNewMoves)(long, Space &,double &),
bool (*pfNewMCMC)(long,Space &,double &))
template <class Space, class Params>
moveset<Space,Params>::moveset(void (*pfInit)(Space &, double &, Params &),long (*pfMoveSelector)(long ,const Space &, const double &, Params &),
long nMoves, void (**pfNewMoves)(long, Space &,double &, Params &),
bool (*pfNewMCMC)(long,Space &,double &, Params &))
{
SetInitialisor(pfInit);
SetMoveSelectionFunction(pfMoveSelector);
Expand All @@ -142,28 +142,28 @@ namespace smc {
SetMCMCFunction(pfNewMCMC);
}

template <class Space>
moveset<Space>::~moveset()
template <class Space, class Params>
moveset<Space,Params>::~moveset()
{ if(pfMoves)
delete [] pfMoves;
}


template <class Space>
void moveset<Space>::DoInit(population<Space> & pFrom, long N) {
template <class Space, class Params>
void moveset<Space,Params>::DoInit(population<Space> & pFrom, long N, Params & params) {
for (long i=0; i<N; i++){
(*pfInitialise)(pFrom.GetValueRefN(i),pFrom.GetLogWeightRefN(i));
(*pfInitialise)(pFrom.GetValueRefN(i),pFrom.GetLogWeightRefN(i),params);
}
}

template <class Space>
bool moveset<Space>::DoMCMC(long lTime, population<Space> & pFrom, long N, int nRepeats, int & nAccepted)
template <class Space, class Params>
bool moveset<Space,Params>::DoMCMC(long lTime, population<Space> & pFrom, long N, int nRepeats, int & nAccepted, Params & params)
{
if(pfMCMC && nRepeats>0) {
nAccepted = 0;
for (int j=0; j<nRepeats; j++){
for (long i=0; i<N; i++){
nAccepted += pfMCMC(lTime,pFrom.GetValueRefN(i),pFrom.GetLogWeightRefN(i));
nAccepted += pfMCMC(lTime,pFrom.GetValueRefN(i),pFrom.GetLogWeightRefN(i),params);
}
}
return TRUE;
Expand All @@ -174,17 +174,17 @@ namespace smc {
}
}

template <class Space>
void moveset<Space>::DoMove(long lTime, population<Space> & pFrom, long N)
template <class Space, class Params>
void moveset<Space,Params>::DoMove(long lTime, population<Space> & pFrom, long N, Params & params)
{
if(number > 1)
for (long i=0; i<N; i++){
pfMoves[pfMoveSelect(lTime,pFrom.GetValueRefN(i),pFrom.GetLogWeightRefN(i))](lTime,pFrom.GetValueRefN(i),pFrom.GetLogWeightRefN(i));
pfMoves[pfMoveSelect(lTime,pFrom.GetValueRefN(i),pFrom.GetLogWeightRefN(i),params)](lTime,pFrom.GetValueRefN(i),pFrom.GetLogWeightRefN(i),params);
}
else

for (long i=0; i<N; i++){
pfMoves[0](lTime,pFrom.GetValueRefN(i),pFrom.GetLogWeightRefN(i));
pfMoves[0](lTime,pFrom.GetValueRefN(i),pFrom.GetLogWeightRefN(i),params);
}
}

Expand All @@ -194,20 +194,20 @@ namespace smc {
/// The move functions accept two arguments, the first of which corresponds to the system evolution time and the
/// second to an initial particle position and the second to a weighted starting position. It returns a new
/// weighted position corresponding to the moved particle.
template <class Space>
void moveset<Space>::SetMoveFunctions(long nMoves,void (**pfNewMoves)(long,Space &,double &))
template <class Space, class Params>
void moveset<Space,Params>::SetMoveFunctions(long nMoves,void (**pfNewMoves)(long,Space &,double &, Params &))
{
number = nMoves;
if(pfMoves)
delete [] pfMoves;
pfMoves = (void (**)(long int, Space &, double &)) new void* [nMoves];
pfMoves = (void (**)(long int, Space &, double &, Params &)) new void* [nMoves];
for(int i = 0; i < nMoves; i++)
pfMoves[i] = pfNewMoves[i];
return;
}

template <class Space>
moveset<Space> & moveset<Space>::operator= (moveset<Space> & pFrom)
template <class Space, class Params>
moveset<Space,Params> & moveset<Space,Params>::operator= (moveset<Space,Params> & pFrom)
{
SetInitialisor(pFrom.pfInitialise);
SetMCMCFunction(pFrom.pfMCMC);
Expand Down