About This Notebook¶

Hello, I am a Computer Engineering student with a strong interest in mathematics and building a deep, first-principles understanding of concepts.
This notebook explores the mathematics behind Singular Value Decomposition (SVD), its implementation from scratch, NumPy in Python, and its application to image compression.

Required Packages:¶

NumPy: For arrays.
Matplotlib: For visualisation and viewing images.
OpenCV Python: For reading and writing images.
Jupyter Lab: For running Jupyter Notebooks.
All these packages are mentioned in requirements.txt

In [1]:
import numpy as np
import matplotlib.pyplot as plt
import cv2 as cv
from IPython.display import display, Math

Singular Value Decomposition (SVD):¶

SVD is a matrix factorisation technique that decomposes a matrix $A\in\mathbb{R}^{m\times n}$ in three parts, $U$, $\Sigma$, $V$ such that $A=U\Sigma V^{T}$
Where:
$U\in\mathbb{R}^{m\times m}$ is an orthogonal matrix whose columns are the normalised eigenvectors of matrix $AA^{T}$
$V\in\mathbb{R}^{n\times n}$ is an orthogonal matrix whose columns are the normalised eigenvectors of matrix $A^{T}A$
$\Sigma\in\mathbb{R}^{m\times n}$ is a diagonal matrix containing singular values of $A$
$\Sigma=\begin{bmatrix}\sigma_{1}&0&\cdots&0\\0&\sigma_{2}&\cdots&0\\\vdots&\vdots&\ddots&\vdots\\0&0&\cdots&\sigma_{\min(m,n)}\end{bmatrix}_{m\times n}$ such that $\sigma_{i}=\sqrt{\lambda_{i}}$ where $\lambda_{i}$ is the ith eigenvalue of $AA^{T}$ or $A^{T}A$. The eigenvalues are arranged in descending order.

Example:¶

Suppose $A=\begin{bmatrix}3&2&2\\2&3&-2\end{bmatrix}$
To find the SVD of $A$ we need to find $AA^{T}$ and $A^{T}A$
$$ AA^{T}=\begin{bmatrix} 3&2&2\\ 2&3&-2 \end{bmatrix} \begin{bmatrix} 3&2\\ 2&3\\ 2&-2 \end{bmatrix}=\begin{bmatrix} 17&8\\ 8&17 \end{bmatrix} $$ Eigenvalues of $AA^{T}$ are $9$ and $25$
Eigenvectors of $AA^{T}$ are $X_{9}=\langle\frac{1}{\sqrt{2}},-\frac{1}{\sqrt{2}}\rangle$ and $X_{25}=\langle\frac{1}{\sqrt{2}},\frac{1}{\sqrt{2}}\rangle$
Note:
$X_{\lambda}$ is the eigenvector corresponding to eigenvalue $\lambda$.
Hence; $$ U=\begin{bmatrix} X_{25}&X_{9} \end{bmatrix}=\begin{bmatrix} \frac{1}{\sqrt{2}}&\frac{1}{\sqrt{2}}\\ \frac{1}{\sqrt{2}}&-\frac{1}{\sqrt{2}} \end{bmatrix} $$ Singular values are:
$\sigma_{1}=\sqrt{25}=5$ and $\sigma_{2}=\sqrt{9}=3$ $$ \Sigma=\text{diag}(\sigma_{1},\sigma_{2})=\begin{bmatrix} 5&0\\ 0&3 \end{bmatrix} $$ $$ A^{T}A=\begin{bmatrix} 3&2\\ 2&3\\ 2&-2 \end{bmatrix} \begin{bmatrix} 3&2&2\\ 2&3&-2 \end{bmatrix}=\begin{bmatrix} 13&12&2\\ 12&13&-2\\ 2&-2&8 \end{bmatrix} $$ Eigenvectors are: $\langle{\frac{1}{\sqrt{2}},\frac{1}{\sqrt{2}},0}\rangle$, $\langle{\frac{1}{\sqrt{18}},-\frac{1}{\sqrt{18}},\frac{4}{\sqrt{18}}}\rangle$, $\langle{-\frac{2}{3},\frac{2}{3},\frac{1}{3}}\rangle$
$$ \therefore V=\begin{bmatrix} \frac{1}{\sqrt{2}}&\frac{1}{\sqrt{18}}&-\frac{2}{3}\\ \frac{1}{\sqrt{2}}&-\frac{1}{\sqrt{18}}&\frac{2}{3}\\ 0&\frac{4}{\sqrt{18}}&\frac{1}{3} \end{bmatrix} $$
$$ \therefore A= \begin{bmatrix} \frac{1}{\sqrt{2}}&\frac{1}{\sqrt{2}}\\ \frac{1}{\sqrt{2}}&-\frac{1}{\sqrt{2}} \end{bmatrix}\begin{bmatrix} 5&0\\ 0&3 \end{bmatrix}\begin{bmatrix} \frac{1}{\sqrt{2}}&\frac{1}{\sqrt{2}}&0\\ \frac{1}{\sqrt{18}}&-\frac{1}{\sqrt{18}}&\frac{4}{\sqrt{18}}\\ -\frac{2}{3}&\frac{2}{3}&\frac{1}{3} \end{bmatrix} $$

Using the function made from scratch:¶

In [2]:
from mathematical_foundation import svd_from_scratch
A = np.array([[3, 2, 2],
              [2, 3, -2]])
U, S, Vt = svd_from_scratch.svd_from_scratch(A)
display(Math(r"U="))
print(U)
display(Math(r"\Sigma="))
print(S)
display(Math(r"V^{T}"))
print(Vt)
$\displaystyle U=$
[[-0.70710678  0.70710678]
 [-0.70710678 -0.70710678]]
$\displaystyle \Sigma=$
[[5. 0. 0.]
 [0. 3. 0.]]
$\displaystyle V^{T}$
[[-7.07106781e-01 -7.07106781e-01 -1.38777878e-17]
 [ 2.35702260e-01 -2.35702260e-01  9.42809042e-01]
 [-6.66666667e-01  6.66666667e-01  3.33333333e-01]]

Using optimised function (NumPy function):¶

This does not construct full matrices, which leads to optimisation.

In [3]:
from optimized_method import optimized_svd
U, S, Vt = optimized_svd.optimized_svd(A)
display(Math(r"U="))
print(U)
display(Math(r"\Sigma="))
print(S)
display(Math(r"V^{T}"))
print(Vt)
$\displaystyle U=$
[[-0.70710678 -0.70710678]
 [-0.70710678  0.70710678]]
$\displaystyle \Sigma=$
[5. 3.]
$\displaystyle V^{T}$
[[-7.07106781e-01 -7.07106781e-01 -6.47932334e-17]
 [-2.35702260e-01  2.35702260e-01 -9.42809042e-01]]

Approximating a matrix:¶

We can use this decomposition to reconstruct an approximation of the original matrix.

The math for approximation using SVD:¶

$$\begin{aligned} A&=U\Sigma V^{T}\\ U&=\begin{bmatrix}u_{1}&u_{2}&\cdots&u_{i}&\cdots&u_{m}\end{bmatrix},\quad U\in\mathbb{R}^{m\times m}\quad u_{i}\in\mathbb{R}^{m\times1}\\ V^{T}&=\begin{bmatrix}v_{1}\\v_{2}\\\vdots\\v_{i}\\\vdots\\v_{n}\end{bmatrix},\quad V\in\mathbb{R}^{n\times n}\quad v_{i}\in\mathbb{R}^{1\times n}\\ \Sigma&=\text{diag}(\sigma_{1},\sigma_{2},\ldots,\sigma_{\min(m,n)})\\ A&=\sum_{i=1}^{\min{(m,n)}}{\sigma_{i}u_{i}v_{i}}\\ &\text{Suppose }k<\min(m, n)\\ A&\approx\sum_{i=1}^{k}{\sigma_{i}u_{i}v_{i}} \end{aligned}$$

In [4]:
from metrics import metrics
double_equal_line = "====================================================================="
rng = np.random.default_rng(seed=42)
A = rng.random((50, 50))
U, S, Vt = optimized_svd.optimized_svd(A)
steps = [5*i for i in range(10, 0, -1)]
print("A:")
print(A)
print(double_equal_line)
for i in steps:
    print(double_equal_line)
    print(f"k={i}:")
    print(optimized_svd.reconstruct(U, S, Vt, i))
    print(f"Reconstruction Error: {metrics.reconstruction_error(S, i)}")
    print(double_equal_line)
A:
[[0.77395605 0.43887844 0.85859792 ... 0.2883281  0.6824955  0.13975248]
 [0.1999082  0.00736227 0.78692438 ... 0.08764992 0.1180059  0.96189766]
 [0.90858069 0.69970713 0.26586996 ... 0.5206724  0.43891146 0.02161208]
 ...
 [0.05628881 0.85363017 0.23470349 ... 0.18949128 0.12392589 0.55525562]
 [0.5964451  0.79227084 0.75307327 ... 0.40427646 0.58297826 0.58519135]
 [0.65733396 0.72021603 0.6415472  ... 0.2117411  0.81945076 0.49860004]]
=====================================================================
=====================================================================
k=50:
[[0.77395605 0.43887844 0.85859792 ... 0.2883281  0.6824955  0.13975248]
 [0.1999082  0.00736227 0.78692438 ... 0.08764992 0.1180059  0.96189766]
 [0.90858069 0.69970713 0.26586996 ... 0.5206724  0.43891146 0.02161208]
 ...
 [0.05628881 0.85363017 0.23470349 ... 0.18949128 0.12392589 0.55525562]
 [0.5964451  0.79227084 0.75307327 ... 0.40427646 0.58297826 0.58519135]
 [0.65733396 0.72021603 0.6415472  ... 0.2117411  0.81945076 0.49860004]]
Reconstruction Error: 0.0
=====================================================================
=====================================================================
k=45:
[[ 0.76246363  0.43189955  0.83867593 ...  0.2666112   0.66119
   0.13321381]
 [ 0.21759429 -0.01080555  0.77590924 ...  0.07432289  0.11254037
   0.9566713 ]
 [ 0.92147975  0.69287097  0.2693694  ...  0.53567002  0.44143431
   0.02166681]
 ...
 [ 0.06579091  0.85207624  0.23599319 ...  0.18406779  0.13033416
   0.5563018 ]
 [ 0.59128546  0.78746045  0.75152602 ...  0.40649403  0.57954252
   0.58213329]
 [ 0.6592782   0.72478211  0.64374869 ...  0.20314578  0.82706457
   0.50078761]]
Reconstruction Error: 0.43681817604114564
=====================================================================
=====================================================================
k=40:
[[ 0.77914445  0.43154519  0.86800254 ...  0.2770824   0.65893636
   0.11642856]
 [ 0.19568832 -0.0211405   0.7631807  ...  0.09996862  0.1390347
   0.93813222]
 [ 0.91743176  0.70873322  0.26255742 ...  0.55434827  0.45407572
   0.02032964]
 ...
 [ 0.07841245  0.84971359  0.26272525 ...  0.20698328  0.1319369
   0.52947423]
 [ 0.58728743  0.79544568  0.75404982 ...  0.43620771  0.5994958
   0.567482  ]
 [ 0.66719921  0.718862    0.63437678 ...  0.21161138  0.80362461
   0.4978805 ]]
Reconstruction Error: 1.154451740075357
=====================================================================
=====================================================================
k=35:
[[0.77694133 0.43773741 0.86087428 ... 0.29740162 0.59023372 0.12386661]
 [0.19090853 0.00390548 0.70495844 ... 0.04169258 0.15981635 0.88411585]
 [0.92665377 0.68967529 0.27757104 ... 0.62206489 0.39708764 0.05162488]
 ...
 [0.07062957 0.84775746 0.29561684 ... 0.19614038 0.14455977 0.53658679]
 [0.61757503 0.78373118 0.75217779 ... 0.50816693 0.56218218 0.6063632 ]
 [0.65706947 0.74195963 0.61745935 ... 0.13275519 0.87091192 0.46787816]]
Reconstruction Error: 2.1103824296413665
=====================================================================
=====================================================================
k=30:
[[0.78251897 0.38005235 0.85748105 ... 0.30598354 0.6041642  0.23826775]
 [0.20083904 0.08026046 0.70310383 ... 0.04434625 0.23322392 0.90823623]
 [0.89855658 0.57957807 0.24054025 ... 0.64620931 0.30983831 0.1022256 ]
 ...
 [0.04799909 0.85658667 0.40033847 ... 0.12800414 0.21889606 0.47279622]
 [0.59502466 0.73817822 0.77252252 ... 0.45789779 0.51914022 0.50333906]
 [0.65191365 0.65430733 0.57360612 ... 0.1781841  0.73924768 0.46554828]]
Reconstruction Error: 3.403421894211069
=====================================================================
=====================================================================
k=25:
[[0.6728902  0.34159925 0.88690843 ... 0.2597378  0.55745265 0.27381008]
 [0.27516757 0.06757712 0.5755006  ... 0.0174259  0.26159473 0.86434522]
 [0.86743985 0.52015796 0.2377472  ... 0.7351008  0.29140948 0.18042703]
 ...
 [0.01718583 0.88052784 0.44545784 ... 0.229741   0.17548115 0.50955312]
 [0.6134489  0.7706738  0.82519995 ... 0.44530478 0.64719963 0.47115084]
 [0.68897941 0.70783342 0.5693888  ... 0.23303359 0.75956407 0.43748992]]
Reconstruction Error: 4.879313828884314
=====================================================================
=====================================================================
k=20:
[[0.56545679 0.38871705 1.02918977 ... 0.15686298 0.57535796 0.34887795]
 [0.2640601  0.08442598 0.55069003 ... 0.04223734 0.27890005 0.8708532 ]
 [0.86683377 0.56748996 0.3216705  ... 0.63328232 0.15754047 0.20213758]
 ...
 [0.09417762 0.80320076 0.38490773 ... 0.29374114 0.19189351 0.43180028]
 [0.56398626 0.78756377 0.78824333 ... 0.49198891 0.67629633 0.51571416]
 [0.82679272 0.80151818 0.54958084 ... 0.1105062  0.69411017 0.38713921]]
Reconstruction Error: 6.4165530030656655
=====================================================================
=====================================================================
k=15:
[[0.55914382 0.38396948 0.85563264 ... 0.23141459 0.67296436 0.36239759]
 [0.31539327 0.1339768  0.67501438 ... 0.14258255 0.24354873 0.74679405]
 [1.01452485 0.53116042 0.48114697 ... 0.51998596 0.2426378  0.0246075 ]
 ...
 [0.18018987 0.71490431 0.35263792 ... 0.30822214 0.32668413 0.38032587]
 [0.636415   0.74280517 0.83205779 ... 0.57660581 0.67412906 0.5666082 ]
 [0.89787772 0.85352677 0.69788615 ... 0.39434983 0.57441715 0.38188592]]
Reconstruction Error: 8.184063246823262
=====================================================================
=====================================================================
k=10:
[[0.37635098 0.35519284 0.86307702 ... 0.17749449 0.64672962 0.51092533]
 [0.33114711 0.14963397 0.64334369 ... 0.26386376 0.35262762 0.63021994]
 [0.82147129 0.44582441 0.50841448 ... 0.61892857 0.33564542 0.17791365]
 ...
 [0.18762034 0.6028922  0.31736134 ... 0.39624113 0.16837416 0.35562548]
 [0.59158462 0.81363973 0.82036406 ... 0.69218288 0.5290473  0.66927607]
 [0.62051426 0.7068069  0.67871574 ... 0.62210548 0.5116244  0.56560119]]
Reconstruction Error: 10.108049225024347
=====================================================================
=====================================================================
k=5:
[[0.38005397 0.27760238 0.68247311 ... 0.40696117 0.57171529 0.53661063]
 [0.23369737 0.10437353 0.59831128 ... 0.2015875  0.46174013 0.47891876]
 [0.76267675 0.63135886 0.65206    ... 0.5956501  0.25800367 0.30566458]
 ...
 [0.37977964 0.29724666 0.58420214 ... 0.39745003 0.45130071 0.42781185]
 [0.5840594  0.58944725 0.58765593 ... 0.52447638 0.58979376 0.55998993]
 [0.61002132 0.61105418 0.52185506 ... 0.62526611 0.54263022 0.43730619]]
Reconstruction Error: 12.380650233070941
=====================================================================

As $k$ decreases the reconstruction error increases.

Application to Image Compression:¶

For better understanding, we will try to compress a grayscale image.

About the image used in this notebook:¶

Image Source: The image used in this notebook is AI-generated using This Person Does Not Exist dot com
The dimension of the originally used image is $(819\times819)$

In [5]:
image_name = "image.png"
import os
path = os.path.join("..", "data", image_name)

What is an Image?¶

An image is a visual representation of the world captured at a moment in time.
In mathematical terms, an image is a matrix whose elements are pixels, each representing a light intensity value.

In [6]:
image = cv.imread(path, cv.IMREAD_GRAYSCALE)
image_dimensions = image.shape
print("Image Matrix:")
print(image)
print("Image Dimensions:", image_dimensions)
plt.title("Original Image")
plt.imshow(image, cmap="gray")
plt.show()
Image Matrix:
[[186 181 178 ... 164 164 159]
 [188 180 179 ... 161 163 168]
 [187 180 179 ... 160 164 165]
 ...
 [134 134 136 ... 132 130 132]
 [137 136 137 ... 130 131 134]
 [146 138 135 ... 130 134 144]]
Image Dimensions: (819, 819)
No description has been provided for this image

Decomposing this image:

How to compress an image?¶

We can use SVD to approximate this matrix — the approximated matrix is the compressed image.

In [7]:
_, S, _ = optimized_svd.optimized_svd(image)
print("First 10 singular values:")
print(S[:10])
print("Last 10 singular values:")
print(S[-10:])
First 10 singular values:
[100957.54465131  21927.03110361  16446.95781641  10535.70882771
   7526.01995845   5964.6588328    4922.29266784   4562.8443728
   4192.81900941   3492.41139879]
Last 10 singular values:
[0.55201224 0.5339179  0.42555796 0.39841393 0.30532788 0.25718209
 0.21223382 0.12631836 0.05953665 0.01881724]

From above we see that the first $10$ singular values are very large where as the last $10$ are very small.
The starting singular values store maximum information of the image where as the last singular values only store minor correction which can be ignored, as human eye does not recognise the difference.

Compression using svd_from_scratch():¶

In [8]:
from application import image_compression
k = [5*i for i in range(2, 11)]
figure, axes = plt.subplots(3, 3, figsize=(16, 13))
_, S, _ = svd_from_scratch.svd_from_scratch(image)
for n, i in enumerate(k):
    compressed_image = image_compression.compress_from_scratch(image, i)
    axes[n // 3, n - 3*(n // 3)].imshow(compressed_image, cmap="gray")
    axes[n // 3, n - 3*(n // 3)].set_title(rf"$k={i}$" + f"\nError: {np.round(metrics.reconstruction_error(S, i), 2)}")
    axes[n // 3, n - 3*(n // 3)].axis("off")
plt.show()
No description has been provided for this image

As $k$ increases, the error decreases and the image quality increases. At $k=45$, the image looks quite similar to the original image.

Error, Energy retained and Compression ratio plots:¶

In [9]:
k = np.arange(1, min(image_dimensions)+1, 1)
reconstruction_error = np.empty(min(image_dimensions))
energy_retained = np.empty(min(image_dimensions))
percent_compression_ratio = np.empty(min(image_dimensions))
for i in k:
    reconstruction_error[i-1] = metrics.reconstruction_error(S, i)
    energy_retained[i-1] = metrics.energy_retained(S, i)
    percent_compression_ratio[i-1] = metrics.percent_compression_ratio(image.shape, i)
figure = plt.figure(figsize=(16, 10))
gs = figure.add_gridspec(2, 2)
ax1 = figure.add_subplot(gs[0, 0])
ax2 = figure.add_subplot(gs[0, 1])
ax3 = figure.add_subplot(gs[1, :])

ax1.plot(k, reconstruction_error)
ax1.set_title(r"Reconstruction Error vs $k$")
ax1.set_xlabel(r"$k$")
ax1.set_ylabel(r"Reconstruction Error")
ax1.axhline(0, linestyle="--", color="black")
ax1.grid(True)

ax2.plot(k, energy_retained)
ax2.set_title(r"Energy Retained vs $k$")
ax2.set_xlabel(r"$k$")
ax2.set_ylabel(r"Energy Retained")
ax2.axhline(1, linestyle="--", color="black")
ax2.grid(True)

ax3.plot(k, percent_compression_ratio)
ax3.set_title(r"Percent Compression Ratio vs $k$")
ax3.set_xlabel(r"$k$")
ax3.set_ylabel(r"Percent Compression Ratio")
k = image_compression.get_k_from_compression_ratio(image_dimensions, 100, percentage=True)
ax3.hlines(100, xmin=0, xmax=k, linestyle="--", color="black")
ax3.vlines(k, ymin=0, ymax=100, linestyle="--", color="black")
ax3.grid(True)

plt.show()
display(Math(r"\text{At }k="+f"{k} "+r"\text{ complete image is stored, using more singular values increases the space taken by the image.}"))
No description has been provided for this image
$\displaystyle \text{At }k=409 \text{ complete image is stored, using more singular values increases the space taken by the image.}$

Compression using optimized_svd():¶

In [10]:
k = [5*i for i in range(2, 11)]
figure, axes = plt.subplots(3, 3, figsize=(16, 13))
_, S, _ = optimized_svd.optimized_svd(image)
for n, i in enumerate(k):
    compressed_image = image_compression.optimized_compress(image, i)
    axes[n // 3, n - 3*(n // 3)].imshow(compressed_image, cmap="gray")
    axes[n // 3, n - 3*(n // 3)].set_title(rf"$k={i}$" + f"\nError: {np.round(metrics.reconstruction_error(S, i), 2)}")
    axes[n // 3, n - 3*(n // 3)].axis("off")
plt.show()
No description has been provided for this image

Again, as we increase $k$, image quality increases.

Error, Energy retained and Compression ratio plots:¶

In [11]:
k = np.arange(1, min(image_dimensions)+1, 1)
reconstruction_error = np.empty(min(image_dimensions))
energy_retained = np.empty(min(image_dimensions))
percent_compression_ratio = np.empty(min(image_dimensions))
for i in k:
    reconstruction_error[i-1] = metrics.reconstruction_error(S, i)
    energy_retained[i-1] = metrics.energy_retained(S, i)
    percent_compression_ratio[i-1] = metrics.percent_compression_ratio(image.shape, i)
figure = plt.figure(figsize=(16, 10))
gs = figure.add_gridspec(2, 2)
ax1 = figure.add_subplot(gs[0, 0])
ax2 = figure.add_subplot(gs[0, 1])
ax3 = figure.add_subplot(gs[1, :])

ax1.plot(k, reconstruction_error)
ax1.set_title(r"Reconstruction Error vs $k$")
ax1.set_xlabel(r"$k$")
ax1.set_ylabel(r"Reconstruction Error")
ax1.axhline(0, linestyle="--", color="black")
ax1.grid(True)

ax2.plot(k, energy_retained)
ax2.set_title(r"Energy Retained vs $k$")
ax2.set_xlabel(r"$k$")
ax2.set_ylabel(r"Energy Retained")
ax2.axhline(1, linestyle="--", color="black")
ax2.grid(True)

ax3.plot(k, percent_compression_ratio)
ax3.set_title(r"Percent Compression Ratio vs $k$")
ax3.set_xlabel(r"$k$")
ax3.set_ylabel(r"Percent Compression Ratio")
k = image_compression.get_k_from_compression_ratio(image_dimensions, 100, percentage=True)
ax3.hlines(100, xmin=0, xmax=k, linestyle="--", color="black")
ax3.vlines(k, ymin=0, ymax=100, linestyle="--", color="black")
ax3.grid(True)

plt.show()
display(Math(r"\text{At }k="+f"{k} "+r"\text{ complete image is stored, using more singular values increases the space taken by the image.}"))
No description has been provided for this image
$\displaystyle \text{At }k=409 \text{ complete image is stored, using more singular values increases the space taken by the image.}$

Comparing original image with compressed image having $\approx100\%$ or $\approx50\%$ compression ratio:¶

In [12]:
k = [k, k // 2]

figure, axes = plt.subplots(2, 2, figsize=(8, 8))

axes[0, 0].imshow(image, cmap="gray")
axes[0, 0].set_title("Original Image")
axes[0, 0].axis("off")

axes[1, 0].imshow(image, cmap="gray")
axes[1, 0].set_title("Original Image")
axes[1, 0].axis("off")

compressed_image = image_compression.optimized_compress(image, k[0])
axes[0, 1].imshow(compressed_image, cmap="gray")
axes[0, 1].set_title(rf"$k={k[0]}$" + f"\nCompression Ratio={np.round(metrics.percent_compression_ratio(image.shape, k[0]), 2)}%")
axes[0, 1].axis("off")

compressed_image = image_compression.optimized_compress(image, k[1])
axes[1, 1].imshow(compressed_image, cmap="gray")
axes[1, 1].set_title(rf"$k={k[1]}$" + f"\nCompression Ratio={np.round(metrics.percent_compression_ratio(image.shape, k[1]), 2)}%")
axes[1, 1].axis("off")

plt.show()
No description has been provided for this image

When $\approx50\%$ of image is stored even then the compressed image looks just like original!

Compression of colored images:¶

Colored images have 3 channels:

  1. Red
  2. Green
  3. Blue We can separate these channels, then compress them like a grayscale image and then combine these channels again — we get the compressed image.
    OpenCV Python reads the image in BGR format, and Matplotlib uses RGB format, hence we need to take note of that and compress them accordingly.
    Since we also have to calculate metrics, we will not use the function from the image_compression module.
In [13]:
image = cv.imread(path)
plt.imshow(cv.cvtColor(image, cv.COLOR_BGR2RGB))
plt.title("Original Image")
plt.axis("off")
plt.show()
No description has been provided for this image
In [14]:
B, G, R = cv.split(image)
k = [5*i for i in range(2, 11)]
figures, axes = plt.subplots(3, 3, figsize=(16, 13))
for n, i in enumerate(k):
    compressed_B = image_compression.optimized_compress(B, i)
    compressed_G = image_compression.optimized_compress(G, i)
    compressed_R = image_compression.optimized_compress(R, i)
    compressed_image = cv.merge([np.clip(compressed_B, 0, 255),
                                 np.clip(compressed_G, 0, 255),
                                 np.clip(compressed_R, 0, 255)]).astype("uint8")
    reconstruction_error = np.sqrt(metrics.frobenius_error_squared(B, compressed_B) + metrics.frobenius_error_squared(G, compressed_G) + metrics.frobenius_error_squared(R, compressed_R))
    axes[n // 3, n - 3*(n // 3)].imshow(cv.cvtColor(compressed_image, cv.COLOR_BGR2RGB))
    axes[n // 3, n - 3*(n // 3)].set_title(rf"$k={i}$" + f"\nReconstruction error: {np.round(reconstruction_error, 2)}")
    axes[n // 3, n - 3*(n // 3)].axis("off")
plt.show()
No description has been provided for this image

Again, as $k$ increases, image quality gets better.

Error, Energy retained and Compression ratio plots:¶

In [15]:
k = np.arange(1, min(image_dimensions)+1, 1)
reconstruction_error = np.empty(min(image_dimensions))
energy_retained = np.empty(min(image_dimensions))
percent_compression_ratio = np.empty(min(image_dimensions))

_, S_blue, _ = optimized_svd.optimized_svd(B)
_, S_green, _ = optimized_svd.optimized_svd(G)
_, S_red, _ = optimized_svd.optimized_svd(R)

for i in k:
    reconstruction_error[i-1] = np.sqrt(metrics.reconstruction_error_squared(S_blue, i) +
                                        metrics.reconstruction_error_squared(S_green, i) +
                                        metrics.reconstruction_error_squared(S_red, i))
    energy_retained[i-1] = np.sum(S_blue[:i]**2 + S_green[:i]**2 + S_red[:i]**2) / np.sum(S_blue**2 + S_green**2 + S_red**2)
    percent_compression_ratio[i-1] = metrics.percent_compression_ratio(image.shape, i)

figure = plt.figure(figsize=(16, 10))
gs = figure.add_gridspec(2, 2)
ax1 = figure.add_subplot(gs[0, 0])
ax2 = figure.add_subplot(gs[0, 1])
ax3 = figure.add_subplot(gs[1, :])

ax1.plot(k, reconstruction_error)
ax1.set_title(r"Reconstruction Error vs $k$")
ax1.set_xlabel(r"$k$")
ax1.set_ylabel(r"Reconstruction Error")
ax1.axhline(0, linestyle="--", color="black")
ax1.grid(True)

ax2.plot(k, energy_retained)
ax2.set_title(r"Energy Retained vs $k$")
ax2.set_xlabel(r"$k$")
ax2.set_ylabel(r"Energy Retained")
ax2.axhline(1, linestyle="--", color="black")
ax2.grid(True)

ax3.plot(k, percent_compression_ratio)
ax3.set_title(r"Percent Compression Ratio vs $k$")
ax3.set_xlabel(r"$k$")
ax3.set_ylabel(r"Percent Compression Ratio")
k = image_compression.get_k_from_compression_ratio(image_dimensions, 1)
ax3.hlines(100, xmin=0, xmax=k, linestyle="--", color="black")
ax3.vlines(k, ymin=0, ymax=100, linestyle="--", color="black")
ax3.grid(True)

plt.show()
display(Math(r"\text{At }k="+f"{k} "+r"\text{ complete image is stored, using more singular values increases the space taken by the image.}"))
No description has been provided for this image
$\displaystyle \text{At }k=409 \text{ complete image is stored, using more singular values increases the space taken by the image.}$

Comparing original image with compressed image having $\approx100\%$ or $\approx50\%$ compression ratio:¶

In [16]:
k = [k, k // 2]

figure, axes = plt.subplots(2, 2, figsize=(8, 8))

axes[0, 0].imshow(cv.cvtColor(image, cv.COLOR_BGR2RGB))
axes[0, 0].set_title("Original Image")
axes[0, 0].axis("off")

axes[1, 0].imshow(cv.cvtColor(image, cv.COLOR_BGR2RGB))
axes[1, 0].set_title("Original Image")
axes[1, 0].axis("off")

compressed_image = image_compression.optimized_compress(image, k[0])
axes[0, 1].imshow(cv.cvtColor(compressed_image, cv.COLOR_BGR2RGB))
axes[0, 1].set_title(rf"$k={k[0]}$" + f"\nCompression Ratio={np.round(metrics.percent_compression_ratio(image.shape, k[0]), 2)}%")
axes[0, 1].axis("off")

compressed_image = image_compression.optimized_compress(image, k[1])
axes[1, 1].imshow(cv.cvtColor(compressed_image, cv.COLOR_BGR2RGB))
axes[1, 1].set_title(rf"$k={k[1]}$" + f"\nCompression Ratio={np.round(metrics.percent_compression_ratio(image.shape, k[1]), 2)}%")
axes[1, 1].axis("off")

plt.show()
No description has been provided for this image

When $\approx50\%$ of image is stored even then the compressed image looks just like original!