Plotting Code Samples
Contents
C. Plotting Code Samples¶
C.1. hymod.ipynb¶
The following are the plotting functions as described in the hymod.ipynb Jupyter notebook tutorial.
The following are the necessary package imports to run these functions:
import numpy as np
import seaborn as sns
import matplotlib.pyplot as plt
from matplotlib.lines import Line2D
C.1.1. plot_observed_vs_simulated_streamflow()¶
def plot_observed_vs_simulated_streamflow(df, hymod_dict, figsize=[12, 6]):
"""Plot observed versus simulated streamflow.
:param df: Dataframe of hymod input data including columns for precip, potential evapotranspiration,
and streamflow
:param hymod_dict: A dictionary of hymod outputs
:type hymod_dict: dict
:param figsize: Matplotlib figure size
:type figsize: list
"""
# set plot style
plt.style.use('seaborn-white')
# set up figure
fig, ax = plt.subplots(figsize=figsize)
# plot observed streamflow
ax.plot(range(0, len(df['Strmflw'])), df['Strmflw'], color='pink')
# plot simulated streamflow
ax.plot(range(0, len(df['Strmflw'])), hymod_dict['Q'], color='black')
# set axis labels
ax.set_ylabel('Streamflow($m^3/s$)')
ax.set_xlabel('Days')
# set plot title
plt.title('Observed vs. Simulated Streamflow')
return ax
C.1.2. plot_observed_vs_sensitivity_streamflow()¶
def plot_observed_vs_sensitivity_streamflow(df_obs, df_sim, figsize=[10, 4]):
"""Plot observed streamflow versus simulations generated from sensitivity analysis.
:param df_obs: Dataframe of mean monthly hymod input data including columns for precip,
potential evapotranspiration, and streamflow
:param df_sim: Dataframe of mean monthly simulation data from sensitivity analysis
:param figsize: Matplotlib figure size
:type figsize: list
"""
month_list = range(len(df_sim))
# set up figure
fig, ax = plt.subplots(figsize=figsize)
# set labels
ax.set_xlabel('Days')
ax.set_ylabel('Flow Discharge (m^3/s)')
# plots all simulated streamflow cases under different sample sets
for i in df_sim.columns:
plt.plot(month_list, df_sim[i], color="pink", alpha=0.2)
# plot observed streamflow
plt.plot(month_list, df_obs['Strmflw'], color="black")
plt.title('Observed vs. Sensitivity Analysis Outputs')
return ax
C.1.3. plot_monthly_heatmap()¶
def plot_monthly_heatmap(arr_sim, df_obs, title='', figsize=[14, 6]):
"""Plot a sensitivity metric overlain by observed flow.
:param arr_sim: Numpy array of simulated metrics
:param df_obs: Dataframe of mean monthly observed data from sensitivity analysis
:param title: Title of plot
:type title: str
:param figsize: Matplotlib figure size
:type figsize: list
"""
# set up figure
fig, ax = plt.subplots(figsize=figsize)
# plot heatmap
sns.heatmap(arr_sim,
ax=ax,
yticklabels=['Kq', 'Ks', 'Alp', 'Huz', 'B'],
cmap=sns.color_palette("ch:s=-.2,r=.6"))
# setup overlay axis
ax2 = ax.twinx()
# plot line
ax2.plot(np.arange(0.5, 12.5), df_obs['Strmflw'], color='slateblue')
# plot points on line
ax2.plot(np.arange(0.5, 12.5), df_obs['Strmflw'], color='slateblue', marker='o')
# set axis limits and labels
ax.set_ylim(0, 5)
ax.set_xlim(0, 12)
ax.set_xticklabels(['jan', 'feb', 'mar', 'apr', 'may', 'jun', 'jul', 'aug', 'sep', 'oct', 'nov', 'dec'])
ax2.set_ylabel('Flow Discharge($m^3/s$)')
plt.title(title)
plt.show()
return ax, ax2
C.1.4. plot_annual_heatmap()¶
def plot_annual_heatmap(arr_sim, df_obs, title='', figsize=[14,5]):
"""Plot a sensitivity metric overlain by observed flow..
:param arr_sim: Numpy array of simulated metrics
:param df_obs: Dataframe of mean monthly observed data from sensitivity analysis
:param title: Title of plot
:type title: str
:param figsize: Matplotlib figure size
:type figsize: list
"""
# set up figure
fig, ax = plt.subplots(figsize=figsize)
# plot heatmap
sns.heatmap(arr_sim, ax=ax, cmap=sns.color_palette("YlOrBr"))
# setup overlay axis
ax2 = ax.twinx()
# plot line
ax2.plot(np.arange(0.5, 10.5), df_obs['Strmflw'], color='slateblue')
# plot points on line
ax2.plot(np.arange(0.5, 10.5), df_obs['Strmflw'], color='slateblue', marker='o')
# set up axis lables and limits
ax.set_ylim(0, 5)
ax.set_xlim(0, 10)
ax.set_yticklabels(['Kq', 'Ks', 'Alp', 'Huz', 'B'])
ax.set_xticklabels(range(2000, 2010))
ax2.set_ylabel('Flow Discharge($m^3/s$)')
plt.title(title)
return ax, ax2
C.1.5. plot_varying_heatmap()¶
def plot_varying_heatmap(arr_sim, df_obs, title='', figsize=[14,5]):
"""Plot a sensitivity metric overlain by observed flow..
:param arr_sim: Numpy array of simulated metrics
:param df_obs: Dataframe of mean monthly observed data from sensitivity analysis
:param title: Title of plot
:type title: str
:param figsize: Matplotlib figure size
:type figsize: list
"""
# set up figure
fig, ax = plt.subplots(figsize=figsize)
# plot heatmap
sns.heatmap(arr_sim,
ax=ax,
yticklabels=['Kq', 'Ks', 'Alp', 'Huz', 'B'],
cmap=sns.light_palette("seagreen", as_cmap=True))
n_years = df_obs.shape[0]
# setup overlay axis
ax2 = ax.twinx()
# plot line
ax2.plot(range(0, n_years), df_obs['Strmflw'], color='slateblue')
# plot points on line
ax2.plot(range(0, n_years), df_obs['Strmflw'], color='slateblue', marker='o')
# set up axis lables and limits
ax.set_ylim(0, 5)
ax.set_xlim(-0.5, 119.5)
ax2.set_ylabel('Flow Discharge')
ax.set_xlabel('Number of Months')
plt.title(title)
return ax, ax2
C.1.6. plot_precalibration_flow()¶
def plot_precalibration_flow(df_sim, df_obs, figsize=[10, 4]):
"""Plot flow discharge provided by the ensemble of parameters sets from Pre-Calibration versus the observed
flow data.
:param df_sim: Dataframe of simulated metrics
:param df_obs: Dataframe of mean monthly observed data from sensitivity analysis
:param figsize: Matplotlib figure size
:type figsize: list
"""
# set up figure
fig, ax = plt.subplots(figsize=figsize)
# set axis labels
ax.set_xlabel('Days')
ax.set_ylabel('Flow Discharge')
# plot pre-calibration results
for i in range(df_sim.shape[1]):
plt.plot(range(len(df_sim)), df_sim.iloc[:, i], color="lightgreen", alpha=0.2)
# plot observed
plt.plot(range(len(df_sim)), df_obs['Strmflw'], color="black")
plt.title('Observed vs. Pre-Calibration Outputs')
# customize legend
custom_lines = [Line2D([0], [0], color="lightgreen", lw=4),
Line2D([0], [0], color="black", lw=4)]
plt.legend(custom_lines, ['Pre-Calibration', 'Observed'])
return ax
C.1.7. plot_precalibration_glue()¶
def plot_precalibration_glue(df_precal, df_glue, df_obs, figsize=[10, 4]):
"""Plot flow discharge provided by the ensemble of parameters sets from Pre-Calibration versus the observed
flow data.
:param df_sim: Dataframe of simulated metrics
:param df_obs: Dataframe of mean monthly observed data from sensitivity analysis
:param figsize: Matplotlib figure size
:type figsize: list
"""
# set up figure
fig, ax = plt.subplots(figsize=figsize)
# set axis labels
ax.set_xlabel('Days')
ax.set_ylabel('Flow Discharge')
# plot pre-calibration results
for i in range(df_precal.shape[1]):
plt.plot(range(len(df_precal)), df_precal.iloc[:, i], color="lightgreen", alpha=0.2)
# plot glue
for i in range(df_glue.shape[1]):
plt.plot(range(len(df_glue)), df_glue.iloc[:, i], color="lightblue", alpha=0.2)
# plot observed
plt.plot(range(len(df_precal)), df_obs['Strmflw'], color="black")
plt.title('Observed vs. Sensitivity Analysis Outputs across GLUE/Pre-Calibration')
# customize legend
custom_lines = [Line2D([0], [0], color="lightgreen", lw=4),
Line2D([0], [0], color="lightblue", lw=4),
Line2D([0], [0], color="black", lw=4)]
plt.legend(custom_lines, ['Pre-Calibration', 'GLUE', 'Observed'])
return ax
C.2. fishery_dynamics.ipynb¶
The following are the plotting functions as described in the fishery_dynamics.ipynb Jupyter notebook tutorial.
The following are the necessary package imports to run these functions:
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import patheffects as pe
C.2.1. plot_objective_performance()¶
def plot_objective_performance(objective_performance, profit_solution, robust_solution, figsize=(18, 9)):
"""Plot the identified solutions with regards to their objective performance
in a parallel axis plot
:param objective_performance: Objective performance array
:param profit_solution: Profitable solutions array
:param robust_solution: Robust solutions array
:param figsize: Figure size
:type figsize: tuple
"""
# create the figure object
fig = plt.figure(figsize=figsize)
# set up subplot axis object
ax = fig.add_subplot(1, 1, 1)
# labels where constraint is always 0
objs_labels = ['Net present\nvalue (NPV)',
'Prey population deficit',
'Longest duration\nof low harvest',
'Worst harvest instance',
'Variance of harvest',
'Duration of predator\npopulation collapse']
# normalization across objectives
mins = objective_performance.min(axis=0)
maxs = objective_performance.max(axis=0)
norm_reference = objective_performance.copy()
for i in range(5):
mm = objective_performance[:, i].min()
mx = objective_performance[:, i].max()
if mm != mx:
norm_reference[:, i] = (objective_performance[:, i] - mm) / (mx - mm)
else:
norm_reference[:, i] = 1
# colormap from matplotlib
cmap = plt.cm.get_cmap("Blues")
# plot all solutions
for i in range(len(norm_reference[:, 0])):
ys = np.append(norm_reference[i, :], 1.0)
xs = range(len(ys))
ax.plot(xs, ys, c=cmap(ys[0]), linewidth=2)
# to highlight robust solutions
ys = np.append(norm_reference[profit_solution, :], 1.0) # Most profitable
xs = range(len(ys))
l1 = ax.plot(xs[0:6],
ys[0:6],
c=cmap(ys[0]),
linewidth=3,
label='Most robust in NPV',
path_effects=[pe.Stroke(linewidth=6, foreground='darkgoldenrod'), pe.Normal()])
ys = np.append(norm_reference[robust_solution, :], 1.0) # Most robust in all criteria
xs = range(len(ys))
l2 = ax.plot(xs[0:6],
ys[0:6],
c=cmap(ys[0]),
linewidth=3,
label='Most robust across criteria',
path_effects=[pe.Stroke(linewidth=6, foreground='gold'), pe.Normal()])
# build colorbar
sm = plt.cm.ScalarMappable(cmap=cmap)
sm.set_array([objective_performance[:, 0].min(), objective_performance[:, 0].max()])
cbar = fig.colorbar(sm)
cbar.ax.set_ylabel("\nNet present value (NPV)")
# tick values
minvalues = ["{0:.3f}".format(mins[0]),
"{0:.3f}".format(-mins[1]),
str(-mins[2]),
"{0:.3f}".format(-mins[3]),
"{0:.2f}".format(-mins[4]),
str(0)]
maxvalues = ["{0:.2f}".format(maxs[0]),
"{0:.3f}".format(-maxs[1]),
str(-maxs[2]),
"{0:.2f}".format(maxs[3]),
"{0:.2f}".format(-maxs[4]),
str(0)]
ax.set_ylabel("Preference ->", size=12)
ax.set_yticks([])
ax.set_xticks([0, 1, 2, 3, 4, 5])
ax.set_xticklabels([minvalues[i] + '\n' + objs_labels[i] for i in range(len(objs_labels))])
# make a twin axis for toplabels
ax1 = ax.twiny()
ax1.set_yticks([])
ax1.set_xticks([0, 1, 2, 3, 4, 5])
ax1.set_xticklabels([maxvalues[i] for i in range(len(maxs) + 1)])
return ax, ax1
C.2.2. plot_factor_performance()¶
def plot_factor_performance(param_values, collapse_days, b, m, a):
"""Visualize the performance of our policies in three-dimensional
parametric space.
:param param_values: Saltelli sample array
:param collapse_days: Simulation array
:param b: b parameter boundary interval
:param m: m parameter boundary interval
:param a: a parameter boundary interval
"""
# set colormap
cmap = plt.cm.get_cmap("RdBu_r")
# build figure object
fig = plt.figure(figsize=plt.figaspect(0.5), dpi=600, constrained_layout=True)
# set up scalable colormap
sm = plt.cm.ScalarMappable(cmap=cmap)
# set up subplot for profit maximizing policy
ax1 = fig.add_subplot(1, 2, 1, projection='3d')
# add point data for profit plot
sows = ax1.scatter(param_values[:,1],
param_values[:,6],
param_values[:,0],
c=collapse_days[:,0],
cmap=cmap,
s=0.5)
# add surface data for boundary separating successful and failed states of the world
pts_ineq = ax1.plot_surface(b, m, a, color='black', alpha=0.25, zorder=1)
# add reference point to plot
pt_ref = ax1.scatter(0.5, 0.7, 0.005, c='black', s=50, zorder=0)
# set up plot aesthetics and labels
ax1.set_xlabel("b")
ax1.set_ylabel("m")
ax1.set_zlabel("a")
ax1.set_zlim([0.0, 2.0])
ax1.set_xlim([0.0, 1.0])
ax1.set_ylim([0.0, 1.5])
ax1.xaxis.set_view_interval(0, 0.5)
ax1.set_facecolor('white')
ax1.view_init(12, -17)
ax1.set_title('Profit maximizing policy')
# set up subplot for robust policy
ax2 = fig.add_subplot(1, 2, 2, projection='3d')
# add point data for robust plot
sows = ax2.scatter(param_values[:,1],
param_values[:,6],
param_values[:,0],
c=collapse_days[:,1],
cmap=cmap,
s=0.5)
# add surface data for boundary separating successful and failed states of the world
pts_ineq = ax2.plot_surface(b, m, a, color='black', alpha=0.25, zorder=1)
# add reference point to plot
pt_ref = ax2.scatter(0.5, 0.7, 0.005, c='black', s=50, zorder=0)
# set up plot aesthetics and labels
ax2.set_xlabel("b")
ax2.set_ylabel("m")
ax2.set_zlabel("a")
ax2.set_zlim([0.0, 2.0])
ax2.set_xlim([0.0, 1.0])
ax2.set_ylim([0.0, 1.5])
ax2.xaxis.set_view_interval(0, 0.5)
ax2.set_facecolor('white')
ax2.view_init(12, -17)
ax2.set_title('Robust policy')
# set up colorbar
sm.set_array([collapse_days.min(), collapse_days.max()])
cbar = fig.colorbar(sm)
cbar.set_label('Days with predator collapse')
return ax1, ax2