Python special.gamma方法代码示例

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def _upper_incomplete_gamma(z, a):
    """
    An implementation of the non-regularised upper incomplete gamma
    function. Computed using the relationship with the regularised 
    lower incomplete gamma function (scipy.special.gammainc). 
    Uses the recurrence relation wherever z<0.
    """
    n = int(-np.floor(z))
    if n > 0:
        z = z + n
        u_gamma = (1. - gammainc(z, a))*gamma(z)
        
        for i in range(n):
            z = z - 1.
            u_gamma = (u_gamma - np.power(a, z)*np.exp(-a))/z
        return u_gamma
    else:
        return (1. - gammainc(z, a))*gamma(z) 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def compute_finite_difference_coefficients(derivative_order, grid_size):

    # from http://www.scientificpython.net/pyblog/uniform-finite-differences-all-orders

    n = 2*grid_size -1
    A = np.tile(np.arange(grid_size), (n,1)).T
    B = np.tile(np.arange(1-grid_size,grid_size), (grid_size,1))
    M = (B**A)/gamma(A+1)

    r = np.zeros(grid_size)
    r[derivative_order] = 1

    D = np.zeros((grid_size, grid_size))
    for k in range(grid_size):
        indexes = k + np.arange(grid_size)
        D[:,k] = solve(M[:,indexes], r)

    return D

#################################################################################################### 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def _pdf_skip(self, x, dfn, dfd, nc):
        # ncf.pdf(x, df1, df2, nc) = exp(nc/2 + nc*df1*x/(2*(df1*x+df2))) *
        #             df1**(df1/2) * df2**(df2/2) * x**(df1/2-1) *
        #             (df2+df1*x)**(-(df1+df2)/2) *
        #             gamma(df1/2)*gamma(1+df2/2) *
        #             L^{v1/2-1}^{v2/2}(-nc*v1*x/(2*(v1*x+v2))) /
        #             (B(v1/2, v2/2) * gamma((v1+v2)/2))
        n1, n2 = dfn, dfd
        term = -nc/2+nc*n1*x/(2*(n2+n1*x)) + sc.gammaln(n1/2.)+sc.gammaln(1+n2/2.)
        term -= sc.gammaln((n1+n2)/2.0)
        Px = np.exp(term)
        Px *= n1**(n1/2) * n2**(n2/2) * x**(n1/2-1)
        Px *= (n2+n1*x)**(-(n1+n2)/2)
        Px *= sc.assoc_laguerre(-nc*n1*x/(2.0*(n2+n1*x)), n2/2, n1/2-1)
        Px /= sc.beta(n1/2, n2/2)
        # This function does not have a return.  Drop it for now, the generic
        # function seems to work OK. 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def _pdf(self, x, df, nc):
        # nct.pdf(x, df, nc) =
        #                               df**(df/2) * gamma(df+1)
        #                ----------------------------------------------------
        #                2**df*exp(nc**2/2) * (df+x**2)**(df/2) * gamma(df/2)
        n = df*1.0
        nc = nc*1.0
        x2 = x*x
        ncx2 = nc*nc*x2
        fac1 = n + x2
        trm1 = n/2.*np.log(n) + sc.gammaln(n+1)
        trm1 -= n*np.log(2)+nc*nc/2.+(n/2.)*np.log(fac1)+sc.gammaln(n/2.)
        Px = np.exp(trm1)
        valF = ncx2 / (2*fac1)
        trm1 = np.sqrt(2)*nc*x*sc.hyp1f1(n/2+1, 1.5, valF)
        trm1 /= np.asarray(fac1*sc.gamma((n+1)/2))
        trm2 = sc.hyp1f1((n+1)/2, 0.5, valF)
        trm2 /= np.asarray(np.sqrt(fac1)*sc.gamma(n/2+1))
        Px *= trm1+trm2
        return Px 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def _munp(self, n, beta, m):
        """
        Returns the n-th non-central moment of the crystalball function.
        """
        N = 1.0 / (m/beta / (m-1) * np.exp(-beta**2 / 2.0) + _norm_pdf_C * _norm_cdf(beta))

        def n_th_moment(n, beta, m):
            """
            Returns n-th moment. Defined only if n+1 < m
            Function cannot broadcast due to the loop over n
            """
            A = (m/beta)**m * np.exp(-beta**2 / 2.0)
            B = m/beta - beta
            rhs = 2**((n-1)/2.0) * sc.gamma((n+1)/2) * (1.0 + (-1)**n * sc.gammainc((n+1)/2, beta**2 / 2))
            lhs = np.zeros(rhs.shape)
            for k in range(n + 1):
                lhs += sc.binom(n, k) * B**(n-k) * (-1)**k / (m - k - 1) * (m/beta)**(-m + k + 1)
            return A * lhs + rhs

        return N * _lazywhere(np.atleast_1d(n + 1 < m),
                              (n, beta, m),
                              np.vectorize(n_th_moment, otypes=[np.float]),
                              np.inf) 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def mean(t, a, b):
        # TODO this is not tested yet.
        # tests:
        #    cemean(0., a, b)==mean(a, b, p)
        #    mean(t, a, 1., p)==mean(0., a, 1., p) == a
        # conditional excess mean
        # E[Y|y>t]
        # (conditional mean age at failure)
        # http://reliabilityanalyticstoolkit.appspot.com/conditional_distribution
        from scipy.special import gamma
        from scipy.special import gammainc
        # Regularized lower gamma
        print('not tested')

        v = 1. + 1. / b
        gv = gamma(v)
        L = np.power((t + .0) / a, b)
        cemean = a * gv * np.exp(L) * (1 - gammainc(v, t / a) / gv)

        return cemean 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def volume_unit_ball(n_dims):
    return (pi ** (.5 * n_dims)) / gamma(.5 * n_dims + 1) 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def ball_volume(ndim, radius):
    """Return the volume of a ball of dimension `ndim` and radius `radius`."""
    n = ndim
    r = radius
    return np.pi ** (n / 2) / special.gamma(n / 2 + 1) * r ** n 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def levy(self, D):
		r"""Levy function.

		Returns:
			float: Next Levy number.
		"""
		sigma = (Gamma(1 + self.beta) * sin(pi * self.beta / 2) / (Gamma((1 + self.beta) / 2) * self.beta * 2 ** ((self.beta - 1) / 2))) ** (1 / self.beta)
		return 0.01 * (self.normal(0, 1, D) * sigma / fabs(self.normal(0, 1, D)) ** (1 / self.beta)) 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def factorial(n):
    return gamma(n + 1) 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def _hs(k, cs, rho, omega):
    c0 = (cs * cs * (1 + rho * rho) / (1 - rho * rho) /
          (1 - 2 * rho * rho * cos(2 * omega) + rho ** 4))
    gamma = (1 - rho * rho) / (1 + rho * rho) / tan(omega)
    ak = abs(k)
    return c0 * rho ** ak * (cos(omega * ak) + gamma * sin(omega * ak)) 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def from_spline(cls, tck, extrapolate=None):
        """
        Construct a piecewise polynomial from a spline

        Parameters
        ----------
        tck
            A spline, as returned by `splrep` or a BSpline object.
        extrapolate : bool or 'periodic', optional
            If bool, determines whether to extrapolate to out-of-bounds points
            based on first and last intervals, or to return NaNs.
            If 'periodic', periodic extrapolation is used. Default is True.
        """
        if isinstance(tck, BSpline):
            t, c, k = tck.tck
            if extrapolate is None:
                extrapolate = tck.extrapolate
        else:
            t, c, k = tck

        cvals = np.empty((k + 1, len(t)-1), dtype=c.dtype)
        for m in xrange(k, -1, -1):
            y = fitpack.splev(t[:-1], tck, der=m)
            cvals[k - m, :] = y/spec.gamma(m+1)

        return cls.construct_fast(cvals, t, extrapolate) 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def fromspline(cls, xk, cvals, order, fill=0.0):
        # Note: this spline representation is incompatible with FITPACK
        N = len(xk)-1
        sivals = np.empty((order+1, N), dtype=float)
        for m in xrange(order, -1, -1):
            fact = spec.gamma(m+1)
            res = _fitpack._bspleval(xk[:-1], xk, cvals, order, m)
            res /= fact
            sivals[order-m, :] = res
        return cls(sivals, xk, fill=fill)


# The 3 private functions below can be called by splmake(). 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def _stats(self, df):
        mu = np.sqrt(2)*sc.gamma(df/2.0+0.5)/sc.gamma(df/2.0)
        mu2 = df - mu*mu
        g1 = (2*mu**3.0 + mu*(1-2*df))/np.asarray(np.power(mu2, 1.5))
        g2 = 2*df*(1.0-df)-6*mu**4 + 4*mu**2 * (2*df-1)
        g2 /= np.asarray(mu2**2.0)
        return mu, mu2, g1, g2 

# 需要导入模块: from scipy import special [as 别名]
# 或者: from scipy.special import gamma [as 别名]
def _rvs(self, a):
        sz, rndm = self._size, self._random_state
        u = rndm.random_sample(size=sz)
        gm = gamma.rvs(a, size=sz, random_state=rndm)
        return gm * np.where(u >= 0.5, 1, -1)