Python triangle.corner函数代码示例

def triangleMAPs(savefilename,basename):
    with open(savefilename,'rb') as savefile:
        bf= numpy.array(pickle.load(savefile))
        samples= numpy.array(pickle.load(savefile))
        bf_g15= numpy.array(pickle.load(savefile))
        samples_g15= numpy.array(pickle.load(savefile))
        bf_zero= numpy.array(pickle.load(savefile))
        samples_zero= numpy.array(pickle.load(savefile))
    labels= []
    for jj in range(samples.shape[2]):
        labels.append(r"$\mathrm{param}\ %i$" % jj)
    maps= define_rcsample.MAPs()
    for ii, map in enumerate(maps.map()):
        if ii >= len(bf): break
        tfeh= numpy.nanmedian(map['FE_H'])
        tafe= numpy.nanmedian(map[define_rcsample._AFETAG])
        for tbf,tsamples,ext in zip([bf,bf_g15,bf_zero],
                                    [samples,samples_g15,samples_zero],
                                    ['fid','g15','zero']):
            try:
                triangle.corner(tsamples[ii,].T,quantiles=[0.16, 0.5, 0.84],
                                labels=labels,
                                show_titles=True,title_args={"fontsize": 12},
                                bins=21)
            except ValueError: pass
            else:
                bovy_plot.bovy_text(r'$[\mathrm{{Fe/H}}] = {feh:.1f},$'\
                                        .format(feh=tfeh)+'\n'
                                    +r'$[\alpha/\mathrm{{Fe}}] = {afe:.2f}$'\
                                        .format(afe=tafe),
                                    top_left=True,size=16.)
                bovy_plot.bovy_end_print(basename+"_%i_%s.png" % (ii,ext))
    return None

def plot_corner(logpost, chain, outdir=None):
    flatchain = chain.reshape((-1, chain.shape[2]))

    tri.corner(flatchain, labels=logpost.pnames, quantiles=[0.05, 0.95])

    if outdir is not None:
        pp.savefig(op.join(outdir, 'corner.pdf'))

 def plot_triangle(self, **kwargs):
     try:
         from triangle import corner
     except ImportError:
         raise ImportError('Plotting requires trianglepy')
     corner(self.posterior_data[self.parameter_names],
            labels=self.parameter_names, **kwargs)

def MCParTemp(nwalkers,steps,ntemps=40,burnin=100,threads=1):
    t1 = time.time()
    # use emcee parallel tempering sampler to find prob. surface
    ndim = 4
    startpos = [[[0,0,0,0] + 5e-3*np.random.randn(ndim) for i in range(nwalkers)] for i in range(ntemps)]

    # time how long it should take
    t0 = time.time()
    testpos = [[[0,0,0,0] + 5e-3*np.random.randn(4) for i in range(8)] for i in range(8)]
    sampler  = emcee.PTSampler(8,8, 4, lnlike, xs.lnprior, threads=threads)
    for ptest,logprobtest, logliketest in sampler.sample(testpos,iterations=10):
        pass
    sampler.reset()
    tdiff = time.time() - t0
    esttime = ntemps/8.*nwalkers/8.*(burnin+steps)/10.*tdiff/60.
    print 'Estimated time = {} minutes'.format(esttime)

    # burn-in for defined number of steps (default=100)
    sampler  = emcee.PTSampler(ntemps,nwalkers, ndim, lnlike, xs.lnprior, threads=threads)
    for p,logprob, loglike in sampler.sample(startpos,iterations=burnin):
        pass
    sampler.reset()

    # run mcmc
    sampler = emcee.PTSampler(ntemps, nwalkers, ndim, lnlike, xs.lnprior, threads=threads)
    for p, logprob, loglike in sampler.sample(p, lnprob0=logprob, lnlike0=loglike,iterations=steps):
        pass

    print 'Time: {} minutes'.format(time.time() - t0)
    samples = sampler.chain[:,:,:,:].reshape((-1, ndim))
    triangle.corner(samples,labels=[r'$\Delta s$',r'$M_A$',r'$F_1^S$',r'$\mu_S$'],quantiles=[0.05,0.5,0.95])
    print sampler.thermodynamic_integration_log_evidence(fburnin=0)
    return samples

def plot_mcmc_results(chain):
    # Pull m and b arrays out of the Markov chain.
    mm = [m for b,m in chain]
    bb = [b for b,m in chain]
    # Scatterplot of m,b posterior samples
    plt.clf()
    plt.contour(bgrid, mgrid, posterior, pdf_contour_levels(posterior))
    plt.plot(bb, mm, 'b.', alpha=0.1)
    plot_mb_setup()
    plt.show()
    
    # Histograms
    import triangle
    triangle.corner(chain, labels=['b','m'], extents=[0.99]*2)
    plt.show()
    
    # Traces
    plt.clf()
    plt.subplot(2,1,1)
    plt.plot(mm, 'k-')
    plt.ylim(mlo,mhi)
    plt.ylabel('m')
    plt.subplot(2,1,2)
    plt.plot(bb, 'k-')
    plt.ylabel('b')
    plt.ylim(blo,bhi)
    plt.show()

    def plotCorner(self):
        model = self.model

        mcmcVParams = self.mcmcRes.vparam_names
        nestVParams = self.nestRes.vparam_names

        mcmcSamples = self.mcmcRes.samples
        nestSamples = self.nestRes.samples

        mcmc_ndim, mcmc_nsamples = len(mcmcVParams), len(mcmcSamples)
        nest_ndim, nest_nsamples = len(nestVParams), len(nestSamples)

        # make figure

        figure_mcmc = triangle.corner(mcmcSamples, labels=[mcmcVParams[0], mcmcVParams[1], mcmcVParams[2], mcmcVParams[3]],
                     truths=[model.get(mcmcVParams[0]), model.get(mcmcVParams[1]),
                             model.get(mcmcVParams[2]), model.get(mcmcVParams[3])],
                     range=mcmc_ndim*[0.9999],
                     show_titles=True, title_args={"fontsize": 12})

        figure_mcmc.gca().annotate("mcmc sampling", xy=(0.5, 1.0), xycoords="figure fraction",
                  xytext=(0, -5), textcoords="offset points",
                  ha="center", va="top")

        figure_nest = triangle.corner(nestSamples, labels=[nestVParams[0], nestVParams[1], nestVParams[2], nestVParams[3]],
                     truths=[model.get(nestVParams[0]), model.get(nestVParams[1]),
                             model.get(nestVParams[2]), model.get(nestVParams[3])],
                     weights=self.nestRes.weights, range=nest_ndim*[0.9999],
                     show_titles=True, title_args={"fontsize": 12})

        figure_nest.gca().annotate("nest sampling", xy=(0.5, 1.0), xycoords="figure fraction",
                  xytext=(0, -5), textcoords="offset points",
                  ha="center", va="top")

        return figure_mcmc, figure_nest

	def plot_PDFtriangle(self,parameterset, labels):		

		if parameterset=='10pars':
			figure = triangle.corner(self.chain.flatchain, labels= labels, plot_contours=True, plot_datapoints = False, show_titles=True, quantiles=[0.16, 0.50, 0.84])
		elif parameterset == 'int_lums':
			figure = triangle.corner(self.int_lums.T, labels= labels,   plot_contours=True, plot_datapoints = False, show_titles=True, quantiles=[0.16, 0.50, 0.84])
		return figure

 def plot_triangle(self, **kwargs):
     try:
         from triangle import corner
     except ImportError:
         raise ImportError('Plotting requires trianglepy')
     data_columns = self.posterior_data.columns[2:]
     corner(self.posterior_data[data_columns],
            weights=self.posterior_data.posterior, **kwargs)

def mcmctriangle(ID):
	filename=('tables/ndim_' + str(ndim) + '_walkers_' + str(nwalkers) + '.fits')
	data = Table.read(filename)
	data_t = np.array([data[key] for key in param_keys]).transpose()
	truths = [np.median(data[key]) for key in param_keys]
	triangle.corner(data_t,labels=param_labels,quantiles=[0.1587,0.5000,0.8413])
	plt.savefig('plots/triangle_ndim_' + str(ndim) + '_walkers_' + str(nwalkers) + '.pdf')
	plt.close()

def convergence_plots(sampler):
    plt.figure()
    plt.plot(sampler.lnprobability.T)
    plt.figure()
    pu.plot_emcee_chains(sampler.chain)
    triangle.corner(sampler.flatchain)

    print 'Autocorrelation lengths: ', ac.emcee_chain_autocorrelation_lengths(sampler.chain)
    print 'Gelman-Rubin R: ', ac.emcee_gelman_rubin_r(sampler.chain)

def PlotTriangle(fileroot,usetruths=True):
	data = np.loadtxt(fileroot+'post_equal_weights.dat', usecols=(0,1,2,3,4,5))
	if (usetruths):
		truths = np.loadtxt(fileroot+'injected_values.txt')
		figure = triangle.corner(data, labels=names, truths=truths, bins=30, quantiles=[0.05, 0.5, 0.95],
								 show_titles=True, title_args={"fontsize": 12})
	else:
		figure = triangle.corner(data, labels=names, bins=30, quantiles=[0.05, 0.5, 0.95],
								 show_titles=True, title_args={"fontsize": 12})
	figure.savefig(fileroot+'posterior_plot.png')

def plot_triangle(samples):
    import triangle
    if samples.shape[1] == 2:
        fig = triangle.corner(samples, labels=["$t_0$", "$P$"])
    if samples.shape[1] == 6:
        fig = triangle.corner(samples, labels=["$t_0$", r"depth", r"duration",
                                               r"$b$", "$q_1$", "$q_2$"])
    elif samples.shape[1] == 8:
        fig = triangle.corner(samples, labels=["$t_0$", r"depth", r"duration",
                                               r"$b$", "$q_1$", "$q_2$", "$q_3$", "$q_4$"])
    plt.show()

def triangle_plot_chain(chain, lnprob, prefix, truth=None, truth_lnprob=0.):
    """
    Make a 7x7 triangle.
    """
    nx, nq = chain.shape
    maxlnp = np.max(lnprob)
    lnpextent = [(maxlnp-14.5, maxlnp+0.5)]
    bar = SixPosition(chain[0]) # temporary variable to get names
    foo = np.concatenate((chain, lnprob.reshape((nx, 1))), axis=1)
    if truth is None:
        truths = None
    else:
        truths = np.append(truth, [truth_lnprob])
    labels = np.append(bar.get_sixpos_names(), [r"$\ln p$"])
    extents = bar.get_sixpos_extents() + lnpextent
    fig = tri.corner(foo, labels=labels, extents=extents,
                     truths=truths, plot_contours=False)
    fn = prefix + "a.png"
    print "triangle_plot_chain(): writing " + fn
    fig.savefig(fn)
    obsfoo = 1. * foo # copy
    for i in range(nx):
        obsfoo[i,:6] = SixPosition(foo[i,:6]).get_observables_array()
    if truth is None:
        trueobs = None
    else:
        trueobs = 1. * truths # copy
        trueobs[:6] = SixPosition(truth).get_observables_array()
    labels = np.append(bar.get_observables_names(), [r"$\ln p$"])
    extents = bar.get_observables_extents() + lnpextent
    fig = tri.corner(obsfoo, labels=labels, extents=extents,
                     truths=trueobs, plot_contours=False)
    fn = prefix + "b.png"
    print "triangle_plot_chain(): writing " + fn
    fig.savefig(fn)
    intfoo = 1. * foo[:,2:] # copy
    for i in range(nx):
        intfoo[i,:4] = SixPosition(foo[i,:6]).get_integrals_of_motion()
    if truth is None:
        trueint = None
    else:
        trueint = 1. * truths[2:] # copy
        trueint[:4] = SixPosition(truth).get_integrals_of_motion()
    labels = np.append(bar.get_integrals_of_motion_names(), [r"$\ln p$"])
    extents = bar.get_integrals_of_motion_extents() + lnpextent
    fig = tri.corner(intfoo, labels=labels, extents=extents,
                     truths=trueint, plot_contours=False)
    fn = prefix + "c.png"
    print "triangle_plot_chain(): writing " + fn
    fig.savefig(fn)
    return None

    def visModel(self, labels=None, fname=None):

        nSamples = 1e6
        samples = self.reg.sample(n_samples=nSamples)
        #quick fix for better plotting
        samples = np.log10(samples)

        if self.pred == None:
            if hasTriangle:
                figure = triangle.corner(samples, labels=labels,
                                         quantiles=[0.16, 0.5, 0.84],
                                         show_titles=True, title_args={"fontsize": 12})
                if fname!=None:
                    plt.savefig('predicted_'+fname)

            elif samples.size[1]==2:
                f, ax = plt.subplots(1)
                ax.hist2d(samples[:,0], samples[:,1])
                if fname!=None:
                    plt.savefig('predicted_'+fname)

            else:
                raise NotImplementedError("Plotting datasets w/ dim > 2 without Triangle not implemented")

        else:
            if hasTriangle:
                figure = triangle.corner(samples, labels=labels,
                                         quantiles=[0.16, 0.5, 0.84],
                                         show_titles=True, title_args={"fontsize": 12})
                if fname!=None:
                    plt.savefig('predicted_'+fname)
                
                figure = triangle.corner(self.X, labels=labels,
                                         quantiles=[0.16, 0.5, 0.84],
                                         show_titles=True, title_args={"fontsize": 12})
                if fname!=None:
                    plt.savefig('original_'+fname)

            elif samples.size[1]==2:
                f, ax = plt.subplots(2)
                ax[0].hist2d(samples[:,0], samples[:,1])
                ax[1].hist2d(self.X[:,0], self.X[:,1])
                ax[1].set_xlabel(labels[0])
                ax[1].set_ylabel(labels[1])
                plt.tight_layout()
                if fname!=None:
                    plt.savefig(fname)

            else:
                raise NotImplementedError("Plotting datasets w/ dim > 2 without Triangle not implemented")

def plot_grid_points(db='phoenix2016.db3',saveplot=False):
    # plot the grid points
    connection = sqlite.connect(db)
    sql_query = 'SELECT teff,logg,mh,alpha FROM parameter_sets'
    tab = pd.read_sql_query(sql_query,connection)

    connection.close()

    arr = np.array([np.array(tab['teff']),np.array(tab['logg']),np.array(tab['mh']),tab['alpha']])
    corner(tab,plot_contours=False,labels=['Teff','log g', '[M/H]', '[alpha/Fe]'])
    if saveplot:
        fname = os.path.basename(db)
        fname = os.path.splitext(fname)
        plt.savefig(fname[0]+'.pdf')