rcpp_distribution.hpp

Go to the documentation of this file.

#ifndef FIMS_INTERFACE_RCPP_RCPP_OBJECTS_RCPP_DISTRIBUTION_HPP
#define FIMS_INTERFACE_RCPP_RCPP_OBJECTS_RCPP_DISTRIBUTION_HPP
 
#include <sstream>
#include <string>
#include <utility>
#include <vector>
 
#include "rcpp_interface_base.hpp"
#include "../../../distributions/distributions.hpp"

class DistributionsInterfaceBase : public FIMSRcppInterfaceBase {
 public:
  static uint32_t id_g;
  uint32_t id_m;
  std::shared_ptr<std::vector<uint32_t>> key_m;
  SharedString input_type_m;
  SharedString use_mean_m = fims::to_string("no");
  static std::map<uint32_t, std::shared_ptr<DistributionsInterfaceBase>>
      live_objects;
  SharedInt interface_observed_data_id_m = -999;

  double lpdf_value = 0;
  DistributionsInterfaceBase() {
    this->key_m = std::make_shared<std::vector<uint32_t>>();
    this->id_m = DistributionsInterfaceBase::id_g++;
    /* Create instance of map: key is id and value is pointer to
    DistributionsInterfaceBase */
    // DistributionsInterfaceBase::live_objects[this->id_m] = this;
  }

  DistributionsInterfaceBase(const DistributionsInterfaceBase &other)
      : id_m(other.id_m),
        key_m(other.key_m),
        input_type_m(other.input_type_m),
        use_mean_m(other.use_mean_m),
        interface_observed_data_id_m(other.interface_observed_data_id_m) {}

  virtual ~DistributionsInterfaceBase() {}

  virtual uint32_t get_id() = 0;

  virtual bool set_distribution_links(std::string input_type,
                                      Rcpp::IntegerVector ids) {
    return false;
  }

  virtual bool set_distribution_mean(double input_value) { return false; }

  virtual bool set_observed_data(int observed_data_id) { return false; }

  virtual double evaluate() = 0;

  std::string parameters_to_json(
      const std::vector<std::pair<std::string, VariableVector *>> &parameters) {
    std::vector<std::string> entries;
    for (const auto &parameter : parameters) {
      VariableVector *values = parameter.second;
      if (values == nullptr || values->size() == 0) {
        continue;
      }
      std::stringstream entry;
      entry << "{\n";
      entry << "   \"name\": \"" << parameter.first << "\",\n";
      entry << "   \"id\": " << values->id_m << ",\n";
      entry << "   \"type\": \"vector\",\n";
      entry << "   \"dimensionality\": {\n";
      entry << "    \"header\": [null],\n";
      entry << "    \"dimensions\": [" << values->size() << "]\n";
      entry << "   },\n";
      entry << "   \"values\": " << *values << "\n";
      entry << "  }";
      entries.push_back(entry.str());
    }
    if (entries.empty()) {
      return "";
    }
    std::stringstream ss;
    ss << " \"parameters\": [\n";
    for (size_t i = 0; i < entries.size(); i++) {
      ss << (i == 0 ? "  " : ",\n  ") << entries[i];
    }
    ss << "\n ],\n";
    return ss.str();
  }

  std::string linked_ids_to_json() {
    std::stringstream ss;
    ss << " \"linked_ids\": [";
    for (size_t i = 0; i < this->key_m->size(); i++) {
      ss << (i == 0 ? "" : ", ") << this->key_m->at(i);
    }
    ss << "],\n";
    return ss.str();
  }
};

class DnormDistributionsInterface : public DistributionsInterfaceBase {
 public:
  VariableVector observed_values;
  VariableVector expected_values;
  VariableVector expected_mean;
  VariableVector log_sd;
  RealVector lpdf_vec; 

  DnormDistributionsInterface() : DistributionsInterfaceBase() {
    DistributionsInterfaceBase::live_objects[this->id_m] =
        std::make_shared<DnormDistributionsInterface>(*this);
    FIMSRcppInterfaceBase::fims_interface_objects.push_back(
        DistributionsInterfaceBase::live_objects[this->id_m]);
  }

  DnormDistributionsInterface(const DnormDistributionsInterface &other)
      : DistributionsInterfaceBase(other),
        observed_values(other.observed_values),
        expected_values(other.expected_values),
        expected_mean(other.expected_mean),
        log_sd(other.log_sd),
        lpdf_vec(other.lpdf_vec) {}

  virtual ~DnormDistributionsInterface() {}

  virtual uint32_t get_id() { return this->id_m; }

  virtual bool set_observed_data(int observed_data_id) {
    this->interface_observed_data_id_m.set(observed_data_id);
    return true;
  }

  virtual bool set_distribution_mean(double input_value) {
    this->expected_mean[0].initial_value_m = input_value;
    this->expected_mean[0].estimation_type_m.set("fixed_effects");
    this->use_mean_m.set(fims::to_string("yes"));
    return true;
  }

  virtual bool set_distribution_links(std::string input_type,
                                      Rcpp::IntegerVector ids) {
    this->input_type_m.set(input_type);
    this->key_m->resize(ids.size());
    for (R_xlen_t i = 0; i < ids.size(); i++) {
      this->key_m->at(i) = ids[i];
    }
    return true;
  }

  virtual double evaluate() {
    fims_distributions::NormalLPDF<double> dnorm;
    dnorm.observed_values.resize(this->observed_values.size());
    dnorm.expected_values.resize(this->expected_values.size());
    dnorm.log_sd.resize(this->log_sd.size());
    dnorm.expected_mean.resize(this->expected_mean.size());
    for (size_t i = 0; i < this->observed_values.size(); i++) {
      dnorm.observed_values[i] = this->observed_values[i].initial_value_m;
    }
    for (size_t i = 0; i < this->expected_values.size(); i++) {
      dnorm.expected_values[i] = this->expected_values[i].initial_value_m;
    }
    for (size_t i = 0; i < this->log_sd.size(); i++) {
      dnorm.log_sd[i] = this->log_sd[i].initial_value_m;
    }
    for (size_t i = 0; i < this->expected_mean.size(); i++) {
      dnorm.expected_mean[i] = this->expected_mean[i].initial_value_m;
    }
    dnorm.use_mean = this->use_mean_m;
    return dnorm.evaluate();
  }

  virtual void finalize() {
    if (this->finalized) {
      // log warning that finalize has been called more than once.
      FIMS_WARNING_LOG("DnormDistribution  " + fims::to_string(this->id_m) +
                       " has been finalized already.");
    }
 
    this->finalized = true;  // indicate this has been called already
 
    std::shared_ptr<fims_info::Information<double>> info =
        fims_info::Information<double>::GetInstance();
 
    fims_info::Information<double>::density_components_iterator it;
 
    // search for density component in Information
    it = info->density_components.find(this->id_m);
    // if not found, just return
    if (it == info->density_components.end()) {
      FIMS_WARNING_LOG("DnormDistribution " + fims::to_string(this->id_m) +
                       " not found in Information.");
      return;
    } else {
      std::shared_ptr<fims_distributions::NormalLPDF<double>> dnorm =
          std::dynamic_pointer_cast<fims_distributions::NormalLPDF<double>>(
              it->second);
 
      this->lpdf_value = dnorm->lpdf;
 
      size_t n_x = dnorm->get_n_x();
 
      // If input log_sd is a scalar, resize to n_x and fill with the scalar
      // value
      if (this->log_sd.size() != n_x) {
        // If log_sd size == 1 (scalar), repeat the entry
        if (this->log_sd.size() == 1) {
          auto tmp = this->log_sd[0];  // copy the one log_sd param
          this->log_sd.resize(n_x);
          for (size_t i = 0; i < n_x; ++i) {
            this->log_sd[i] = tmp;  // copies all fields in Param
          }
        } else {
          // Handle error
          FIMS_WARNING_LOG(
              "log_sd size does not match number of observations and is not "
              "scalar.");
        }
      }
      for (size_t i = 0; i < n_x; i++) {
        if (this->log_sd[i].estimation_type_m.get() == "constant") {
          this->log_sd[i].final_value_m = this->log_sd[i].initial_value_m;
        } else {
          this->log_sd[i].final_value_m = dnorm->log_sd.get_force_scalar(i);
        }
      }
 
      for (size_t i = 0; i < this->expected_mean.size(); i++) {
        if (this->expected_mean[i].estimation_type_m.get() == "constant") {
          this->expected_mean[i].final_value_m =
              this->expected_mean[i].initial_value_m;
        } else {
          this->expected_mean[i].final_value_m = dnorm->expected_mean[i];
        }
      }
 
      this->lpdf_vec = RealVector(n_x);
      if (this->expected_values.size() == 1) {
        this->expected_values.resize(n_x);
      }
      if (this->observed_values.size() == 1) {
        this->observed_values.resize(n_x);
      }
 
      for (size_t i = 0; i < this->lpdf_vec.size(); i++) {
        this->lpdf_vec[i] = dnorm->lpdf_vec[i];
        this->expected_values[i].final_value_m = dnorm->get_expected(i);
        this->observed_values[i].final_value_m = dnorm->get_observed(i);
      }
    }
  }

  virtual std::string to_json() {
    std::stringstream ss;
 
    ss << "{\n";
    ss << " \"module_name\": \"density\",\n";
    ss << " \"module_id\": " << this->id_m << ",\n";
    ss << " \"module_type\": \"normal\",\n";
    ss << " \"observed_data_id\" : " << this->interface_observed_data_id_m
       << ",\n";
    ss << " \"input_type\" : \"" << this->input_type_m << "\",\n";
    ss << this->linked_ids_to_json();
    ss << this->parameters_to_json({{"log_sd", &this->log_sd}});
    ss << " \"density_component\": {\n";
    ss << "  \"lpdf_value\": " << sanitize_val(this->lpdf_value) << ",\n";
    ss << "  \"value\":[";
    if (this->lpdf_vec.size() == 0) {
      ss << "],\n";
    } else {
      for (size_t i = 0; i < this->lpdf_vec.size() - 1; i++) {
        ss << this->value_to_string(this->lpdf_vec[i]);
        ss << ", ";
      }
      ss << this->value_to_string(this->lpdf_vec[this->lpdf_vec.size() - 1]);
 
      ss << "],\n";
    }
    ss << "  \"expected_values\":[";
    if (this->expected_values.size() == 0) {
      ss << "],\n";
    } else {
      for (size_t i = 0; i < this->expected_values.size() - 1; i++) {
        ss << this->value_to_string(this->expected_values[i].final_value_m)
           << ", ";
      }
      ss << this->value_to_string(
          this->expected_values[this->expected_values.size() - 1]
              .final_value_m);
      ss << "],\n";
    }
    ss << "  \"log_sd_values\":[";
    if (this->log_sd.size() == 0) {
      ss << "],\n";
    } else {
      for (R_xlen_t i = 0; i < static_cast<R_xlen_t>(this->log_sd.size()) - 1;
           i++) {
        ss << this->value_to_string(this->log_sd[i].final_value_m) << ", ";
      }
      ss << this->value_to_string(
                this->log_sd[this->log_sd.size() - 1].final_value_m)
         << "],\n";
    }
    ss << "  \"observed_values\":[";
    if (this->observed_values.size() == 0) {
      ss << "]\n";
    } else {
      for (size_t i = 0; i < this->observed_values.size() - 1; i++) {
        ss << this->observed_values[i].final_value_m << ", ";
      }
      ss << this->observed_values[this->observed_values.size() - 1]
                .final_value_m
         << "]\n";
    }
    ss << " }}\n";
    return ss.str();
  }
 
#ifdef TMB_MODEL
 
  template <typename Type>
  bool add_to_fims_tmb_internal() {
    std::shared_ptr<fims_info::Information<Type>> info =
        fims_info::Information<Type>::GetInstance();
 
    std::shared_ptr<fims_distributions::NormalLPDF<Type>> distribution =
        std::make_shared<fims_distributions::NormalLPDF<Type>>();
 
    // interface to data/parameter value
 
    distribution->observed_data_id_m = interface_observed_data_id_m;
    std::stringstream ss;
    distribution->input_type = this->input_type_m;
    distribution->key.resize(this->key_m->size());
    for (size_t i = 0; i < this->key_m->size(); i++) {
      distribution->key[i] = this->key_m->at(i);
    }
    distribution->id = this->id_m;
    distribution->observed_values.resize(this->observed_values.size());
    for (size_t i = 0; i < this->observed_values.size(); i++) {
      distribution->observed_values[i] =
          this->observed_values[i].initial_value_m;
    }
    // set relative info
    distribution->expected_values.resize(this->expected_values.size());
    for (size_t i = 0; i < this->expected_values.size(); i++) {
      distribution->expected_values[i] =
          this->expected_values[i].initial_value_m;
    }
    distribution->log_sd.resize(this->log_sd.size());
    for (size_t i = 0; i < this->log_sd.size(); i++) {
      distribution->log_sd[i] = this->log_sd[i].initial_value_m;
      if (this->log_sd[i].estimation_type_m.get() == "fixed_effects") {
        ss.str("");
        ss << "dnorm." << this->id_m << ".log_sd." << this->log_sd[i].id_m;
        info->RegisterParameterName(ss.str());
        info->RegisterParameter(distribution->log_sd[i]);
      }
      if (this->log_sd[i].estimation_type_m.get() == "random_effects") {
        FIMS_ERROR_LOG("standard deviations cannot be set to random effects");
      }
    }
    info->variable_map[this->log_sd.id_m] = &(distribution)->log_sd;
 
    distribution->use_mean = this->use_mean_m.get();
    distribution->expected_mean.resize(this->expected_mean.size());
    for (size_t i = 0; i < this->expected_mean.size(); i++) {
      distribution->expected_mean[i] = this->expected_mean[i].initial_value_m;
      if (this->expected_mean[i].estimation_type_m.get() == "fixed_effects") {
        ss.str("");
        ss << "dnorm." << this->id_m << ".expected_mean."
           << this->expected_mean[i].id_m;
        info->RegisterParameterName(ss.str());
        info->RegisterParameter(distribution->expected_mean[i]);
      }
      if (this->expected_mean[i].estimation_type_m.get() == "random_effects") {
        FIMS_ERROR_LOG("expected_mean cannot be set to random effects");
      }
    }
    info->variable_map[this->expected_mean.id_m] =
        &(distribution)->expected_mean;
 
    info->density_components[distribution->id] = distribution;
 
    return true;
  }

  virtual bool add_to_fims_tmb() {
    this->add_to_fims_tmb_internal<TMB_FIMS_REAL_TYPE>();
    this->add_to_fims_tmb_internal<TMBAD_FIMS_TYPE>();
 
    return true;
  }
 
#endif
};

class DlnormDistributionsInterface : public DistributionsInterfaceBase {
 public:
  VariableVector observed_values;
  VariableVector expected_values;
  VariableVector log_sd;
  RealVector lpdf_vec; 

  DlnormDistributionsInterface() : DistributionsInterfaceBase() {
    DistributionsInterfaceBase::live_objects[this->id_m] =
        std::make_shared<DlnormDistributionsInterface>(*this);
    FIMSRcppInterfaceBase::fims_interface_objects.push_back(
        DistributionsInterfaceBase::live_objects[this->id_m]);
  }

  DlnormDistributionsInterface(const DlnormDistributionsInterface &other)
      : DistributionsInterfaceBase(other),
        observed_values(other.observed_values),
        expected_values(other.expected_values),
        log_sd(other.log_sd),
        lpdf_vec(other.lpdf_vec) {}

  virtual ~DlnormDistributionsInterface() {}

  virtual uint32_t get_id() { return this->id_m; }

  virtual bool set_observed_data(int observed_data_id) {
    this->interface_observed_data_id_m.set(observed_data_id);
    return true;
  }

  virtual bool set_distribution_links(std::string input_type,
                                      Rcpp::IntegerVector ids) {
    this->input_type_m.set(input_type);
    this->key_m->resize(ids.size());
    for (R_xlen_t i = 0; i < ids.size(); i++) {
      this->key_m->at(i) = ids[i];
    }
    return true;
  }

  virtual double evaluate() {
    fims_distributions::LogNormalLPDF<double> dlnorm;
    dlnorm.observed_values.resize(this->observed_values.size());
    dlnorm.expected_values.resize(this->expected_values.size());
    dlnorm.log_sd.resize(this->log_sd.size());
    // dlnorm.input_type = "prior";
    for (size_t i = 0; i < this->observed_values.size(); i++) {
      dlnorm.observed_values[i] = this->observed_values[i].initial_value_m;
    }
    for (size_t i = 0; i < this->expected_values.size(); i++) {
      dlnorm.expected_values[i] = this->expected_values[i].initial_value_m;
    }
    for (size_t i = 0; i < this->log_sd.size(); i++) {
      dlnorm.log_sd[i] = this->log_sd[i].initial_value_m;
    }
    return dlnorm.evaluate();
  }

  virtual void finalize() {
    if (this->finalized) {
      // log warning that finalize has been called more than once.
      FIMS_WARNING_LOG("LogNormalLPDF  " + fims::to_string(this->id_m) +
                       " has been finalized already.");
    }
 
    this->finalized = true;  // indicate this has been called already
 
    std::shared_ptr<fims_info::Information<double>> info =
        fims_info::Information<double>::GetInstance();
 
    fims_info::Information<double>::density_components_iterator it;
 
    // search for density component in Information
    it = info->density_components.find(this->id_m);
    // if not found, just return
    if (it == info->density_components.end()) {
      FIMS_WARNING_LOG("LogNormalLPDF " + fims::to_string(this->id_m) +
                       " not found in Information.");
      return;
    } else {
      std::shared_ptr<fims_distributions::LogNormalLPDF<double>> dlnorm =
          std::dynamic_pointer_cast<fims_distributions::LogNormalLPDF<double>>(
              it->second);
 
      this->lpdf_value = dlnorm->lpdf;
 
      size_t n_x = dlnorm->get_n_x();
 
      if (this->log_sd.size() != n_x) {
        // If log_sd size == 1 (scalar), repeat the entry
        if (this->log_sd.size() == 1) {
          auto tmp = this->log_sd[0];  // copy the one log_sd param
          this->log_sd.resize(n_x);
          for (size_t i = 0; i < n_x; ++i) {
            this->log_sd[i] = tmp;  // copies all fields in Param
          }
        } else {
          // Handle error
          FIMS_WARNING_LOG(
              "log_sd size does not match number of observations and is not "
              "scalar.");
        }
      }
 
      for (size_t i = 0; i < n_x; i++) {
        if (this->log_sd[i].estimation_type_m.get() == "constant") {
          this->log_sd[i].final_value_m = this->log_sd[i].initial_value_m;
        } else {
          this->log_sd[i].final_value_m = dlnorm->log_sd.get_force_scalar(i);
        }
      }
 
      this->lpdf_vec = RealVector(n_x);
      if (this->expected_values.size() == 1) {
        this->expected_values.resize(n_x);
      }
      if (this->observed_values.size() == 1) {
        this->observed_values.resize(n_x);
      }
      for (size_t i = 0; i < this->lpdf_vec.size(); i++) {
        this->lpdf_vec[i] = dlnorm->lpdf_vec[i];
        this->expected_values[i].final_value_m = dlnorm->get_expected(i);
        this->observed_values[i].final_value_m = dlnorm->get_observed(i);
      }
    }
  }

  virtual std::string to_json() {
    std::stringstream ss;
 
    ss << "{\n";
    ss << " \"module_name\": \"density\",\n";
    ss << " \"module_id\": " << this->id_m << ",\n";
    ss << " \"module_type\": \"log_normal\",\n";
    ss << " \"observed_data_id\" : " << this->interface_observed_data_id_m
       << ",\n";
    ss << " \"input_type\" : \"" << this->input_type_m << "\",\n";
    ss << this->linked_ids_to_json();
    ss << this->parameters_to_json({{"log_sd", &this->log_sd}});
    ss << " \"density_component\": {\n";
    ss << "  \"lpdf_value\": " << sanitize_val(this->lpdf_value) << ",\n";
    ss << "  \"value\":[";
    if (this->lpdf_vec.size() == 0) {
      ss << "],\n";
    } else {
      for (size_t i = 0; i < this->lpdf_vec.size() - 1; i++) {
        ss << this->value_to_string(this->lpdf_vec[i]) << ", ";
      }
      ss << this->value_to_string(this->lpdf_vec[this->lpdf_vec.size() - 1]);
 
      ss << "],\n";
    }
    ss << "  \"expected_values\":[";
    if (this->expected_values.size() == 0) {
      ss << "],\n";
    } else {
      for (size_t i = 0; i < this->expected_values.size() - 1; i++) {
        ss << this->value_to_string(this->expected_values[i].final_value_m)
           << ", ";
      }
      ss << this->value_to_string(
          this->expected_values[this->expected_values.size() - 1]
              .final_value_m);
 
      ss << "],\n";
    }
    ss << "  \"log_sd_values\":[";
    if (this->log_sd.size() == 0) {
      ss << "],\n";
    } else {
      for (R_xlen_t i = 0; i < static_cast<R_xlen_t>(this->log_sd.size()) - 1;
           i++) {
        ss << this->value_to_string(this->log_sd[i].final_value_m) << ", ";
      }
      ss << this->value_to_string(
                this->log_sd[this->log_sd.size() - 1].final_value_m)
         << "],\n";
    }
    ss << "  \"observed_values\":[";
    if (this->observed_values.size() == 0) {
      ss << "]\n";
    } else {
      for (size_t i = 0; i < this->observed_values.size() - 1; i++) {
        ss << this->observed_values[i].final_value_m << ", ";
      }
      ss << this->observed_values[this->observed_values.size() - 1]
                .final_value_m
         << "]\n";
    }
    ss << " }}\n";
    return ss.str();
  }
 
#ifdef TMB_MODEL
 
  template <typename Type>
  bool add_to_fims_tmb_internal() {
    std::shared_ptr<fims_info::Information<Type>> info =
        fims_info::Information<Type>::GetInstance();
 
    std::shared_ptr<fims_distributions::LogNormalLPDF<Type>> distribution =
        std::make_shared<fims_distributions::LogNormalLPDF<Type>>();
 
    // set relative info
    distribution->id = this->id_m;
    std::stringstream ss;
    distribution->observed_data_id_m = interface_observed_data_id_m;
    distribution->input_type = this->input_type_m;
    distribution->key.resize(this->key_m->size());
    for (size_t i = 0; i < this->key_m->size(); i++) {
      distribution->key[i] = this->key_m->at(i);
    }
    distribution->observed_values.resize(this->observed_values.size());
    for (size_t i = 0; i < this->observed_values.size(); i++) {
      distribution->observed_values[i] =
          this->observed_values[i].initial_value_m;
    }
    // set relative info
    distribution->expected_values.resize(this->expected_values.size());
    for (size_t i = 0; i < this->expected_values.size(); i++) {
      distribution->expected_values[i] =
          this->expected_values[i].initial_value_m;
    }
    distribution->log_sd.resize(this->log_sd.size());
    for (size_t i = 0; i < this->log_sd.size(); i++) {
      distribution->log_sd[i] = this->log_sd[i].initial_value_m;
      if (this->log_sd[i].estimation_type_m.get() == "fixed_effects") {
        ss.str("");
        ss << "dlnorm." << this->id_m << ".log_sd." << this->log_sd[i].id_m;
        info->RegisterParameterName(ss.str());
        info->RegisterParameter(distribution->log_sd[i]);
      }
      if (this->log_sd[i].estimation_type_m.get() == "random_effects") {
        FIMS_ERROR_LOG("standard deviations cannot be set to random effects");
      }
    }
    info->variable_map[this->log_sd.id_m] = &(distribution)->log_sd;
 
    info->density_components[distribution->id] = distribution;
 
    return true;
  }

  virtual bool add_to_fims_tmb() {
    this->add_to_fims_tmb_internal<TMB_FIMS_REAL_TYPE>();
    this->add_to_fims_tmb_internal<TMBAD_FIMS_TYPE>();
 
    return true;
  }
 
#endif
};

class DmultinomDistributionsInterface : public DistributionsInterfaceBase {
 public:
  VariableVector observed_values;
  VariableVector expected_values;
  RealVector dims;
  RealVector lpdf_vec; 

  SharedString notes;

  DmultinomDistributionsInterface() : DistributionsInterfaceBase() {
    DistributionsInterfaceBase::live_objects[this->id_m] =
        std::make_shared<DmultinomDistributionsInterface>(*this);
    FIMSRcppInterfaceBase::fims_interface_objects.push_back(
        DistributionsInterfaceBase::live_objects[this->id_m]);
  }

  DmultinomDistributionsInterface(const DmultinomDistributionsInterface &other)
      : DistributionsInterfaceBase(other),
        observed_values(other.observed_values),
        expected_values(other.expected_values),
        dims(other.dims),
        lpdf_vec(other.lpdf_vec),
        notes(other.notes) {}

  virtual ~DmultinomDistributionsInterface() {}
  virtual uint32_t get_id() { return this->id_m; }

  virtual bool set_observed_data(int observed_data_id) {
    this->interface_observed_data_id_m.set(observed_data_id);
    return true;
  }

  virtual bool set_distribution_links(std::string input_type,
                                      Rcpp::IntegerVector ids) {
    this->input_type_m.set(input_type);
    this->key_m->resize(ids.size());
    for (R_xlen_t i = 0; i < ids.size(); i++) {
      this->key_m->at(i) = ids[i];
    }
    return true;
  }

  void set_note(std::string note) { this->notes.set(note); }

  virtual double evaluate() {
    fims_distributions::MultinomialLPMF<double> dmultinom;
    // Declare TMBVector in this scope
    dmultinom.observed_values.resize(this->observed_values.size());
    dmultinom.expected_values.resize(this->expected_values.size());
    for (size_t i = 0; i < observed_values.size(); i++) {
      dmultinom.observed_values[i] = this->observed_values[i].initial_value_m;
    }
    for (size_t i = 0; i < expected_values.size(); i++) {
      dmultinom.expected_values[i] = this->expected_values[i].initial_value_m;
    }
    dmultinom.dims.resize(2);
    dmultinom.dims[0] = this->dims[0];
    dmultinom.dims[1] = this->dims[1];
    return dmultinom.evaluate();
  }
 
  void finalize() {
    if (this->finalized) {
      // log warning that finalize has been called more than once.
      FIMS_WARNING_LOG("DmultinomDistributions  " +
                       fims::to_string(this->id_m) +
                       " has been finalized already.");
    }
 
    this->finalized = true;  // indicate this has been called already
 
    std::shared_ptr<fims_info::Information<double>> info =
        fims_info::Information<double>::GetInstance();
 
    fims_info::Information<double>::density_components_iterator it;
 
    // search for density component in Information
    it = info->density_components.find(this->id_m);
    // if not found, just return
    if (it == info->density_components.end()) {
      FIMS_WARNING_LOG("DmultinomDistributions " + fims::to_string(this->id_m) +
                       " not found in Information.");
      return;
    } else {
      std::shared_ptr<fims_distributions::MultinomialLPMF<double>> dmultinom =
          std::dynamic_pointer_cast<
              fims_distributions::MultinomialLPMF<double>>(it->second);
 
      this->lpdf_value = dmultinom->lpdf;
 
      size_t n_x = dmultinom->lpdf_vec.size();
      this->lpdf_vec = Rcpp::NumericVector(n_x);
      if (this->expected_values.size() != n_x) {
        this->expected_values.resize(n_x);
      }
      if (this->observed_values.size() != n_x) {
        this->observed_values.resize(n_x);
      }
      for (size_t i = 0; i < this->lpdf_vec.size(); i++) {
        this->lpdf_vec[i] = dmultinom->lpdf_vec[i];
        this->expected_values[i].final_value_m = dmultinom->get_expected(i);
        if (dmultinom->input_type != "data") {
          this->observed_values[i].final_value_m = dmultinom->get_observed(i);
        }
      }
      if (dmultinom->input_type == "data") {
        dims.resize(2);
        dims[0] = dmultinom->dims[0];
        dims[1] = dmultinom->dims[1];
        for (size_t i = 0; i < dims[0]; i++) {
          for (size_t j = 0; j < dims[1]; j++) {
            size_t idx = (i * dims[1]) + j;
            this->observed_values[idx].final_value_m = dmultinom->get_observed(
                static_cast<size_t>(i), static_cast<size_t>(j));
          }
        }
      }
    }
  }

  virtual std::string to_json() {
    std::stringstream ss;
 
    ss << "{\n";
    ss << " \"module_name\": \"density\",\n";
    ss << " \"module_id\": " << this->id_m << ",\n";
    ss << " \"module_type\": \"multinomial\",\n";
    ss << "\"observed_data_id\" : " << this->interface_observed_data_id_m
       << ",\n";
    ss << " \"input_type\" : \"" << this->input_type_m << "\",\n";
    ss << " \"density_component\": {\n";
    ss << "  \"lpdf_value\": " << sanitize_val(this->lpdf_value) << ",\n";
    ss << "  \"value\":[";
    if (this->lpdf_vec.size() == 0) {
      ss << "],\n";
    } else {
      for (size_t i = 0; i < this->lpdf_vec.size() - 1; i++) {
        ss << this->value_to_string(this->lpdf_vec[i]);
        ss << ", ";
      }
      ss << this->value_to_string(this->lpdf_vec[this->lpdf_vec.size() - 1]);
 
      ss << "],\n";
    }
    ss << "  \"expected_values\":[";
    if (this->expected_values.size() == 0) {
      ss << "],\n";
    } else {
      for (size_t i = 0; i < this->expected_values.size() - 1; i++) {
        ss << this->value_to_string(this->expected_values[i].final_value_m)
           << ", ";
      }
      ss << this->value_to_string(
          this->expected_values[this->expected_values.size() - 1]
              .final_value_m);
 
      ss << "],\n";
    }
    // no log_sd_values for multinomial
    ss << "  \"observed_values\":[";
    if (this->observed_values.size() == 0) {
      ss << "]\n";
    } else {
      for (size_t i = 0; i < this->observed_values.size() - 1; i++) {
        ss << this->observed_values[i].final_value_m << ", ";
      }
      ss << this->observed_values[this->observed_values.size() - 1]
                .final_value_m
         << "]\n";
    }
    ss << " }}\n";
    return ss.str();
  }
 
#ifdef TMB_MODEL
 
  template <typename Type>
  bool add_to_fims_tmb_internal() {
    std::shared_ptr<fims_info::Information<Type>> info =
        fims_info::Information<Type>::GetInstance();
 
    std::shared_ptr<fims_distributions::MultinomialLPMF<Type>> distribution =
        std::make_shared<fims_distributions::MultinomialLPMF<Type>>();
 
    distribution->id = this->id_m;
    distribution->observed_data_id_m = interface_observed_data_id_m;
    distribution->input_type = this->input_type_m;
    distribution->key.resize(this->key_m->size());
    for (size_t i = 0; i < this->key_m->size(); i++) {
      distribution->key[i] = this->key_m->at(i);
    }
    distribution->observed_values.resize(this->observed_values.size());
    for (size_t i = 0; i < this->observed_values.size(); i++) {
      distribution->observed_values[i] =
          this->observed_values[i].initial_value_m;
    }
    // set relative info
    distribution->expected_values.resize(this->expected_values.size());
    for (size_t i = 0; i < this->expected_values.size(); i++) {
      distribution->expected_values[i] =
          this->expected_values[i].initial_value_m;
    }
 
    info->density_components[distribution->id] = distribution;
    return true;
  }
 
  virtual bool add_to_fims_tmb() {
    this->add_to_fims_tmb_internal<TMB_FIMS_REAL_TYPE>();
    this->add_to_fims_tmb_internal<TMBAD_FIMS_TYPE>();
 
    return true;
  }
 
#endif
};
 
#endif