python类svd()的实例源码

def test_svd_no_uv(self):
        # gh-4733
        for shape in (3, 4), (4, 4), (4, 3):
            for t in float, complex:
                a = np.ones(shape, dtype=t)
                w = linalg.svd(a, compute_uv=False)
                c = np.count_nonzero(np.absolute(w) > 0.5)
                assert_equal(c, 1)
                assert_equal(np.linalg.matrix_rank(a), 1)
                assert_array_less(1, np.linalg.norm(a, ord=2))

def do(self, a, b):
        u, s, vt = linalg.svd(a, 0)
        assert_allclose(a, dot_generalized(np.asarray(u) * np.asarray(s)[..., None, :],
                                           np.asarray(vt)),
                        rtol=get_rtol(u.dtype))
        assert_(imply(isinstance(a, matrix), isinstance(u, matrix)))
        assert_(imply(isinstance(a, matrix), isinstance(vt, matrix)))

def test_types(self):
        def check(dtype):
            x = np.array([[1, 0.5], [0.5, 1]], dtype=dtype)
            u, s, vh = linalg.svd(x)
            assert_equal(u.dtype, dtype)
            assert_equal(s.dtype, get_real_dtype(dtype))
            assert_equal(vh.dtype, dtype)
            s = linalg.svd(x, compute_uv=False)
            assert_equal(s.dtype, get_real_dtype(dtype))

        for dtype in [single, double, csingle, cdouble]:
            yield check, dtype

def do(self, a, b):
        c = asarray(a)  # a might be a matrix
        s = linalg.svd(c, compute_uv=False)
        old_assert_almost_equal(
            s[..., 0] / s[..., -1], linalg.cond(a, 2), decimal=5)

def _get_orthogonal_init_weights(weights):
    fan_out = weights.size(0)
    fan_in = weights.size(1) * weights.size(2)*weights.size(3)*weights.size(4)

    u, _, v = svd(normal(0.0, 0.01, (fan_out, fan_in)), full_matrices=False)

    if u.shape == (fan_out, fan_in):
        return torch.Tensor(u.reshape(weights.size()))
    else:
        return torch.Tensor(v.reshape(weights.size()))

def PerspectiveFromPointsOld(source, dest, new_size):
    '''
    Python/Scipy implementation implementation which finds a perspective 
    transform between points.

    Most users should use PerspectiveFromPoints instead.  This method
    may be eliminated in the future.
    '''
    assert len(source) == len(dest)

    src_nrm = pv.AffineNormalizePoints(source)
    source = src_nrm.transformPoints(source)
    dst_nrm = pv.AffineNormalizePoints(dest)
    dest   = dst_nrm.transformPoints(dest)

    A = []
    for i in range(len(source)):
        src = source[i]
        dst = dest[i]

        # See Hartley and Zisserman Ch. 4.1, 4.1.1, 4.4.4
        row1 = [0.0,0.0,0.0,-dst.w*src.x,-dst.w*src.y,-dst.w*src.w,dst.y*src.x,dst.y*src.y,dst.y*src.w]
        row2 = [dst.w*src.x,dst.w*src.y,dst.w*src.w,0.0,0.0,0.0,-dst.x*src.x,-dst.x*src.y,-dst.x*src.w]
        #row3 = [-dst.y*src.x,-dst.y*src.y,-dst.y*src.w,dst.x*src.x,dst.x*src.y,dst.x*src.w,0.0,0.0,0.0]
        A.append(row1)
        A.append(row2)
        #A.append(row3)
    A = np.array(A)
    U,D,Vt = la.svd(A)
    H = Vt[8,:].reshape(3,3)

    matrix = np.dot(dst_nrm.inverse,np.dot(H,src_nrm.matrix))

    return PerspectiveTransform(matrix,new_size)

def rmsd(X, Y):
    """
    Calculate the root mean squared deviation (RMSD) using Kabsch' formula.

    @param X: (n, d) input vector
    @type X: numpy array

    @param Y: (n, d) input vector
    @type Y: numpy array

    @return: rmsd value between the input vectors
    @rtype: float
    """

    from numpy import sum, dot, sqrt, clip, average
    from numpy.linalg import svd, det

    X = X - X.mean(0)
    Y = Y - Y.mean(0)

    R_x = sum(X ** 2)
    R_y = sum(Y ** 2)

    V, L, U = svd(dot(Y.T, X))

    if det(dot(V, U)) < 0.:
        L[-1] *= -1

    return sqrt(clip(R_x + R_y - 2 * sum(L), 0., 1e300) / len(X))

def wrmsd(X, Y, w):
    """
    Calculate the weighted root mean squared deviation (wRMSD) using Kabsch'
    formula.

    @param X: (n, d) input vector
    @type X: numpy array

    @param Y: (n, d) input vector
    @type Y: numpy array

    @param w: input weights
    @type w: numpy array

    @return: rmsd value between the input vectors
    @rtype: float
    """

    from numpy import sum, dot, sqrt, clip, average
    from numpy.linalg import svd

    ## normalize weights

    w = w / w.sum()

    X = X - dot(w, X)
    Y = Y - dot(w, Y)

    R_x = sum(X.T ** 2 * w)
    R_y = sum(Y.T ** 2 * w)

    L = svd(dot(Y.T * w, X))[1]

    return sqrt(clip(R_x + R_y - 2 * sum(L), 0., 1e300))

def main():
    """
    """

    A = np.array([[1,-2],[3,5]])

    print A

    res = svd(A)
    u = res[0]
    s = res[1]
    v = res[2]

    print s

def pca(m, k):
    from numpy.linalg import svd
    from numpy.linalg import eig
    from numpy.linalg import det
    u,s,v = svd(m)
    rs = np.sqrt(np.diag(s[:k]))
    x=np.dot(u[:,:k], rs)
    y=np.dot(rs, v[:k])
    mhat=np.dot(x, y)
    return s, x, y, mhat

def test_svd_build(self, level=rlevel):
        # Ticket 627.
        a = array([[0., 1.], [1., 1.], [2., 1.], [3., 1.]])
        m, n = a.shape
        u, s, vh = linalg.svd(a)

        b = dot(transpose(u[:, n:]), a)

        assert_array_almost_equal(b, np.zeros((2, 2)))

def test_large_svd_32bit(self):
        # See gh-4442, 64bit would require very large/slow matrices.
        x = np.eye(1000, 66)
        np.linalg.svd(x)

def test_svd_no_uv(self):
        # gh-4733
        for shape in (3, 4), (4, 4), (4, 3):
            for t in float, complex:
                a = np.ones(shape, dtype=t)
                w = linalg.svd(a, compute_uv=False)
                c = np.count_nonzero(np.absolute(w) > 0.5)
                assert_equal(c, 1)
                assert_equal(np.linalg.matrix_rank(a), 1)
                assert_array_less(1, np.linalg.norm(a, ord=2))

def do(self, a, b):
        u, s, vt = linalg.svd(a, 0)
        assert_allclose(a, dot_generalized(np.asarray(u) * np.asarray(s)[..., None, :],
                                           np.asarray(vt)),
                        rtol=get_rtol(u.dtype))
        assert_(imply(isinstance(a, matrix), isinstance(u, matrix)))
        assert_(imply(isinstance(a, matrix), isinstance(vt, matrix)))

def test_types(self):
        def check(dtype):
            x = np.array([[1, 0.5], [0.5, 1]], dtype=dtype)
            u, s, vh = linalg.svd(x)
            assert_equal(u.dtype, dtype)
            assert_equal(s.dtype, get_real_dtype(dtype))
            assert_equal(vh.dtype, dtype)
            s = linalg.svd(x, compute_uv=False)
            assert_equal(s.dtype, get_real_dtype(dtype))

        for dtype in [single, double, csingle, cdouble]:
            yield check, dtype

def do(self, a, b):
        c = asarray(a)  # a might be a matrix
        s = linalg.svd(c, compute_uv=False)
        old_assert_almost_equal(
            s[..., 0] / s[..., -1], linalg.cond(a, 2), decimal=5)

def normRotation(self, tmpSkel = None, refPtIdx = None):
        '''tmpSkel: normalize every palm pose to the tmpSkel pose (template skeleton)
           refPtIdx: indexes of joints on the palm
        '''
        if tmpSkel is None:
            tmpSkel = self.frmList[0].norm_skel

        if refPtIdx is None:
            refPtIdx = self.refPtIdx
        refIdx = []
        for idx in refPtIdx:
            refIdx += [idx*3, idx*3+1, idx*3+2]

        keep_list = set(range(3*self.skel_num)).\
            difference(set(refIdx+range(self.centerPtIdx, self.centerPtIdx+3)))
        keep_list = list(keep_list)

        temp = tmpSkel[refIdx].copy()
        temp.shape = (-1,3)

        for frm in self.frmList:
           model = frm.norm_skel[refIdx] 
           model.shape = (-1,3)

           R = np.zeros((3,3), np.float32)
           for vt, vm in zip(temp, model):
               R = R + np.dot(vm.reshape(3,1), vt.reshape(1,3))

           U,s,V = svd(R, full_matrices=True) 
           R = np.dot(V.transpose(), U.transpose())
           frm.quad = Quaternion(R)
           frm.norm_skel.shape = (-1,3)
           frm.norm_skel = np.dot(R,frm.norm_skel.transpose())
           frm.norm_skel = frm.norm_skel.flatten('F')
           # frm.norm_skel = frm.norm_skel[keep_list]

def _calc(self, data, ret_obj):
        """
        Calculate SVD (wrap numpy SVD).
        """
        try:
            if self.rng is None:
                self._U, self._s, self._Vh = _svd(data, full_matrices=False)
            else:
                self._U, self._s, self._Vh = _svd(data[..., self.rng],
                                               full_matrices=False)
        except:
            return False
        else:
            return True

def truncated_svd(W, k):
    """ Given input filters, return a set of basis and the linear combination
        required to approximate the original input filters
        Input: 
            W: [dxc] matrix, where c is the input dimension, 
                d is the output dimension
        Output:
            B: [kxc] matrix, where c is the input dimension, 
                k is the maximum rank of output filters
            L: [dxk] matrix, where k is the maximum rank of the
                output filters, d is the output dimension

        Note that k <= min(c,d). It is an error if that is encountered.
    """
    d, c = W.shape
    assert k <= min(c,d), 'k={} is too large for c={}, d={}'.format(k,c,d)
    # S in this case is a vector with len=K=min(c,d), and U is [d x K], V is [K x c]
    u, s, v = svd(W, full_matrices=False)
    # compute square of s -> s_sqrt
    s_sqrt = np.sqrt(s[:k])
    # extract L from u
    B = v[:k, :] * s_sqrt[:, np.newaxis]
    # extract B from v
    L = u[:, :k] * s_sqrt

    return B, L

def test_svd_build(self, level=rlevel):
        # Ticket 627.
        a = array([[0., 1.], [1., 1.], [2., 1.], [3., 1.]])
        m, n = a.shape
        u, s, vh = linalg.svd(a)

        b = dot(transpose(u[:, n:]), a)

        assert_array_almost_equal(b, np.zeros((2, 2)))