# Gram-Schmidt process

The Gram-Schmidt procedure takes a list of (column) vectors and forms an orthonormal basis from this set.
As a corollary, the procedure allows us to determine the dimension of the space spanned by the basis vectors, which is equal to or less than the space which the vectors sit.


```python
# !pip install sympy
```


```python
import numpy as np
import numpy.linalg as la

verySmallNumber = 1e-14 # That's 1×10⁻¹⁴ = 0.00000000000001
```


```python
# sympy converts numpy as laTex
# from IPython.display import display, Math, Latex
import sympy as sp
from sympy import Matrix as spm
sp.init_printing(use_unicode=True)
```

### Gram-Schmidt procedure
Take 4 basis vectors as a list of vectors as the columns of a matrix, A.
Go through the vectors one at a time and set them to be orthogonal to all the vectors that came before it and normalise it.


```python
def gsBasis4(A):
    """Gram-Schmidt for 4 vectors"""
    
    # Work with a copy - vectors are mutable
    B = np.array(A, dtype=np.float_) 
    
    # Column(vector) 0: Normalise(divide by its modulus or norm)
    B[:, 0] = B[:, 0] / la.norm(B[:, 0])
    
    # Column 1: - subtract any overlap with the zeroth vector
    B[:, 1] = B[:, 1] - B[:, 1] @ B[:, 0] * B[:, 0]
    
    # If there's anything left after that subtraction, then B[:, 1] is linearly independant of B[:, 0]
    # Normalise - norm(indepent)=1, norm(dependent)=0
    if la.norm(B[:, 1]) > verySmallNumber :
        B[:, 1] = B[:, 1] / la.norm(B[:, 1])
    else :
        B[:, 1] = np.zeros_like(B[:, 1])
        
    # Column 2: - subtract the overlap with the zeroth vector
    #           - subtract the overlap with the first
    B[:, 2] = B[:, 2] - B[:, 2] @ B[:, 0] * B[:, 0] - B[:, 2] @ B[:, 1] * B[:, 1]
    
    # Normalise - norm(indepent)=1, norm(dependent)=0
    if la.norm(B[:, 2]) > verySmallNumber :
        B[:, 2] = B[:, 2] / la.norm(B[:, 2])
    else :
        B[:, 2] = np.zeros_like(B[:, 2])
    
    # Column 2: - subtract the overlap with the zeroth vector
    #           - subtract the overlap with the first
    #           - subtract the overlap with the second
    B[:, 3] = B[:, 3] - (B[:, 3] @ B[:, 0] * B[:, 0]) - (B[:, 3] @ B[:, 1] * B[:, 1]) - (B[:, 3] @ B[:, 2] * B[:, 2])
    
    # Normalise - norm(indepent)=1, norm(dependent)=0
    if la.norm(B[:, 3]) > verySmallNumber :
        B[:, 3] = B[:, 3] / la.norm(B[:, 3])
    else :
        B[:, 3] = np.zeros_like(B[:, 3])
    
    return B
```


```python
def gsBasis(A):
    """Gram-Schmidt for n vectors"""
    
    # Work with a copy - vectors are mutable
    B = np.array(A, dtype=np.float_)
    
    # Loop over all vectors
    for i in range(B.shape[1]):
        # Loop over all previous vectors, j, to subtract.
        for j in range(i):
            # Subtract the overlap with previous vectors
            # you'll need the current vector B[:, i] and a previous vector B[:, j]
            B[:, i] -= (B[:, i] @ B[:, j] * B[:, j])
            
        # Normalise - norm(indepent)=1, norm(dependent)=0
        if la.norm(B[:, i]) > verySmallNumber :
            B[:, i] = B[:, i] / la.norm(B[:, i])
        else :
            B[:, i] = np.zeros_like(B[:, i])

    return B

def dimensions(A):
    """Gram-schmidt process to calculate the dimension spanned by a list of vectors.
    Independent vectors are normalised to one and dependent vectors to zero
    Thus the sum of all the norms will be the dimension"""
    return np.sum(la.norm(gsBasis(A), axis=0))
```

## Test your code before submission
To test the code you've written above, run the cell (select the cell above, then press the play button [ ▶| ] or press shift-enter).
You can then use the code below to test out your function.
You don't need to submit this cell; you can edit and run it as much as you like.

Try out your code on tricky test cases!

### 4 vector function


```python
V = np.array([[1,0,2,6],
              [0,1,8,2],
              [2,8,3,1],
              [1,-6,2,3]], dtype=np.float_)
spm(gsBasis4(V).round(2))
```


```python
# Once you've done Gram-Schmidt once, doing it again should give you the same result.
U = gsBasis4(V)
spm(U.round(2)), spm(gsBasis4(U).round(2))
np.testing.assert_almost_equal(gsBasis4(U), gsBasis4(V))
```

### Generic function


```python
spm(gsBasis(V).round(2))
```

### Non-square matrices


```python
A = np.array([[3,2,3],
              [2,5,-1],
              [2,4,8],
              [12,2,1]], dtype=np.float_)
spm(gsBasis(A).round(2)), dimensions(A)
```

### Dependent vectors - a linear combination of the others


```python
B = np.array([[6,2,1,7,5],
              [2,8,5,-4,1],
              [1,-6,3,2,8]], dtype=np.float_)
spm(gsBasis(B).round(2)), dimensions(B)
```


```python
C = np.array([[1,0,2],
              [0,1,-3],
              [1,0,2]], dtype=np.float_)
spm(gsBasis(C).round(2)), dimensions(C)
```


```python

```

45f6f8fb4ba9289c8eed7a12617e742b3022a0e2
ipynb
Jupyter Notebook
_math/Imperial College London - Math for ML/GramSchmidtProcess.ipynb
aixpact/data-science
f04a54595fbc2d797918d450b979fd4c2eabac15
MIT
2020-07-22T23:12:39.000Z
2020-07-25T02:30:48.000Z
_math/Imperial College London - Math for ML/GramSchmidtProcess.ipynb
aixpact/data-science
f04a54595fbc2d797918d450b979fd4c2eabac15
MIT
_math/Imperial College London - Math for ML/GramSchmidtProcess.ipynb
aixpact/data-science
f04a54595fbc2d797918d450b979fd4c2eabac15
MIT
Qwen/Qwen-72B
1. YES
2. YES