How to use the pyriemann.utils.base.invsqrtm function in pyriemann

To help you get started, we’ve selected a few pyriemann examples, based on popular ways it is used in public projects.

Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.

github alexandrebarachant / pyRiemann / tests / test_utils_base.py View on Github external
def test_invsqrtm():
    """Test matrix inverse square root"""
    C = 2*np.eye(3)
    Ctrue = (1.0/np.sqrt(2))*np.eye(3)
    assert_array_almost_equal(invsqrtm(C), Ctrue)
github alexandrebarachant / pyRiemann / pyriemann / utils / tangentspace.py View on Github external
def transport(Covs, Cref, metric='riemann'):
    """Parallel transport of two set of covariance matrix.

    """
    C = mean_covariance(Covs, metric=metric)
    iC = invsqrtm(C)
    E = sqrtm(numpy.dot(numpy.dot(iC, Cref), iC))
    out = numpy.array([numpy.dot(numpy.dot(E, c), E.T) for c in Covs])
    return out
github alexandrebarachant / pyRiemann / pyriemann / utils / mean.py View on Github external
Nt, Ne, Ne = covmats.shape
    crit = numpy.inf
    k = 0

    # init with AJD
    B, _ = ajd_pham(covmats)
    while (crit > tol) and (k < maxiter):
        k += 1
        J = numpy.zeros((Ne, Ne))

        for index, Ci in enumerate(covmats):
            tmp = logm(numpy.dot(numpy.dot(B.T, Ci), B))
            J += sample_weight[index] * tmp

        update = numpy.diag(numpy.diag(expm(J)))
        B = numpy.dot(B, invsqrtm(update))

        crit = distance_riemann(numpy.eye(Ne), update)

    A = numpy.linalg.inv(B)

    J = numpy.zeros((Ne, Ne))
    for index, Ci in enumerate(covmats):
        tmp = logm(numpy.dot(numpy.dot(B.T, Ci), B))
        J += sample_weight[index] * tmp

    C = numpy.dot(numpy.dot(A.T, expm(J)), A)
    return C
github alexandrebarachant / pyRiemann / pyriemann / utils / mean.py View on Github external
# init
    sample_weight = _get_sample_weight(sample_weight, covmats)
    Nt, Ne, Ne = covmats.shape
    if init is None:
        C = numpy.mean(covmats, axis=0)
    else:
        C = init
    k = 0
    nu = 1.0
    tau = numpy.finfo(numpy.float64).max
    crit = numpy.finfo(numpy.float64).max
    # stop when J<10^-9 or max iteration = 50
    while (crit > tol) and (k < maxiter) and (nu > tol):
        k = k + 1
        C12 = sqrtm(C)
        Cm12 = invsqrtm(C)
        J = numpy.zeros((Ne, Ne))

        for index in range(Nt):
            tmp = numpy.dot(numpy.dot(Cm12, covmats[index, :, :]), Cm12)
            J += sample_weight[index] * logm(tmp)

        crit = numpy.linalg.norm(J, ord='fro')
        h = nu * crit
        C = numpy.dot(numpy.dot(C12, expm(nu * J)), C12)
        if h < tau:
            nu = 0.95 * nu
            tau = h
        else:
            nu = 0.5 * nu

    return C
github alexandrebarachant / pyRiemann / pyriemann / utils / tangentspace.py View on Github external
def tangent_space(covmats, Cref):
    """Project a set of covariance matrices in the tangent space. according to
    the reference point Cref

    :param covmats: np.ndarray
        Covariance matrices set, Ntrials X Nchannels X Nchannels
    :param Cref: np.ndarray
        The reference covariance matrix
    :returns: np.ndarray
        the Tangent space , a matrix of Ntrials X (Nchannels*(Nchannels+1)/2)

    """
    Nt, Ne, Ne = covmats.shape
    Cm12 = invsqrtm(Cref)
    idx = numpy.triu_indices_from(Cref)
    Nf = int(Ne * (Ne + 1) / 2)
    T = numpy.empty((Nt, Nf))
    coeffs = (numpy.sqrt(2) * numpy.triu(numpy.ones((Ne, Ne)), 1) +
              numpy.eye(Ne))[idx]
    for index in range(Nt):
        tmp = numpy.dot(numpy.dot(Cm12, covmats[index, :, :]), Cm12)
        tmp = logm(tmp)
        T[index, :] = numpy.multiply(coeffs, tmp[idx])
    return T
github alexandrebarachant / pyRiemann / pyriemann / utils / geodesic.py View on Github external
def geodesic_riemann(A, B, alpha=0.5):
    """Return the matrix at the position alpha on the riemannian geodesic between A and B  :

    .. math::
            \mathbf{C} = \mathbf{A}^{1/2} \left( \mathbf{A}^{-1/2} \mathbf{B} \mathbf{A}^{-1/2} \\right)^\\alpha \mathbf{A}^{1/2}

    C is equal to A if alpha = 0 and B if alpha = 1

    :param A: the first coavriance matrix
    :param B: the second coavriance matrix
    :param alpha: the position on the geodesic
    :returns: the covariance matrix

    """
    sA = sqrtm(A)
    isA = invsqrtm(A)
    C = numpy.dot(numpy.dot(isA, B), isA)
    D = powm(C, alpha)
    E = numpy.dot(numpy.dot(sA, D), sA)
    return E