netutils.py 文件源码

def sample_gamma_posterior (data, shape0, scale0):
    S = numpy.sum(data)
    Slog = numpy.sum(map(numpy.log, data))
    N = len(data)

    logP = 1.0+Slog
    q = 1.0 + S
    r = 1.0 + N
    s = 1.0 + N
#    lnf_shape = lambda al: (al-1)*logP + special.gammaln(s*al+1) - (1-al*s)*numpy.log(q) -r*special.gammaln(al)

    lnf_shape = lambda al: (al-1)*logP - q/scale0 - r*special.gammaln(al) - (al*s)*numpy.log(scale0)
    shape_new = sampling.slice (lnf_shape, shape0, thin=5, N=1, lower=0, upper=numpy.inf)[0]
#    print "Sfoo", shape0, lnf_shape(shape0)
#    print "....", shape_new, lnf_shape(shape_new)

#    lnf_rate = lambda beta: -beta*q + s*shape_new*numpy.log(beta)
#    rate_new = sampling.slice (lnf_rate, 1.0/scale0, thin=5, N=1, lower=0, upper=numpy.inf)[0]
    lnf_scale = lambda beta: (shape_new-1)*logP - q/beta - r*special.gammaln(shape_new) - (shape_new*s)*numpy.log(beta)
    scale_new = sampling.slice (lnf_scale, scale0, thin=5, N=1, lower=0, upper=numpy.inf)[0]
#    print "Rfoo", 1.0/scale0, lnf_scale(scale0)
#    print "...", scale_new, lnf_scale(scale_new)

    return [shape_new, scale_new]