python类erfc()的实例源码

def runs_test(bits):
    n = len(bits)
    zeroes,ones = count_ones_zeroes(bits)

    prop = float(ones)/float(n)
    print "  prop ",prop

    tau = 2.0/math.sqrt(n)
    print "  tau ",tau

    if abs(prop-0.5) > tau:
        return (False,0.0,None)

    vobs = 1.0
    for i in xrange(n-1):
        if bits[i] != bits[i+1]:
            vobs += 1.0

    print "  vobs ",vobs

    p = math.erfc(abs(vobs - (2.0*n*prop*(1.0-prop)))/(2.0*math.sqrt(2.0*n)*prop*(1-prop) ))
    success = (p >= 0.01)
    return (success,p,None)

def monobit_test(bits):
    n = len(bits)

    zeroes,ones = count_ones_zeroes(bits)
    s = abs(ones-zeroes)
    print "  Ones count   = %d" % ones
    print "  Zeroes count = %d" % zeroes

    p = math.erfc(float(s)/(math.sqrt(float(n)) * math.sqrt(2.0)))

    success = (p >= 0.1)
    return (success,p,None)

def dft_test(bits):
    n = len(bits)
    if (n % 2) == 1:        # Make it an even number
        bits = bits[:-1]

    ts = list()             # Convert to +1,-1
    for bit in bits:
        ts.append((bit*2)-1)

    ts_np = numpy.array(ts)
    fs = numpy.fft.fft(ts_np)  # Compute DFT

    mags = abs(fs)[:n/2] # Compute magnitudes of first half of sequence

    T = math.sqrt(math.log(1.0/0.05)*n) # Compute upper threshold
    N0 = 0.95*n/2.0
    print "  N0 = %f" % N0

    N1 = 0.0   # Count the peaks above the upper theshold
    for mag in mags:
        if mag < T:
            N1 += 1.0
    print "  N1 = %f" % N1
    d = (N1 - N0)/math.sqrt((n*0.95*0.05)/4) # Compute the P value
    p = math.erfc(abs(d)/math.sqrt(2))

    success = (p >= 0.01)
    return (success,p,None)

def normcdf(n):
    return 0.5 * math.erfc(-n * math.sqrt(0.5))

def pvalue(x,sigma):

    return 0.5*erfc(x/(sigma*np.sqrt(2)))

def test_one():
    from math import sin, cos, tan, asin, acos, atan
    from math import sinh, cosh, tanh, asinh, acosh, atanh
    from math import exp, expm1, log, log10, log1p, sqrt, lgamma
    from math import fabs, ceil, floor, trunc, erf, erfc
    try:
        from math import log2
    except ImportError:
        def log2(x):
            return log(x) / log(2)

    def wrapper(f, v):
        try:
            return f(v)
        except ValueError:
            if f == sqrt:
                return float('nan')
            if v >= 0:
                return float('inf')
            else:
                return -float('inf')

    def compare(a, b):
        if isfinite(a) and isfinite(b):
            return assert_almost_equals(a, b)
        return str(a) == str(b)

    for f in [sin, cos, tan, asin, acos, atan,
              sinh, cosh, tanh, asinh, acosh, atanh,
              exp, expm1, log, log2, log10, log1p, sqrt,
              lgamma,
              fabs, ceil, floor, trunc,
              erf, erfc]:
        for p in [0.5, 1]:
            a = random_lst(p=p)
            b = SparseArray.fromlist(a)
            c = getattr(b, f.__name__)()
            res = [wrapper(f, x) for x in a]
            index = [k for k, v in enumerate(res) if v != 0]
            res = [x for x in res if x != 0]
            print(f, p, c.non_zero, len(res))
            assert c.non_zero == len(res)
            [assert_almost_equals(v, w) for v, w in zip(index,
                                                        c.index)]
            [compare(v, w) for v, w in zip(res,
                                           c.data)]