<jupyter_start><jupyter_text>The beauty of the Gershgorin disk theoremIn this post I'll talk about one of the most beatiful theorems I've encountered while studying linear algebra. I bumped into it while taking the ACM104 Applied Linear Algebra course at Caltech, taught by Prof. <NAME>. Btw he is amazing - check out his math videos on Youtube  !  I find this theorem very aesthetically pleasing because it has a visual representation, and also a remarkable number of applications. Moreover, I think that the wit with which it came about is to be praised for. It is a clear example of how mathematics are *just right there*, waiting to be discovered. We follow the derivation of the book *Applied Linear Algebra* by Shakiban.  The statement Before stating the theorem, we'll need some definitions. Also, recall that for a complex number $z = a + ib$ its magnitude is defined by $|z| = \sqrt{a^2 + b^2}$. A nice equality is also that if we define the conjugate $\bar{z} = a - ib$, then $|z|^2 = z\bar{z}$. **Definition:** *Gershgorin disk*. Let $A \in \mathbb{M}_{n \times n}$ be a real matrix over $\mathbb{F}$. For each $1 \le i \le n$, define the $i$-th Gershgorin disk as : $$D_i = \{ |A_{ii} - z_i | \le r_i : z \in \mathbb{C} \}$$where $r_i = \sum_{j = 1, j \neq i} |A_ij| $**Definition** *Gershgorin domain*. The union of the $n$ Gershgorin disks is the Gershgorin domain. $$\mathfrak{D}_A = \{ \cup_{i = 1}^n D_i \} \subseteq \mathbb{C}$$**Definition** *Spectrum of a matrix*. We call the spectrum of a matrix A to the set of eigenvalues associated with A. We denote it as $\mathrm{spec} A$.  **Theorem** $$\mathrm{spec A} \subseteq \mathfrak{D}_A \subseteq \mathbb{C}$$*Proof* The constructive proof is surprisingly straightforward. Let $\lambda$ be an eigenvalue of A, $\vec{v}$ be its associated eigenvector. Let $\vec{u} = \frac{v}{ || v ||_{\infty} }$ be the corresponding unit eigenvector with respect to (w.r.t.) the  $\infty$ norm, i.e. $$|| u ||_\infty  = \mathrm{max} \{ |u_1|, ..., |u_n| \} = 1$$Let $u_i$ be an entry of $\vec{u}$ that achievs the maximum $|u_i| = 1$. Writing out the $i$-th component of the equation $\mathbf{A}\vec{v} = \lambda \vec{v}$ we obtain: $$\sum_{j = 1}^ n \mathbf{A}_{ij} u_j = \lambda u_i$$which we can rewrite as: $$\sum_{j \neq i } \mathbf{A}_{ij} u_j = \lambda u_i - \mathbf{A}_{ii} u_i = (\lambda - \mathbf{A}_{ii} ) u_i$$Note: in the last step we just subtracted $\mathbf{A}_{ii} u_i $ from both sides. Thus, since all $|u_j| \le 1 $ while $|u_i|$ = 1 we have that \begin{align}|\mathbf{A}_{ii} - \lambda| &= |\lambda - \mathbf{A}_{ii}|\\[1em]                          &= |\lambda - \mathbf{A}_{ii}| |u_i|\\[1em]                          &= | (\lambda - \mathbf{A}_{ii}) u_i | \\[.7em]                          &= | \sum_{j \neq i} \mathbf{A}_{ij} u_j | \\[.7em]                          &\le \sum_{j \neq i } |\mathbf{A}_{ij}||u_j| \\[.7em]                          &\le \sum_{j \neq i} |\mathbf{A}_{ij}|\\[1em]                          &= r_i.\end{align}This implies that $|\lambda - \mathbf{A}_{ii}| \le r_i$, i.e. $\lambda \in D_i$. We have equality in the third step as can be checked for all cases ( ++, +-, -+, -- ). The fourth step is just substituting the equation above. The fifth step holds by the triangle inequality $ ||x+y|| \le ||x|| + ||y|| $. The sixth line holds as $|u_j| \le 1 \forall j \neq i$ by construction of $\vec{u}$.  **Definition.** A square matrix is called diagonally dominant if $$|a_{ii}| > \sum_{j \neq i} |a_{ij}| \forall i = 1, ..., n.$$**Corollary.** A strictly diagonally dominant matrix is **nonsingular**.   ExamplesHere are some visualizations in a jupyter notebook if you want to get a feel of the theorem.<jupyter_code>import numpy as np 
import matplotlib.pyplot as plt 
from numpy import linalg as la
from matplotlib import rc
import panel as pn
pn.extension()

#rc('font',**{'family':'sans-serif','sans-serif':['Helvetica']})
## for Palatino and other serif fonts use:
#rc('font',**{'family':'serif','serif':['Palatino']})
#rc('text', usetex=True)
%config InlineBackend.figure_format = "svg"<jupyter_output><empty_output><jupyter_text>Let's test the hypothesis that a random matrix is invertible/non-singular with high probability.<jupyter_code>rng = np.random.default_rng(seed=47)
A = rng.random((3,3))
A<jupyter_output><empty_output><jupyter_text>We can quickly check if the matrix is non-singular by checking that the determinant is not zero.<jupyter_code>la.det(A)<jupyter_output><empty_output><jupyter_text>Indeed it is not singular, but close to being singular. Let's plot its Gershgorin disks.<jupyter_code>def get_gershgorin_disk_radius(A, idx): 
    "Returns Gershgorin disk radius for square matrix."
    r_i = np.delete(A[idx, :], idx).sum()
    return r_i
def plot_gershgorin_domain(A):
    n = A.shape[0]
    fig, ax = plt.subplots(figsize=(4,4))
    # Extract diagonal (A_ii's)
    diag = np.diag(A)
    diag_x,diag_y = np.real(diag),np.imag(diag)
    radius = np.zeros(n)

    for i in range(n):
        a_ii = diag[i]
        r_i = get_gershgorin_disk_radius(A, i)
        radius[i] = r_i # Save radius
        
        real_aii, imag_aii = diag_x[i], diag_y[i]

        plt.scatter(
            real_aii, imag_aii, marker = "D", color = "lightgrey", 
            s= 5
        )
        ax.add_patch(
            plt.Circle(
                (real_aii, imag_aii), r_i, alpha = 0.3, facecolor="lightblue", 
                edgecolor = "lightgrey"
            )
        )

    plt.scatter(0, 0, marker = "X", color = "k", label = "origin",s= 10)
    
    # Make symmetric limits 
    max_radius = max(radius) + 0.5
    maxmin_centers = max(diag_x), min(diag_x), max(diag_y), min(diag_y)
    maxmin_centers = [abs(x) for x in maxmin_centers]
    lim = np.max(maxmin_centers)
    
    plt.xlim(-lim - max_radius , lim + max_radius)
    plt.ylim(-lim - max_radius , lim + max_radius)
    plt.legend()
    
    plt.xlabel("$\mathrm{Re}(z)$", fontsize=18)
    plt.ylabel("$\mathrm{Im}(z)$", fontsize =18)
    plt.title("$\mathcal{D}_A$", fontsize = 20)
    
    plt.close(fig);
    return fig
fig=plot_gershgorin_domain(A)
fig<jupyter_output><empty_output><jupyter_text>Indeed, we can see that the origin is in the Gershgorin domain $ \mathcal{D}_A$. We can use our corollary by increasing the values along the diagonal.<jupyter_code>def get_diag_dominant(A_, eps):
    """Scales values along the diagonal of a square matrix."""
    return A_ + np.diag(np.diag(A_)*eps)
def get_diag_dominant_iden(A_, eps):
    """Regularization, applies A = A_ + ϵI."""
    n = A.shape[0]
    
    return A_ + eps*np.eye(n)
eps = 20 # multiplier

A_diag_dom = get_diag_dominant(A,eps)
fig = plot_gershgorin_domain(A_diag_dom)
fig<jupyter_output><empty_output><jupyter_text>By heavily increasing the diagonal we can help making a matrix invertible. This process is usually called regularization. Let's visualize the process of how both the Gershgorin disks and the determinant changes as we increase the value of our multiplier.<jupyter_code>data_dict = {
    "A":A, 
    # add another matrix here 
    # "name_of_mat": np.array
}
eps_slider = pn.widgets.FloatSlider(start= 0.0, end = 20, step = 0.5, value =1)
data_widget = pn.widgets.Select(
    name='data', 
    options=list(data_dict.keys())
)
@pn.depends(eps_slider.param.value, data_widget.param.value)
def diag_dom_gershgorin(eps, data_name): 
    """
    Visualize "diagonally dominant-alizing" a matrix.
    """
    A_ = data_dict[data_name]
    
    A = get_diag_dominant(A_, eps)
    
    n = A.shape[0]
    
    #fig, ax1 = plt.subplots(figsize=(4,4))
    fig = plt.figure(figsize = (10, 5))
    ax1 = fig.add_subplot(121)
    
    diag = np.diag(A) # Extract diagonal (A_ii's)
    diag_x,diag_y = np.real(diag),np.imag(diag)
    radius = np.zeros(n)

    for i in range(n):
        a_ii = diag[i]
        r_i = get_gershgorin_disk_radius(A, i)
        radius[i] = r_i # Save radius
        real_aii, imag_aii = diag_x[i], diag_y[i]
        
        plt.scatter(
            real_aii, imag_aii, marker = "D", color = "lightgrey", 
            s= 5
        )
        ax1.add_patch(
            plt.Circle(
                (real_aii, imag_aii), r_i, alpha = 0.3, facecolor="lightblue", 
                edgecolor = "lightgrey"
            )
        )


    plt.scatter(0, 0, marker = "X", color = "k", label = "origin",s= 5)
    
    # Make symmetric limits 
    max_radius = max(radius) + 0.5
    maxmin_centers = max(diag_x), min(diag_x), max(diag_y), min(diag_y)
    maxmin_centers = [abs(x) for x in maxmin_centers]
    lim = np.max(maxmin_centers)
    
    plt.xlim(-lim - max_radius , lim + max_radius)
    plt.ylim(-lim - max_radius , lim + max_radius)
    plt.legend()
    
    plt.xlabel("$\mathrm{Re}(z)$", fontsize=18)
    plt.ylabel("$\mathrm{Im}(z)$", fontsize =18)
    plt.title("$\mathcal{D}_A$", fontsize = 20)
    
    ax2 = fig.add_subplot(122)
    evals = np.linspace(0, eps, 500)
    dets = [la.det(get_diag_dominant(A_, eps)) for eps in evals]
    #plt.scatter(evals, dets, s = 5, color = "lightgrey")
    plt.semilogy(evals, np.abs(dets), color = "lightblue", lw=1)
    
    plt.xlim(0, 23)
    plt.ylim(0,500)
    plt.xlabel("$\epsilon$", fontsize=18)
    plt.ylabel("$\mathrm{det}(A)$", fontsize=18)
    #plt.yscale("log")
    plt.tight_layout()
    plt.close(fig)
    return fig
pn.Column(
    eps_slider, 
    data_widget, 
    diag_dom_gershgorin
)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.14581791265729088, 'fraction_numerical': 0.010785661414824997, 'mean_word_length': 3.1193379790940767, 'pattern_counts': {'":': 2, '<': 26, '<?xml version=': 0, '>': 27, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '6', 'dir': 'jupyter-structured-clean-dedup', 'id': '252777', 'license': 'MIT', 'n_tokens_mistral': 3486, 'n_tokens_neox': 3303, 'n_words': 979, 'path': 'notebooks/gershgorin.ipynb', 'repo_name': 'manuflores/sandbox', 'size': '9462'}