28 const std::vector<double>* parameters,
const std::vector<double>* fixed) {
51 if(!adapt_manager[
"AdaptionOptions"][
"Covariance"][matrix_name_str]) {
52 MACH3LOG_WARN(
"Adaptive Settings not found for {}, this is fine if you don't want adaptive MCMC", matrix_name_str);
57 if(!GetFromManager<bool>(adapt_manager[
"AdaptionOptions"][
"Covariance"][matrix_name_str][
"DoAdaption"],
false)) {
62 if(!
CheckNodeExists(adapt_manager,
"AdaptionOptions",
"Settings",
"OutputFileName")) {
63 MACH3LOG_ERROR(
"No OutputFileName specified in AdaptionOptions::Settings into your config file");
68 start_adaptive_throw = GetFromManager<int>(adapt_manager[
"AdaptionOptions"][
"Settings"][
"StartThrow"], 10);
69 start_adaptive_update = GetFromManager<int>(adapt_manager[
"AdaptionOptions"][
"Settings"][
"StartUpdate"], 0);
70 end_adaptive_update = GetFromManager<int>(adapt_manager[
"AdaptionOptions"][
"Settings"][
"EndUpdate"], 10000);
71 adaptive_update_step = GetFromManager<int>(adapt_manager[
"AdaptionOptions"][
"Settings"][
"UpdateStep"], 100);
73 output_file_name = GetFromManager<std::string>(adapt_manager[
"AdaptionOptions"][
"Settings"][
"OutputFileName"],
"");
84 auto matrix_blocks = GetFromManager<std::vector<std::vector<int>>>(adapt_manager[
"AdaptionOptions"][
"Covariance"][matrix_name_str][
"MatrixBlocks"], {{}});
114 if(block_indices.size()==0 || block_indices[0].size()==0)
return;
116 int block_size =
static_cast<int>(block_indices.size());
118 for(
int iblock=0; iblock < block_size; iblock++){
120 int sub_block_size =
static_cast<int>(block_indices.size()-1);
121 for(
int isubblock=0; isubblock < sub_block_size ; isubblock+=2){
122 int block_lb = block_indices[iblock][isubblock];
123 int block_ub = block_indices[iblock][isubblock+1];
126 MACH3LOG_ERROR(
"Cannot set matrix block with edges {}, {} for matrix of size {}",
130 for(
int ipar = block_lb; ipar < block_ub; ipar++){
142 const std::string &systematicName,
bool is_final) {
146 TFile *outFile =
new TFile(outFileName.c_str(),
"UPDATE");
147 if (outFile->IsZombie()) {
152 TVectorD *outMeanVec =
new TVectorD(
int(
par_means.size()));
153 for (
int i = 0; i < int(
par_means.size()); i++) {
157 std::string adaptive_cov_name = systematicName +
"_posfit_matrix";
158 std::string mean_vec_name = systematicName +
"_mean_vec";
161 std::string total_steps_str =
163 std::string syst_name_str = systematicName;
165 total_steps_str +
'_' + syst_name_str + std::string(
"_throw_matrix");
167 total_steps_str +
'_' + syst_name_str + std::string(
"_mean_vec");
172 outMeanVec->Write(mean_vec_name.c_str());
183 const std::string& matrix_name,
184 const std::string& means_name,
185 bool& use_adaptive) {
189 std::unique_ptr<TFile>matrix_file(
new TFile(matrix_file_name.c_str()));
192 if(matrix_file->IsZombie()){
200 MACH3LOG_ERROR(
"Couldn't find {} in {}", matrix_name, matrix_file_name);
205 TVectorD* means_vector =
static_cast<TVectorD*
>(matrix_file->Get(means_name.c_str()));
211 if(means_vector->GetNrows()){
220 MACH3LOG_INFO(
"Found Means in External File, Will be able to adapt");
225 MACH3LOG_WARN(
"Cannot find means vector in {}, therefore I will not be able to adapt!", matrix_file_name);
226 use_adaptive =
false;
229 matrix_file->Close();
236 std::vector<double> par_means_prev =
par_means;
243 par_means[i] = (
CurrVal(i) + par_means_prev[i]*steps_post_burn)/(steps_post_burn+1);
249 #pragma omp parallel for
253 (*adaptive_covariance)(i, i) = 1.0;
259 for (
int j = 0; j <= i; ++j) {
261 (*adaptive_covariance)(i, j) = 0.0;
262 (*adaptive_covariance)(j, i) = 0.0;
266 double cov_prev = (*adaptive_covariance)(i, j);
267 double cov_updated = 0.0;
269 if (steps_post_burn > 0) {
271 double cov_t = cov_prev * (steps_post_burn - 1) / steps_post_burn;
272 double prev_means_t = steps_post_burn * par_means_prev[i] * par_means_prev[j];
276 cov_updated = cov_t +
adaption_scale * (prev_means_t - curr_means_t + curr_step_t) / steps_post_burn;
279 (*adaptive_covariance)(i, j) = cov_updated;
280 (*adaptive_covariance)(j, i) = cov_updated;
bool CheckNodeExists(const YAML::Node &node, Args... args)
KS: Wrapper function to call the recursive helper.
Custom exception class for MaCh3 errors.
std::vector< double > par_means
Mean values for all parameters.
bool IsFixed(const int ipar) const
Check if a parameter is fixed.
bool UpdateMatrixAdapt()
Tell whether matrix should be updated.
std::vector< int > adapt_block_matrix_indices
Indices for block-matrix adaption.
double CurrVal(const int par_index)
Get Current value of parameter.
bool InitFromConfig(const YAML::Node &adapt_manager, const std::string &matrix_name_str, const std::vector< double > *parameters, const std::vector< double > *fixed)
Read initial values from config file.
std::string output_file_name
Name of the file to save the adaptive matrices into.
bool AdaptionUpdate()
To be fair not a clue...
double adaption_scale
Scaling factor.
void SaveAdaptiveToFile(const std::string &outFileName, const std::string &systematicName, const bool is_final=false)
HW: Save adaptive throw matrix to file.
void SetAdaptiveBlocks(const std::vector< std::vector< int > > &block_indices)
HW: sets adaptive block matrix.
const std::vector< double > * _fCurrVal
Current values of parameters.
void Print()
Print all class members.
int start_adaptive_update
When do we stop update the adaptive matrix.
std::vector< int > adapt_block_sizes
Size of blocks for adaption.
virtual ~AdaptiveMCMCHandler()
Destructor.
bool IndivStepScaleAdapt()
Tell whether we want reset step scale or not.
int adaptive_save_n_iterations
void SetFixed(const std::vector< double > *fix)
Set the fixed parameters.
void SetParams(const std::vector< double > *params)
Set the current values of the parameters.
AdaptiveMCMCHandler()
Constructor.
void CreateNewAdaptiveCovariance()
If we don't have a covariance matrix to start from for adaptive tune we need to make one!
int end_adaptive_update
Steps between changing throw matrix.
int total_steps
Total number of MCMC steps.
void UpdateAdaptiveCovariance()
Method to update adaptive MCMC .
int adaptive_update_step
Steps between changing throw matrix.
void SetThrowMatrixFromFile(const std::string &matrix_file_name, const std::string &matrix_name, const std::string &means_name, bool &use_adaptive)
sets throw matrix from a file
int GetNumParams() const
Get the current values of the parameters.
bool SkipAdaption()
Tell if we are Skipping Adaption.
TMatrixDSym * adaptive_covariance
Full adaptive covariance matrix.