Modelling.h

Go to the documentation of this file.

/********************************************************************
 *  Copyright (C) 2016 by Federico Marulli and Alfonso Veropalumbo  *
 *  federico.marulli3@unibo.it                                      *
 *                                                                  *
 *  This program is free software; you can redistribute it and/or   *
 *  modify it under the terms of the GNU General Public License as  *
 *  published by the Free Software Foundation; either version 2 of  * 
 *  the License, or (at your option) any later version.             *
 *                                                                  *
 *  This program is distributed in the hope that it will be useful, *
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of  *
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the   *
 *  GNU General Public License for more details.                    *
 *                                                                  *
 *  You should have received a copy of the GNU General Public       *
 *  License along with this program; if not, write to the Free      *
 *  Software Foundation, Inc.,                                      *
 *  59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.       *
 ********************************************************************/
 
#ifndef __MODELLING__
#define __MODELLING__
 
 
#include "Posterior.h"
 
 
// ===================================================================================================
 
 
namespace cbl {
 
  namespace modelling {
 
    class Modelling {
 
    protected:
 
      std::shared_ptr<data::Data> m_data = NULL;
 
      bool m_fit_range = false;
 
      std::shared_ptr<data::Data> m_data_fit;
 
      std::shared_ptr<statistics::Model> m_model = NULL;
      
      std::shared_ptr<statistics::Model> m_response_func = NULL;
 
      std::shared_ptr<statistics::Likelihood> m_likelihood = NULL;
 
      std::vector<std::shared_ptr<statistics::PriorDistribution>> m_parameter_priors;
 
      std::shared_ptr<statistics::Posterior> m_posterior = NULL;
 
      void m_set_prior (std::vector<statistics::PriorDistribution> prior_distribution);
      
      void m_isSet_response () {if (m_response_func == NULL) ErrorCBL("the response function for the SSC is not set!", "m_isSet_response", "Modelling.h");};
 
    public:
 
      
      Modelling () = default;
 
      virtual ~Modelling () = default;
 
      
      void m_set_posterior (const int seed);
      
      
      std::shared_ptr<data::Data> data () { return m_data; }
 
      std::shared_ptr<data::Data> data_fit () 
      { 
    if (!m_fit_range) 
      ErrorCBL("no fit range has been set!", "data_fit", "Modelling.h");
 
    return m_data_fit; 
      }
      
      std::shared_ptr<statistics::Likelihood> likelihood ();
 
      std::shared_ptr<statistics::Posterior> posterior ();
 
      std::shared_ptr<statistics::ModelParameters> likelihood_parameters ();
 
      std::shared_ptr<statistics::ModelParameters> posterior_parameters ();
      
      virtual void set_parameter_from_string (const std::string parameter, const double value)
      { (void)parameter; (void)value; ErrorCBL("", "set_parameter_from_string", "Modelling.h"); }
      
      virtual double get_parameter_from_string (const std::string parameter) const
      { return ErrorCBL("", "get_parameter_from_string ("+parameter+")", "Modelling.h"); }
      
      std::shared_ptr<statistics::PriorDistribution> get_prior (const int i) { return m_parameter_priors[i]; }
 
      std::shared_ptr<statistics::Model> get_response_function () {m_isSet_response(); return m_response_func;};
 
      
 
      
      void reset_fit_range () { m_fit_range = false; }
 
      void set_fit_range (const double xmin, const double xmax);
      
      void set_fit_range (const double xmin, const double xmax, const double ymin, const double ymax);
 
      void set_data (const std::shared_ptr<data::Data> dataset) { m_data = move(dataset); }
 
      void set_likelihood (const statistics::LikelihoodType likelihood_type, const std::vector<size_t> x_index={0, 2}, const int w_index=-1, const double prec=1.e-10, const int Nres=-1);
      
      void set_likelihood (const cbl::statistics::Likelihood_function log_likelihood_function);
 
      
 
      
      void maximize_likelihood (const std::vector<double> start, const std::vector<std::vector<double>> parameter_limits, const unsigned int max_iter=10000, const double tol=1.e-6, const double epsilon=1.e-3);
 
      void maximize_posterior (const std::vector<double> start, const unsigned int max_iter=10000, const double tol=1.e-6, const double epsilon=1.e-3, const int seed=666);
 
      void sample_posterior (const int chain_size, const int nwalkers, const int seed=666, const double aa=2, const bool parallel=true);
 
      void sample_posterior (const int chain_size, const int nwalkers, const double radius, const std::vector<double> start, const unsigned int max_iter=10000, const double tol=1.e-6, const double epsilon=1.e-3, const int seed=666, const double aa=2, const bool parallel=true);
 
      void sample_posterior (const int chain_size, const int nwalkers, std::vector<double> &value, const double radius, const int seed=666, const double aa=2, const bool parallel=true);
 
      void sample_posterior (const int chain_size, const std::vector<std::vector<double>> chain_value, const int seed=666, const double aa=2, const bool parallel=true);
 
      void sample_posterior (const int chain_size, const int nwalkers, const std::string input_dir, const std::string input_file, const int seed=666, const double aa=2, const bool parallel=true);
 
      void importance_sampling (const std::string input_dir, const std::string input_file, const int seed=666, const std::vector<size_t> column={}, const int header_lines_to_skip=1, const bool is_FITS_format=false, const bool apply_to_likelihood=false);
 
 
      
      
      void write_chain (const std::string output_dir, const std::string output_file, const int start=0, const int thin=1, const bool is_FITS_format=false, const int prec=5, const int ww=14);
 
      void read_chain (const std::string input_dir, const std::string input_file, const int nwalkers, const std::vector<size_t> columns={}, const int skip_header=1, const bool fits=false);
 
      void show_results (const int start=0, const int thin=1, const int nbins=50, const bool show_mode=false, const int ns=-1);
 
      void write_results (const std::string output_dir, const std::string root_file, const int start=0, const int thin=1, const int nbins=50, const bool fits=false, const bool compute_mode=false, const int ns=-1);
 
      virtual void write_model (const std::string output_dir, const std::string output_file, const std::vector<double> xx, const std::vector<double> parameters);
 
      virtual void write_model (const std::string output_dir, const std::string output_file, const std::vector<double> xx, const std::vector<double> yy, const std::vector<double> parameters);
 
      virtual void write_model_at_bestfit (const std::string output_dir, const std::string output_file, const std::vector<double> xx);
 
      virtual void write_model_at_bestfit (const std::string output_dir, const std::string output_file, const std::vector<double> xx, const std::vector<double> yy);
 
      virtual void write_model_from_chains (const std::string output_dir, const std::string output_file, const std::vector<double> xx, const int start=0, const int thin=1);
 
      virtual void write_model_from_chains (const std::string output_dir, const std::string output_file, const std::vector<double> xx, const std::vector<double> yy, const int start=0, const int thin=1);
 
      double reduced_chi2 (const std::vector<double> parameter={});
      
      
    };
  }
}
 
#endif