Skip to content

Commit

Permalink
Convert ReporterValue to include confidence intervals
Browse files Browse the repository at this point in the history 
(refs idaholab#16584)
  • Loading branch information
@aeslaughter
aeslaughter committed Jun 2, 2021
1 parent e9171b1 commit 7b9c52c
Show file tree
Hide file tree
Showing 8 changed files with 308 additions and 260 deletions.
59 changes: 25 additions & 34 deletions modules/stochastic_tools/include/reporters/StatisticsReporter.h
View file
Expand Up @@ -19,17 +19,17 @@
* ReporterContext that utilizes a Calculator object to compute its value and confidence levels
*/
template <typename InType, typename OutType>
class ReporterStatisticsContext : public ReporterContext<OutType>
class ReporterStatisticsContext : public ReporterContext<std::pair<OutType, std::vector<OutType>>>
{
public:
ReporterStatisticsContext(const libMesh::ParallelObject & other,
ReporterState<OutType> & state,
ReporterState<std::pair<OutType, std::vector<OutType>>> & state,
const InType & data,
const ReporterProducerEnum & mode,
const MooseEnumItem & stat);

ReporterStatisticsContext(const libMesh::ParallelObject & other,
ReporterState<OutType> & state,
ReporterState<std::pair<OutType, std::vector<OutType>>> & state,
const InType & data,
const ReporterProducerEnum & mode,
const MooseEnumItem & stat,
Expand All @@ -54,19 +54,16 @@ class ReporterStatisticsContext : public ReporterContext<OutType>
/// Storage for the BootstrapCalculator for the desired confidence interval calculations (optional)
std::unique_ptr<const StochasticTools::BootstrapCalculator<InType, OutType>> _ci_calc_ptr =
nullptr;

/// The results
std::vector<OutType> _ci_results;
};

template <typename InType, typename OutType>
ReporterStatisticsContext<InType, OutType>::ReporterStatisticsContext(
const libMesh::ParallelObject & other,
ReporterState<OutType> & state,
ReporterState<std::pair<OutType, std::vector<OutType>>> & state,
const InType & data,
const ReporterProducerEnum & mode,
const MooseEnumItem & stat)
: ReporterContext<OutType>(other, state),
: ReporterContext<std::pair<OutType, std::vector<OutType>>>(other, state),
_data(data),
_data_mode(mode),
_calc_ptr(StochasticTools::makeCalculator<InType, OutType>(stat, other))
Expand All @@ -76,7 +73,7 @@ ReporterStatisticsContext<InType, OutType>::ReporterStatisticsContext(
template <typename InType, typename OutType>
ReporterStatisticsContext<InType, OutType>::ReporterStatisticsContext(
const libMesh::ParallelObject & other,
ReporterState<OutType> & state,
ReporterState<std::pair<OutType, std::vector<OutType>>> & state,
const InType & data,
const ReporterProducerEnum & mode,
const MooseEnumItem & stat,
Expand All @@ -94,22 +91,22 @@ template <typename InType, typename OutType>
void
ReporterStatisticsContext<InType, OutType>::finalize()
{
this->_state.value() = _calc_ptr->compute(_data, _data_mode == REPORTER_MODE_DISTRIBUTED);
ReporterContext<OutType>::finalize();
this->_state.value().first = _calc_ptr->compute(_data, _data_mode == REPORTER_MODE_DISTRIBUTED);
ReporterContext<std::pair<OutType, std::vector<OutType>>>::finalize();

if (_ci_calc_ptr)
_ci_results = _ci_calc_ptr->compute(_data, _data_mode == REPORTER_MODE_DISTRIBUTED);
this->_state.value().second =
_ci_calc_ptr->compute(_data, _data_mode == REPORTER_MODE_DISTRIBUTED);
}

template <typename InType, typename OutType>
void
ReporterStatisticsContext<InType, OutType>::store(nlohmann::json & json) const
{
ReporterContext<OutType>::store(json);
ReporterContext<std::pair<OutType, std::vector<OutType>>>::store(json);
json["stat"] = _calc_ptr->name();
if (_ci_calc_ptr)
json["confidence_intervals"] = {{"method", _ci_calc_ptr->name()},
{"values", _ci_results},
{"levels", _ci_calc_ptr->levels()},
{"replicates", _ci_calc_ptr->replicates()},
{"seed", _ci_calc_ptr->seed()}};
Expand Down Expand Up @@ -140,23 +137,15 @@ class StatisticsReporter : public GeneralReporter
// CI Method to be computed (optional)
const MooseEnum & _ci_method;

// CI levels to be computed (not a ref. by design since these are computed from input parameter)
const std::vector<Real> _ci_levels;
// CI levels to be computed
const std::vector<Real> & _ci_levels;

// Number of CI replicates to use in Bootstrap methods
const unsigned int & _ci_replicates;

// Random seed for producing CI replicates
const unsigned int & _ci_seed;

/**
* Helper function for converting confidence levels given in (0, 0.5] into levels in (0, 1).
* For example, levels_in = {0.05, 0.1, 0.5} converts to {0.05 0.1, 0.5, 0.9, 0.95}.
*
* This also performs error checking on the supplied "ci_levels".
*/
std::vector<Real> computeLevels(const std::vector<Real> & levels_in) const;

/**
* Helper for adding statistic reporter values
*
Expand All @@ -176,17 +165,19 @@ StatisticsReporter::declareValueHelper(const ReporterName & r_name)
{
const std::string s_name = r_name.getCombinedName() + "_" + item.name();
if (_ci_method.isValid())
declareValueByName<StatType, ReporterStatisticsContext<InType, OutType>>(s_name,
REPORTER_MODE_ROOT,
data,
mode,
item,
_ci_method,
_ci_levels,
_ci_replicates,
_ci_seed);
declareValueByName<std::pair<OutType, std::vector<OutType>>,
ReporterStatisticsContext<InType, OutType>>(s_name,
REPORTER_MODE_ROOT,
data,
mode,
item,
_ci_method,
_ci_levels,
_ci_replicates,
_ci_seed);
else
declareValueByName<StatType, ReporterStatisticsContext<InType, OutType>>(
declareValueByName<std::pair<OutType, std::vector<OutType>>,
ReporterStatisticsContext<InType, OutType>>(
s_name, REPORTER_MODE_ROOT, data, mode, item);
}
}
51 changes: 19 additions & 32 deletions modules/stochastic_tools/src/reporters/StatisticsReporter.C
View file
Expand Up @@ -46,7 +46,9 @@ StatisticsReporter::validParams()
"ci_method", ci, "The method to use for computing confidence level intervals.");

params.addParam<std::vector<Real>>(
"ci_levels", "A vector of confidence levels to consider, values must be in (0, 0.5].");
"ci_levels",
std::vector<Real>({0.1, 0.9}),
"A vector of confidence levels to consider, values must be in (0, 1).");
params.addParam<unsigned int>(
"ci_replicates",
10000,
Expand All @@ -62,11 +64,25 @@ StatisticsReporter::StatisticsReporter(const InputParameters & parameters)
: GeneralReporter(parameters),
_compute_stats(getParam<MultiMooseEnum>("compute")),
_ci_method(getParam<MooseEnum>("ci_method")),
_ci_levels(_ci_method.isValid() ? computeLevels(getParam<std::vector<Real>>("ci_levels"))
: std::vector<Real>()),
_ci_levels(getParam<std::vector<Real>>("ci_levels")),
_ci_replicates(getParam<unsigned int>("ci_replicates")),
_ci_seed(getParam<unsigned int>("ci_seed"))
{
// CI levels error checking
if (_ci_method.isValid())
{
if (_ci_levels.empty())
paramError("ci_levels",
"If the 'ci_method' parameter is supplied then the 'ci_levels' must also be "
"supplied with values in (0, 1).");

else if (*std::min_element(_ci_levels.begin(), _ci_levels.end()) <= 0)
paramError("ci_levels", "The supplied levels must be greater than zero.");

else if (*std::max_element(_ci_levels.begin(), _ci_levels.end()) >= 1)
paramError("ci_levels", "The supplied levels must be less than 1.0");
}

// Stats for Reporters
if (isParamValid("reporters"))
{
Expand Down Expand Up @@ -109,32 +125,3 @@ StatisticsReporter::StatisticsReporter(const InputParameters & parameters)
else
mooseError("The 'vectorpostprocessors' and/or 'reporters' parameters must be defined.");
}

std::vector<Real>
StatisticsReporter::computeLevels(const std::vector<Real> & levels_in) const
{
if (levels_in.empty())
paramError("ci_levels",
"If the 'ci_method' parameter is supplied then the 'ci_levels' must also be "
"supplied with values in (0, 0.5].");

else if (*std::min_element(levels_in.begin(), levels_in.end()) <= 0)
paramError("ci_levels", "The supplied levels must be greater than zero.");

else if (*std::max_element(levels_in.begin(), levels_in.end()) > 0.5)
paramError("ci_levels", "The supplied levels must be less than or equal to 0.5");

std::list<Real> levels_out;
for (auto it = levels_in.rbegin(); it != levels_in.rend(); ++it)
{
if (*it == 0.5)
levels_out.push_back(*it);

else
{
levels_out.push_front(*it);
levels_out.push_back(1 - *it);
}
}
return std::vector<Real>(levels_out.begin(), levels_out.end());
}

0 comments on commit 7b9c52c

Please sign in to comment.