<jupyter_start><jupyter_text>Singular Value Decomposition If a square matrix, X, is diagonalizable (ie. an n x n matrix has n linearly independent eigenvectors), then it has what is known as an eigendecomposition:\begin{equation*}\mathbf{X}_{nxn} = \mathbf{P}_{nxn} \times \mathbf{D}_{nxn} \times \mathbf{P}^{-1}_{nxn}\end{equation*}Here, the columns of P are the eigenvectors of X, with their corresponding eigenvalues as the diagonal of D. Since P is an orthonormal matrix, it can be written:\begin{equation*}\mathbf{X}_{nxn} = \mathbf{P}_{nxn} \times \mathbf{D}_{nxn} \times \mathbf{P}^{T}_{nxn}\end{equation*}There exist several other forms of matrix decomposition such as QR decomposition or LU factorization. The singular value decomposition of a matrix, X, is a generalization of matrix decomposition that can be applied to any (including non-square) matricies, written in the form:\begin{equation*}\mathbf{X}_{nxm} = \mathbf{U}_{nxn} \times \mathbf{\Sigma}_{nxm} \times \mathbf{V}^{T}_{mxm}\end{equation*}<jupyter_code>import numpy as np
import numpy.linalg as la
np.random.seed(0)<jupyter_output><empty_output><jupyter_text>Let's start by generating a random matrix<jupyter_code>X = np.random.rand(10, 5)*10
X<jupyter_output><empty_output><jupyter_text>From here we can easily use a linear algebra library to compute the SVD decomposition<jupyter_code>u, s, vt = la.svd(X)
print(f"Shape of U: {u.shape}\n{u}")
print(f"Shape of S: {s.shape}\n{s}")
print(f"Shape of VT: {vt.shape}\n{vt}")<jupyter_output>Shape of U: (10, 10)
[[-0.31808548  0.07010396 -0.32968818 -0.20815326  0.03691277 -0.02946223
  -0.25638887 -0.46442258 -0.58817451 -0.33508442]
 [-0.39389315 -0.05336169  0.31509607 -0.00527192 -0.21680262 -0.75122741
   0.31358241  0.10104013 -0.09260469 -0.11925904]
 [-0.34811653  0.3079775   0.34400872  0.11394361 -0.10649223  0.32986776
   0.05509008 -0.32520012  0.46238299 -0.45814648]
 [-0.29454003 -0.66054805  0.30903705 -0.16418078  0.24752302  0.04648673
  -0.2394105  -0.32836927  0.16263075  0.31534429]
 [-0.36601696  0.4769409  -0.08382718  0.15530076  0.04631442 -0.20468167
  -0.58210152  0.13210523  0.16860938  0.42267119]
 [-0.3135048  -0.23949071  0.04253854  0.38383416 -0.59753613  0.40323796
  -0.02800053  0.22059419 -0.32276126  0.14815472]
 [-0.2430838   0.27000639 -0.05134591 -0.64772875 -0.15503673  0.23960138
   0.41862828 -0.02566453 -0.00955963  0.43401955]
 [-0.3932126  -0.06265564  0.01596145 -0.10592092  0.53700403  0.22321976
   0.04292737  0.63470004 -0.1[...]<jupyter_text>Let's transform s (a single vector) into the diagonal matrix we're looking for.<jupyter_code># Eigenvectors of covariance matrix
sigma = np.diag(s)

add_rows = max(u.shape[0] - vt.shape[0], 0)
add_columns = max(vt.shape[0] - u.shape[0], 0)

# Fix n
sigma = np.concatenate(
    (sigma, np.zeros((add_rows, sigma.shape[1]))),
    axis=0
)
# Fix m
sigma = np.concatenate(
    (sigma, np.zeros((sigma.shape[0], add_columns))),
    axis=1
)
print(f"Shape of sigma: {sigma.shape}\n{sigma}")<jupyter_output>Shape of sigma: (10, 5)
[[39.45197377  0.          0.          0.          0.        ]
 [ 0.         12.11383468  0.          0.          0.        ]
 [ 0.          0.          6.97052559  0.          0.        ]
 [ 0.          0.          0.          6.30369226  0.        ]
 [ 0.          0.          0.          0.          5.11348667]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]]<jupyter_text>Note that $\Sigma$ is a diagonal matrix of the form:\begin{equation*}\Sigma = \begin{bmatrix}     \sigma_{0} & 0          & \dots  & \dots & 0 & \dots \\    0          & \sigma_{1} & 0      & \dots & 0 & \dots \\    \vdots     & \ddots     & \ddots & \ddots & 0 & \dots \\    \vdots      & \ddots      & \ddots  & \ddots & \sigma_{k-1} & \dots\\    \vdots      & \ddots      & \ddots  & \ddots & \ddots & \ddots \\    \end{bmatrix}\end{equation*}Such that:\begin{equation*}\sigma_{0} > \sigma_{1} > ... > \sigma_{k-1}\end{equation*}These diagonal entries are known as the __singular values__ of X.Now that we've reconstructed $\Sigma$, we can easily reconstruct the original matrix.<jupyter_code>def check_two_matricies_almost_equal(x, y, eps=1e-10):
    return np.sum((np.abs(x - y) < eps)) == np.multiply(*list(x.shape))
X_reconstructed = u @ sigma @ vt
check_two_matricies_almost_equal(X, X_reconstructed)<jupyter_output><empty_output><jupyter_text>Note that U and V are orthonormal matricies.<jupyter_code>def check_number_almost_equal(x, y, eps=1e-10):
    return (np.abs(x-y) < eps)

def show_orthonormality(x):
    for col in range(x.shape[1]):
        self_dot = x[:, col] @ x[:, col].T
        equal_one = check_number_almost_equal(self_dot, 1.0)
        print(f"Column {col}: Dot Product With Itself = 1 Is {equal_one}")
        for other_col in range(col+1, x.shape[1]):
            other_dot = x[:, col] @ x[:, other_col].T
            equal_zero = check_number_almost_equal(other_dot, 0.0)
            print(f"Column {col}: Dot Product With {other_col} = 0 Is {equal_zero}")
print(f"MATRIX U: {u}\n")
show_orthonormality(u)
print(f"MATRIX V: {vt.T}\n")
show_orthonormality(vt.T)<jupyter_output>MATRIX V: [[-0.457698    0.4019739  -0.1829935   0.7710481  -0.03045884]
 [-0.37968881  0.47891457 -0.51517818 -0.59365513  0.09296378]
 [-0.50493055 -0.41420274 -0.0135855  -0.11656918 -0.74813655]
 [-0.54138922  0.09006533  0.74727999 -0.17794927  0.32968565]
 [-0.31351523 -0.65527258 -0.37748228  0.08834084  0.56747606]]

Column 0: Dot Product With Itself = 1 Is True
Column 0: Dot Product With 1 = 0 Is True
Column 0: Dot Product With 2 = 0 Is True
Column 0: Dot Product With 3 = 0 Is True
Column 0: Dot Product With 4 = 0 Is True
Column 1: Dot Product With Itself = 1 Is True
Column 1: Dot Product With 2 = 0 Is True
Column 1: Dot Product With 3 = 0 Is True
Column 1: Dot Product With 4 = 0 Is True
Column 2: Dot Product With Itself = 1 Is True
Column 2: Dot Product With 3 = 0 Is True
Column 2: Dot Product With 4 = 0 Is True
Column 3: Dot Product With Itself = 1 Is True
Column 3: Dot Product With 4 = 0 Is True
Column 4: Dot Product With Itself = 1 Is True<jupyter_text>It follows naturally that $U^{-1} = U^{T}$ and $V^{-1} = V^{T}$<jupyter_code>check_two_matricies_almost_equal(u.T @ u, np.identity(u.shape[0]))
check_two_matricies_almost_equal(u.T @ u, u @ u.T)
check_two_matricies_almost_equal(vt.T @ vt, np.identity(vt.shape[0]))
check_two_matricies_almost_equal(vt.T @ vt, vt @ vt.T)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1381051381051381, 'fraction_numerical': 0.16394416394416395, 'mean_word_length': 2.7604690117252932, 'pattern_counts': {'":': 0, '<': 28, '<?xml version=': 0, '>': 29, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 9, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '7', 'dir': 'jupyter-structured-clean-dedup', 'id': '396798', 'license': 'MIT', 'n_tokens_mistral': 3221, 'n_tokens_neox': 2604, 'n_words': 671, 'path': 'p0_svd_intro/Singular Value Decomposition.ipynb', 'repo_name': 'ehelenowski/learn_by_notebooks', 'size': '6734'}