Simplify _Q_opt in parmest with block scenario structure.#3789
Simplify _Q_opt in parmest with block scenario structure.#3789sscini wants to merge 138 commits intoPyomo:mainfrom
Conversation
|
@djlaky @adowling2 Here are my initial thoughts on the redesign. Please provide feedback on code layout to this point. |
|
Dan and I had a good conversation today about the functional things needed in the redesign, and how to go about adding components to the overall block and sub-models. Working to implement changes, closely related to what is done currently in doe's _generate_scenario_blocks. |
|
@adowling2 @djlaky Added case for bootstrap, now works with theta_est, theta_est_bootstrap, and theta_est_leaveNout if I replace _Q_opt with _Q_opt_blocks. Please review. I am aware of the duplicated code, but when I attempted to unify them using n_scenarios = len(self.exp_list) or len(bootlist), I am getting an error for bootstrap I am still working out. Would we want to fully transition to using theta_est for obj and theta, and cov_est for covariance and add an error/warning? Or should I make it work to calculate cov if calc_cov=True? Thanks! |
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## main #3789 +/- ##
==========================================
- Coverage 89.40% 89.31% -0.10%
==========================================
Files 909 909
Lines 105541 105646 +105
==========================================
- Hits 94364 94360 -4
- Misses 11177 11286 +109
Flags with carried forward coverage won't be shown. Click here to find out more. ☔ View full report in Codecov by Sentry. 🚀 New features to boost your workflow:
|
|
@adowling2 @djlaky Should I tag larger team on this? Please review when available. Thanks! |
…ocks" This reverts commit 1aea99f.
|
From development meeting conversation, I moved all additional tests to test_parmest.py, and will work on removing deprecated theta_est arguments for covariance. Some tests may need to be adjusted. Review still would be greatly appreciated, and recommendation on whether the full MPI-SPPY dependent interface can be cleanly removed here, or address in a follow-up PR. Thank you! |
blnicho
left a comment
There was a problem hiding this comment.
@sscini I'm still going through parmest.py but here are my review comments so far on the rest of the files.
In response to your specific review questions
- How to handle the import structure that only applies to the deprecated interface
Could you leave a comment pointing to the "import structure" you're referring to. I'm not following. By deprecated interface do you mean the ParmEst interface before the Experiment class or the old interface for calculating the covariance?
- Are there any changes needed in the new implementation I did not consider?
Not sure yet.
- Now that this will be Pyomo 6.10.1, should I remove deprecation warnings for covariance structural change? (separate theta_est and cov_est)
Yes, I think you can remove the deprecated covariance functionality.
- Can elements from the new interface be updated into the deprecated and the old structural utility files be removed?
Are create_ef.py and mpi_utils.py the structural utility files you're referring to? I think we should leave these files for now and let's address the rest of the mpisppy infrastructure in a separate PR. Could you be more specific with which "elements" and which deprecated interface you're referring to?
|
@blnicho Thank you very much for your comments so far! I will add a specific comment in the code for the import clarification. Can elements from the new interface be updated into the deprecated and the old structural utility files be removed? Are create_ef.py and mpi_utils.py the structural utility files you're referring to? I think we should leave these files for now and let's address the rest of the mpisppy infrastructure in a separate PR. Could you be more specific with which "elements" and which deprecated interface you're referring to? Leaving create_ef.py and mpi_utils makes sense. To clarify this: I was referring to the _DeprecatedEstimator class as the deprecated interface. I updated _Q_opt in the Estimator class to use blocks and combined it with _Q_at_theta. The _DeprecatedEstimator still uses both _Q_opt and _Q_at_theta, which rely on the _experiment_instance_creation_callback function (line 95), and have separate booleans to connect to either mpisppy, or create_ef files. If I updated the _DeprecatedEstimator to use my _Q_opt implementation, since we previously discussed that the main point of the deprecated estimator is to preserve old user facing functions, then _experiment_instance_creation_callback, and the mentioned files could be removed. The main deprecated function to preserve now would be the use of theta_est for covariance calculations. Otherwise, I tried my best to avoid the user interface in this PR. Hopefully this clarifies, but please let me know if more information is still needed. As you mentioned, this will most likely be done in its own PR. |
|
@blnicho I removed the covariance support from theta_est and adjusted testing. Please review changes. Thank you! |
|
@blnicho Please share remaining comments mentioned for PR and parmest.py when available. Thank you! |
| total_data_points = 0 | ||
| for experiment in experiment_list: | ||
| total_number_data += len(experiment.get_labeled_model().experiment_outputs) | ||
| # 1. Identify the first parent component of the experiment outputs | ||
| output_vars = experiment.get_labeled_model().experiment_outputs | ||
|
|
||
| # 1. Identify the first parent component | ||
| first_var_key = list(output_vars.keys())[0] | ||
| first_parent = first_var_key.parent_component() | ||
| # 2. Count only the keys that belong to this specific parent | ||
| first_parent_indices = [ | ||
| v for v in output_vars.keys() if v.parent_component() is first_parent | ||
| ] | ||
| total_data_points += len(first_parent_indices) |
There was a problem hiding this comment.
This implementation seems fragile to me. What if the experiment outputs are in a different order? What if the first experiment output component is indexed by data point (i.e. time) and another index? What if you're trying to combine data from experiments with different outputs?
Is this calculating the n value needed for doing the covariance calculation or what is this value being used for?
| cov_n=NOTSET, | ||
| ): | ||
| def _create_scenario_blocks(self, bootlist=None, theta_vals=None, fix_theta=False): | ||
| # Create scenario block structure |
There was a problem hiding this comment.
| # Create scenario block structure |
This comment can be removed now that there is a docstring
| # 1. Identify the first parent component of the experiment outputs | ||
| output_vars = experiment.get_labeled_model().experiment_outputs | ||
|
|
||
| # 1. Identify the first parent component |
There was a problem hiding this comment.
You have two "# 1" comments here
| expanded_theta_names = self._expand_indexed_unknowns(template_model) | ||
| model._parmest_theta_names = tuple(expanded_theta_names) | ||
| model.parmest_theta = pyo.Var(model._parmest_theta_names) | ||
|
|
||
| for name in expanded_theta_names: | ||
| template_theta_var = template_model.find_component(name) | ||
| parent_theta_var = model.parmest_theta[name] |
There was a problem hiding this comment.
Using names to find equivalent theta components across the different experiment models kind of defeats the purpose of the unknown_parameters Suffix. I think a better solution might be to enforce that only unindexed components, or "ComponentData", are allowed in the unknown_parameters Suffix. Then the Suffix keys can be used directly and all these extra steps expanding potentially indexed things goes away. If this is too restrictive then in _expand_indexed_unknowns we could return a ComponentList or list of cuids instead of names.
| raise ValueError("At least one scenario is required to build the EF model.") | ||
|
|
||
| # Create indexed block for holding scenario models | ||
| model.exp_scenarios = pyo.Block(range(self.obj_probability_constant)) |
There was a problem hiding this comment.
Why not index this block by the scenario_numbers list you created above? Maybe add scenario_numbers to the model as a Pyomo Set so that you have easy access to it later?
| if self.diagnostic_mode: | ||
| print("Parmest _Q_opt model with scenario blocks:") | ||
| model.pprint() |
There was a problem hiding this comment.
You should use a logger here instead of printing directly to stdout
| # Initialize worst_status to optimal, update if not optimal | ||
| worst_status = pyo.TerminationCondition.optimal | ||
| # Get termination condition from solve result | ||
| status = solve_result.solver.termination_condition | ||
|
|
||
| # In case of fixing theta, just log a warning if not optimal | ||
| if status != pyo.TerminationCondition.optimal: | ||
| # logger.warning( | ||
| # "Solver did not terminate optimally when thetas were fixed. " | ||
| # "Termination condition: %s", | ||
| # str(status), | ||
| # ) | ||
| # Unless infeasible, update worst_status | ||
| if worst_status != pyo.TerminationCondition.infeasible: | ||
| worst_status = status |
There was a problem hiding this comment.
Why do you obscure an "infeasible" termination condition here? Why not just return the termination condition directly?
| if return_values is not None and len(return_values) > 0: | ||
| var_values = [] | ||
| # In the scenario blocks structure, exp_scenarios is an IndexedBlock | ||
| exp_blocks = self.ef_instance.exp_scenarios.values() |
There was a problem hiding this comment.
This is the only place you use self.ef_instance in this method instead of model so I suggest making this change for consistency
| exp_blocks = self.ef_instance.exp_scenarios.values() | |
| exp_blocks = model.exp_scenarios.values() |
| # retrieve the independent variables (i.e., estimated parameters) | ||
| ind_vars = [] | ||
| for nd_name, Var, sol_val in ef_nonants(self.ef_instance): | ||
| ind_vars.append(Var) | ||
| # calculate the reduced hessian | ||
| for name in self.ef_instance._parmest_theta_names: | ||
| var = self.ef_instance.parmest_theta[name] | ||
| ind_vars.append(var) | ||
|
|
There was a problem hiding this comment.
I suspect you can simplify this to:
| # retrieve the independent variables (i.e., estimated parameters) | |
| ind_vars = [] | |
| for nd_name, Var, sol_val in ef_nonants(self.ef_instance): | |
| ind_vars.append(Var) | |
| # calculate the reduced hessian | |
| for name in self.ef_instance._parmest_theta_names: | |
| var = self.ef_instance.parmest_theta[name] | |
| ind_vars.append(var) | |
| # retrieve the independent variables (i.e., estimated parameters) | |
| ind_vars = list(self.ef_instance.parmest_theta.values()) |
| deprecation_warning( | ||
| "The `initialize_parmest_model` option in `objective_at_theta()` is " | ||
| "deprecated and will be removed in future releases. Please ensure the" | ||
| "model is initialized within the Experiment class definition.", | ||
| version="6.9.5", |
There was a problem hiding this comment.
If this is a new deprecation warning introduced in this PR then our standard is to set the version to the current dev version of Pyomo. We'll update it to the next official Pyomo version as part of our release process.
| deprecation_warning( | |
| "The `initialize_parmest_model` option in `objective_at_theta()` is " | |
| "deprecated and will be removed in future releases. Please ensure the" | |
| "model is initialized within the Experiment class definition.", | |
| version="6.9.5", | |
| deprecation_warning( | |
| "The `initialize_parmest_model` option in `objective_at_theta()` is " | |
| "deprecated and will be removed in future releases. Please ensure the" | |
| "model is initialized within the Experiment class definition.", | |
| version="6.10.1.dev0", |
Fixes # .
Summary/Motivation:
_Q_opt, the main model building and solving function within parmest, is still heavily dependent on MPISPPY, and the code is becoming outdated and difficult to interpret. The goal of this PR is to redesign _Q_opt to work without this dependence, retaining all the current functionality
Changes proposed in this PR:
-_Q_opt redesign, using a block structure, with minimal changes needed to functions that utilize it.
Legal Acknowledgement
By contributing to this software project, I have read the contribution guide and agree to the following terms and conditions for my contribution: