Skip to content

plotting

plot_3d(fs, vmin=None, vmax=None, max_sum=None, max_index=None, show=True, colorbar=False)

Parameters:

Name Type Description Default
fs TLSpectrum

the TLSpectrum to plot

 required
vmin float

minimum value to plot

 None
vmax float

maximum value to plot

 None
max_sum int

slice spectrum on diagonal plane, plotting points closer to the origin than the plane

 None
max_index int

only show entries this distance from one of the axes planes

 None
show bool

True to show the plot

 True
colorbar bool

True to show colorbar

 False
Source code in dadi/TwoLocus/plotting.py 
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 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
 95
 96
 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
def plot_3d(fs, vmin=None, vmax=None, max_sum=None, max_index=None, show=True, colorbar=False):
    """
    Args:
        fs (TLSpectrum): the TLSpectrum to plot
        vmin (float, optional): minimum value to plot
        vmax (float, optional): maximum value to plot
        max_sum (int, optional): slice spectrum on diagonal plane, plotting points closer to the origin than the plane
        max_index (int, optional): only show entries this distance from one of the axes planes
        show (bool, optional): True to show the plot
        colorbar (bool, optional): True to show colorbar
    """
    fig = plt.figure()
    plt.clf()
    ax = mplot3d.Axes3D(fig)

    if vmin==None:
        vmin = np.min(fs[fs>0])
    if vmax==None:
        vmax = np.max(fs[fs>0])

    fs = np.swapaxes(fs,0,2)

    toplot = np.logical_not(fs.mask)
    toplot = np.logical_and(toplot, fs.data >= vmin)
    if max_sum != None:
        iis = np.where(toplot)[0]
        jjs = np.where(toplot)[1]
        kks = np.where(toplot)[2]
        for ll in range(len(iis)):
            ii = iis[ll]
            jj = jjs[ll]
            kk = kks[ll]
            if ii+jj+kk > max_sum:
                toplot[ii,jj,kk] = False
    if max_index != None:
        iis = np.where(toplot)[0]
        jjs = np.where(toplot)[1]
        kks = np.where(toplot)[2]
        for ll in range(len(iis)):
            ii = iis[ll]
            jj = jjs[ll]
            kk = kks[ll]
            if ii > max_index and jj > max_index and kk > max_index:
                toplot[ii,jj,kk] = False

    normalized = (np.log(fs)-np.log(vmin))\
            /(np.log(vmax)-np.log(vmin))
    normalized = np.minimum(normalized, 1)
    # scrunch by a factor 
    # XXX: this is really hacky
    factor = .1
    normalized = (1-2*factor)*normalized + .2
    colors = plt.cm.cubehelix_r(normalized)

    # We draw by calculating which faces are visible and including each as a
    # polygon.
    polys, polycolors = [],[]
    for ii in range(fs.shape[0]):
        for jj in range(fs.shape[1]):
            for kk in range(fs.shape[2]):
                if not toplot[ii,jj,kk]:
                    continue
                if kk < fs.shape[2]-1 and toplot[ii,jj,kk+1]:
                    pass
                else:
                    polys.append([[ii-0.5,jj+0.5,kk+0.5],[ii+0.5,jj+0.5,kk+0.5],
                                  [ii+0.5,jj-0.5,kk+0.5],[ii-0.5,jj-0.5,kk+0.5]]
                                 )
                    polycolors.append(colors[ii,jj,kk])
                if kk > 0 and toplot[ii,jj,kk-1]:
                    pass
                else:
                    polys.append([[ii-0.5,jj+0.5,kk-0.5],[ii+0.5,jj+0.5,kk-0.5],
                                  [ii+0.5,jj-0.5,kk-0.5],[ii-0.5,jj-0.5,kk-0.5]]
                                 )
                    polycolors.append(colors[ii,jj,kk])
                if jj < fs.shape[1]-1 and toplot[ii,jj+1,kk]:
                    pass
                else:
                    polys.append([[ii-0.5,jj+0.5,kk+0.5],[ii+0.5,jj+0.5,kk+0.5],
                                  [ii+0.5,jj+0.5,kk-0.5],[ii-0.5,jj+0.5,kk-0.5]]
                                 )
                    polycolors.append(colors[ii,jj,kk])
                if jj > 0 and toplot[ii,jj-1,kk]:
                    pass
                else:
                    polys.append([[ii-0.5,jj-0.5,kk+0.5],[ii+0.5,jj-0.5,kk+0.5],
                                  [ii+0.5,jj-0.5,kk-0.5],[ii-0.5,jj-0.5,kk-0.5]]
                                 )
                    polycolors.append(colors[ii,jj,kk])
                if ii < fs.shape[0]-1 and toplot[ii+1,jj,kk]:
                    pass
                else:
                    polys.append([[ii+0.5,jj-0.5,kk+0.5],[ii+0.5,jj+0.5,kk+0.5],
                                  [ii+0.5,jj+0.5,kk-0.5],[ii+0.5,jj-0.5,kk-0.5]]
                                 )
                    polycolors.append(colors[ii,jj,kk])
                if ii > 0 and toplot[ii-1,jj,kk]:
                    pass
                else:
                    polys.append([[ii-0.5,jj-0.5,kk+0.5],[ii-0.5,jj+0.5,kk+0.5],
                                  [ii-0.5,jj+0.5,kk-0.5],[ii-0.5,jj-0.5,kk-0.5]]
                                 )
                    polycolors.append(colors[ii,jj,kk])

    polycoll = mplot3d.art3d.Poly3DCollection(polys, facecolor=polycolors, 
                                              edgecolor='k', linewidths=0.5)
    ax.add_collection(polycoll)

    # Set the limits
    ax.set_xlim3d(-0.5,fs.shape[0]-0.5)
    ax.set_ylim3d(-0.5,fs.shape[1]-0.5)
    ax.set_zlim3d(-0.5,fs.shape[2]-0.5)

    ax.set_xlabel('aB', horizontalalignment='left')
    ax.set_ylabel('Ab', verticalalignment='bottom')
    ax.set_zlabel('AB', verticalalignment='bottom')

    ax.view_init(elev=30., azim=45)
    ax.colorbar()
    if show == True:
        plt.show()

plot_3d_comp(model, data, resid_range=1, max_sum=None, max_index=None, show=True)

plots (model - data)/data

Parameters:

Name Type Description Default
model TLSpectrum

the model TLSpectrum

 required
data TLSpectrum

the data TLSpectrum

 required
resid_range float

range of residuals (+/-)

 1
max_sum int

slice spectrum on diagonal plane, plotting points closer to the origin than the plane

 None
max_index int

only show entries this distance from one of the axes planes

 None
show bool

True to show the plot

 True
Source code in dadi/TwoLocus/plotting.py 
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
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
def plot_3d_comp(model,data, resid_range=1, max_sum=None, max_index=None, show=True):
    """
    plots (model - data)/data

    Args:
        model (TLSpectrum): the model TLSpectrum
        data (TLSpectrum): the data TLSpectrum
        resid_range (float, optional): range of residuals (+/-)
        max_sum (int, optional): slice spectrum on diagonal plane, plotting points closer to the origin than the plane
        max_index (int, optional): only show entries this distance from one of the axes planes
        show (bool, optional): True to show the plot
    """
    model *= np.sum(data)/np.sum(model)
    resids = (model-data)/np.sqrt(data)
    for ii in range(len(resids)):
        for jj in range(len(resids)):
            for kk in range(len(resids)):
                if resids.mask[ii,jj,kk] == True:
                    continue
                elif np.isnan(resids[ii,jj,kk]) == True:
                    resids.mask[ii,jj,kk] = True
                elif np.isinf(resids[ii,jj,kk]) == True:
                    resids.mask[ii,jj,kk] = True

    fig = plt.figure()
    plt.clf()
    ax = mplot3d.Axes3D(fig)

    vmin = -resid_range
    vmax = resid_range

    resids = np.swapaxes(resids,0,2)

    toplot = np.logical_not(resids.mask)
    #toplot = np.logical_and(toplot, resids.data >= vmin)
    if max_sum != None:
        iis = np.where(toplot)[0]
        jjs = np.where(toplot)[1]
        kks = np.where(toplot)[2]
        for ll in range(len(iis)):
            ii = iis[ll]
            jj = jjs[ll]
            kk = kks[ll]
            if ii+jj+kk > max_sum:
                toplot[ii,jj,kk] = False
    if max_index != None:
        iis = np.where(toplot)[0]
        jjs = np.where(toplot)[1]
        kks = np.where(toplot)[2]
        for ll in range(len(iis)):
            ii = iis[ll]
            jj = jjs[ll]
            kk = kks[ll]
            if ii > max_index and jj > max_index and kk > max_index:
                toplot[ii,jj,kk] = False


    normalized = resids/(vmax-vmin) + .5
    #normalized = np.minimum(normalized, 1)
    # scrunch by a factor 
    # XXX: this is really hacky
    colors = plt.cm.RdBu(normalized)

    # We draw by calculating which faces are visible and including each as a
    # polygon.
    polys, polycolors = [],[]
    for ii in range(resids.shape[0]):
        for jj in range(resids.shape[1]):
            for kk in range(resids.shape[2]):
                if not toplot[ii,jj,kk]:
                    continue
                if kk < resids.shape[2]-1 and toplot[ii,jj,kk+1]:
                    pass
                else:
                    polys.append([[ii-0.5,jj+0.5,kk+0.5],[ii+0.5,jj+0.5,kk+0.5],
                                  [ii+0.5,jj-0.5,kk+0.5],[ii-0.5,jj-0.5,kk+0.5]]
                                 )
                    polycolors.append(colors[ii,jj,kk])
                if kk > 0 and toplot[ii,jj,kk-1]:
                    pass
                else:
                    polys.append([[ii-0.5,jj+0.5,kk-0.5],[ii+0.5,jj+0.5,kk-0.5],
                                  [ii+0.5,jj-0.5,kk-0.5],[ii-0.5,jj-0.5,kk-0.5]]
                                 )
                    polycolors.append(colors[ii,jj,kk])
                if jj < resids.shape[1]-1 and toplot[ii,jj+1,kk]:
                    pass
                else:
                    polys.append([[ii-0.5,jj+0.5,kk+0.5],[ii+0.5,jj+0.5,kk+0.5],
                                  [ii+0.5,jj+0.5,kk-0.5],[ii-0.5,jj+0.5,kk-0.5]]
                                 )
                    polycolors.append(colors[ii,jj,kk])
                if jj > 0 and toplot[ii,jj-1,kk]:
                    pass
                else:
                    polys.append([[ii-0.5,jj-0.5,kk+0.5],[ii+0.5,jj-0.5,kk+0.5],
                                  [ii+0.5,jj-0.5,kk-0.5],[ii-0.5,jj-0.5,kk-0.5]]
                                 )
                    polycolors.append(colors[ii,jj,kk])
                if ii < resids.shape[0]-1 and toplot[ii+1,jj,kk]:
                    pass
                else:
                    polys.append([[ii+0.5,jj-0.5,kk+0.5],[ii+0.5,jj+0.5,kk+0.5],
                                  [ii+0.5,jj+0.5,kk-0.5],[ii+0.5,jj-0.5,kk-0.5]]
                                 )
                    polycolors.append(colors[ii,jj,kk])
                if ii > 0 and toplot[ii-1,jj,kk]:
                    pass
                else:
                    polys.append([[ii-0.5,jj-0.5,kk+0.5],[ii-0.5,jj+0.5,kk+0.5],
                                  [ii-0.5,jj+0.5,kk-0.5],[ii-0.5,jj-0.5,kk-0.5]]
                                 )
                    polycolors.append(colors[ii,jj,kk])

    polycoll = mplot3d.art3d.Poly3DCollection(polys, facecolor=polycolors, 
                                              edgecolor='k', linewidths=0.5)
    ax.add_collection(polycoll)

    # Set the limits
    ax.set_xlim3d(-0.5,resids.shape[0]-0.5)
    ax.set_ylim3d(-0.5,resids.shape[1]-0.5)
    ax.set_zlim3d(-0.5,resids.shape[2]-0.5)

    ax.set_xlabel('aB', horizontalalignment='left')
    ax.set_ylabel('Ab', verticalalignment='bottom')
    ax.set_zlabel('AB', verticalalignment='bottom')

    ax.view_init(elev=30., azim=45)
    if show == True:
        plt.show()