Source code for TransportMaps.Diagnostics.Plotting

#
# This file is part of TransportMaps.
#
# TransportMaps is free software: you can redistribute it and/or modify
# it under the terms of the GNU Lesser General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# TransportMaps is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU Lesser General Public License for more details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with TransportMaps.  If not, see <http://www.gnu.org/licenses/>.
#
# Transport Maps Library
# Copyright (C) 2015-2018 Massachusetts Institute of Technology
# Uncertainty Quantification group
# Department of Aeronautics and Astronautics
#
# Author: Transport Map Team
# Website: transportmaps.mit.edu
# Support: transportmaps.mit.edu/qa/
#

import numpy as np
import scipy.stats as scistat

from TransportMaps.External import PLOT_SUPPORT
from TransportMaps.MPI import mpi_map
from TransportMaps.ObjectBase import TMO
from TransportMaps.Distributions.FrozenDistributions import StandardNormalDistribution

if PLOT_SUPPORT:
    import matplotlib as mpl
    import matplotlib.pyplot as plt

__all__ = ['AlignedConditionalsObject', 'computeAlignedConditionals',
           'plotAlignedConditionals',
           'RandomConditionalsObject', 'computeRandomConditionals',
           'plotRandomConditionals',
           'plotAlignedSliceMap', 'plotAlignedMarginals',
           'plotAlignedScatters',
           'plotGradXMap', 'plotLinearityMap', 'niceSpy']

def frexp10(x,numdecimals=1):  #Returns mantissa and exponent of a given number
    exp = int( np.floor( np.log10(x) ) )
    factor = 10.**(numdecimals)
    mant =   x / 10.**exp
    return mant, exp

[docs]class AlignedConditionalsObject(TMO): def __init__(self, nplots, X_list, Y_list, pdfEval_list): self.nplots = nplots self.X_list = X_list self.Y_list = Y_list self.pdfEval_list = pdfEval_list
[docs]def computeAlignedConditionals(distribution, dimensions_vec=0, range_vec=[-3,3], numPointsXax=30, do_diag=True, mpi_pool=None): r""" Compute the conditionals aligned with the axis Args: distribution (:class:`Distribution<Distribution>`): distribution :math:`\pi` dimensions_vec (list of int): list of dimensions to be displayed. Default 0: display 10 dimensions at most. range_vec (:class:`list<list>`): range to be displayed. Either a :class:`list<list>` [2] of integers, or a :class:`list<list>` [d] of :class:`list<list>` [2] of integers. numPointsXax (int): number of points for each axis. do_diag (bool): whether to include the one dimensional conditionals on the diagonal mpi_pool (:class:`mpi_map.MPI_Pool<mpi_map.MPI_Pool>`): pool of processes Returns: (:class:`AlignedConditionalsObject<AlignedConditionalsObject>`) -- object storing all the necessary evaluated values """ dim = distribution.dim #Assign optional parameters if not defined if dimensions_vec == 0: dimensions_vec = list(range( min(dim,30) )) else: dimensions_vec = list(dimensions_vec) nplots = len(dimensions_vec) # number of subplots per dimension # Prepare range_vec if not isinstance(range_vec[0], list): range_vec = [range_vec] * distribution.dim cntr_vec = [ (a+b)/2 for a, b in range_vec ] # Create list of points for every 2d plot X_list = [] Y_list = [] tempMat_list = [] for ii in range(nplots): for jj in range(ii+1): if jj<ii: # 2d plot #The x variable is jj, whereas the y-variable is ii. #Create evaluation grid xlin = np.linspace( range_vec[jj][0] , range_vec[jj][1], num=numPointsXax) ylin = np.linspace( range_vec[ii][0] , range_vec[ii][1], num=numPointsXax) X, Y = np.meshgrid(xlin, ylin) X_eval = X.flat #Flat iterator over the array Y_eval = Y.flat #Slice 2d pdf npointEval =len(X_eval) #Number of points where to evaluate the slice tempMat = np.tile(cntr_vec,(npointEval,1)) tempMat[:,dimensions_vec[ii]] = Y_eval tempMat[:,dimensions_vec[jj]] = X_eval # Append X_list.append(X) Y_list.append(Y) tempMat_list.append(tempMat) # Create list of points for every 1d plot if do_diag: for jj in range(nplots): #1D plot xlin = np.linspace( range_vec[jj][0] , range_vec[jj][1], num=numPointsXax) #Slice 1d pdf npointEval =len(xlin) #Number of points where to evaluate the slice tempMat = np.tile(cntr_vec,(npointEval,1)) tempMat[:,dimensions_vec[jj]] = xlin # Append X_list.append(xlin) tempMat_list.append(tempMat) # Evaluate npoints_list = [0] + [ tmat.shape[0] for tmat in tempMat_list ] npoints_pos = np.cumsum(npoints_list) tempMat_mat = np.concatenate(tempMat_list, axis=0) scatter_tuple = (['x'], [tempMat_mat]) pdfEval_mat = mpi_map("pdf", obj=distribution, scatter_tuple=scatter_tuple, mpi_pool=mpi_pool) pdfEval_list = [ pdfEval_mat[npoints_pos[i]:npoints_pos[i+1]] for i in range(len(tempMat_list)) ] data = AlignedConditionalsObject(nplots, X_list, Y_list, pdfEval_list) return data
[docs]def plotAlignedConditionals(distribution=None, data=None, dimensions_vec=0, range_vec=[-3,3], numPointsXax = 30, numCont = 15, figname=None, show_flag=True, do_diag=True, show_title=False, show_axis=True, title='Aligned conditionals', vartitles=None, fig=None, ret_handles=False, mpi_pool=None): r""" Plot the conditionals aligned with the axis Args: distribution (:class:`Distribution<Distribution>`): distribution :math:`\pi` data (:class:`AlignedConditionalsObject`): output of :func:`computeAlignedConditionals` dimensions_vec (list of int): list of dimensions to be displayed. Default 0: display 10 dimensions at most. range_vec (:class:`list<list>`): range to be displayed. Either a :class:`list<list>` [2] of integers, or a :class:`list<list>` [d] of :class:`list<list>` [2] of integers. numPointsXax (int): number of points for each axis. numCont (int): number of contours in the contour plots. figname (str): if defined, store the figure in the provided path. show_flag (bool): whether to show the plot before returning do_diag (bool): whether to include the one dimensional conditionals on the diagonal show_title (bool): whether to show the title show_axis (bool): whether to show the axis title (str): title for the figure vartitles (list): list of titles for each variable fig (figure): matplotlib figure object if one wants to re-useit. mpi_pool (:class:`mpi_map.MPI_Pool<mpi_map.MPI_Pool>`): pool of processes ret_handles (bool): whether to return the axes handles """ if distribution is None and data is None: raise ValueError("Either distribution or data need to be provided") if vartitles is not None and len(vartitles) != distribution.dim: raise ValueError("You need to provide the right number of titles.") if not PLOT_SUPPORT: raise RuntimeError("Plotting is not suported") if data is None: data = computeAlignedConditionals( distribution, dimensions_vec=dimensions_vec, range_vec=range_vec, numPointsXax=numPointsXax, mpi_pool=mpi_pool) nplots = data.nplots X_list = data.X_list Y_list = data.Y_list pdfEval_list = data.pdfEval_list if fig is None: fig = plt.figure(figsize=(8,8)); else: plt.figure(fig.number); fig.clf(); naxes = nplots if do_diag else nplots - 1 axarr = np.empty((naxes,naxes), dtype=object) for i in range(naxes): for j in range(1,naxes+1): if i == 0 and j == 1: axarr[0,0] = plt.subplot( naxes, naxes, (i*naxes)+j); else: axarr[i,j-1] = plt.subplot( naxes, naxes, (i*naxes)+j, sharex=axarr[i-1,j-1]); if show_title: fig.suptitle(title,fontsize=16); plt.set_cmap( 'jet' ); #Set colormap #plt.subplots_adjust(wspace=.15) #Reduce horizontal distance between subplots colorboxes = 'black' lw_spines = .5 # Plot idx = 0 for ii in range(nplots): for jj in range(ii+1): if jj<ii: # 2d plot X = X_list[idx] Y = Y_list[idx] pdfEval = pdfEval_list[idx] #Reshape results for plotting pdfEval_mat = pdfEval.reshape( X.shape[0], X.shape[1] ) #Set current axes if do_diag: ax = axarr[ii,jj] else: ax = axarr[ii-1,jj] cont_h = ax.contourf(X, Y, pdfEval_mat, numCont); cmin, cmax = cont_h.get_clim(); cmin = cont_h.levels[2] # cont_h.cmap.set_under('white'); #Set color below levels cont_h.set_clim(cmin, cmax); ax.spines['top'].set_color(colorboxes); ax.spines['right'].set_color(colorboxes); ax.spines['left'].set_color(colorboxes); ax.spines['bottom'].set_color(colorboxes); ax.spines['top'].set_linewidth(lw_spines); ax.spines['right'].set_linewidth(lw_spines); ax.spines['left'].set_linewidth(lw_spines); ax.spines['bottom'].set_linewidth(lw_spines); ax.get_xaxis().set_visible(False); ax.get_yaxis().set_visible(False); ax.set_aspect('equal'); if not show_axis: ax.axis('off') idx += 1 if do_diag: for jj in range(nplots): #1D plot xlin = X_list[idx] pdfEval = pdfEval_list[idx] maxPdf = np.max(pdfEval) #Set current axes ax = axarr[jj,jj] ax.plot( xlin , pdfEval ); ax.plot( xlin[0] , maxPdf*1.1 , color ='white' ); ax.get_xaxis().set_visible(False); ax.get_yaxis().set_visible(False); ax.spines['top'].set_color(colorboxes); ax.spines['right'].set_color(colorboxes); ax.spines['left'].set_color(colorboxes); ax.spines['bottom'].set_color(colorboxes); ax.spines['top'].set_linewidth(lw_spines); ax.spines['right'].set_linewidth(lw_spines); ax.spines['left'].set_linewidth(lw_spines); ax.spines['bottom'].set_linewidth(lw_spines); if vartitles is not None: tit = vartitles[jj] ax.set_title(tit); idx += 1 for ii in range(naxes): for jj in range(ii+1,naxes): axarr[ii,jj].set_visible(False); # fig.canvas.draw(); if not figname is None: plt.savefig(figname, bbox_inches = 'tight') else: if show_flag: plt.show(False); if ret_handles: axhandles = {'axarr': axarr} return fig, axhandles else: return fig
[docs]def plotAlignedSliceMap(tr_map, dimensions_vec=0, pointEval=0, range_vec=[-4,4], numPointsXax = 30, numCont = 30, figname=None, show_flag=True, tickslabelsize = 6, show_title=False, fig=None, mpi_pool=None): r""" Plot the conditionals aligned with the axis Args: tr_map (:class:`TriangularTransportMap`): Triangular transport map dimensions_vec (list of int): list of dimensions to be displayed. Default 0: display 10 dimensions at most. pointEval (:class:`ndarray`[:math:`d`]): anchor point. Default is zero. range_vec (:class:`list<list>`): range to be displayed. Either a :class:`list<list>` [2] of integers, or a :class:`list<list>` [d] of :class:`list<list>` [2] of integers. numPointsXax (int): number of points for each axis. numCont (int): number of contours in the contour plots. figname (str): if defined, store the figure in the provided path. show_flag (bool): whether to show the plot before returning show_title (bool): whether to show a title on top of the figure fig (figure): matplotlib figure object if one wants to re-use it. mpi_pool (:class:`mpi_map.MPI_Pool<mpi_map.MPI_Pool>`): pool of processes """ if not PLOT_SUPPORT: raise RuntimeError("Plotting is not suported") dim = tr_map.dim #Assign optional parameters if not defined if dimensions_vec == 0: dimensions_vec = range( min(dim,10) ) if pointEval ==0: pointEval=np.zeros(dim) nplots = len(dimensions_vec) # number of subplots per dimension if fig is None: fig = plt.figure(figsize=(8,8)) else: plt.figure(fig.number) fig.clf() axarr = np.empty((nplots,nplots), dtype=object) for i in range(nplots): for j in range(1,nplots+1): if i != 0 or j != 0: axarr[i,j-1] = plt.subplot( nplots, nplots, (i*nplots)+j, sharex=axarr[0,0]) else: axarr[i,j-1] = plt.subplot( nplots, nplots, (i*nplots)+j) # Prepare range_vec if isinstance(range_vec[0], int): range_vec = [range_vec] * tr_map.dim if show_title: fig.suptitle('Slices map along coordinate axes',fontsize=16) #plt.set_cmap( 'jet' ) #Set colormap plt.set_cmap( 'Blues' ) #plt.subplots_adjust(wspace=.15) #Reduce horizontal distance between subplots colorboxes = 'black' lw_spines = .5 # Create list of points (for every plot) X_list = [] Y_list = [] tempMat_list = [] for ii in range(nplots): for jj in range(ii+1): if jj<ii: # 2d plot #The x variable is jj, whereas the y-variable is ii. #Create evaluation grid xlin = np.linspace( range_vec[jj][0] , range_vec[jj][1], num=numPointsXax) + \ pointEval[jj] #Notice x -- jj ylin = np.linspace( range_vec[ii][0] , range_vec[ii][1], num=numPointsXax) + \ pointEval[ii] X, Y = np.meshgrid(xlin, ylin) X_eval = X.flat #Flat iterator over the array Y_eval = Y.flat #Slice 2d pdf npointEval =len(X_eval) #Number of points where to evaluate the slice tempMat = np.tile(pointEval,(npointEval,1)) tempMat[:,dimensions_vec[ii]] = Y_eval tempMat[:,dimensions_vec[jj]] = X_eval X_list.append(X) Y_list.append(Y) tempMat_list.append(tempMat) # Evaluate npoints_list = [0] + [ tmat.shape[0] for tmat in tempMat_list ] npoints_pos = np.cumsum(npoints_list) tempMat_mat = np.concatenate(tempMat_list, axis=0) scatter_tuple = (['x'], [tempMat_mat]) mapEval_mat = mpi_map("evaluate", obj=tr_map, scatter_tuple=scatter_tuple, mpi_pool=mpi_pool) mapEval_list = [ mapEval_mat[npoints_pos[i]:npoints_pos[i+1]] for i in range(len(tempMat_list)) ] # Plot idx = 0 for ii in range(nplots): for jj in range(ii+1): if jj<ii: # 2d plot X = X_list[idx] Y = Y_list[idx] tempMat = tempMat_list[idx] mapEval = mapEval_list[idx] #Reshape results for plotting mapEval_mat = mapEval[:,ii].reshape( X.shape[0], X.shape[1] ) #Set current axes ax = axarr[ii,jj] cont_h = ax.contourf(X, Y, mapEval_mat, numCont) # cmin, cmax = cont_h.get_clim() #cmin = cont_h.levels[1] #cont_h.cmap.set_under('white') #Set color below levels #cont_h.set_clim(cmin, cmax+10) ax.spines['top'].set_color(colorboxes) ax.spines['right'].set_color(colorboxes) ax.spines['left'].set_color(colorboxes) ax.spines['bottom'].set_color(colorboxes) ax.spines['top'].set_linewidth(lw_spines) ax.spines['right'].set_linewidth(lw_spines) ax.spines['left'].set_linewidth(lw_spines) ax.spines['bottom'].set_linewidth(lw_spines) ax.get_xaxis().set_visible(False) ax.get_yaxis().set_visible(False) ax.set_aspect('equal') idx += 1 max1dplot = 0 min1dplot = 0 for jj in range(nplots): #1D plot xlin = np.linspace( range_vec[jj][0] , range_vec[jj][1], num=numPointsXax) + pointEval[jj] #Slice 1d pdf npointEval =len(xlin) #Number of points where to evaluate the slice tempMat = np.tile(pointEval,(npointEval,1)) tempMat[:,dimensions_vec[jj]] = xlin #Evaluate distribution on the slice mapEval = tr_map.evaluate(tempMat)[:,jj] max1dplot = max(max1dplot , max(mapEval) ) min1dplot = min(min1dplot , min(mapEval) ) #maxPdf = np.max(pdfEval) #Set current axes if nplots == 1: ax = axarr else: ax = axarr[jj,jj] ax.plot( xlin , mapEval ) #ax.plot( xlin[0] , maxPdf*1.1 , color ='white' ) ax.get_xaxis().set_visible(False) #ax.get_yaxis().set_visible(False) ax.yaxis.tick_right() ax.yaxis.set_ticks_position('right') # Set the y-ticks to only the right ax.tick_params(direction='out') ax.tick_params(axis='y', colors=colorboxes, labelsize = tickslabelsize) ax.spines['top'].set_color(colorboxes) ax.spines['right'].set_color(colorboxes) ax.spines['left'].set_color(colorboxes) ax.spines['bottom'].set_color(colorboxes) ax.spines['top'].set_linewidth(lw_spines) ax.spines['right'].set_linewidth(lw_spines) ax.spines['left'].set_linewidth(lw_spines) ax.spines['bottom'].set_linewidth(lw_spines) max1dplot = int( np.ceil(max1dplot*1.1) ) min1dplot = int ( np.floor(min1dplot*1.1) ) delta = max1dplot - min1dplot quart_delta = delta/4. tick_quarter = int(min1dplot+quart_delta) tick_triquarter = int(min1dplot+3.*quart_delta) for jj in range(nplots): if nplots == 1: ax = axarr else: ax = axarr[jj,jj] ax.set_ylim([min1dplot, max1dplot]) plt.setp(ax, yticks=[min1dplot, tick_quarter, 0, tick_triquarter, max1dplot]) for ii in range(nplots): for jj in range(ii+1,nplots): axarr[ii,jj].set_visible(False) fig.canvas.draw() if not figname is None: plt.savefig(figname, bbox_inches = 'tight') else: if show_flag: plt.show(False) return fig
[docs]class RandomConditionalsObject(TMO): def __init__(self, nplots, X_list, Y_list, pdfEval_list, Q_rand_list): self.nplots = nplots self.X_list = X_list self.Y_list = Y_list self.pdfEval_list = pdfEval_list self.Q_rand_list = Q_rand_list
[docs]def computeRandomConditionals( distribution, num_conditionalsXax=0, pointEval=None, range_vec=[-3,3], numPointsXax = 30, Q_rand_list = None, mpi_pool=None): r""" Compute the random conditionals Args: distribution (:class:`Distribution`): distribution :math:`\pi` num_conditionalsXax (int): number of random conditionals per axis pointEval (:class:`ndarray`[:math:`d`]): anchor point. Default is zero. range_vec (:class:`tuple`[2]): range to be displayed. numPointsXax (int): number of points for each axis. Q_rand_list (list): list of random directions. mpi_pool (:class:`mpi_map.MPI_Pool<mpi_map.MPI_Pool>`): pool of processes """ dim = distribution.dim #Assign optional parameters if not defined if num_conditionalsXax == 0: num_conditionalsXax = 4 # number of subplots per dimension if pointEval is None: pointEval=np.zeros(dim) # Create list of points (for every plot) X_list = [] Y_list = [] tempMat_list = [] if Q_rand_list is None: Q_rand_list = [None for i in range(num_conditionalsXax**2)] for ii in range(num_conditionalsXax): for jj in range(num_conditionalsXax): #The x variable is jj, whereas the y-variable is ii. #Create evaluation grid xlin = np.linspace( range_vec[0] , range_vec[1], num=numPointsXax) #Notice x -- jj ylin = np.linspace( range_vec[0] , range_vec[1], num=numPointsXax) X, Y = np.meshgrid(xlin, ylin) X_eval = X.flat #Flat iterator over the array Y_eval = Y.flat #Generate a pair of random orthogonal directions if Q_rand_list[ii*num_conditionalsXax + jj] is None: Q_rand = np.random.randn( dim , 2 ) Q_rand[:,0] = Q_rand[:,0]/np.linalg.norm(Q_rand[:,0]) Q_rand[:,1] = Q_rand[:,1] - np.dot( Q_rand[:,0], Q_rand[:,1])*Q_rand[:,0] Q_rand[:,1] = Q_rand[:,1]/np.linalg.norm(Q_rand[:,1]) Q_rand_list[ii*num_conditionalsXax + jj] = Q_rand else: Q_rand = Q_rand_list[ii*num_conditionalsXax + jj] #Slice 2d pdf tempMat = np.dot( np.vstack((X_eval, Y_eval)).T, Q_rand.T ) tempMat += pointEval # Append X_list.append(X) Y_list.append(Y) tempMat_list.append(tempMat) # Evaluate npoints_list = [0] + [ tmat.shape[0] for tmat in tempMat_list ] npoints_pos = np.cumsum(npoints_list) tempMat_mat = np.concatenate(tempMat_list, axis=0) scatter_tuple = (['x'], [tempMat_mat]) pdfEval_mat = mpi_map("pdf", obj=distribution, scatter_tuple=scatter_tuple, mpi_pool=mpi_pool) pdfEval_list = [ pdfEval_mat[npoints_pos[i]:npoints_pos[i+1]] for i in range(len(tempMat_list)) ] data = RandomConditionalsObject(num_conditionalsXax, X_list, Y_list, pdfEval_list, Q_rand_list) return data
[docs]def plotRandomConditionals(distribution=None, data=None, num_conditionalsXax=0, pointEval=None, range_vec=[-3,3], numPointsXax = 30, numCont = 15, Q_rand_list=None, figname=None, show_flag=True, show_title=False, title="Random conditionals", fig=None, mpi_pool=None): r""" Plot the random conditionals Args: distribution (:class:`Distribution`): distribution :math:`\pi` data (:class:`RandomConditionalsObject`): output of :func:`computeRandomConditionals` num_conditionalsXax (int): number of random conditionals per axis pointEval (:class:`ndarray`[:math:`d`]): anchor point. Default is zero. range_vec (:class:`tuple`[2]): range to be displayed. numPointsXax (int): number of points for each axis. numCont (int): number of contours in the contour plots. figname (str): if defined, store the figure in the provided path. show_flag (bool): whether to show the plot before returning show_title (bool): whether to show the title fig (figure): matplotlib figure object if one wants to re-use it. mpi_pool (:class:`mpi_map.MPI_Pool<mpi_map.MPI_Pool>`): pool of processes """ if distribution is None and data is None: raise ValueError("Either distribution or data need to be provided") if not PLOT_SUPPORT: raise RuntimeError("Plotting is not suported") if data is None: data = computeRandomConditionals( distribution, num_conditionalsXax=num_conditionalsXax, pointEval=pointEval, range_vec=range_vec, numPointsXax=numPointsXax, Q_rand_list=Q_rand_list, mpi_pool=mpi_pool) num_conditionalsXax = data.nplots X_list = data.X_list Y_list = data.Y_list pdfEval_list = data.pdfEval_list if fig is None: fig = plt.figure(figsize=(8,8)); else: plt.figure(fig.number); fig.clf() axarr = np.empty((num_conditionalsXax,num_conditionalsXax), dtype=object) for i in range(num_conditionalsXax): for j in range(1,num_conditionalsXax+1): if i != 0 or j != 0: axarr[i,j-1] = plt.subplot(num_conditionalsXax, num_conditionalsXax, (i*num_conditionalsXax)+j, sharex=axarr[0,0]); else: axarr[i,j-1] = plt.subplot( num_conditionalsXax, num_conditionalsXax, (i*num_conditionalsXax)+j); if show_title: fig.suptitle(title,fontsize=16); plt.set_cmap( 'jet' ); #Set colormap #plt.subplots_adjust(wspace=.15) #Reduce horizontal distance between subplots colorboxes = 'black' lw_spines = .5 # Plot idx = 0 for ii in range(num_conditionalsXax): for jj in range(num_conditionalsXax): X = X_list[idx] Y = Y_list[idx] pdfEval = pdfEval_list[idx] #Reshape results for plotting pdfEval_mat = pdfEval.reshape( X.shape[0], X.shape[1] ) #Set current axes ax = axarr[ii,jj] cont_h = ax.contourf(X, Y, pdfEval_mat, numCont); cmin, cmax = cont_h.get_clim(); cmin = cont_h.levels[1] cont_h.cmap.set_under('white'); #Set color below levels cont_h.set_clim(cmin, cmax); ax.spines['top'].set_color(colorboxes); ax.spines['right'].set_color(colorboxes); ax.spines['left'].set_color(colorboxes); ax.spines['bottom'].set_color(colorboxes); ax.spines['top'].set_linewidth(lw_spines); ax.spines['right'].set_linewidth(lw_spines); ax.spines['left'].set_linewidth(lw_spines); ax.spines['bottom'].set_linewidth(lw_spines); ax.get_xaxis().set_visible(False); ax.get_yaxis().set_visible(False); ax.set_aspect('equal'); idx += 1 fig.canvas.draw(); if not figname is None: plt.savefig(figname, bbox_inches = 'tight'); else: if show_flag: plt.show(False); return fig
# def niceBoxplot(data, box_loc=0, widhts_val = .5, # xlabel="", ylabel="", title="", # xtickslabel = 0, logyscale = 0): # if not PLOT_SUPPORT: # raise RuntimeError("Plotting is not suported") # nboxes = data.shape[1] # #Assign unassigned optional arguments # if box_loc == 0: # box_loc = np.linspace(1,3*nboxes,nboxes) # fig = plt.figure() # ax = fig.add_subplot(111) # beta = .25 # bplot = ax.boxplot( data, notch=False, positions = box_loc, # widths = widhts_val, patch_artist=True, # showmeans=False, showfliers=False) # ax.set_xlabel(xlabel) # ax.set_ylabel(ylabel) # ax.set_title(title) # ax.set_xlim(-1, max(box_loc)+2) # if logyscale==1: # ax.set_yscale('log') # if xtickslabel != 0: # ax.set_xticklabels( [ "$%d\\times10^{%d}$"%frexp10( xtickslabel[ii] ) for ii in range(len(box_loc)) ] ) #List comprehension # ylab = ax.get_yticklabels() # def mjrFormatter(x, pos): # return "$%.1f\\times10^{%d}$"%frexp10(x,numdecimals=2) # ax.yaxis.set_major_formatter(mpl.ticker.FuncFormatter(mjrFormatter)) # for patch in bplot['boxes']: # patch.set_facecolor('lightblue') # patch.set_edgecolor('darkblue') # patch.set_alpha(.5) # plt.setp(bplot['whiskers'],linestyle = '-',color='darkblue',alpha=.5) # plt.setp(bplot['caps'],color='darkblue',alpha=.5) # coloraxes = 'gray' # ax.spines['right'].set_visible(False) # Remove the right axis boundary # ax.spines['top'].set_visible(False) # Remove the top axis boundary # ax.xaxis.set_ticks_position('bottom') # Set the x-ticks to only the bottom # ax.yaxis.set_ticks_position('left') # Set the y-ticks to only the left # ax.spines['bottom'].set_position(('axes',-0.04)) # Offset the bottom scale from the axis # ax.spines['left'].set_position(('axes',-0.04)) # Offset the left scale from the axis # ax.spines['left'].set_linewidth(.5) # ax.spines['bottom'].set_linewidth(.5) # ax.spines['bottom'].set_color(coloraxes) # ax.spines['left'].set_color(coloraxes) # ax.xaxis.label.set_color(coloraxes) # ax.yaxis.label.set_color(coloraxes) # ax.tick_params(axis='x', colors=coloraxes) # ax.tick_params(axis='y', colors=coloraxes) # ax.yaxis.label.set_size(18) # ax.xaxis.label.set_size(18) # ax.title.set_color(coloraxes) # ax.title.set_size(16) # return fig # def accuracyIntegratedExponentialDiagnostic( distribution1, distribution2, # numMCsamples = 1e4, # order_increase = 2, # figname=None, show_flag=True ): # r""" Variance diagnositc # Statistical analysis of the variance diagnostic # .. math:: # \mathcal{D}_{KL}(\pi_1 \Vert \pi_2) \approx \frac{1}{2} \mathbb{V}_{\pi_1} \left( \log \frac{\pi_1}{\pi_2}\right) # Args: # distribution1 (:class:`DIST.Distribution`): distribution :math:`\pi_1`, with which one must be able to approximate # integrals (i.e. must be able to sample from :math:`\pi_1`) # distribution2 (:class:`DIST.Distribution`): distribution :math:`\pi_2` # numMCsamples (int): number of samples # order_increase (int): number of increases of the order for the # integral approximation. # figname (str): if defined, store the figure in the provided path. # show_flag (bool): whether to show the plot before returning # """ # if not PLOT_SUPPORT: # raise RuntimeError("Plotting is not suported") # d2 = copy.deepcopy(distribution2) # tm = d2.transport_map # if not isinstance(tm, MAPS.IntegratedExponentialTriangularTransportMap): # raise ValueError("The transport map must be of type " + # "IntegratedExponentialTriangularTransportMap") # base_ord_list = [ t.integ_ord_mult for t in tm.approx_list ] # dim = distribution1.dim # #Will compute the variance diagnostic for an increasing number of samples # allVar = np.zeros( order_increase + 1 ) # samplesMat =distribution1.rvs( int( numMCsamples ) ) # for ii in range(order_increase+1): # logger.info("accuracyIntegratedExponentialDiagnostic: Iteration %d" % ii) # for t,bo in zip(tm.approx_list, base_ord_list): # t.integ_ord_mult = bo + ii # intEval = distribution1.log_pdf(samplesMat)-d2.log_pdf(samplesMat) # allVar[ii] = np.var( intEval ) # allVar = allVar/2. # #fig = plt.semilogy(range(order_increase+1), allVar) # fig = nicePlot( range(order_increase+1), allVar , logyscale=True , title="Accuracy integrated exponential", xlabel="Integration order", ylabel="Variance diagnostic" ) # if not figname is None: # fig.savefig(figname, bbox_inches = 'tight' ) # else: # if show_flag: # plt.show(False) # def nicePlot( x , y, xlabel="", ylabel="", title="", xtickslabel = 0, # linewidth=2, logyscale = False, logxscale = False, linestyle='b', # figname=None, show_flag=True): # if not PLOT_SUPPORT: # raise RuntimeError("Plotting is not suported") # fig = plt.figure() # ax = fig.add_subplot(111) # if logyscale and logxscale: # plot1 = ax.loglog( x, y, linestyle ) # elif logyscale: # plot1 = ax.semilogy( x, y, linestyle ) # elif logxscale: # plot1 = ax.semilogx( x, y, linestyle ) # else: # plot1 = ax.plot( x, y, linestyle ) # plt.setp(plot1, linewidth=linewidth) # ax.set_xlabel(xlabel) # ax.set_ylabel(ylabel) # ax.set_title(title) # coloraxes = 'gray' # ax.spines['right'].set_visible(False) # Remove the right axis boundary # ax.spines['top'].set_visible(False) # Remove the top axis boundary # ax.xaxis.set_ticks_position('bottom') # Set the x-ticks to only the bottom # ax.yaxis.set_ticks_position('left') # Set the y-ticks to only the left # ax.spines['bottom'].set_position(('axes',-0.04)) # Offset the bottom scale from the axis # ax.spines['left'].set_position(('axes',-0.04)) # Offset the left scale from the axis # ax.spines['left'].set_linewidth(.5) # ax.spines['bottom'].set_linewidth(.5) # ax.spines['bottom'].set_color(coloraxes) # ax.spines['left'].set_color(coloraxes) # ax.xaxis.label.set_color(coloraxes) # ax.yaxis.label.set_color(coloraxes) # ax.tick_params(axis='x', colors=coloraxes) # ax.tick_params(axis='y', colors=coloraxes) # ax.yaxis.label.set_size(18) # ax.xaxis.label.set_size(18) # ax.title.set_color(coloraxes) # ax.title.set_size(16) # fig.subplots_adjust(bottom=0.14) # if not figname is None: # fig.savefig(figname, bbox_inches = 'tight' ) # else: # if show_flag: # plt.show(False) # return fig # def niceContour( x , y, z, xlabel="", ylabel="", title="", # numCont = 30, xlim=[-3,3], ylim=[-3,3], # colormap='Blues', contourf = False, colorbar=False, # figname=None, show_flag=True): # if not PLOT_SUPPORT: # raise RuntimeError("Plotting is not suported") # plt.set_cmap( colormap ) # fig = plt.figure() # ax = fig.add_subplot(111) # if contourf: # cont = ax.contourf(x, y, z, numCont) # else: # cont = ax.contour(x, y, z, numCont) # coloraxes = 'gray' # #Set axis limits # ax.set_ylim(ylim[0],ylim[1]) # ax.set_xlim(xlim[0],xlim[1]) # #Remove last level # cmin, cmax = cont.get_clim() # cont.cmap.set_under('white') #Set color below level # cmin = cont.levels[1] # cont.set_clim(cmin, cmax) #Set limit levels # plt.draw() # #Colorbar # if colorbar: # cbar = plt.colorbar(cont) # prop = plt.getp(cbar.ax.axes, 'yticklabels') # plt.setp(prop, color=coloraxes) # cbar.outline.set_edgecolor(coloraxes) # cbar.ax.tick_params(color=coloraxes) # #Assign labels and title # ax.set_xlabel(xlabel) # ax.set_ylabel(ylabel) # ax.set_title(title) # #Axis aspect # ax.spines['right'].set_visible(False) # Remove the right axis boundary # ax.spines['top'].set_visible(False) # Remove the top axis boundary # ax.xaxis.set_ticks_position('bottom') # Set the x-ticks to only the bottom # ax.yaxis.set_ticks_position('left') # Set the y-ticks to only the left # ax.spines['bottom'].set_position(('axes',-0.04)) # Offset the bottom scale from the axis # ax.spines['left'].set_position(('axes',-0.04)) # Offset the left scale from the axis # ax.spines['left'].set_linewidth(.5) # ax.spines['bottom'].set_linewidth(.5) # ax.spines['bottom'].set_color(coloraxes) # ax.spines['left'].set_color(coloraxes) # ax.xaxis.label.set_color(coloraxes) # ax.yaxis.label.set_color(coloraxes) # ax.tick_params(axis='x', colors=coloraxes) # ax.tick_params(axis='y', colors=coloraxes) # ax.yaxis.label.set_size(18) # ax.xaxis.label.set_size(18) # ax.title.set_color(coloraxes) # ax.title.set_size(16) # fig.subplots_adjust(bottom=0.14) # ax.set_aspect('equal') # if not figname is None: # fig.savefig(figname, bbox_inches = 'tight' ) # else: # if show_flag: # plt.show(False) # return fig # def plot2distribution( distribution, xlabel="", ylabel="", title="", # numCont = 30, xlim=[-3,3], ylim=[-3,3], numPointsXax = 150, # colormap='Blues', contourf = False, colorbar = False, # figname=None, show_flag=True): # if not PLOT_SUPPORT: # raise RuntimeError("Plotting is not suported") # #Create evaluation grid # xlin = np.linspace( xlim[0] , xlim[1], num=numPointsXax) # ylin = np.linspace( ylim[0] , ylim[1], num=numPointsXax) # X, Y = np.meshgrid(xlin, ylin) # X_eval = X.flatten() #Flat iterator over the array # Y_eval = Y.flatten() # npointEval =len(X_eval) #Number of points where to evaluate the distribution # tempMat = np.vstack( ( X_eval , Y_eval ) ) # tempMat = tempMat.transpose() # #Evaluate the distribution # pdfEval = distribution.pdf(tempMat) # #Reshape results for plotting # pdfEval_mat = pdfEval.reshape( X.shape[0], X.shape[1] ) # niceContour( X , Y, pdfEval_mat, xlabel=xlabel, ylabel=ylabel, # title=title, numCont = numCont, xlim=xlim, ylim=ylim, # colormap=colormap, contourf = contourf, colorbar=colorbar, # figname=figname, show_flag=show_flag)
[docs]def plotAlignedMarginals(mat_points, mat_points2=None, dimensions_vec=0, range_vec=None, scatter=False, colormap = 'jet', white_background=True, levels=10, do_diag=True, figname=None, show_flag=True, show_axis=False, title='Marginals along coordinate axes', vartitles=None, fig=None, ret_handles=False, mpi_pool=None): r""" Plot the marginals aligned with the axis Args: mat_points (ndarray): first dataset mat_points2 (ndarray): second dataset (optional) dimensions_vec (list of int): list of dimensions to be displayed. Default 0: display 10 dimensions at most. range_vec (list of :class:`tuple<tuple>` [2]): range to be displayed. scatter (bool): whether to plot a scatter plots instead of densities colormap (str): colormap to be used white_background (bool): whether to have a white background of to use the last layer of the colormap levels (int or list): number of levels to be displayed or list of values defining the levels. do_diag (bool): whether to include the one dimensional marginals on the diagonal numPointsXax (int): number of points for each axis. numCont (int): number of contours in the contour plots. figname (str): if defined, store the figure in the provided path. show_flag (bool): whether to show the plot before returning show_axis (bool): whether to show the axis of the plot vartitles (list): list of titles for each variable fig (figure): matplotlib figure object if one wants to re-use it. mpi_pool (:class:`mpi_map.MPI_Pool<mpi_map.MPI_Pool>`): pool of processes ret_handles (bool): whether to return the axes handles Returns: (fig, [list]) -- figure handle and dictionary of handles """ if not PLOT_SUPPORT: raise RuntimeError("Plotting is not suported") dim = mat_points.shape[1] #Assign optional parameters if not defined if dimensions_vec == 0: dimensions_vec = list(range( min(dim,20) )) else: dimensions_vec = list(dimensions_vec) if (np.max(dimensions_vec)+1)>dim: raise RuntimeError("Desired dimensions not available") if vartitles is not None and len(vartitles) != len(dimensions_vec): raise ValueError("You need to provide the right number of titles.") if fig is None: fig = plt.figure(figsize=(8,8)) else: plt.figure(fig.number) fig.clf() nplots = len(dimensions_vec) # number of subplots per dimension naxes = nplots if do_diag else nplots - 1 axarr = np.empty((naxes,naxes), dtype=object) for i in range(naxes): for j in range(1,naxes+1): if i == 0 and j == 1: axarr[0,0] = plt.subplot( naxes, naxes, (i*naxes)+j ) else: axarr[i,j-1] = plt.subplot( naxes, naxes, (i*naxes)+j, sharex=axarr[i-1,j-1] ) if title is not None: fig.suptitle(title,fontsize=16) # plt.set_cmap( 'jet' ) #Set colormap # plt.subplots_adjust(wspace=.15) #Reduce horizontal distance between subplots colorboxes = 'black' lw_spines = .5 xx_list = [] yy_list = [] zz_list = [] handles = [] # Plot idx = 0 for ii in range(nplots): for jj in range(ii+1): if jj<ii: # 2d plot index_xx = dimensions_vec[jj] index_yy = dimensions_vec[ii] xx = mat_points[:,index_xx] yy = mat_points[:,index_yy] if mat_points2 is not None: xx2 = mat_points2[:,index_xx] yy2 = mat_points2[:,index_yy] if range_vec is not None: xx_min = range_vec[jj][0] xx_max = range_vec[jj][1] yy_min = range_vec[ii][0] yy_max = range_vec[ii][1] else: xx_min = xx.min() xx_max = xx.max() yy_min = yy.min() yy_max = yy.max() if mat_points2 is not None: xx_min = min( xx_min, xx2.min() ) xx_max = max( xx_min, xx2.max() ) yy_min = min( yy_min, yy2.min() ) yy_min = min( yy_min, yy2.min() ) if do_diag: ax = axarr[ii,jj] else: ax = axarr[ii-1, jj] if scatter: cont_h = ax.scatter(xx, yy, marker='.', s=.5) else: XX, YY = np.mgrid[xx_min:xx_max:100j, yy_min:yy_max:100j] positions = np.vstack([XX.ravel(), YY.ravel()]) values = np.vstack([xx, yy]) kernel = scistat.gaussian_kde(values) if mat_points2 is not None: values2 = np.vstack([xx2, yy2]) kernel2 = scistat.gaussian_kde(values2) ZZ = np.reshape(kernel(positions).T, XX.shape) zz_max = ZZ.max() if mat_points2 is not None: ZZ2 = np.reshape(kernel2(positions).T, XX.shape) zz_max = max( zz_max, ZZ.max() ) lvls = np.linspace(0, zz_max*0.95, levels) cont_h = ax.contour(XX,YY,ZZ, levels=lvls, cmap=plt.get_cmap(colormap)) if mat_points2 is not None: cont_h2 = ax.contour(XX,YY,ZZ2, linestyles='dashed', levels=lvls, cmap=plt.get_cmap(colormap)) # if white_background: # cmin, cmax = cont_h.get_clim() # cmin = cont_h.levels[1] # cont_h.cmap.set_under('white') #Set color below levels # cont_h.set_clim(cmin, cmax) xx_list.append(XX) yy_list.append(YY) zz_list.append(ZZ) handles.append(cont_h) ax.set_xlim([xx_min, xx_max]) ax.set_ylim([yy_min, yy_max]) ax.spines['top'].set_color(colorboxes) ax.spines['right'].set_color(colorboxes) ax.spines['left'].set_color(colorboxes) ax.spines['bottom'].set_color(colorboxes) ax.spines['top'].set_linewidth(lw_spines) ax.spines['right'].set_linewidth(lw_spines) ax.spines['left'].set_linewidth(lw_spines) ax.spines['bottom'].set_linewidth(lw_spines) ax.get_xaxis().set_visible(show_axis) ax.get_yaxis().set_visible(show_axis) #ax.set_aspect('equal') if do_diag: for jj in range(nplots): #1D plot index_xx = dimensions_vec[jj] xx = mat_points[:,index_xx] if mat_points2 is not None: xx2 = mat_points2[:,index_xx] if range_vec is not None: xx_min = range_vec[jj][0] xx_max = range_vec[jj][1] else: xx_min = xx.min() xx_max = xx.max() if mat_points2 is not None: xx_min = min( xx_min, xx2.min() ) xx_max = max( xx_max, xx2.max() ) XX = np.linspace(xx_min, xx_max,100) kernel = scistat.gaussian_kde(xx) if mat_points2 is not None: kernel2 = scistat.gaussian_kde(xx2) pdfEval = kernel.pdf(XX) if mat_points2 is not None: pdfEval2 = kernel2.pdf(XX) ax = axarr[jj,jj] ph = ax.plot( XX , pdfEval ) if mat_points2 is not None: ax.plot(XX, pdfEval2, '--b') xx_list.append(XX) zz_list.append(pdfEval) handles.append( ph ) ax.set_xlim([xx_min, xx_max]) ax.get_xaxis().set_visible(show_axis) ax.get_yaxis().set_visible(False) ax.spines['top'].set_color(colorboxes) ax.spines['right'].set_color(colorboxes) ax.spines['left'].set_color(colorboxes) ax.spines['bottom'].set_color(colorboxes) ax.spines['top'].set_linewidth(lw_spines) ax.spines['right'].set_linewidth(lw_spines) ax.spines['left'].set_linewidth(lw_spines) ax.spines['bottom'].set_linewidth(lw_spines) # Variables titles if vartitles is not None: if do_diag: for jj, tit in enumerate(vartitles): ax = axarr[jj,jj] ax.set_title(tit) else: for jj, tit in enumerate(vartitles): if jj < naxes: ax = axarr[naxes-1,jj] ax.set_xlabel(tit) if jj > 0: ax = axarr[jj-1,0] ax.set_ylabel(tit) for ii in range(naxes): for jj in range(ii+1,naxes): axarr[ii,jj].set_visible(False) # fig.canvas.draw() if not figname is None: plt.savefig(figname, bbox_inches = 'tight') else: if show_flag: plt.show(False) if ret_handles: axhandles = {'axarr': axarr, 'handles': handles, 'xx_list': xx_list, 'yy_list': yy_list, 'zz_list': zz_list} return (fig, axhandles) else: return fig
[docs]def plotAlignedScatters(mat_points, dimensions_vec=0, do_diag=True, s=5, bins=10, show_axis=True, axis_fmt=None, figname=None, show_flag=True, show_title=False, title='Marginals along coordinate axes', vartitles=None, fig=None): r""" Plot the marginals aligned with the axis Args: mat_points (:class:`ndarray<numpy.ndarray>` [:math:`m,d`]): samples dimensions_vec (list of int): list of dimensions to be displayed. Default 0: display 10 dimensions at most. do_diag (bool): whether to include the one dimensional marginals on the diagonal s (int): size of scatter points bins (int): number of bins for one dimensional plots show_axis (bool): whether to show the axis axis_fmt (list): list of matplotlib formatters figname (str): if defined, store the figure in the provided path. show_flag (bool): whether to show the plot before returning show_title (bool): whether to show the title title (str): title for the figure vartitles (list): list of titles for each variable fig (figure): matplotlib figure object if one wants to re-use it. """ if not PLOT_SUPPORT: raise RuntimeError("Plotting is not suported") dim = mat_points.shape[1] #Assign optional parameters if not defined if dimensions_vec == 0: dimensions_vec = list(range( min(dim,20) )) else: dimensions_vec = list(dimensions_vec) if (np.max(dimensions_vec)+1)>dim: raise RuntimeError("Desired dimensions not available") if vartitles is not None and len(vartitles) != dim: raise ValueError("You need to provide the right number of titles.") if fig is None: fig = plt.figure(figsize=(8,8)); else: plt.figure(fig.number); fig.clf(); nplots = len(dimensions_vec) # number of subplots per dimension naxes = nplots if do_diag else nplots - 1 axarr = np.empty((naxes,naxes), dtype=object) for i in range(naxes): for j in range(1,naxes+1): if i == 0 and j == 1: axarr[0,0] = plt.subplot( naxes, naxes, (i*naxes)+j ); else: axarr[i,j-1] = plt.subplot( naxes, naxes, (i*naxes)+j, sharex=axarr[i-1,j-1] ); if show_title: fig.suptitle(title,fontsize=16); # plt.set_cmap( 'jet' ); #Set colormap # plt.subplots_adjust(wspace=.15) #Reduce horizontal distance between subplots colorboxes = 'black' lw_spines = .5 # Plot idx = 0 for ii in range(nplots): for jj in range(ii+1): if jj<ii: # 2d plot index_xx = dimensions_vec[jj] index_yy = dimensions_vec[ii] xx = mat_points[:,index_xx] yy = mat_points[:,index_yy] xx_min = xx.min() xx_max = xx.max() yy_min = yy.min() yy_max = yy.max() if do_diag: ax = axarr[ii,jj] else: ax = axarr[ii-1, jj] ax.scatter(xx, yy, s=s) ax.set_xlim([xx_min, xx_max]); ax.set_ylim([yy_min, yy_max]); ax.spines['top'].set_color(colorboxes); ax.spines['right'].set_color(colorboxes); ax.spines['left'].set_color(colorboxes); ax.spines['bottom'].set_color(colorboxes); ax.spines['top'].set_linewidth(lw_spines); ax.spines['right'].set_linewidth(lw_spines); ax.spines['left'].set_linewidth(lw_spines); ax.spines['bottom'].set_linewidth(lw_spines); ax.get_xaxis().set_visible(ii==(nplots-1) and show_axis); ax.get_yaxis().set_visible(jj==0 and show_axis); ax.get_xaxis().set_major_locator(mpl.ticker.MaxNLocator(1)) ax.get_yaxis().set_major_locator(mpl.ticker.MaxNLocator(1)) if axis_fmt is not None: ax.get_xaxis().set_major_formatter(axis_fmt[jj]) ax.get_yaxis().set_major_formatter(axis_fmt[ii]) if do_diag: for jj in range(nplots): #1D plot index_xx = dimensions_vec[jj] xx = mat_points[:,index_xx] xx_min = xx.min() xx_max = xx.max() #Set current axes ax = axarr[jj,jj] ax.hist(xx, bins=bins) ax.set_xlim([xx_min, xx_max]) ax.spines['top'].set_color(colorboxes); ax.spines['right'].set_color(colorboxes); ax.spines['left'].set_color(colorboxes); ax.spines['bottom'].set_color(colorboxes); ax.spines['top'].set_linewidth(lw_spines); ax.spines['right'].set_linewidth(lw_spines); ax.spines['left'].set_linewidth(lw_spines); ax.spines['bottom'].set_linewidth(lw_spines); ax.get_xaxis().set_visible(jj==(nplots-1) and show_axis); ax.get_yaxis().set_visible(False); ax.get_xaxis().set_major_locator(mpl.ticker.MaxNLocator(1)) if axis_fmt is not None: ax.get_xaxis().set_major_formatter(axis_fmt[jj]) if vartitles is not None: tit = vartitles[jj] ax.set_title(tit); for ii in range(naxes): for jj in range(ii+1,naxes): axarr[ii,jj].set_visible(False); fig.canvas.draw(); if not figname is None: plt.savefig(figname, bbox_inches = 'tight'); else: if show_flag: plt.show(False); return fig
[docs]def plotGradXMap(tmap, base_distribution = None, nsamples = 1000, show_cbar=True, show_ticks=True, title = "Intensity coefficients map", cmap='Blues', mpi_pool=None, show_flag=True): dim = tmap.dim if base_distribution == None: base_distribution = StandardNormalDistribution(dim) X = base_distribution.rvs(nsamples) #Parallel evaluation scatter_tuple = (['x'], [X]) gradient_map = mpi_map("grad_x", obj=tmap, scatter_tuple=scatter_tuple, mpi_pool=mpi_pool) intensityCoeff = np.sum( np.abs( gradient_map),axis = 0)/nsamples intensityCoeff_normalized = intensityCoeff / np.max(intensityCoeff) return niceSpy( intensityCoeff_normalized, title=title, cmap=cmap, show_flag=show_flag, show_cbar=show_cbar, show_ticks=show_ticks)
[docs]def plotLinearityMap(tmap, nlevels=2, threshold=1e-2, title = "Linearity pattern", cmap='Blues', show_flag=True): dim = tmap.dim lin_mat = np.zeros((dim,dim)) lvs = np.linspace(0,1,nlevels+1) for d, (tm_comp, tm_comp_avar) in enumerate(zip(tmap.approx_list, tmap.active_vars)): if len(tm_comp.c.multi_idxs) == 0: const_midxs_mat = np.zeros((0,len(tm_comp_avar)), dtype=int) else: const_midxs_mat = np.asarray( tm_comp.c.multi_idxs ) exp_midxs_mat = np.asarray( tm_comp.h.multi_idxs ) # Adjust multi-index of the exponential part to map to the linear span approx exp_midxs_mat[:,-1] += 1 # Threshold parameters const_threshold_mask = np.abs(tm_comp.c.coeffs) > threshold const_midxs_mat = const_midxs_mat[const_threshold_mask,:] exp_threshold_mask = np.abs(tm_comp.h.coeffs) > threshold exp_midxs_mat = exp_midxs_mat[exp_threshold_mask,:] if exp_midxs_mat.shape[0] == 0: exp_midxs_mat = np.zeros((1,len(tm_comp_avar)), dtype=int) exp_midxs_mat[0,-1] = 1 # Merge multi-index matrices midxs_mat = np.vstack( (const_midxs_mat, exp_midxs_mat) ) # Compute mixed order mixed_ord = np.sum(midxs_mat, axis=1) for nl in range(nlevels): # Find dimensions of mixed order nl+1 or >=nlevels if nl+1=nlevels if nl+1 == nlevels: order_idxs = np.where( mixed_ord > nl )[0] else: order_idxs = np.where( mixed_ord == nl+1 )[0] for idx in order_idxs: # Find active variables for selected multi-indices avar_idxs = np.where( midxs_mat[idx,:] > 0 )[0] avar = [tm_comp_avar[i] for i in avar_idxs] # Set value on matrix lin_mat[d,avar] = lvs[nl+1] return niceSpy(lin_mat, title=title, cmap=cmap, show_flag=show_flag, vrange=[0,1])
[docs]def niceSpy(mat, title = None, cmap = 'Blues', show_flag=True, vrange=None, show_cbar=True, show_ticks=True): fig = plt.figure() ax = fig.add_subplot(111) #ims = ax.imshow( intensityCoeff, interpolation='none', title = title ) plt.set_cmap( cmap ) if vrange is None: vmin = np.min(mat) vmax = np.max(mat) else: vmin = vrange[0] vmax = vrange[1] ims = ax.matshow(mat, vmin=vmin, vmax=vmax) coloraxes = 'gray' #Colorbar if show_cbar: cbar = plt.colorbar(ims) prop = plt.getp(cbar.ax.axes, 'yticklabels') plt.setp(prop, color=coloraxes) cbar.outline.set_edgecolor(coloraxes) cbar.ax.tick_params(color=coloraxes) #Assign title if title is not None: ax.set_title(title) if not show_ticks: ax.set_xticks([]) ax.set_yticks([]) else: ax.spines['left'].set_linewidth(.5) ax.spines['bottom'].set_linewidth(.5) ax.spines['top'].set_linewidth(.5) ax.spines['right'].set_linewidth(.5) ax.spines['bottom'].set_color(coloraxes) ax.spines['left'].set_color(coloraxes) ax.spines['right'].set_color(coloraxes) ax.spines['top'].set_color(coloraxes) ax.tick_params(axis='x', colors=coloraxes) ax.tick_params(axis='y', colors=coloraxes) ax.xaxis.label.set_color(coloraxes) ax.yaxis.label.set_color(coloraxes) ax.yaxis.label.set_size(18) ax.xaxis.label.set_size(18) ax.title.set_color(coloraxes) ax.title.set_size(16) ax.set_aspect('equal') if show_flag: plt.show(False) return fig