Skip to content

plotting

plot_single_trispectrum(sfs, folded=False, cmap=_new_cmap, vmin=None, vmax=None, colorbar=False, fraction=0.046)

Plots a single triallelic spectrum (sfs)

Note

~~~~8/13: to check - does it need to be given an unfolded spectrum? what if we want to plot a folded spectrum...

Parameters:

Name Type Description Default
sfs Spectrum

the triallelic spectrum to plot

 required
folded bool

True if we want to fold the spectrum before plotting, False if to remain unfolded

 False
cmap

Define the colormap to use. Default is same colormap as used in figures from Ragsdale et al (2016)

 _new_cmap
vmin float

lower limit of colormap and smallest value to be plotted

 None
vmax float

upper limit of colormap

 None
colorbar bool

True to show colorbar for Count

 False
fraction float

fraction of the axis to use for colorbar

 0.046
Source code in dadi/Triallele/plotting.py 
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
def plot_single_trispectrum(sfs, folded=False, cmap=_new_cmap, vmin=None, vmax=None, colorbar=False, fraction=.046):
    """
    Plots a single triallelic spectrum (sfs)

    Note:
        ~~~~8/13: to check - does it need to be given an unfolded spectrum? what if we want to plot a folded spectrum...

    Args:
        sfs (Spectrum): the triallelic spectrum to plot
        folded (bool): True if we want to fold the spectrum before plotting, False if to remain unfolded
        cmap (): Define the colormap to use. Default is same colormap as used in figures from Ragsdale et al (2016)
        vmin (float): lower limit of colormap and smallest value to be plotted
        vmax (float): upper limit of colormap
        colorbar (bool): True to show colorbar for Count
        fraction (float): fraction of the axis to use for colorbar
    """
    ax = plt.gca()
    if vmin is None:
        vmin = sfs.min()
    if vmax is None:
        vmax = sfs.max()

    if folded == True:
        sfs = _fold(sfs)
        mappable = ax.pcolor(np.transpose(np.ma.masked_where(sfs < vmin, sfs)), 
                    norm=LogNorm(vmin=vmin, vmax=vmax), 
                    vmin = vmin, vmax = vmax, cmap=_cmap)
    else:
        mappable = ax.pcolor(np.ma.masked_where(sfs < vmin, sfs), 
                    norm=LogNorm(vmin=vmin, vmax=vmax), 
                    vmin = vmin, vmax = vmax, cmap=_cmap)

    ax.set_xlim(0, sfs.shape[1])
    if folded == True:
        ax.set_ylim(0, sfs.shape[0]/2)
        ax.set_xlabel('Major derived allele frequency')
        ax.set_ylabel('Minor derived allele frequency')
    else:
        ax.set_ylim(0, sfs.shape[0])
        ax.set_xlabel('Frequency')
        ax.set_ylabel('Frequency')

    if colorbar == True:
        cbar = ax.figure.colorbar(mappable,fraction=fraction, pad=0.04,
                                  shrink=0.85)
        #cbar.ax.tick_params(labelsize=6)
        cbar.set_label(label='Count', fontsize=10)
        return cbar

plot_trispectrum_comp(sfs1, sfs2, folded=False, cmap=_new_cmap, vmin=None, vmax=None, resid_range=None, colorbar=False, title1='sfs1', title2='sfs2', title3='Residual', fraction=0.046)

Plots the two spectra and the residual ( (sfs1-sfs2)/sqrt(sfs1) )

Parameters:

Name Type Description Default
sfs1 Spectrum

the first triallelic spectrum to plot

 required
sfs2 Spectrum

the second triallelic spectrum to plot

 required
folded bool

True if we want to fold the spectrum before plotting, False if to remain unfolded

 False
cmap

Define the colormap to use. Default is same colormap as used in figures from Ragsdale et al (2016)

 _new_cmap
vmin float

lower limit of colormap and smallest value to be plotted for frequency spectra

 None
vmax float

upper limit of colormap for frequency spectra

 None
resid_range float

residual colormap ranges from -resid_range to +resid_range

 None
colorbar bool

True to show colorbars for Count and Residual

 False
title1 str

title for first spectrum

 'sfs1'
title2 str

title for second spectrum

 'sfs2'
title3 str

title for residuals

 'Residual'
fraction float

fraction of the axis to use for colorbars

 0.046
Source code in dadi/Triallele/plotting.py 
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
def plot_trispectrum_comp(sfs1, sfs2, folded=False, cmap=_new_cmap, vmin=None, vmax=None, resid_range=None, colorbar=False, title1="sfs1", title2="sfs2", title3="Residual", fraction=.046):
    """
    Plots the two spectra and the residual ( (sfs1-sfs2)/sqrt(sfs1) )

    Args:
        sfs1 (Spectrum): the first triallelic spectrum to plot
        sfs2 (Spectrum): the second triallelic spectrum to plot
        folded (bool): True if we want to fold the spectrum before plotting, False if to remain unfolded
        cmap (): Define the colormap to use. Default is same colormap as used in figures from Ragsdale et al (2016)
        vmin (float): lower limit of colormap and smallest value to be plotted for frequency spectra
        vmax (float): upper limit of colormap for frequency spectra
        resid_range (float): residual colormap ranges from -resid_range to +resid_range
        colorbar (bool): True to show colorbars for Count and Residual
        title1 (str): title for first spectrum
        title2 (str): title for second spectrum
        title3 (str): title for residuals
        fraction (float): fraction of the axis to use for colorbars
    """
    fig = plt.figure(np.random.randint(1000))

    if vmin is None:
        vmin = np.min((sfs1.min(),sfs2.min()))
    if vmax is None:
        vmax = np.max((sfs1.max(),sfs2.max()))

    ax1 = plt.subplot(3,1,1)
    plot_single_trispectrum(sfs1, folded=folded, cmap=cmap, vmin=vmin, vmax=vmax, colorbar=colorbar, fraction=.046)
    ax1.set_title(title1)

    ax2 = plt.subplot(3,1,2)
    plot_single_trispectrum(sfs2, folded=folded, cmap=cmap, vmin=vmin, vmax=vmax, colorbar=colorbar, fraction=.046)
    ax2.set_title(title2)

    ax3 = plt.subplot(3,1,3)
    ax3.set_title(title3)

    residuals = (sfs1 - sfs2)/np.sqrt(sfs1)
    if folded == True:
        residuals = _fold(residuals)
        if resid_range == None:
            mappable = ax3.pcolor(np.transpose(residuals), vmin = -np.max(abs(residuals)), vmax = np.max(abs(residuals)), cmap='RdBu')
        else:
            mappable = ax3.pcolor(np.transpose(residuals), vmin = -resid_range, vmax = resid_range, cmap='RdBu')
    else:
        if resid_range == None:
            mappable = ax3.pcolor(residuals, vmin = -np.max(abs(residuals)), vmax = np.max(abs(residuals)), cmap='RdBu')
        else:
            mappable = ax3.pcolor(residuals, vmin = -resid_range, vmax = resid_range, cmap='RdBu')

    ax3.set_xlim(0, sfs1.shape[1])
    if folded == True:
        ax3.set_ylim(0, sfs1.shape[0]/2)
    else:
        ax3.set_ylim(0, sfs1.shape[0])

    if colorbar == True:
        cbar = ax3.figure.colorbar(mappable,fraction=fraction, pad=0.04,
                                  shrink=0.85)
        #cbar.ax.tick_params(labelsize=6)
        cbar.set_label(label='Residual', fontsize=10)
        return cbar