+ This notebook is part of lecture 17 *Orthogonal matrices and Gram-Scmidt* in the OCW MIT course 18.06 by Prof Gilbert Strang [1]
+ Created by me, Dr Juan H Klopper
    + Head of Acute Care Surgery
    + Groote Schuur Hospital
    + University Cape Town
    + <a href="mailto:juan.klopper@uct.ac.za">Email me with your thoughts, comments, suggestions and corrections</a> 
<a rel="license" href="http://creativecommons.org/licenses/by-nc/4.0/"></a><br /><span xmlns:dct="http://purl.org/dc/terms/" href="http://purl.org/dc/dcmitype/InteractiveResource" property="dct:title" rel="dct:type">Linear Algebra OCW MIT18.06</span> <span xmlns:cc="http://creativecommons.org/ns#" property="cc:attributionName">IPython notebook [2] study notes by Dr Juan H Klopper</span> is licensed under a <a rel="license" href="http://creativecommons.org/licenses/by-nc/4.0/">Creative Commons Attribution-NonCommercial 4.0 International License</a>.

+ [1] <a href="http://ocw.mit.edu/courses/mathematics/18-06sc-linear-algebra-fall-2011/index.htm">OCW MIT 18.06</a>
+ [2] Fernando Pérez, Brian E. Granger, IPython: A System for Interactive Scientific Computing, Computing in Science and Engineering, vol. 9, no. 3, pp. 21-29, May/June 2007, doi:10.1109/MCSE.2007.53. URL: http://ipython.org


```python
from IPython.core.display import HTML, Image
css_file = 'style.css'
HTML(open(css_file, 'r').read())
```


```python
from sympy import init_printing, symbols, Matrix, sin, cos, sqrt, Rational, GramSchmidt
from warnings import filterwarnings
```


```python
init_printing(use_latex = 'mathjax')
filterwarnings('ignore')
```


```python
theta = symbols('theta')
```

# Orthogonal basis
# Orthogonal matrix
# Gram-Schmidt

## Orthogonal basis

* Here we mean vectors *q*<sub>1</sub>,*q*<sub>2</sub>,...,*q*<sub>n</sub>

* We actually mean orthonormal vectors (for orthogonal or perpendicular and of unit length / normalized)
* Vectors that are orthogonal have a dot product equal to zero
    * If they are orthogonal
$$ {q}_{i}^{T}{q}_{j}={0} $$
    * If they are not
$$ {q}_{i}^{T}{q}_{j}\neq{0} $$

## Orthogonal matrix

* We can now put these (column) basis vectors into a matrix Q

* This brings about
$$ {Q}^{T}{Q}={I} $$

* In the case of the matrix Q being square the word *orthogonal matrix* is used
* When it is square we can calculate the inverse making
$$ {Q}^{T}={Q}^{-1} $$

* Consider the following permutation matrix with orthonormal column vectors


```python
Q = Matrix([[0, 0, 1], [1, 0, 0], [0, 1, 0]])
Q, Q.transpose()
```

* In this example the transpose also contains orthonormal column vectors
* Multiplication gives the identity matrix


```python
Q.transpose() * Q
```

* Consider this example


```python
Q = Matrix([[cos(theta), -sin(theta)], [sin(theta), cos(theta)]])
Q, Q.transpose()
```

* The two column vectors are orthogonal and the length of each column vector is 1
* It is thus an orthogonal matrix


```python
Q.transpose() * Q
```

* The example below certainly has orthogonal column vectors, but they are not of unit length
$$ {Q}=\begin{bmatrix} 1 & 1 \\ 1 & {-1} \end{bmatrix} $$

* Well, we can change them into unit vectors by dividing each component by the length of that vector
$$ \sqrt { { \left( 1 \right)  }^{ 2 }+{ \left( 1 \right)  }^{ 2 } } =\sqrt { 2 } \\ \sqrt { { \left( 1 \right)  }^{ 2 }+{ \left( -1 \right)  }^{ 2 } } =\sqrt { 2 }  $$
$$ Q=\frac { 1 }{ \sqrt { 2 }  } \begin{bmatrix} 1 & 1 \\ 1 & -1 \end{bmatrix} $$

* As it stands Q<sup>T</sup>Q is not the identity matrix


```python
Q = Matrix([[1, 1], [1, -1]])
Q.transpose() * Q
```

* Turning it into an orthogonal matrix


```python
Q = (1 / sqrt(2)) * Matrix([[1, 1], [1, -1]])
Q 
```


```python
Q.transpose() * Q
```

* Consider this example with orthogonal (but not orthonormal) column vectors


```python
Q = Matrix([[1, 1, 1, 1], [1, -1, 1, -1], [1, 1, -1, -1], [1, -1, -1, 1]])
Q
```

* Again, as it stands Q<sup>T</sup>Q is not the identity matrix


```python
Q.transpose() * Q
```

* But turning it into an orthogonal matrix works


```python
Q = Rational(1, 2) * Matrix([[1, 1, 1, 1], [1, -1, 1, -1], [1, 1, -1, -1], [1, -1, -1, 1]])
# Rational() creates a mathematical fraction instead of a decimal
Q
```


```python
Q.transpose() * Q
```

* Consider this matrix Q with orthogonal column vectors, but that is not square


```python
Q = Rational(1, 3) * Matrix([[1, -2], [2, -1], [2, 2]])
Q
```

* We now have a matrix with two column vectors that are normalized and orthogonal to each other and they form a basis for a plane (subspace) in &#8477;<sup>3</sup>

* There must be a third column matrix of unit length, orthogonal to the other two so we end up with an orthogonal matrix 


```python
Q = Rational(1, 3) * Matrix([[1, -2, 2], [2, -1, -2], [2, 2, 1]])
Q
```


```python
Q.transpose() * Q
```

* Let's make use of these matrices with orthonormal columns (which we will always denote with a letter Q) and project them onto their columnspace
* What would the projection matrix be?
$$ Q\underline{x}={b} \\ {P}={Q}{\left({Q}^{T}{Q}\right)}^{-1}{Q}^{T} $$

* Remember, though that for matrices with orthonormal column vectors we have Q<sup>T</sup>Q is the identity matrix and we have
$$ {P}={Q}{Q}^{T} $$

* If additionally, Q is square, then we have independent columns and the columnspace contain the whole space &#8477;<sup>n</sup> and the projection matrix is the identity matrix in *n*
    * Remember Q<sup>T</sup> = Q<sup>-1</sup> in these cases making it easy to see that we get the identity matrix
    * Remember also that the projection matrix is symmetric
    * Lastly the projection matrix has the property of squaring it leaves us in the same spot, so here we will have (QQ<sup>T</sup>)<sup>2</sup>=QQ<sup>T</sup>

* All of this has the final consequence that
$$ {Q}^{T}Q\hat{x}={Q}^{T}\underline{b} \\ \hat{x}={Q}^{T}\underline{b} \\ \hat{x}_{i}={q}_{i}^{T}{b} $$

## Gram-Schmidt

* All of the above makes things quite easy, so we should try and create orthogonal matrices

* Good, let's start with two independent vectors **a** and **b** and try and create two orthogonal vectors **A** and **B** and then create two orthonormal vectors
$$ { q }_{ 1 }=\frac { A }{ \left\| A \right\|  } \\ { q }_{ 2 }=\frac { B }{ \left\| B \right\|  }  $$

* We can choose one of them as our initial vector, say **a** = **A**, so we have to get an orthogonal projection (to **a**) for **B**
* This is what we previously called the error vector **e**
$$ \underline{e}=\underline{b}-\underline{p} $$
* Remembering how to get **p** we have the following
$$ B=\underline { b } -\frac { { A }^{ T }\underline { b }  }{ { A }^{ T }A } A $$

* Let's do an example


```python
a = Matrix([1, 1, 1])
b = Matrix([1, 0, 2])
a, b
```


```python
A = a
A
```


```python
A.transpose() * b
```


```python
A.transpose() * A
```


```python
B = b - A
B
```

* Checking that they are perpendicular


```python
A.transpose() * B
```

* Now we have to create Q by turning **A** and **B** into unit vectors and place them in the same matrix


```python
A.normalized() # Easy way no normalize a matrix
```


```python
B.normalized()
```


```python
Q = Matrix([[sqrt(3) / 3, 0], [sqrt(3) / 3, -sqrt(2) / 2], [sqrt(3) / 3, sqrt(2) / 2]])
Q
```

* The columnspace of the original matrix (of two column vectors) and Q are the same

* In python&#8482; we can use the following code


```python
# The column matrices (independant orthogonal column vectors) are entered indivisually inside square bracket []
A = [Matrix([1, 1, 1]), Matrix([1, 0, 2])]
A
```


```python
Q = GramSchmidt(A, True) # The True argument normalizes the columns
Q
```

## Example problems

### Example problem 1

* Create an orthogonal matrix from the following matrix
$$ \begin{bmatrix} 1 & 2 & 4 \\ 0 & 0 & 5 \\ 0 & 3 & 6 \end{bmatrix} $$

#### Solution


```python
A = [Matrix([1, 0, 0]), Matrix([2, 0, 3]), Matrix([4, 5, 6])]
A
```


```python
Q = GramSchmidt(A, True)
Q
```

* We can also consider QR-factorization


```python
from sympy.mpmath import matrix, qr
```


```python
A = matrix([[1, 2, 4], [0, 0, 5], [0, 3, 6]])
print(A)
```


```python
Q, R = qr(A)
```


```python
print(Q)
```


```python
print(R)
```


```python

```

f67ff3740ba70da6a7238cd6768a1e227753c7dc
ipynb
Jupyter Notebook
_math/MIT_OCW_18_06_Linear_algebra/II_04_Orthogonal_matrices_Gram_Schmidt.ipynb
aixpact/data-science
f04a54595fbc2d797918d450b979fd4c2eabac15
MIT
2020-07-22T23:12:39.000Z
2020-07-25T02:30:48.000Z
_math/MIT_OCW_18_06_Linear_algebra/II_04_Orthogonal_matrices_Gram_Schmidt.ipynb
aixpact/data-science
f04a54595fbc2d797918d450b979fd4c2eabac15
MIT
_math/MIT_OCW_18_06_Linear_algebra/II_04_Orthogonal_matrices_Gram_Schmidt.ipynb
aixpact/data-science
f04a54595fbc2d797918d450b979fd4c2eabac15
MIT
Qwen/Qwen-72B
1. YES
2. YES