Python linalg.inv方法代码示例

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def test_scipy_gmres_linop_parameter(self):
    """ This is a test on gmres method with a parameter-dependent linear operator """
    for omega in np.linspace(-10.0, 10.0, 10):
      for eps in np.linspace(-10.0, 10.0, 10):
        
        linop_param = linalg.aslinearoperator(vext2veff_c(omega, eps, n))
        
        Aparam = np.zeros((n,n), np.complex64)
        for i in range(n):
          uv = np.zeros(n, np.complex64); uv[i] = 1.0
          Aparam[:,i] = linop_param.matvec(uv)
        x_ref = np.dot(inv(Aparam), b)
    
        x_itr,info = linalg.lgmres(linop_param, b)
        derr = abs(x_ref-x_itr).sum()/x_ref.size
        self.assertLess(derr, 1e-6) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def test_predict(self):
        X_test = np.random.rand(10, 2)

        m, v = self.model.predict(X_test)

        assert len(m.shape) == 1
        assert m.shape[0] == X_test.shape[0]
        assert len(v.shape) == 1
        assert v.shape[0] == X_test.shape[0]

        m, v = self.model.predict(X_test, full_cov=True)

        assert len(m.shape) == 1
        assert m.shape[0] == X_test.shape[0]
        assert len(v.shape) == 2
        assert v.shape[0] == X_test.shape[0]
        assert v.shape[1] == X_test.shape[0]

        K_zz = self.kernel.get_value(X_test)
        K_zx = self.kernel.get_value(X_test, self.X)
        K_nz = self.kernel.get_value(self.X) + self.model.noise * np.eye(self.X.shape[0])
        inv = spla.inv(K_nz)
        K_zz_x = K_zz - np.dot(K_zx, np.inner(inv, K_zx))
        assert np.mean((K_zz_x - v) ** 2) < 10e-5 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def _b_orthonormalize(B, blockVectorV, blockVectorBV=None, retInvR=False):
    if blockVectorBV is None:
        if B is not None:
            blockVectorBV = B(blockVectorV)
        else:
            blockVectorBV = blockVectorV  # Shared data!!!
    gramVBV = np.dot(blockVectorV.T, blockVectorBV)
    gramVBV = cholesky(gramVBV)
    gramVBV = inv(gramVBV, overwrite_a=True)
    # gramVBV is now R^{-1}.
    blockVectorV = np.dot(blockVectorV, gramVBV)
    if B is not None:
        blockVectorBV = np.dot(blockVectorBV, gramVBV)

    if retInvR:
        return blockVectorV, blockVectorBV, gramVBV
    else:
        return blockVectorV, blockVectorBV 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def H(self):
        k = self.neqs
        Lk = tsa.elimination_matrix(k)
        Kkk = tsa.commutation_matrix(k, k)
        Ik = np.eye(k)

        # B = chain_dot(Lk, np.eye(k**2) + commutation_matrix(k, k),
        #               np.kron(self.P, np.eye(k)), Lk.T)

        # return np.dot(Lk.T, L.inv(B))

        B = chain_dot(Lk,
                      np.dot(np.kron(Ik, self.P), Kkk) + np.kron(self.P, Ik),
                      Lk.T)

        return np.dot(Lk.T, L.inv(B)) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def mvn_loglike(x, sigma):
    '''loglike multivariate normal

    assumes x is 1d, (nobs,) and sigma is 2d (nobs, nobs)

    brute force from formula
    no checking of correct inputs
    use of inv and log-det should be replace with something more efficient
    '''
    #see numpy thread
    #Sturla: sqmahal = (cx*cho_solve(cho_factor(S),cx.T).T).sum(axis=1)
    sigmainv = linalg.inv(sigma)
    logdetsigma = np.log(np.linalg.det(sigma))
    nobs = len(x)

    llf = - np.dot(x, np.dot(sigmainv, x))
    llf -= nobs * np.log(2 * np.pi)
    llf -= logdetsigma
    llf *= 0.5
    return llf 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def kalman_filter(self, z): 
        '''
        Implement the Kalman Filter, including the predict and the update stages,
        with the measurement z
        '''
        x = self.x_state
        # Predict
        x = dot(self.F, x)
        self.P = dot(self.F, self.P).dot(self.F.T) + self.Q

        #Update
        S = dot(self.H, self.P).dot(self.H.T) + self.R
        K = dot(self.P, self.H.T).dot(inv(S)) # Kalman gain
        y = z - dot(self.H, x) # residual
        x += dot(K, y)
        self.P = self.P - dot(K, self.H).dot(self.P)
        self.x_state = x.astype(int) # convert to integer coordinates 
                                     #(pixel values) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def getElemBezierData( elemCoords , C , order=4 , method="Gauss" , elemType = 'default' ):

  elemData = elemShapeData()
    
  if elemType == 'default':
    elemType = getElemType( elemCoords )
    
  (intCrds,intWghts) = getIntegrationPoints( "Line3" , order , method )
    
  for xi,intWeight in zip( real(intCrds) , intWghts ):    
    try:
      sData = eval( 'getBezier'+elemType+'(xi,C)' )
    except:
      raise NotImplementedError('Unknown type :'+elemType)

    jac = dot ( sData.dhdxi.transpose() , elemCoords )

    if jac.shape[0] is jac.shape[1]:
      sData.dhdx   = (dot ( inv( jac ) , sData.dhdxi.transpose() )).transpose()
   
    sData.weight = calcWeight( jac ) * intWeight

    elemData.sData.append(sData)

  return elemData 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def _mri_landmarks_to_mri_voxels(mri_landmarks, t1_mgh):
    """Convert landmarks from MRI RAS space to MRI voxel space.

    Parameters
    ----------
    mri_landmarks : array, shape (3, 3)
        The MRI RAS landmark data: rows LPA, NAS, RPA, columns x, y, z.
    t1_mgh : nib.MGHImage
        The image data in MGH format.

    Returns
    -------
    mri_landmarks : array, shape (3, 3)
        The MRI voxel-space landmark data.
    """
    # Get landmarks in voxel space, using the T1 data
    vox2ras_tkr = t1_mgh.header.get_vox2ras_tkr()
    ras2vox_tkr = linalg.inv(vox2ras_tkr)
    mri_landmarks = apply_trans(ras2vox_tkr, mri_landmarks)  # in vox
    return mri_landmarks 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def test_twodiags(self):
        A = spdiags([[1, 2, 3, 4, 5], [6, 5, 8, 9, 10]], [0, 1], 5, 5)
        b = array([1, 2, 3, 4, 5])

        # condition number of A
        cond_A = norm(A.todense(),2) * norm(inv(A.todense()),2)

        for t in ['f','d','F','D']:
            eps = finfo(t).eps  # floating point epsilon
            b = b.astype(t)

            for format in ['csc','csr']:
                Asp = A.astype(t).asformat(format)

                x = spsolve(Asp,b)

                assert_(norm(b - Asp*x) < 10 * cond_A * eps) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def b_orthonormalize(B, blockVectorV,
                      blockVectorBV=None, retInvR=False):
    """Internal."""
    import scipy.linalg as sla
    if blockVectorBV is None:
        if B is not None:
            blockVectorBV = B(blockVectorV)
        else:
            blockVectorBV = blockVectorV  # Shared data!!!
    gramVBV = sp.dot(blockVectorV.T, blockVectorBV)
    gramVBV = sla.cholesky(gramVBV)
    gramVBV = sla.inv(gramVBV, overwrite_a=True)
    # gramVBV is now R^{-1}.
    blockVectorV = sp.dot(blockVectorV, gramVBV)
    if B is not None:
        blockVectorBV = sp.dot(blockVectorBV, gramVBV)

    if retInvR:
        return blockVectorV, blockVectorBV, gramVBV
    else:
        return blockVectorV, blockVectorBV 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def test_graphical_lasso_cv(random_state=1):
    # Sample data from a sparse multivariate normal
    dim = 5
    n_samples = 6
    random_state = check_random_state(random_state)
    prec = make_sparse_spd_matrix(dim, alpha=.96,
                                  random_state=random_state)
    cov = linalg.inv(prec)
    X = random_state.multivariate_normal(np.zeros(dim), cov, size=n_samples)
    # Capture stdout, to smoke test the verbose mode
    orig_stdout = sys.stdout
    try:
        sys.stdout = StringIO()
        # We need verbose very high so that Parallel prints on stdout
        GraphicalLassoCV(verbose=100, alphas=5, tol=1e-1).fit(X)
    finally:
        sys.stdout = orig_stdout

    # Smoke test with specified alphas
    GraphicalLassoCV(alphas=[0.8, 0.5], tol=1e-1, n_jobs=1).fit(X) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def test_graph_lasso_cv(random_state=1):
    # Sample data from a sparse multivariate normal
    dim = 5
    n_samples = 6
    random_state = check_random_state(random_state)
    prec = make_sparse_spd_matrix(dim, alpha=.96,
                                  random_state=random_state)
    cov = linalg.inv(prec)
    X = random_state.multivariate_normal(np.zeros(dim), cov, size=n_samples)
    # Capture stdout, to smoke test the verbose mode
    orig_stdout = sys.stdout
    try:
        sys.stdout = StringIO()
        # We need verbose very high so that Parallel prints on stdout
        GraphLassoCV(verbose=100, alphas=5, tol=1e-1).fit(X)
    finally:
        sys.stdout = orig_stdout

    # Smoke test with specified alphas
    GraphLassoCV(alphas=[0.8, 0.5], tol=1e-1, n_jobs=1).fit(X) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def test_property():
    rng = np.random.RandomState(0)
    rand_data = RandomData(rng, scale=7)
    n_components = rand_data.n_components

    for covar_type in COVARIANCE_TYPE:
        X = rand_data.X[covar_type]
        gmm = GaussianMixture(n_components=n_components,
                              covariance_type=covar_type, random_state=rng,
                              n_init=5)
        gmm.fit(X)
        if covar_type == 'full':
            for prec, covar in zip(gmm.precisions_, gmm.covariances_):

                assert_array_almost_equal(linalg.inv(prec), covar)
        elif covar_type == 'tied':
            assert_array_almost_equal(linalg.inv(gmm.precisions_),
                                      gmm.covariances_)
        else:
            assert_array_almost_equal(gmm.precisions_, 1. / gmm.covariances_) 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def predict(self, a_hist, t):
        """
        This function implements the prediction formula discussed is section 6 (1.59)
        It takes a realization for a^N, and the period in which the prediciton is formed

        Output:  E[abar | a_t, a_{t-1}, ..., a_1, a_0]
        """

        N = np.asarray(a_hist).shape[0] - 1
        a_hist = np.asarray(a_hist).reshape(N + 1, 1)
        V = self.construct_V(N + 1)

        aux_matrix = np.zeros((N + 1, N + 1))
        aux_matrix[:(t + 1), :(t + 1)] = np.eye(t + 1)
        L = la.cholesky(V).T
        Ea_hist = la.inv(L) @ aux_matrix @ L @ a_hist

        return Ea_hist 

# 需要导入模块: from scipy import linalg [as 别名]
# 或者: from scipy.linalg import inv [as 别名]
def plot4():
    # Density 1
    Z = gen_gaussian_plot_vals(x_hat, Σ)
    cs1 = ax.contour(X, Y, Z, 6, colors="black")
    ax.clabel(cs1, inline=1, fontsize=10)
    # Density 2
    M = Σ * G.T * linalg.inv(G * Σ * G.T + R)
    x_hat_F = x_hat + M * (y - G * x_hat)
    Σ_F = Σ - M * G * Σ
    Z_F = gen_gaussian_plot_vals(x_hat_F, Σ_F)
    cs2 = ax.contour(X, Y, Z_F, 6, colors="black")
    ax.clabel(cs2, inline=1, fontsize=10)
    # Density 3
    new_x_hat = A * x_hat_F
    new_Σ = A * Σ_F * A.T + Q
    new_Z = gen_gaussian_plot_vals(new_x_hat, new_Σ)
    cs3 = ax.contour(X, Y, new_Z, 6, colors="black")
    ax.clabel(cs3, inline=1, fontsize=10)
    ax.contourf(X, Y, new_Z, 6, alpha=0.6, cmap=cm.jet)
    ax.text(float(y[0]), float(y[1]), r"$y$", fontsize=20, color="black")

# == Choose a plot to generate == #