Skip to content

Commit

Permalink
refactor nll cleanup
Browse files Browse the repository at this point in the history 
*issue #644
*refactor nll to positive log in density functions
*sum negative log densities in model.hpp
*rename observed_values to x

*only add standardization to dlnorm if data
*add lpdf_type to dlnorm interface

*fix do_log in density evaluations to true
  • Loading branch information
@Andrea-Havron-NOAA
Andrea-Havron-NOAA committed Jul 26, 2024
1 parent 7af5dc7 commit 66050f6
Show file tree
Hide file tree
Showing 15 changed files with 228 additions and 233 deletions.
6 changes: 3 additions & 3 deletions inst/include/common/model.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -113,7 +113,7 @@ class Model { // may need singleton
// Evaluate population
(*it).second->Evaluate();
// Recrtuiment negative log-likelihood
rec_nll += (*it).second->recruitment->evaluate_nll();
rec_nll -= (*it).second->recruitment->evaluate_lpdf();
MODEL_LOG << "Recruitment negative log-likelihood is: " << rec_nll
<< std::endl;
}
Expand All @@ -134,10 +134,10 @@ class Model { // may need singleton
#endif
MODEL_LOG << "Setting up pointer to fleet " << (*jt).second->id << "."
<< std::endl;
age_comp_nll += (*jt).second->evaluate_age_comp_nll();
age_comp_nll -= (*jt).second->evaluate_age_comp_lpmf();
MODEL_LOG << "Sum of survey and age comp negative log-likelihood is: "
<< age_comp_nll << std::endl;
index_nll += (*jt).second->evaluate_index_nll();
index_nll -= (*jt).second->evaluate_index_lpdf();
}
MODEL_LOG << "All fleets successfully evaluated." << std::endl;
// Loop over populations and fleets/surveys and fill in reporting
Expand Down
9 changes: 4 additions & 5 deletions inst/include/distributions/functors/density_components_base.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -34,10 +34,10 @@ struct DensityComponentBase : public fims_model_object::FIMSObject<Type> {
// Assigning each one its own ID is a way to keep track of
// all the instances of the DensityComponentBase class.
static uint32_t id_g; /**< global unique identifier for distribution modules */
fims::Vector<Type> observed_values; /**< input value of distribution function */
fims::Vector<Type> x; /**< input value of distribution function */
fims::Vector<Type> expected_values; /**< expected value of distribution function */
fims::Vector<Type> nll_vec; /**< vector to record observation level negative log-likelihood values */
std::string nll_type; /**< string classifies the type of the negative log-likelihood; options are: prior, re, data */
fims::Vector<Type> lpdf_vec; /**< vector to record observation level negative log-likelihood values */
std::string lpdf_type; /**< string classifies the type of the negative log-likelihood; options are: prior, re, data */
bool osa_flag = false; /**< Boolean; if true, osa residuals are calculated */
bool simulate_flag = false; /**< Boolean; if true, data are simulated from the distribution */

Expand All @@ -49,9 +49,8 @@ struct DensityComponentBase : public fims_model_object::FIMSObject<Type> {
/**
* @brief Generic probability density function. Calculates the pdf at the
* independent variable value.
* @param do_log Boolean; if true, log densities are returned
*/
virtual const Type evaluate(const bool& do_log) = 0;
virtual const Type evaluate() = 0;
};

/** @brief Default id of the singleton distribution class
Expand Down
63 changes: 31 additions & 32 deletions inst/include/distributions/functors/lognormal_lpdf.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -20,19 +20,19 @@
namespace fims_distributions
{
/**
* LogNormal Negative Log-Likelihood
* LogNormal Log Probability Density Function
*/
template <typename Type>
struct LogNormalLPDF : public DensityComponentBase<Type>
{
fims::Vector<Type> log_sd; /**< log of the standard deviation of the distribution on the log scale; can be a vector or scalar */
fims::Vector<Type> mu; /**< mean of the distribution on the log scale; can be a vector or scalar */
fims::Vector<Type> sd; /**< standard deviation of the distribution on the log scale; can be a vector or scalar */
fims::Vector<Type> log_logsd; /**< log of the standard deviation of the distribution on the log scale; can be a vector or scalar */
fims::Vector<Type> logmu; /**< mean of the distribution on the log scale; can be a vector or scalar */
fims::Vector<Type> logsd; /**< standard deviation of the distribution on the log scale; can be a vector or scalar */
std::vector<bool> is_na; /**< Boolean; if true, data observation is NA and the likelihood contribution is skipped */
#ifdef TMB_MODEL
::objective_function<Type> *of; /**< Pointer to the TMB objective function */
#endif
Type nll = 0.0; /**< total negative log-likelihood contribution of the distribution */
Type lpdf = 0.0; /**< total log probability density contribution of the distribution */
// data_indicator<tmbutils::vector<Type> , Type> keep; /**< Indicator used in TMB one-step-ahead residual calculations */

/** @brief Constructor.
Expand All @@ -46,74 +46,73 @@ namespace fims_distributions
virtual ~LogNormalLPDF() {}

/**
* @brief Evaluates the negative log-likelihood of the lognormal probability density function
* @param do_log Boolean; if true, log densities are returned
* @brief Evaluates the lognormal probability density function
*/
virtual const Type evaluate(const bool &do_log)
virtual const Type evaluate()
{
this->mu.resize(this->observed_values.size());
this->sd.resize(this->observed_values.size());
is_na.resize(this->observed_values.size());
this->nll_vec.resize(this->observed_values.size());
for (size_t i = 0; i < this->observed_values.size(); i++)
this->logmu.resize(this->x.size());
this->logsd.resize(this->x.size());
is_na.resize(this->x.size());
this->lpdf_vec.resize(this->x.size());
for (size_t i = 0; i < this->x.size(); i++)
{
if (this->expected_values.size() == 1)
{
this->mu[i] = this->expected_values[0];
this->logmu[i] = this->expected_values[0];
} else {
if(this->observed_values.size() != this->expected_values.size()){
if(this->x.size() != this->expected_values.size()){
/* move error handling to CreateModel in information so not to crash R
Rcpp::stop("the dimensions of the observed and expected values from lognormal negative log likelihood do not match");
*/
} else {
this->mu[i] = this->expected_values[i];
this->logmu[i] = this->expected_values[i];
}
}
if (log_sd.size() == 1)
if (log_logsd.size() == 1)
{
sd[i] = fims_math::exp(log_sd[0]);
logsd[i] = fims_math::exp(log_logsd[0]);
} else {
if(this->observed_values.size() != this->log_sd.size()){
if(this->x.size() != this->log_logsd.size()){
/* move error handling to CreateModel in information so not to crash R
Rcpp::stop("the dimensions of the observed and log sd values from lognormal negative log likelihood do not match");
Rcpp::stop("the dimensions of the observed and log logsd values from lognormal negative log likelihood do not match");
*/
} else {
sd[i] = fims_math::exp(log_sd[i]);
logsd[i] = fims_math::exp(log_logsd[i]);
}
}

if(!is_na[i])
{
#ifdef TMB_MODEL
// this->nll_vec[i] = this->keep[i] * -dnorm(this->observed_values[i], mu[i], sd[i], do_log);
this->nll_vec[i] = -dnorm(log(this->observed_values[i]), mu[i], sd[i], true) + log(this->observed_values[i]);
if(!do_log){
this->nll_vec[i] = -exp(-this->nll_vec[i]);
// this->lpdf_vec[i] = this->keep[i] * dnorm(this->x[i], logmu[i], logsd[i], true);
this->lpdf_vec[i] = dnorm(log(this->x[i]), logmu[i], logsd[i], true);
if(this->lpdf_type == "data"){
this->lpdf_vec[i] -= log(this->x[i]);
}
nll += this->nll_vec[i];
lpdf += this->lpdf_vec[i];
if (this->simulate_flag)
{
FIMS_SIMULATE_F(this->of)
{ // preprocessor definition in interface.hpp
// this simulates data that is mean biased
this->observed_values[i] = fims_math::exp(rnorm(mu[i], sd[i]));
this->x[i] = fims_math::exp(rnorm(logmu[i], logsd[i]));
}
}
#endif

/* osa not working yet
if(osa_flag){//data observation type implements osa residuals
//code for osa cdf method
this->nll_vec[i] = this->keep.cdf_lower[i] * -log( pnorm(this->observed_values[i], mu[i], sd[i]) );
this->nll_vec[i] = this->keep.cdf_upper[i] * -log( 1.0 - pnorm(this->observed_values[i], mu[i], sd[i]) );
this->lpdf_vec[i] = this->keep.cdf_lower[i] * log( pnorm(this->x[i], logmu[i], logsd[i]) );
this->lpdf_vec[i] = this->keep.cdf_upper[i] * log( 1.0 - pnorm(this->x[i], logmu[i], logsd[i]) );
} */
}
}
#ifdef TMB_MODEL
vector<Type> lognormal_observed_values = this->observed_values;
// FIMS_REPORT_F(lognormal_observed_values, this->of);
vector<Type> lognormal_x = this->x;
// FIMS_REPORT_F(lognormal_x, this->of);
#endif
return (nll);
return (lpdf);
}
};
} // namespace fims_distributions
Expand Down
40 changes: 20 additions & 20 deletions inst/include/distributions/functors/multinomial_lpmf.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -21,12 +21,12 @@
namespace fims_distributions
{
/**
* Multinomial Negative Log-Likelihood
* Multinomial Log Probability Mass Function
*/
template <typename Type>
struct MultinomialLPMF : public DensityComponentBase<Type>
{
Type nll = 0.0; /**< total negative log-likelihood contribution of the distribution */
Type lpmf = 0.0; /**< total negative log-likelihood contribution of the distribution */
fims::Vector<size_t> dims; /**< Dimensions of the number of rows and columns of the multivariate dataset */
std::vector<bool> is_na; /**< Boolean; if true, data observation is NA and the likelihood contribution for the entire row is skipped */
#ifdef TMB_MODEL
Expand All @@ -45,42 +45,42 @@ namespace fims_distributions
virtual ~MultinomialLPMF() {}

/**
* @brief Evaluates the negative log-likelihood of the multinomial probability mass function
* @param do_log Boolean; if true, log densities are returned
* @brief Evaluates the multinomial probability mass function
*/
virtual const Type evaluate(const bool& do_log)
virtual const Type evaluate()
{
this->nll_vec.resize(dims[0]);
fims::Vector<Type> observed_vector;
fims::Vector<Type> expected_vector;
observed_vector.resize(dims[1]);
expected_vector.resize(dims[1]);
this->lpdf_vec.resize(dims[0]);
fims::Vector<Type> x_vector;
fims::Vector<Type> prob_vector;
x_vector.resize(dims[1]);
prob_vector.resize(dims[1]);
Type lpdf = 0.0; /**< total log probability mass contribution of the distribution */
for (size_t i = 0; i < dims[0]; i++)
{
if(!is_na[i]){
#ifdef TMB_MODEL
for (size_t j = 0; j < dims[1]; j++)
{
size_t idx = (i * dims[1]) + j;
observed_vector[j] = this->observed_values[idx];
expected_vector[j] = this->expected_values[idx];
x_vector[j] = this->x[idx];
prob_vector[j] = this->expected_values[idx];
}

this->nll_vec[i] = -dmultinom((vector<Type>)observed_vector, (vector<Type>)expected_vector, do_log);
nll += this->nll_vec[i];
this->lpdf_vec[i] = dmultinom((vector<Type>)x_vector, (vector<Type>)prob_vector, true);
lpdf += this->lpdf_vec[i];
/*
if (this->simulate_flag)
{
FIMS_SIMULATE_F(this->of)
{
fims::Vector<Type> sim_observed;
sim_observed.resize(dims[1]);
sim_observed = rmultinom(expected_vector);
sim_observed.resize(this->observed_values);
sim_observed = rmultinom(prob_vector);
sim_observed.resize(this->x);
for (size_t j = 0; j < dims[1]; j++)
{
idx = (i * dims[1]) + j;
this->observed_values[idx] = sim_observed[j];
this->x[idx] = sim_observed[j];
}
}
}
Expand All @@ -89,10 +89,10 @@ namespace fims_distributions
}
}
#ifdef TMB_MODEL
vector<Type> observed_values = this->observed_values;
// FIMS_REPORT_F(observed_values, this->of);
vector<Type> x = this->x;
// FIMS_REPORT_F(x, this->of);
#endif
return (nll);
return (lpdf);
}

};
Expand Down
39 changes: 19 additions & 20 deletions inst/include/distributions/functors/normal_lpdf.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -20,14 +20,14 @@

namespace fims_distributions {
/**
* Normal Negative Log-Likelihood
* Normal Log Probability Density Function
*/
template<typename Type>
struct NormalLPDF : public DensityComponentBase<Type> {
fims::Vector<Type> log_sd; /**< log of the standard deviation of the distribution; can be a vector or scalar */
fims::Vector<Type> mu; /**< mean of the distribution; can be a vector or scalar */
fims::Vector<Type> sd; /**< standard deviation of the distribution; can be a vector or scalar */
Type nll = 0.0; /**< total negative log-likelihood contribution of the distribution */
Type lpdf = 0.0; /**< total log probability density contribution of the distribution */
std::vector<bool> is_na; /**< Boolean; if true, data observation is NA and the likelihood contribution is skipped */
#ifdef TMB_MODEL
::objective_function<Type> *of; /**< Pointer to the TMB objective function */
Expand All @@ -45,18 +45,17 @@ struct NormalLPDF : public DensityComponentBase<Type> {
virtual ~NormalLPDF() {}

/**
* @brief Evaluates the negative log-likelihood of the normal probability density function
* @param do_log Boolean; if true, log densities are returned
* @brief Evaluates the normal probability density function
*/
virtual const Type evaluate(const bool& do_log){
this->mu.resize(this->observed_values.size());
this->sd.resize(this->observed_values.size());
this->nll_vec.resize(this->observed_values.size());
for(size_t i=0; i<this->observed_values.size(); i++){
virtual const Type evaluate(){
this->mu.resize(this->x.size());
this->sd.resize(this->x.size());
this->lpdf_vec.resize(this->x.size());
for(size_t i=0; i<this->x.size(); i++){
if(this->expected_values.size() == 1){
this->mu[i] = this->expected_values[0];
} else {
if(this->observed_values.size() != this->expected_values.size()){
if(this->x.size() != this->expected_values.size()){
/* move error handling to CreateModel in information so not to crash R
Rcpp::stop("the dimensions of the observed and expected values from normal negative log likelihood do not match");
*/
Expand All @@ -67,7 +66,7 @@ struct NormalLPDF : public DensityComponentBase<Type> {
if(log_sd.size() == 1){
sd[i] = fims_math::exp(log_sd[0]);
} else {
if(this->observed_values.size() != this->log_sd.size()){
if(this->x.size() != this->log_sd.size()){
/* move error handling to CreateModel in information so not to crash R
Rcpp::stop("the dimensions of the observed and log sd values from normal negative log likelihood do not match");
*/
Expand All @@ -77,30 +76,30 @@ struct NormalLPDF : public DensityComponentBase<Type> {
}
if(!is_na[i])
{
// this->nll_vec[i] = this->keep[i] * -dnorm(this->observed_values[i], mu[i], sd[i], do_log);
// this->lpdf_vec[i] = this->keep[i] * -dnorm(this->x[i], mu[i], sd[i], true);
#ifdef TMB_MODEL
this->nll_vec[i] = -dnorm(this->observed_values[i], mu[i], sd[i], do_log);
this->lpdf_vec[i] = dnorm(this->x[i], mu[i], sd[i], true);

nll += this->nll_vec[i];
lpdf += this->lpdf_vec[i];
if(this->simulate_flag){
FIMS_SIMULATE_F(this->of){
this->observed_values[i] = rnorm(mu[i], sd[i]);
this->x[i] = rnorm(mu[i], sd[i]);
}
}
#endif
/* osa not working yet
if(osa_flag){//data observation type implements osa residuals
//code for osa cdf method
this->nll_vec[i] = this->keep.cdf_lower[i] * -log( pnorm(this->observed_values[i], mu[i], sd[i]) );
this->nll_vec[i] = this->keep.cdf_upper[i] * -log( 1.0 - pnorm(this->observed_values[i], mu[i], sd[i]) );
this->lpdf_vec[i] = this->keep.cdf_lower[i] * log( pnorm(this->x[i], mu[i], sd[i]) );
this->lpdf_vec[i] = this->keep.cdf_upper[i] * log( 1.0 - pnorm(this->x[i], mu[i], sd[i]) );
} */
}
}
#ifdef TMB_MODEL
vector<Type> normal_observed_values = this->observed_values;
//FIMS_REPORT_F(normal_observed_values, this->of);
vector<Type> normal_x = this->x;
//FIMS_REPORT_F(normal_x, this->of);
#endif
return(nll);
return(lpdf);
}

};
Expand Down
Loading

0 comments on commit 66050f6

Please sign in to comment.