def MakePlans(Data, model, LP, Q=None, **kwargs):
''' Create list of Plans
'''
newTopics, Info = makeCandidateTopics(Data, Q, model, LP, **kwargs)
if 'doVizBirth' in kwargs and kwargs['doVizBirth']:
from matplotlib import pylab
pylab.imshow(newTopics, vmin=0, vmax=0.01,
aspect=Data.vocab_size / newTopics.shape[0],
interpolation='nearest')
Plans = list()
for kk in xrange(newTopics.shape[0]):
Plan = dict(ktarget=None, Data=None, targetWordIDs=None,
targetWordFreq=newTopics[kk])
# Add material to the log
topWords = np.argsort(-1 * Plan['targetWordFreq'])[:10]
if hasattr(Data, 'vocab_dict'):
Vocab = getVocab(Data)
topWordStr = ' '.join([Vocab[w] for w in topWords])
else:
topWordStr = ' '.join([str(w) for w in topWords])
Plan['log'] = list()
Plan['topWordIDs'] = topWords
Plan['log'].append(topWordStr)
if 'anchorID' in Info:
anchorWordStr = 'Anchor: ' + str(Info['anchorID'][kk])
Plan['anchorWordID'] = anchorWordStr
Plan['log'].append(anchorWordStr)
Plans.append(Plan)
return Plans
python类imshow()的实例源码
def roll_2Dprofile(profile,position,padvalue=0.0,showprofiles=False):
"""
Move 2D profile to given position in array by rolling it in x and y.
Note that the roll does not handle sub-pixel moves.
tu.shift_2Dprofile() does this using interpolation
--- INPUT ---
profile profile to shift
position position to move center of image (profile) to: [ypos,xpos]
padvalue the values to padd the images with when shifting profile
showprofiles Show profile when shifted?
--- EXAMPLE OF USE ---
tu.roll_2Dprofile(gauss2D,)
"""
profile_dim = profile.shape
yroll = np.int(np.round(position[0]-profile_dim[0]/2.))
xroll = np.int(np.round(position[1]-profile_dim[1]/2.))
profile_shifted = np.roll(np.roll(profile,yroll,axis=0),xroll,axis=1)
if showprofiles:
vmaxval = np.max(profile_shifted)
plt.imshow(profile_shifted,interpolation='none',vmin=-vmaxval, vmax=vmaxval)
plt.title('Positioned Source')
plt.show()
if yroll < 0:
profile_shifted[yroll:,:] = padvalue
else:
profile_shifted[:yroll,:] = padvalue
if xroll < 0:
profile_shifted[:,xroll:] = padvalue
else:
profile_shifted[:,:xroll] = padvalue
if showprofiles:
plt.imshow(profile_shifted,interpolation='none',vmin=-vmaxval, vmax=vmaxval)
plt.title('Positioned Source with 0s inserted')
plt.show()
return profile_shifted
# = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
def optimize_img_scale(img_data,img_std,img_model,optimizer='curve_fit',show_residualimg=False,verbose=True):
"""
optimize the (flux) scaling of an image with respect to a (noisy) data image
--- INPUT ---
img_data The (noisy) data image to scale model image provide in img_model to
img_std Standard deviation image for data to use in optimization
img_model Model image to (flux) scale to match img_data
optimizer The optimizer to use when scaling the layers
show_residualimg To show the residual image (data - model) for the optimize layers, set this to true
verbose Toggle verbosity
--- EXAMPLE OF USE ---
import tdose_model_cube as tmc
scale, cov = tmc.optimize_img_scale()
"""
if verbose: print ' - Optimize residual between model (multiple Gaussians) and data with least squares in 2D'
if verbose: print ' ----------- Started on '+tu.get_now_string()+' ----------- '
# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
if optimizer == 'leastsq':
sys.exit('optimizer = "leastsq" no enabled')
# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
elif optimizer == 'curve_fit':
imgsize = img_data.shape
xgrid, ygrid = tu.gen_gridcomponents(imgsize)
scale_best, scale_cov = opt.curve_fit(lambda (xgrid, ygrid), scale:
tmc.curve_fit_fct_wrapper_imgscale((xgrid, ygrid), scale, img_model),
(xgrid, ygrid),
img_data.ravel(), p0 = [1.0], sigma=img_std.ravel() )
output = scale_best, scale_cov
else:
sys.exit(' ---> Invalid optimizer ('+optimizer+') chosen in optimize_img_scale()')
# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
if verbose: print ' ----------- Finished on '+tu.get_now_string()+' ----------- '
# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
if show_residualimg:
if verbose: print ' - Displaying the residual image between data and scaled model image '
res_img = img_model-img_data
plt.imshow(res_img,interpolation='none', vmin=1e-5, vmax=np.max(res_img), norm=mpl.colors.LogNorm())
plt.title('Initial Residual = Initial Model Image - Data Image')
plt.show()
res_img = img_model*scale_best-img_data
plt.imshow(res_img,interpolation='none', vmin=1e-5, vmax=np.max(res_img), norm=mpl.colors.LogNorm())
plt.title('Best Residual = Scaled (by '+str(scale_best)+') Model Image - Data Image')
plt.show()
# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
return output
# = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
def run_mnist():
np.random.seed(42)
# import dataset
f = gzip.open('./tmp/data/mnist.pkl.gz', 'rb')
(x_train, t_train), (x_valid, t_valid), (x_test, t_test) = cPickle.load(f)
f.close()
Y = x_train[:100, :]
labels = t_train[:100]
Y[Y < 0.5] = -1
Y[Y > 0.5] = 1
# inference
print "inference ..."
M = 30
D = 2
# lvm = vfe.SGPLVM(Y, D, M, lik='Gaussian')
lvm = vfe.SGPLVM(Y, D, M, lik='Probit')
# lvm.train(alpha=0.5, no_epochs=10, n_per_mb=100, lrate=0.1, fixed_params=['sn'])
lvm.optimise(method='L-BFGS-B')
plt.figure()
mx, vx = lvm.get_posterior_x()
zu = lvm.sgp_layer.zu
plt.scatter(mx[:, 0], mx[:, 1], c=labels)
plt.plot(zu[:, 0], zu[:, 1], 'ko')
nx = ny = 30
x_values = np.linspace(-5, 5, nx)
y_values = np.linspace(-5, 5, ny)
sx = 28
sy = 28
canvas = np.empty((sx * ny, sy * nx))
for i, yi in enumerate(x_values):
for j, xi in enumerate(y_values):
z_mu = np.array([[xi, yi]])
x_mean, x_var = lvm.predict_f(z_mu)
t = x_mean / np.sqrt(1 + x_var)
Z = 0.5 * (1 + special.erf(t / np.sqrt(2)))
canvas[(nx - i - 1) * sx:(nx - i) * sx, j *
sy:(j + 1) * sy] = Z.reshape(sx, sy)
plt.figure(figsize=(8, 10))
Xi, Yi = np.meshgrid(x_values, y_values)
plt.imshow(canvas, origin="upper", cmap="gray")
plt.tight_layout()
plt.show()
def plot_model_no_control(model, plot_title='', name_suffix=''):
# plot function
mx, vx = model.get_posterior_x()
mins = np.min(mx, axis=0) - 0.5
maxs = np.max(mx, axis=0) + 0.5
nGrid = 50
xspaced = np.linspace(mins[0], maxs[0], nGrid)
yspaced = np.linspace(mins[1], maxs[1], nGrid)
xx, yy = np.meshgrid(xspaced, yspaced)
Xplot = np.vstack((xx.flatten(), yy.flatten())).T
mf, vf = model.predict_f(Xplot)
fig = plt.figure()
plt.imshow((mf[:, 0]).reshape(*xx.shape),
vmin=mf.min(), vmax=mf.max(), origin='lower',
extent=[mins[0], maxs[0], mins[1], maxs[1]], aspect='auto')
plt.colorbar()
plt.contour(
xx, yy, (mf[:, 0]).reshape(*xx.shape),
colors='k', linewidths=2, zorder=100)
zu = model.dyn_layer.zu
plt.plot(zu[:, 0], zu[:, 1], 'wo', mew=0, ms=4)
for i in range(mx.shape[0] - 1):
plt.plot(mx[i:i + 2, 0], mx[i:i + 2, 1],
'-bo', ms=3, linewidth=2, zorder=101)
plt.xlabel(r'$x_{t, 1}$')
plt.ylabel(r'$x_{t, 2}$')
plt.xlim([mins[0], maxs[0]])
plt.ylim([mins[1], maxs[1]])
plt.title(plot_title)
plt.savefig('/tmp/hh_gpssm_dim_0' + name_suffix + '.pdf')
fig = plt.figure()
plt.imshow((mf[:, 1]).reshape(*xx.shape),
vmin=mf.min(), vmax=mf.max(), origin='lower',
extent=[mins[0], maxs[0], mins[1], maxs[1]], aspect='auto')
plt.colorbar()
plt.contour(
xx, yy, (mf[:, 1]).reshape(*xx.shape),
colors='k', linewidths=2, zorder=100)
zu = model.dyn_layer.zu
plt.plot(zu[:, 0], zu[:, 1], 'wo', mew=0, ms=4)
for i in range(mx.shape[0] - 1):
plt.plot(mx[i:i + 2, 0], mx[i:i + 2, 1],
'-bo', ms=3, linewidth=2, zorder=101)
plt.xlabel(r'$x_{t, 1}$')
plt.ylabel(r'$x_{t, 2}$')
plt.xlim([mins[0], maxs[0]])
plt.ylim([mins[1], maxs[1]])
plt.title(plot_title)
plt.savefig('/tmp/hh_gpssm_dim_1' + name_suffix + '.pdf')
def run_mnist():
np.random.seed(42)
# import dataset
f = gzip.open('./tmp/data/mnist.pkl.gz', 'rb')
(x_train, t_train), (x_valid, t_valid), (x_test, t_test) = cPickle.load(f)
f.close()
Y = x_train[:100, :]
labels = t_train[:100]
Y[Y < 0.5] = -1
Y[Y > 0.5] = 1
# inference
print "inference ..."
M = 30
D = 2
# lvm = aep.SGPLVM(Y, D, M, lik='Gaussian')
lvm = aep.SGPLVM(Y, D, M, lik='Probit')
# lvm.train(alpha=0.5, no_epochs=10, n_per_mb=100, lrate=0.1, fixed_params=['sn'])
lvm.optimise(method='L-BFGS-B', alpha=0.1)
plt.figure()
mx, vx = lvm.get_posterior_x()
zu = lvm.sgp_layer.zu
plt.scatter(mx[:, 0], mx[:, 1], c=labels)
plt.plot(zu[:, 0], zu[:, 1], 'ko')
nx = ny = 30
x_values = np.linspace(-5, 5, nx)
y_values = np.linspace(-5, 5, ny)
sx = 28
sy = 28
canvas = np.empty((sx * ny, sy * nx))
for i, yi in enumerate(x_values):
for j, xi in enumerate(y_values):
z_mu = np.array([[xi, yi]])
x_mean, x_var = lvm.predict_f(z_mu)
t = x_mean / np.sqrt(1 + x_var)
Z = 0.5 * (1 + special.erf(t / np.sqrt(2)))
canvas[(nx - i - 1) * sx:(nx - i) * sx, j *
sy:(j + 1) * sy] = Z.reshape(sx, sy)
plt.figure(figsize=(8, 10))
Xi, Yi = np.meshgrid(x_values, y_values)
plt.imshow(canvas, origin="upper", cmap="gray")
plt.tight_layout()
plt.show()
def plot_waterfall(fil, f_start=None, f_stop=None, if_id=0, logged=True,cb=False,freq_label=False,MJD_time=False, **kwargs):
""" Plot waterfall of data
Args:
f_start (float): start frequency, in MHz
f_stop (float): stop frequency, in MHz
logged (bool): Plot in linear (False) or dB units (True),
cb (bool): for plotting the colorbar
kwargs: keyword args to be passed to matplotlib imshow()
"""
matplotlib.rc('font', **font)
plot_f, plot_data = fil.grab_data(f_start, f_stop, if_id)
# Make sure waterfall plot is under 4k*4k
dec_fac_x, dec_fac_y = 1, 1
if plot_data.shape[0] > MAX_IMSHOW_POINTS[0]:
dec_fac_x = plot_data.shape[0] / MAX_IMSHOW_POINTS[0]
if plot_data.shape[1] > MAX_IMSHOW_POINTS[1]:
dec_fac_y = plot_data.shape[1] / MAX_IMSHOW_POINTS[1]
plot_data = rebin(plot_data, dec_fac_x, dec_fac_y)
if MJD_time:
extent=(plot_f[0], plot_f[-1], fil.timestamps[-1], fil.timestamps[0])
else:
extent=(plot_f[0], plot_f[-1], (fil.timestamps[-1]-fil.timestamps[0])*24.*60.*60, 0.0)
this_plot = plt.imshow(plot_data,
aspect='auto',
rasterized=True,
interpolation='nearest',
extent=extent,
cmap='viridis_r',
**kwargs
)
if cb:
plt.colorbar()
if freq_label:
plt.xlabel("Frequency [Hz]",fontdict=font)
if MJD_time:
plt.ylabel("Time [MJD]",fontdict=font)
else:
plt.ylabel("Time [s]",fontdict=font)
return this_plot
