QR factorisation using modified Gram Schmidt

Question

The question: For this problem, you are given a list of matrices called As, and your job is to find the QR factorization for each of them.

Implement qr_by_gram_schmidt: This function takes as input a matrix A and computes a QR decomposition, returning two variables, Q and R where A=QR, with Q orthogonal and R zero below the diagonal.

A is an n×m matrix with n≥m (i.e. more rows than columns).

You should implement this function using the modified Gram-Schmidt procedure.

INPUT:

As: List of arrays

OUTPUT:

Qs: List of the Q matrices output by qr_by_gram_schmidt, in the same order as As. For a matrix A of shape n×m, Q should have shape n×m.
Rs: List of the R matrices output by qr_by_gram_schmidt, in the same order as As. For a matrix A of shape n×m, R should have shape m×m

I have written the code for the QR factorization which I believe is correct:

import numpy as np
def qr_by_gram_schmidt(A):
m = np.shape(A)[0]
n = np.shape(A)[1]
Q =  np.zeros((m, m))
R =  np.zeros((n, n)) 
for j in xrange(n):
    v = A[:,j]
    for i in xrange(j):
        R[i,j] = Q[:,i].T * A[:,j]
        v = v.squeeze() - (R[i,j] * Q[:,i])
    R[j,j] =  np.linalg.norm(v)
    Q[:,j] = (v / R[j,j]).squeeze()
return Q, R

How do I write the loop to calculate the the QR factorization of each of the matrices in As and storing them in that order?

edit: The code has some error too. I will appreciate it if you can help me in debugging it.

Thanks

score 0 · Answer 1 · answered Apr 15 '16 at 02:44

I didn't check your GS code, but had to make a change (may not be correct!) to make it compile. You just have to set up a list of your matrices, I made 2 of them and then loop through that list and apply your function.

import numpy as np

def gs(A):
    m = np.shape(A)[0]
    n = np.shape(A)[1]
    Q =  np.zeros((m, m))
    R =  np.zeros((n, n)) 
    print m,n,Q,R
    for j in xrange(n):
        v = A[:,j]
        for i in xrange(j):
            R[i,j] =  np.dot(Q[:,i].T , A[:,j])   # I made an arbitrary change here!!!
            v = v.squeeze() - (R[i,j] * Q[:,i])
        R[j,j] =  np.linalg.norm(v)
        Q[:,j] = (v / R[j,j]).squeeze()
    return Q, R

As= np.random.rand(2,3,3)  # list of 2 (3x3) matrices
print As

for A in As:
    print gs(A)

Output:

[[[ 0.9599614   0.02213113  0.43343881]
  [ 0.44202415  0.6816688   0.88321052]
  [ 0.93098107  0.80528361  0.88473308]]

 [[ 0.41794678  0.10762796  0.42110659]
  [ 0.89598082  0.81225543  0.52947205]
  [ 0.0621515   0.59826789  0.14021332]]]
(array([[ 0.68158915, -0.67980134,  0.27075149],
       [ 0.31384477,  0.60583989,  0.73106736],
       [ 0.66101262,  0.41331364, -0.626286  ]]), array([[ 1.40841649,  0.76132516,  1.15743793],
       [ 0.        ,  0.73077208,  0.60610414],
       [ 0.        ,  0.        ,  0.20894464]]))
(array([[ 0.42190511, -0.39510208,  0.81602109],
       [ 0.90446656,  0.121136  , -0.40898205],
       [ 0.06274013,  0.91061541,  0.40846452]]), array([[ 0.99061796,  0.81760207,  0.66535379],
       [ 0.        ,  0.6006613 ,  0.02543844],
       [ 0.        ,  0.        ,  0.18435946]]))

Thank you for your help. But this code does not give the correct output and I cant understand why. The shape of Q and R at each iteration is not correct — amateur_programmer, Apr 15 '16 at 06:35

QR factorisation using modified Gram Schmidt

1 Answers1