dhmc_sampler.py 文件源码

def run_sampler(self, theta0, dt_range, nstep_range, n_burnin, n_sample, seed=None, n_update=10):
        """Run DHMC and return samples and some additional info."""

        np.random.seed(seed)

        # Run HMC.
        theta = theta0
        n_per_update = math.ceil((n_burnin + n_sample) / n_update)
        pathlen_ave = 0
        samples = np.zeros((n_sample + n_burnin, len(theta)))
        logp_samples = np.zeros(n_sample + n_burnin)
        accept_prob = np.zeros(n_sample + n_burnin)

        tic = time.process_time()  # Start clock
        logp, grad, aux = self.f(theta)
        for i in range(n_sample + n_burnin):
            dt = np.random.uniform(dt_range[0], dt_range[1])
            nstep = np.random.randint(nstep_range[0], nstep_range[1] + 1)
            theta, logp, grad, aux, accept_prob[i], pathlen \
                = self.hmc(dt, nstep, theta, logp, grad, aux)
            pathlen_ave = i / (i + 1) * pathlen_ave + 1 / (i + 1) * pathlen
            samples[i, :] = theta
            logp_samples[i] = logp
            if (i + 1) % n_per_update == 0:
                print('{:d} iterations have been completed.'.format(i + 1))

        toc = time.process_time()
        time_elapsed = toc - tic
        print(('The average path length of each DHMC iteration was '
               '{:.2f}.'.format(pathlen_ave)))

        return samples, logp_samples, accept_prob, pathlen_ave, time_elapsed