Python sklearn.utils.extmath 模块，cartesian() 实例源码

def test_cartesian():
    # Check if cartesian product delivers the right results

    axes = (np.array([1, 2, 3]), np.array([4, 5]), np.array([6, 7]))

    true_out = np.array([[1, 4, 6],
                         [1, 4, 7],
                         [1, 5, 6],
                         [1, 5, 7],
                         [2, 4, 6],
                         [2, 4, 7],
                         [2, 5, 6],
                         [2, 5, 7],
                         [3, 4, 6],
                         [3, 4, 7],
                         [3, 5, 6],
                         [3, 5, 7]])

    out = cartesian(axes)
    assert_array_equal(true_out, out)

    # check single axis
    x = np.arange(3)
    assert_array_equal(x[:, np.newaxis], cartesian((x,)))

def test_cgauss_likelihood():
    mu = np.array([0], dtype='float')
    sigma = np.array([2], dtype='float')
    x = np.linspace(-1, 2, 2)
    lapse = np.array([0], dtype='float')
    parameters = cartesian((mu, sigma, lapse, x))
    proportionMethod = PsiMarginal.pf(parameters, psyfun='cGauss')
    samples = np.random.normal(mu, sigma, (200000, 1))
    proportionSamples = np.empty([2, ])
    proportionSamples[0] = np.mean(samples <= x[0])  # cdf is p(X<=x), compute this through sampling to check likelihood
    proportionSamples[1] = np.mean(samples <= x[1])
    np.testing.assert_almost_equal(proportionSamples, proportionMethod, decimal=2) == 1

def compute_overlap(mat1, mat2):
    s1 = mat1.shape[0]
    s2 = mat2.shape[0]
    area1 = (mat1[:, 2] - mat1[:, 0]) * (mat1[:, 3] - mat1[:, 1])
    if mat2.shape[1] == 5:
        area2 = mat2[:, 4]
    else:
        area2 = (mat2[:, 2] - mat2[:, 0]) * (mat2[:, 3] - mat2[:, 1])

    x1 = cartesian([mat1[:, 0], mat2[:, 0]])

    x1 = np.amax(x1, axis=1)
    x2 = cartesian([mat1[:, 2], mat2[:, 2]])
    x2 = np.amin(x2, axis=1)
    com_zero = np.zeros(x2.shape[0])
    w = x2 - x1
    w = w - 1

    w = np.maximum(com_zero, w)

    y1 = cartesian([mat1[:, 1], mat2[:, 1]])
    y1 = np.amax(y1, axis=1)
    y2 = cartesian([mat1[:, 3], mat2[:, 3]])
    y2 = np.amin(y2, axis=1)
    h = y2 - y1
    h = h - 1
    h = np.maximum(com_zero, h)

    oo = w * h

    aa = cartesian([area1[:], area2[:]])
    aa = np.sum(aa, axis=1)

    ooo = oo / (aa - oo)

    overlap = np.transpose(ooo.reshape(s1, s2), (1, 0))

    return overlap

def weights(dim, degree):
        # 1D sigma-points (x) and weights (w)
        x, w = hermegauss(degree)
        # hermegauss() provides weights that cause posdef errors
        w = factorial(degree) / (degree ** 2 * hermeval(x, [0] * (degree - 1) + [1]) ** 2)
        return np.prod(cartesian([w] * dim), axis=1)

def unit_sigma_points(dim, degree):
        # 1D sigma-points (x) and weights (w)
        x, w = hermegauss(degree)
        # nD sigma-points by cartesian product
        return cartesian([x] * dim).T  # column/sigma-point