了解 Gram-Schmidt Process 是如何翻译成这段代码作为实现的

Question

试图从这个解释中理解 Gram-Schmidt 过程：

http://mlwiki.org/index.php/Gram-Schmidt_Process

计算的步骤对我来说很有意义。但是，同一篇文章中包含的 Python 实现似乎并不一致。

def normalize(v):
    return v / np.sqrt(v.dot(v))

n = len(A)

A[:, 0] = normalize(A[:, 0])

for i in range(1, n):
    Ai = A[:, i]
    for j in range(0, i):
        Aj = A[:, j]
        t = Ai.dot(Aj)
        Ai = Ai - t * Aj
    A[:, i] = normalize(Ai)

从上面的代码中，我们看到它对 V1 和 b 做了 点积，但是 (V1,V1) 部分没有作为分母（参考下面的等式）。我想知道下面的等式是如何转换为驻留在 for 循环中的代码的？

Answer 1

这正是代码的作用

基本上，它将前一个向量（A 中的列）归一化并将当前向量投影到它并减去当前向量。

为了简洁的计算，每个向量都会进行归一化。

上面的 V2 方程没有对前一个向量进行归一化，因此存在差异。

Answer 2

试试这个矢量化实现。

我还建议阅读 David C 的理论书。

def replace_zero(array): 
    
    for i in range(len(array)) :
        if array[i] == 0 : 
            array[i] = 1
    return array

def gram_schmidt(self,A, norm=True, row_vect=False):
        """Orthonormalizes vectors by gram-schmidt process

        Parameters
        -----------
        A : ndarray,
        Matrix having vectors in its columns

        norm : bool,
        Do you need Normalized vectors?

        row_vect: bool,
        Does Matrix A has vectors in its rows?

        Returns
        -------
        G : ndarray,
        Matrix of orthogonal vectors
        
        Gram-Schmidt Process
        --------------------
        The Gram–Schmidt process is a simple algorithm for 
        producing an orthogonal or orthonormal basis for any 
        nonzero subspace of Rn.

        Given a basis {x1,....,xp} for a nonzero subspace W of Rn,
        define

        v1 = x1
        v2 = x2 - (x2.v1/v1.v1) * v1
        v3 = x3 - (x3.v1/v1.v1) * v1  - (x3.v2/v2.v2) * v2
         .
         .
         .
        vp = xp - (xp.v1/v1.v1) * v1  - (xp.v2/v2.v2) * v2 - ....... 
             .... - (xp.v(p-1) / v(p-1).v(p-1) ) * v(p-1)


        Then {v1,.....,vp} is an orthogonal basis for W .
        In addition,
        Span {v1,.....,vp} = Span {x1,.....,xp}    for 1 <= k <= p

        References
        ----------
        Linear Algebra and Its Applications - By David.C.Lay

        """
        if row_vect :
            # if true, transpose it to make column vector matrix
            A = A.T

        no_of_vectors = A.shape[1]
        G = A[:,0:1].copy() # copy the first vector in matrix
        # 0:1 is done to to be consistent with dimensions - [[1,2,3]]

        # iterate from 2nd vector to number of vectors
        for i in range(1,no_of_vectors):

            # calculates weights(coefficents) for every vector in G
            numerator = A[:,i].dot(G)
            denominator = np.diag(np.dot(G.T,G)) #to get elements in diagonal
            weights = np.squeeze(numerator/denominator)

            # projected vector onto subspace G
            projected_vector = np.sum(weights * G,
                                      axis=1,
                                      keepdims=True)

            # orthogonal vector to subspace G
            orthogonalized_vector = A[:,i:i+1] - projected_vector

            # now add the orthogonal vector to our set
            G = np.hstack((G,orthogonalized_vector))

        if norm :
            # to get orthoNormal vectors (unit orthogonal vectors)
            # replace zero to 1 to deal with division by 0 if matrix has 0 vector
            # or normazalization value comes out to be zero
            G = G/self.replace_zero(np.linalg.norm(G,axis=0))

        if row_vect:
            return G.T

        return G
        

G = np.array([[1,0,0],[1,1,0],[1,1,1],[1,1,1]])
gram_schmidt(G)
>
array([[ 0.5       , -0.8660254 ,  0.        ],
       [ 0.5       ,  0.28867513, -0.81649658],
       [ 0.5       ,  0.28867513,  0.40824829],
       [ 0.5       ,  0.28867513,  0.40824829]])