From 97ca231afc4521a0a846db23bc133c34a1d2eaab Mon Sep 17 00:00:00 2001
From: titusquah <46580668+titusquah@users.noreply.github.com>
Date: Tue, 16 Jun 2020 09:27:19 -0600
Subject: [PATCH] Edited parity plot to allow color to represent 3rd dimension.
Still need to improve colorbar axis name.
---
.idea/dictionaries/Titus.xml | 2 +
.idea/workspace.xml | 126 ++++++++++++++++++-----------------
reeps/reeps.py | 106 ++++++++++++++++++++---------
tests/test_multi_reeps.py | 6 +-
4 files changed, 147 insertions(+), 93 deletions(-)
diff --git a/.idea/dictionaries/Titus.xml b/.idea/dictionaries/Titus.xml
index 43b55a4..c7fc181 100644
--- a/.idea/dictionaries/Titus.xml
+++ b/.idea/dictionaries/Titus.xml
@@ -19,11 +19,13 @@
ndarray
pred
quah
+ rbfopt
reeps
scipy
seaborn
slsqp
thermo
+ viridis
xmls
diff --git a/.idea/workspace.xml b/.idea/workspace.xml
index d476773..e28ee0d 100644
--- a/.idea/workspace.xml
+++ b/.idea/workspace.xml
@@ -4,13 +4,8 @@
-
-
-
-
-
+
-
@@ -47,14 +42,11 @@
-
+
-
-
-
@@ -69,7 +61,7 @@
-
+
@@ -77,20 +69,20 @@
-
+
-
+
-
+
-
+
@@ -98,32 +90,29 @@
-
+
-
+
-
+
-
+
-
-
-
-
-
+
+
-
+
@@ -132,27 +121,34 @@
-
+
+
+
+
-
-
+
+
-
-
-
+
+
+
+
+
+
+
-
@@ -245,7 +241,14 @@
1591388459379
-
+
+ 1591972483195
+
+
+
+ 1591972483196
+
+
@@ -266,7 +269,8 @@
-
+
+
@@ -286,43 +290,43 @@
-
+
-
+
-
+
-
-
-
+
+
+
-
+
-
-
-
+
+
+
-
+
-
-
-
+
+
+
-
+
-
-
-
+
+
+
-
+
-
+
@@ -335,12 +339,12 @@
-
+
-
+
@@ -351,10 +355,10 @@
-
+
-
+
diff --git a/reeps/reeps.py b/reeps/reeps.py
index 8f7d339..8c43006 100644
--- a/reeps/reeps.py
+++ b/reeps/reeps.py
@@ -64,7 +64,8 @@ class REEPS:
The ordering of the columns needs to be:
- [h_i, h_eq, z_i, z_eq, {RE_1}_aq_i, {RE_1}_aq_eq, {RE_1}_d_eq,
+ [h_i, h_eq, z_i, z_eq,
+ {RE_1}_aq_i, {RE_1}_aq_eq, {RE_1}_d_eq,
{RE_2}_aq_i, {RE_2}_aq_eq, {RE_2}_d_eq,...
{RE_N}_aq_i, {RE_N}_aq_eq, {RE_N}_d_eq]
@@ -347,15 +348,17 @@ class REEPS:
self.update_predicted_dict()
@staticmethod
- def slsqp_optimizer(objective, x_guess):
+ def scipy_minimize(objective, x_guess, optimizer_kwargs=None):
""" The default optimizer for REEPS
- Uses scipy.minimize with options
+ Uses scipy.minimize
+
+ By default, options are
.. code-block:: python
default_kwargs= {"method": 'SLSQP',
- "bounds": [(1e-1, 1e1)*len(x_guess)],
+ "bounds": [(1e-1, 1e1)]*len(x_guess),
"constraints": (),
"options": {'disp': True,
'maxiter': 1000,
@@ -363,14 +366,16 @@ class REEPS:
:param objective: (func) the objective function
:param x_guess: (np.ndarray) the initial guess (always 1)
+ :param optimizer_kwargs: (dict) dictionary of options for minimize
:returns: (np.ndarray) Optimized parameters
"""
- optimizer_kwargs = {"method": 'SLSQP',
- "bounds": [(1e-1, 1e1)] * len(x_guess),
- "constraints": (),
- "options": {'disp': True,
- 'maxiter': 1000,
- 'ftol': 1e-6}}
+ if optimizer_kwargs is None:
+ optimizer_kwargs = {"method": 'SLSQP',
+ "bounds": [(1e-1, 1e1)] * len(x_guess),
+ "constraints": (),
+ "options": {'disp': True,
+ 'maxiter': 1000,
+ 'ftol': 1e-6}}
res = minimize(objective, x_guess, **optimizer_kwargs)
est_parameters = res.x
return est_parameters
@@ -450,7 +455,7 @@ class REEPS:
"""Change list of Cantera solutions by inputting
new xml file name and phase names
- Also runs set_in_moles to set initial molality to 1 g/L
+ Also runs set_in_moles to set feed volume to 1 L
:param phases_xml_filename: (str) xml file with parameters
for equilibrium calc
@@ -656,7 +661,7 @@ class REEPS:
This function also calls update_predicted_dict
- :param feed_vol: (float) feed volume of mixture (g/L)
+ :param feed_vol: (float) feed volume of mixture (L)
"""
phases_copy = self._phases.copy()
exp_df = self._exp_df.copy()
@@ -809,7 +814,7 @@ class REEPS:
"with at least 2 arguments: "
"f(objective_func,x_guess, kwargs)")
if optimizer == 'SLSQP':
- optimizer = self.slsqp_optimizer
+ optimizer = self.scipy_minimize
self._optimizer = optimizer
return None
@@ -939,7 +944,9 @@ class REEPS:
i = 0
for species_name in opt_dict.keys():
for thermo_prop in opt_dict[species_name].keys():
- opt_dict[species_name][thermo_prop] *= x[i]
+ if not np.isnan(
+ x[i]): # if nan, do not update xml with nan
+ opt_dict[species_name][thermo_prop] *= x[i]
i += 1
self.update_xml(opt_dict, temp_xml_file_path)
@@ -969,9 +976,12 @@ class REEPS:
with optimizer. Returns dictionary with opt_dict structure
:param objective_function: (function) function to compute objective
+ If 'None', last set objective or default function is used
:param optimizer: (function) function to perform optimization
- :param optimizer_kwargs: (dict) arguments for optimizer
- :param objective_kwargs: (dict) arguments for objective function
+ If 'None', last set optimizer or default is used
+ :param optimizer_kwargs: (dict) optional arguments for optimizer
+ :param objective_kwargs: (dict) optional arguments
+ for objective function
:returns opt_dict: (dict) optimized opt_dict. Has identical structure
as opt_dict
"""
@@ -1043,6 +1053,8 @@ class REEPS:
def parity_plot(self,
compared_value=None,
+ color_axis=None,
+ plot_title=None,
save_path=None,
print_r_squared=False):
"""
@@ -1052,6 +1064,16 @@ class REEPS:
:param compared_value: (str) Quantity to compare predicted and
experimental data. Can be any column containing "eq" in exp_df i.e.
h_eq, z_eq, {RE}_d_eq, etc.
+ :param plot_title: (str or boolean)
+
+ If None (default): Plot title will be generated from compared_value
+ Recommend to just explore. If h_eq, plot_title is
+ "H^+ eq conc".
+
+ If str: Plot title will be plot_title string
+
+ If "False": No plot title
+ :param color_axis: (dict)
:param save_path: (str) save path for parity plot
:param print_r_squared: (boolean) To plot or not to plot r-squared
value. Prints 2 places past decimal
@@ -1071,33 +1093,57 @@ class REEPS:
name_breakdown = re.findall('[^_\W]+', compared_value)
compared_species = name_breakdown[0]
if compared_species == 'h':
- species_name = '$H^+$'
+ default_title = '$H^+$ eq. conc. (mol/L)'
elif compared_species == 'z':
- species_name = extractant_name
+ default_title = '{0} eq. conc. (mol/L)'.format(extractant_name)
else:
phase = name_breakdown[1]
if phase == 'aq':
re_charge = re_charges[re_species_list.index(compared_species)]
- species_name = '$%s^{%d+}$' % (compared_species, re_charge)
+ default_title = '$%s^{%d+}$ eq. conc. (mol/L)' \
+ % (compared_species, re_charge)
elif phase == 'd':
- species_name = '{0} distribution ratio'.format(
+ default_title = '{0} distribution ratio'.format(
compared_species)
else:
- species_name = '{0} complex'.format(compared_species)
+ default_title = '{0} complex eq. conc. (mol/L)'.format(
+ compared_species)
fig, ax = plt.subplots()
- p1 = sns.scatterplot(meas, pred, color="r",
- label="{0} eq. conc. (mol/L)".format(
- species_name),
- legend=False)
+ if color_axis is None:
+ sns.scatterplot(meas, pred, color="r",
+ legend=False)
+ else:
+ key = list(color_axis.keys())[0]
+ value = list(color_axis.values())[0]
+ if key == 'predicted':
+ y = self.get_predicted_dict()[value]
+ elif key == 'measured':
+ y = self.get_exp_df()[value].values
+ else:
+ raise Exception('color_axis must be a dictionary with key'
+ '"predicted" or "measured"')
+ y = np.array(y)
+ meas = np.array(meas)
+ pred = np.array(pred)
+
+ p1 = ax.scatter(meas, pred, c=y, alpha=1, cmap='viridis')
+ c_bar = fig.colorbar(p1, format='%.2f')
+ # Fix next line. value is just the dictionary value.
+ c_bar.set_label(value, rotation=270, labelpad=20)
sns.lineplot(min_max_data, min_max_data, color="b", label="")
+
if print_r_squared:
- p1.text(min_max_data[0],
+ ax.text(min_max_data[0],
min_max_data[1] * 0.9,
'$R^2$={0:.2f}'.format(self.r_squared(compared_value)))
- plt.legend(loc='lower right')
- else:
- plt.legend()
+ # plt.legend(loc='lower right')
+ # else:
+ # plt.legend()
ax.set(xlabel='Measured', ylabel='Predicted')
+ if plot_title is None:
+ ax.set_title(default_title)
+ elif isinstance(plot_title, str):
+ ax.set_title(plot_title)
plt.show()
if save_path is not None:
plt.savefig(save_path, bbox_inches='tight')
@@ -1110,7 +1156,7 @@ class REEPS:
:param compared_value: (str) Quantity to compare predicted and
experimental data. Can be any column containing "eq" in exp_df i.e.
- h_eq, z_eq, {RE}_d_eq, etc.
+ h_eq, z_eq, {RE}_d_eq, etc. default is {RE}_aq_eq
"""
exp_df = self.get_exp_df()
predicted_dict = self.get_predicted_dict()
diff --git a/tests/test_multi_reeps.py b/tests/test_multi_reeps.py
index e84c7d4..e2017bc 100644
--- a/tests/test_multi_reeps.py
+++ b/tests/test_multi_reeps.py
@@ -15,8 +15,10 @@ def optimizer(func, x_guess):
ub = np.array([1e1])
bounds = (lb, ub)
options = {'c1': 1e-3, 'c2': 1e-3, 'w': 0.9}
- mini_optimizer = ps.single.global_best.GlobalBestPSO(n_particles=100, dimensions=1,
- options=options, bounds=bounds)
+ mini_optimizer = ps.single.global_best.GlobalBestPSO(n_particles=100,
+ dimensions=1,
+ options=options,
+ bounds=bounds)
f_opt, x_opt = mini_optimizer.optimize(func, iters=100)
return x_opt