Source code for pysb.tools.sensitivity_analysis

import os
from itertools import product
import matplotlib.pyplot as plt
from matplotlib import gridspec
import numpy as np
from pysb.bng import generate_equations
import pysb.simulator.base
from pysb.logging import get_logger
import warnings


[docs]class PairwiseSensitivity(object): """ Pairwise sensitivity analysis of model parameters This class calculates the sensitivity of a specified model Observable to changes in pairs of initial species concentrations. The results are stored in matrices described in Attributes. .. warning:: The interface for this class is considered experimental and may change without warning as PySB is updated. Parameters ---------- solver : pysb.simulator.Simulator Simulator instance used to perform models. Must be initialized with tspan argument set. values_to_sample : vector_like Values to sample for each initial concentration of the model.parameters values. objective_function : function A function that returns a scalar value. Used to calculate fraction of changed that is used for calculating sensitivity. See Example. observable : str Observable name used in the objective_function. sens_type: {'params', 'initials', 'all'} Type of sensitivity analysis to perform. sample_list: list List of model pysb.Parameters names to be used. Attributes ---------- b_matrix: numpy.ndarray Matrix of 2-tuples containing (perturbation, species index) b_prime_matrix: numpy.ndarray Same as b_matrix, only where one of the species concentrations is unchanged (i.e. with the single variable contribution removed) index : list List of model parameter names that will be used in analysis index_of_param : dict Dictionary that maps parameters name to index in orig_values array objective_function : Identical to Parameters (see above). orig_vals : numpy.array Original values of the model.Parameters. p_matrix: numpy.ndarray Pairwise sensitivity matrix p_prime_matrix: numpy.ndarray Normalized pairwise sensitivity matrix (in the sense that it contains changes from the baseline, unperturbed case) params_to_run : np.array Parameter sets to be passed to simulator References ---------- 1. Harris, L.A., Nobile, M.S., Pino, J.C., Lubbock, A.L.R., Besozzi, D., Mauri, G., Cazzaniga, P., and Lopez, C.F. 2017. GPU-powered model analysis with PySB/cupSODA. Bioinformatics 33, pp.3492-3494. https://academic.oup.com/bioinformatics/article/33/21/3492/3896987 Examples -------- Sensitivity analysis on the Tyson cell cycle model >>> from pysb.examples.tyson_oscillator import model >>> import numpy as np >>> from pysb.simulator.scipyode import ScipyOdeSimulator >>> np.set_printoptions(precision=4, suppress=True) >>> tspan=np.linspace(0, 200, 201) >>> observable = 'Y3' >>> values_to_sample = [.8, 1.2] >>> def obj_func_cell_cycle(out): ... timestep = tspan[:-1] ... y = out[:-1] - out[1:] ... freq = 0 ... local_times = [] ... prev = y[0] ... for n in range(1, len(y)): ... if y[n] > 0 > prev: ... local_times.append(timestep[n]) ... freq += 1 ... prev = y[n] ... local_times = np.array(local_times) ... local_freq = np.average(local_times)/len(local_times)*2 ... return local_freq >>> solver = ScipyOdeSimulator(model, tspan, integrator='lsoda',\ integrator_options={'atol' : 1e-8,\ 'rtol' : 1e-8,\ 'mxstep' :20000}) >>> sens = PairwiseSensitivity(\ values_to_sample=values_to_sample,\ observable=observable,\ objective_function=obj_func_cell_cycle,\ solver=solver\ ) >>> print(sens.b_matrix) [[((0.8, 'cdc0'), (0.8, 'cdc0')) ((0.8, 'cdc0'), (1.2, 'cdc0')) ((0.8, 'cdc0'), (0.8, 'cyc0')) ((0.8, 'cdc0'), (1.2, 'cyc0'))] [((1.2, 'cdc0'), (0.8, 'cdc0')) ((1.2, 'cdc0'), (1.2, 'cdc0')) ((1.2, 'cdc0'), (0.8, 'cyc0')) ((1.2, 'cdc0'), (1.2, 'cyc0'))] [((0.8, 'cyc0'), (0.8, 'cdc0')) ((0.8, 'cyc0'), (1.2, 'cdc0')) ((0.8, 'cyc0'), (0.8, 'cyc0')) ((0.8, 'cyc0'), (1.2, 'cyc0'))] [((1.2, 'cyc0'), (0.8, 'cdc0')) ((1.2, 'cyc0'), (1.2, 'cdc0')) ((1.2, 'cyc0'), (0.8, 'cyc0')) ((1.2, 'cyc0'), (1.2, 'cyc0'))]] >>> sens.run() >>> print(sens.p_matrix)#doctest: +NORMALIZE_WHITESPACE [[ 0. 0. 5.0243 -4.5381] [ 0. 0. 5.0243 -4.5381] [ 5.0243 5.0243 0. 0. ] [-4.5381 -4.5381 0. 0. ]] >>> print(sens.p_prime_matrix) #doctest: +NORMALIZE_WHITESPACE [[ 0. 0. 5.0243 -4.5381] [ 0. 0. 5.0243 -4.5381] [ 0. 0. 0. 0. ] [ 0. 0. 0. 0. ]] >>> print(sens.p_matrix - sens.p_prime_matrix) \ #doctest: +NORMALIZE_WHITESPACE [[ 0. 0. 0. 0. ] [ 0. 0. 0. 0. ] [ 5.0243 5.0243 0. 0. ] [-4.5381 -4.5381 0. 0. ]] >>> sens.create_boxplot_and_heatplot() #doctest: +SKIP >>> values_to_sample = [.9, 1.1] >>> sens = PairwiseSensitivity(\ values_to_sample=values_to_sample,\ observable=observable,\ objective_function=obj_func_cell_cycle,\ solver=solver,\ sens_type='params'\ ) >>> print(sens.b_matrix.shape == (14, 14)) True >>> sens.run() >>> print(sens.p_matrix)#doctest: +NORMALIZE_WHITESPACE [[ 0. 0. 13.6596 13.6596 24.3955 4.7909 16.4603 11.3258 0.1621 31.2804 13.6596 13.6596 13.6596 13.6596] [ 0. 0. -10.3728 -10.3728 -3.7277 -14.9803 -7.2934 -12.2416 -18.3144 0. -10.3728 -10.3728 -10.3728 -10.3728] [ 13.6596 -10.3728 0. 0. 7.3582 -6.483 3.0794 -2.269 -10.6969 12.7261 0. 0. 0. 0. ] [ 13.6596 -10.3728 0. 0. 7.3582 -6.483 3.0794 -2.269 -10.6969 12.7261 0. 0. 0. 0. ] [ 24.3955 -3.7277 7.3582 7.3582 0. 0. 10.859 5.2577 -4.376 23.2285 7.3582 7.3582 7.3582 7.3582] [ 4.7909 -14.9803 -6.483 -6.483 0. 0. -3.4036 -9.0762 -15.2185 3.8574 -6.483 -6.483 -6.483 -6.483 ] [ 16.4603 -7.2934 3.0794 3.0794 10.859 -3.4036 0. 0. -7.9417 15.5267 3.0794 3.0794 3.0794 3.0794] [ 11.3258 -12.2416 -2.269 -2.269 5.2577 -9.0762 0. 0. -13.128 10.859 -2.269 -2.269 -2.269 -2.269 ] [ 0.1621 -18.3144 -10.6969 -10.6969 -4.376 -15.2185 -7.9417 -13.128 0. 0. -10.6969 -10.6969 -10.6969 -10.6969] [ 31.2804 0. 12.7261 12.7261 23.2285 3.8574 15.5267 10.859 0. 0. 12.7261 12.7261 12.7261 12.7261] [ 13.6596 -10.3728 0. 0. 7.3582 -6.483 3.0794 -2.269 -10.6969 12.7261 0. 0. 0. 0. ] [ 13.6596 -10.3728 0. 0. 7.3582 -6.483 3.0794 -2.269 -10.6969 12.7261 0. 0. 0. 0. ] [ 13.6596 -10.3728 0. 0. 7.3582 -6.483 3.0794 -2.269 -10.6969 12.7261 0. 0. 0. 0. ] [ 13.6596 -10.3728 0. 0. 7.3582 -6.483 3.0794 -2.269 -10.6969 12.7261 0. 0. 0. 0. ]] >>> print(sens.p_matrix - sens.p_prime_matrix) \ #doctest: +NORMALIZE_WHITESPACE [[ 0. 0. 13.6596 13.6596 17.0373 11.2739 13.3809 13.5948 10.859 18.5543 13.6596 13.6596 13.6596 13.6596] [ 0. 0. -10.3728 -10.3728 -11.0859 -8.4973 -10.3728 -9.9725 -7.6175 -12.7261 -10.3728 -10.3728 -10.3728 -10.3728] [ 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ] [ 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ] [ 10.7358 6.6451 7.3582 7.3582 0. 0. 7.7796 7.5267 6.3209 10.5024 7.3582 7.3582 7.3582 7.3582] [ -8.8687 -4.6075 -6.483 -6.483 0. 0. -6.483 -6.8071 -4.5215 -8.8687 -6.483 -6.483 -6.483 -6.483 ] [ 2.8006 3.0794 3.0794 3.0794 3.5008 3.0794 0. 0. 2.7553 2.8006 3.0794 3.0794 3.0794 3.0794] [ -2.3339 -1.8688 -2.269 -2.269 -2.1005 -2.5932 0. 0. -2.4311 -1.8671 -2.269 -2.269 -2.269 -2.269 ] [-13.4976 -7.9417 -10.6969 -10.6969 -11.7342 -8.7355 -11.0211 -10.859 0. 0. -10.6969 -10.6969 -10.6969 -10.6969] [ 17.6207 10.3728 12.7261 12.7261 15.8703 10.3404 12.4473 13.128 0. 0. 12.7261 12.7261 12.7261 12.7261] [ 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ] [ 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ] [ 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ] [ 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. ]] >>> sens.create_boxplot_and_heatplot() #doctest: +SKIP >>> sens = PairwiseSensitivity(\ values_to_sample=values_to_sample,\ observable=observable,\ objective_function=obj_func_cell_cycle,\ solver=solver,\ sample_list=['k1', 'cdc0']\ ) >>> print(sens.b_matrix) [[((0.9, 'k1'), (0.9, 'k1')) ((0.9, 'k1'), (1.1, 'k1')) ((0.9, 'k1'), (0.9, 'cdc0')) ((0.9, 'k1'), (1.1, 'cdc0'))] [((1.1, 'k1'), (0.9, 'k1')) ((1.1, 'k1'), (1.1, 'k1')) ((1.1, 'k1'), (0.9, 'cdc0')) ((1.1, 'k1'), (1.1, 'cdc0'))] [((0.9, 'cdc0'), (0.9, 'k1')) ((0.9, 'cdc0'), (1.1, 'k1')) ((0.9, 'cdc0'), (0.9, 'cdc0')) ((0.9, 'cdc0'), (1.1, 'cdc0'))] [((1.1, 'cdc0'), (0.9, 'k1')) ((1.1, 'cdc0'), (1.1, 'k1')) ((1.1, 'cdc0'), (0.9, 'cdc0')) ((1.1, 'cdc0'), (1.1, 'cdc0'))]] """ def __init__(self, solver, values_to_sample, objective_function, observable, sens_type='initials', sample_list=None): if not isinstance(solver, pysb.simulator.base.Simulator): raise TypeError("solver must be a pysb.simulator object") self._model = solver.model self._logger = get_logger(__name__, model=self._model) self._logger.info('%s created for observable %s' % ( self.__class__.__name__, observable)) generate_equations(self._model) self._values_to_sample = values_to_sample self._solver = solver self.objective_function = objective_function self._observable = observable self._sens_type = sens_type if self._sens_type not in ('params', 'initials', 'all'): if sample_list is None: raise ValueError("Please provide 'sens_type' or 'sample_list'") if sample_list is not None: _valid_options = [i.name for i in self._model.parameters] for i in sample_list: if i not in _valid_options: raise ValueError("{} not in model.parameters".format(i)) self.index = sample_list elif self._sens_type == 'params': self.index = [i.name for i in self._model.parameters_rules()] elif self._sens_type == 'initials': self.index = [i[1].name for i in self._model.initial_conditions] elif self._sens_type == 'all': self.index = [i.name for i in self._model.parameters] self.orig_vals = [i.value for i in self._model.parameters] self.index_of_param = {i.name: n for n, i in enumerate(self._model.parameters)} self._n_sam = len(self._values_to_sample) self._n_species = len(self.index) self._nm = self._n_species * self._n_sam self._size_of_matrix = self._nm ** 2 self._shape_of_matrix = self._nm, self._nm # Outputs self.b_matrix = [] self.b_prime_matrix = [] self.params_to_run = self._setup_simulations() self.p_prime_matrix = np.zeros(self._size_of_matrix) self.p_matrix = np.zeros(self._size_of_matrix) # Stores the objective function value for the original unperturbed # model self._objective_fn_standard = None @property def sensitivity_multiset(self): """ Sensitivity analysis multiset (also called "Q" matrix) Returns ------- list List of lists containing the sensitivity analysis multiset """ sens_ij_nm = [] sens_matrix = self.p_matrix - self.p_prime_matrix # separate each species sensitivity for j in range(0, self._nm, self._n_sam): per_species1 = [] for i in range(0, self._nm, self._n_sam): if i != j: tmp = sens_matrix[j:j + self._n_sam, i:i + self._n_sam].copy() per_species1.append(tmp) sens_ij_nm.append(per_species1) return sens_ij_nm def _calculate_objective(self, function_value): """ Calculate fraction of change for obj value and standard Parameters ---------- function_value : scalar scalar value provided by objective function """ return (self.objective_function(function_value) - self._objective_fn_standard) / \ self._objective_fn_standard * 100.
[docs] def run(self, save_name=None, out_dir=None): """ Run sensitivity analysis Parameters ---------- save_name : str, optional prefix of saved files out_dir : str, optional location to save output if required """ self._solver.initials = None self._solver.param_values = None sim_results = self._solver.run(param_values=None, initials=None) self._objective_fn_standard = self.objective_function( np.array(sim_results.observables[self._observable]) ) traj = self._solver.run(param_values=self.params_to_run) output = np.array(traj.observables)[self._observable].T p_matrix = np.zeros(self._size_of_matrix) p_prime_matrix = np.zeros(self._size_of_matrix) counter = 0 # places values in p matrix that are unique for i in range(len(p_matrix)): if i in self._b_index: p_matrix[i] = self._calculate_objective(output[:, counter]) counter += 1 # places values in p matrix that are duplicated for i in range(len(p_matrix)): if i in self._b_prime_not_in_b: p_prime_matrix[i] = self._calculate_objective( output[:, self._b_prime_not_in_b[i]] ) elif i in self._b_prime_in_b: p_prime_matrix[i] = p_matrix[self._b_prime_in_b[i]] p_matrix = p_matrix.reshape(self._shape_of_matrix) # Project the mirrored image self.p_matrix = p_matrix + p_matrix.T self.p_prime_matrix = p_prime_matrix.reshape(self._shape_of_matrix) # save output if desired if save_name is not None: if out_dir is None: out_dir = '.' elif not os.path.exists(out_dir): os.mkdir(out_dir) p_name = os.path.join(out_dir, '{}_p_matrix.csv'.format(save_name)) p_prime_name = os.path.join( out_dir, '{}_p_prime_matrix.csv'.format(save_name)) self._logger.debug("Saving p matrix and p' matrix to {} and {}". format(p_name, p_prime_name)) np.savetxt(p_name, self.p_matrix) np.savetxt(p_prime_name, self.p_prime_matrix)
def _g_function(self): """ Create sample matrix, index of samples values, and shows bij """ counter = -1 sampled_values_index = set() bij_unique = dict() sampling_matrix = np.zeros( (self._size_of_matrix, len(self.orig_vals))) sampling_matrix[:, :] = self.orig_vals matrix = self.b_matrix for j in range(len(matrix)): for i in matrix[j, :]: sigma_i, index_i = i[0] sigma_j, index_j = i[1] s_1 = str(i[0]) + str(i[1]) s_2 = str(i[1]) + str(i[0]) counter += 1 if index_i == index_j \ or s_1 in bij_unique or s_2 in bij_unique: continue else: x = self.index_of_param[index_i] y = self.index_of_param[index_j] sampling_matrix[counter, x] *= sigma_i sampling_matrix[counter, y] *= sigma_j bij_unique[s_1] = counter bij_unique[s_2] = counter sampled_values_index.add(counter) return sampling_matrix, sampled_values_index, bij_unique def _setup_simulations(self): """ Create initial conditions matrix for sensitivity analysis Returns ------- numpy.ndarray Matrix of initial conditions """ # create matrix (cartesian product of sample vals vs index of species a_matrix = cartesian_product(self._values_to_sample, self.index).T.flatten() # creates matrix b self.b_matrix = cartesian_product(a_matrix, a_matrix) # create matrix a' a_prime = cartesian_product(np.ones(self._n_sam), self.index).T.flatten() # creates matrix b prime self.b_prime_matrix = cartesian_product(a_prime, a_matrix) b_to_run, self._b_index, in_b = self._g_function() n_b_index = len(self._b_index) b_prime = np.zeros((self._size_of_matrix, len(self.orig_vals))) b_prime[:, :] = self.orig_vals counter = -1 bp_not_in_b_raw = set() bp_dict = dict() bp_not_in_b_dict = dict() bp_not_in_b_visited = dict() new_sim_counter = -1 # checks for and removes duplicates of simulations initial conditions for j in range(len(self.b_prime_matrix)): for i in self.b_prime_matrix[j, :]: sigma_i, index_i = i[0] sigma_j, index_j = i[1] s_1 = str(i[0]) + str(i[1]) counter += 1 # no need for doing if same index if index_i == index_j: continue # pointing to same simulation if already in b elif s_1 in in_b: bp_dict[counter] = in_b[s_1] elif s_1 in bp_not_in_b_visited: bp_not_in_b_dict[counter] = bp_not_in_b_visited[s_1] else: new_sim_counter += 1 x = self.index_of_param[index_i] y = self.index_of_param[index_j] b_prime[new_sim_counter, x] *= sigma_i b_prime[new_sim_counter, y] *= sigma_j bp_not_in_b_visited[s_1] = new_sim_counter + n_b_index bp_not_in_b_dict[counter] = new_sim_counter + n_b_index bp_not_in_b_raw.add(new_sim_counter) self._b_prime_in_b = bp_dict self._b_prime_not_in_b = bp_not_in_b_dict x = b_to_run[sorted(self._b_index)] y = b_prime[sorted(bp_not_in_b_raw)] simulations = np.vstack((x, y)) self._logger.debug("{} simulations to run".format(len(simulations))) return simulations
[docs] def create_plot_p_h_pprime(self, save_name=None, out_dir=None, show=False): """ Plot of P, H(B), and P' See :class:`PairwiseSensitivity` attributes for descriptions of these matrices Parameters ---------- save_name : str, optional name to save figure as out_dir : str, optional location to save figure show : bool show the plot if True Returns ------- matplotlib.figure.Figure The matplotlib figure object for further adjustments, if required """ colors = 'seismic' sens_matrix = self.p_matrix - self.p_prime_matrix v_max = max(np.abs(self.p_matrix.min()), self.p_matrix.max()) v_min = -1 * v_max fig = plt.figure(figsize=(12, 6)) ax1 = fig.add_subplot(131) ax2 = fig.add_subplot(132) ax3 = fig.add_subplot(133) ax1.imshow(self.p_matrix, interpolation='nearest', origin='upper', cmap=plt.get_cmap(colors), vmin=v_min, vmax=v_max, extent=[0, self._nm, 0, self._nm]) ax2.imshow(self.p_prime_matrix, interpolation='nearest', origin='upper', cmap=plt.get_cmap(colors), vmin=v_min, vmax=v_max, extent=[0, self._nm, 0, self._nm]) ax3.imshow(sens_matrix, interpolation='nearest', origin='upper', cmap=plt.get_cmap(colors), vmin=v_min, vmax=v_max) for ax in (ax1, ax2, ax3): ax.set_xticks([]) ax.set_yticks([]) ax1.set_title('P', fontsize=22) ax2.set_title('H(B\')', fontsize=22) ax3.set_title('P\' = P - H(B\')', fontsize=22) fig.subplots_adjust(wspace=0, hspace=0.0) if save_name is not None: if out_dir is None: out_dir = '.' if not os.path.exists(out_dir): os.mkdir(out_dir) fig.savefig(os.path.join(out_dir, '{}_P_H_P_prime.png'.format(save_name)), bbox_inches='tight') if show: plt.show() plt.close() return fig
[docs] def create_individual_pairwise_plots(self, save_name=None, out_dir=None, show=False): """ Single plot containing heat plot of each specie pair Parameters ---------- save_name : str, optional name ot save figure as out_dir : str, optional output directory show : bool show figure Returns ------- matplotlib.figure.Figure The matplotlib figure object for further adjustments, if required """ colors = 'seismic' sens_matrix = self.p_matrix - self.p_prime_matrix v_max = max(np.abs(self.p_matrix.min()), self.p_matrix.max()) v_min = -1 * v_max fig = plt.figure(figsize=(self._n_species + 6, self._n_species + 6)) gs = gridspec.GridSpec(self._n_species, self._n_species) # creates a plot of each species vs each species # adds space between plots so you can zoom in on output pairs for n, j in enumerate(range(0, self._nm, self._n_sam)): for m, i in enumerate(range(0, self._nm, self._n_sam)): ax2 = plt.subplot(gs[n, m]) if n == 0: ax2.set_xlabel(self.index[m], fontsize=20) ax2.xaxis.set_label_position('top') if m == 0: ax2.set_ylabel(self.index[n], fontsize=20) plt.xticks([]) plt.yticks([]) if i != j: tmp = sens_matrix[j:j + self._n_sam, i:i + self._n_sam].copy() ax2.imshow(tmp, interpolation='nearest', origin='upper', cmap=plt.get_cmap(colors), vmin=v_min, vmax=v_max) else: ax2.imshow(np.zeros((self._n_sam, self._n_sam)), interpolation='nearest', origin='upper', cmap=plt.get_cmap(colors), vmin=v_min, vmax=v_max) plt.tight_layout() if save_name is not None: if out_dir is None: out_dir = '.' if not os.path.exists(out_dir): os.mkdir(out_dir) plt.savefig(os.path.join(out_dir, '{}_subplots.png'.format(save_name)), bbox_inches='tight') if show: plt.show() plt.close() return fig
[docs] def create_boxplot_and_heatplot(self, x_axis_label=None, save_name=None, out_dir=None, show=False): """ Heat map and box plot of sensitivities Parameters ---------- x_axis_label : str, optional label for x asis save_name : str, optional name of figure to save out_dir : str, option output directory to save figures show : bool Show plot if True Returns ------- matplotlib.figure.Figure The matplotlib figure object for further adjustments, if required """ colors = 'seismic' sens_ij_nm = self.sensitivity_multiset # Create heatmap and boxplot of data fig = plt.figure(figsize=(14, 10)) plt.subplots_adjust(hspace=0.1) # use gridspec to scale colorbar nicely outer = gridspec.GridSpec(2, 1, width_ratios=[1.], height_ratios=[0.03, 1]) gs1 = gridspec.GridSpecFromSubplotSpec(1, 1, subplot_spec=outer[0]) gs2 = gridspec.GridSpecFromSubplotSpec(2, 1, subplot_spec=outer[1], hspace=.35) ax0 = plt.subplot(gs1[0]) ax1 = plt.subplot(gs2[0]) # scale the colors to minimum or maximum of p matrix v_max = max(np.abs(self.p_matrix.min()), self.p_matrix.max()) v_min = -1 * v_max # create heatmap of sensitivities im = ax1.imshow(self.p_matrix, interpolation='nearest', origin='upper', cmap=plt.get_cmap(colors), vmin=v_min, vmax=v_max, extent=[0, self._nm, 0, self._nm]) shape_label = np.arange(self._n_sam / 2, self._nm, self._n_sam) plt.xticks(shape_label, self.index, rotation='vertical', fontsize=12) plt.yticks(shape_label, reversed(self.index), fontsize=12) x_ticks = list(range(0, self._nm, self._n_sam)) ax1.set_xticks(x_ticks, minor=True) ax1.set_yticks(x_ticks, minor=True) plt.grid(True, which='minor', linestyle='--') color_bar = plt.colorbar(im, cax=ax0, orientation='horizontal', use_gridspec=True) color_bar.set_label('% change', y=1, labelpad=5) color_bar.ax.xaxis.set_label_position('top') ticks = np.linspace(v_min, v_max, 5, dtype=int) color_bar.set_ticks(ticks) color_bar.ax.set_xticklabels(ticks) # create boxplot of single parameter sensitivities ax2 = plt.subplot(gs2[1]) x = [np.array(mat).flatten() for mat in sens_ij_nm[::-1]] ax2.boxplot(x, vert=False, labels=None, showfliers=True, whis=(0, 100)) ax2.set_xlim(v_min - 2, v_max + 2) if x_axis_label is not None: ax2.set_xlabel(x_axis_label, fontsize=12) plt.setp(ax2, yticklabels=reversed(self.index)) ax2.yaxis.tick_left() ax2.set_aspect(1. / ax2.get_data_ratio(), adjustable='box') if save_name is not None: if out_dir is None: out_dir = '.' if not os.path.exists(out_dir): os.mkdir(out_dir) plt.savefig(os.path.join(out_dir, save_name + '.png'), bbox_inches='tight') plt.savefig(os.path.join(out_dir, save_name + '.eps'), bbox_inches='tight') plt.savefig(os.path.join(out_dir, save_name + '.svg'), bbox_inches='tight') if show: plt.show() plt.close() return fig
[docs]def cartesian_product(array_1, array_2): """ Cartesian product between two lists Parameters ---------- array_1 : list_like array_2 : list_like Returns ------- np.array array of shape (len(array_1), len(array_2)) """ a = list(product(array_1, array_2)) a = np.asarray(a, dtype=','.join(['object'] * len(a[0]))) return a.reshape(len(array_1), len(array_2))
[docs]class InitialsSensitivity(PairwiseSensitivity): """ Deprecated; use :class:`PairwiseSensitivity` instead. """ def __init__(self, *args, **kwargs): warnings.warn("InitialsSensitivity will be removed in a future " "version of PySB. Use PairwiseSensitivity instead.", DeprecationWarning, stacklevel=2) super(InitialsSensitivity, self).__init__(*args, **kwargs)