eigen_II.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "a8d2740f",
   "metadata": {},
   "source": [
    "# The Perron-Frobenius Theorem\n",
    "\n",
    "\n",
    "<a id='index-0'></a>\n",
    "In addition to what’s in Anaconda, this lecture will need the following libraries:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "103a682c",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "!pip install quantecon"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ac590d87",
   "metadata": {},
   "source": [
    "In this lecture we will begin with the foundational concepts in spectral theory.\n",
    "\n",
    "Then we will explore the Perron-Frobenius theorem and connect it to applications in Markov chains and networks.\n",
    "\n",
    "We will use the following imports:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9ffd91a9",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from numpy.linalg import eig\n",
    "import scipy as sp\n",
    "import quantecon as qe"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e6a3b9ab",
   "metadata": {},
   "source": [
    "## Nonnegative matrices\n",
    "\n",
    "Often, in economics, the matrix that we are dealing with is nonnegative.\n",
    "\n",
    "Nonnegative matrices have several special and useful properties.\n",
    "\n",
    "In this section we will discuss some of them — in particular, the connection\n",
    "between nonnegativity and eigenvalues.\n",
    "\n",
    "An $ n \\times m $ matrix $ A $ is called **nonnegative** if every element of $ A $\n",
    "is nonnegative, i.e., $ a_{ij} \\geq 0 $ for every $ i,j $.\n",
    "\n",
    "We denote this as $ A \\geq 0 $.\n",
    "\n",
    "\n",
    "<a id='irreducible'></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1ec526e7",
   "metadata": {},
   "source": [
    "### Irreducible matrices\n",
    "\n",
    "We introduced irreducible matrices in the [Markov chain lecture](https://intro.quantecon.org/markov_chains_II.html#mc-irreducible).\n",
    "\n",
    "Here we generalize this concept:\n",
    "\n",
    "Let $ a^{k}_{ij} $ be element $ (i,j) $ of $ A^k $.\n",
    "\n",
    "An $ n \\times n $ nonnegative matrix $ A $ is called irreducible if $ A + A^2 + A^3 + \\cdots \\gg 0 $, where $ \\gg 0 $ indicates that every element in $ A $ is strictly positive.\n",
    "\n",
    "In other words, for each $ i,j $ with $ 1 \\leq i, j \\leq n $, there exists a $ k \\geq 0 $ such that $ a^{k}_{ij} > 0 $."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6ea55225",
   "metadata": {},
   "source": [
    "### \n",
    "\n",
    "Here are some examples to illustrate this further:\n",
    "\n",
    "$$\n",
    "A = \\begin{bmatrix} 0.5 & 0.1 \\\\ \n",
    "                    0.2 & 0.2 \n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "$ A $ is irreducible since $ a_{ij}>0 $ for all $ (i,j) $.\n",
    "\n",
    "$$\n",
    "B = \\begin{bmatrix} 0 & 1 \\\\ \n",
    "                    1 & 0 \n",
    "\\end{bmatrix}\n",
    ", \\quad\n",
    "B^2 = \\begin{bmatrix} 1 & 0 \\\\ \n",
    "                      0 & 1\n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "$ B $ is irreducible since $ B + B^2 $ is a matrix of ones.\n",
    "\n",
    "$$\n",
    "C = \\begin{bmatrix} 1 & 0 \\\\ \n",
    "                    0 & 1 \n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "$ C $ is not irreducible since $ C^k = C $ for all $ k \\geq 0 $ and thus\n",
    "$ c^{k}_{12},c^{k}_{21} = 0 $ for all $ k \\geq 0 $."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a25b60e5",
   "metadata": {},
   "source": [
    "### Left eigenvectors\n",
    "\n",
    "Recall that we previously discussed eigenvectors in [Eigenvalues and Eigenvectors](https://intro.quantecon.org/eigen_I.html#la-eigenvalues).\n",
    "\n",
    "In particular, $ \\lambda $ is an eigenvalue of $ A $ and $ v $ is an eigenvector of $ A $ if $ v $ is nonzero and satisfy\n",
    "\n",
    "$$\n",
    "Av = \\lambda v.\n",
    "$$\n",
    "\n",
    "In this section we introduce left eigenvectors.\n",
    "\n",
    "To avoid confusion, what we previously referred to as “eigenvectors” will be called “right eigenvectors”.\n",
    "\n",
    "Left eigenvectors will play important roles in what follows, including that of stochastic steady states for dynamic models under a Markov assumption.\n",
    "\n",
    "A vector $ w $ is called a left eigenvector of $ A $ if $ w $ is a right eigenvector of $ A^\\top $.\n",
    "\n",
    "In other words, if $ w $ is a left eigenvector of matrix $ A $, then $ A^\\top w = \\lambda w $, where $ \\lambda $ is the eigenvalue associated with the left eigenvector $ v $.\n",
    "\n",
    "This hints at how to compute left eigenvectors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d52d2eef",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "A = np.array([[3, 2],\n",
    "              [1, 4]])\n",
    "\n",
    "# Compute eigenvalues and right eigenvectors\n",
    "λ, v = eig(A)\n",
    "\n",
    "# Compute eigenvalues and left eigenvectors\n",
    "λ, w = eig(A.T)\n",
    "\n",
    "# Keep 5 decimals\n",
    "np.set_printoptions(precision=5)\n",
    "\n",
    "print(f\"The eigenvalues of A are:\\n {λ}\\n\")\n",
    "print(f\"The corresponding right eigenvectors are: \\n {v[:,0]} and {-v[:,1]}\\n\")\n",
    "print(f\"The corresponding left eigenvectors are: \\n {w[:,0]} and {-w[:,1]}\\n\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e6fa55e1",
   "metadata": {},
   "source": [
    "We can also use `scipy.linalg.eig` with argument `left=True` to find left eigenvectors directly"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5143059f",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "eigenvals, ε, e = sp.linalg.eig(A, left=True)\n",
    "\n",
    "print(f\"The eigenvalues of A are:\\n {eigenvals.real}\\n\")\n",
    "print(f\"The corresponding right eigenvectors are: \\n {e[:,0]} and {-e[:,1]}\\n\")\n",
    "print(f\"The corresponding left eigenvectors are: \\n {ε[:,0]} and {-ε[:,1]}\\n\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "948bd38c",
   "metadata": {},
   "source": [
    "The eigenvalues are the same while the eigenvectors themselves are different.\n",
    "\n",
    "(Also note that we are taking the nonnegative value of the eigenvector of [dominant eigenvalue](#perron-frobe), this is because `eig` automatically normalizes the eigenvectors.)\n",
    "\n",
    "We can then take transpose to obtain $ A^\\top w = \\lambda w $ and obtain $ w^\\top A= \\lambda w^\\top $.\n",
    "\n",
    "This is a more common expression and where the name left eigenvectors originates.\n",
    "\n",
    "\n",
    "<a id='perron-frobe'></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e9715d48",
   "metadata": {},
   "source": [
    "### The Perron-Frobenius theorem\n",
    "\n",
    "For a square nonnegative matrix $ A $, the behavior of $ A^k $ as $ k \\to \\infty $ is controlled by the eigenvalue with the largest\n",
    "absolute value, often called the **dominant eigenvalue**.\n",
    "\n",
    "For any such matrix $ A $, the Perron-Frobenius theorem characterizes certain\n",
    "properties of the dominant eigenvalue and its corresponding eigenvector."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cf1492dd",
   "metadata": {},
   "source": [
    "###  (Perron-Frobenius Theorem)\n",
    "\n",
    "If a matrix $ A \\geq 0 $ then,\n",
    "\n",
    "1. the dominant eigenvalue of $ A $, $ r(A) $, is real-valued and nonnegative.  \n",
    "1. for any other eigenvalue (possibly complex) $ \\lambda $ of $ A $, $ |\\lambda| \\leq r(A) $.  \n",
    "1. we can find a nonnegative and nonzero eigenvector $ v $ such that $ Av = r(A)v $.  \n",
    "\n",
    "\n",
    "Moreover if $ A $ is also irreducible then,\n",
    "\n",
    "1. the eigenvector $ v $ associated with the eigenvalue $ r(A) $ is strictly positive.  \n",
    "1. there exists no other positive eigenvector $ v $ (except scalar multiples of $ v $) associated with $ r(A) $.  \n",
    "\n",
    "\n",
    "(More of the Perron-Frobenius theorem about primitive matrices will be introduced [below](#prim-matrices).)\n",
    "\n",
    "(This is a relatively simple version of the theorem — for more details see\n",
    "[here](https://en.wikipedia.org/wiki/Perron%E2%80%93Frobenius_theorem)).\n",
    "\n",
    "We will see applications of the theorem below.\n",
    "\n",
    "Let’s build our intuition for the theorem using a simple example we have seen [before](https://intro.quantecon.org/markov_chains_I.html#mc-eg1).\n",
    "\n",
    "Now let’s consider examples for each case."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "195da87e",
   "metadata": {},
   "source": [
    "#### Example: irreducible matrix\n",
    "\n",
    "Consider the following irreducible matrix $ A $:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "04e958bf",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "A = np.array([[0, 1, 0],\n",
    "              [.5, 0, .5],\n",
    "              [0, 1, 0]])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5f0ed823",
   "metadata": {},
   "source": [
    "We can compute the dominant eigenvalue and the corresponding eigenvector"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d1ecfe52",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "eig(A)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b9ab2073",
   "metadata": {},
   "source": [
    "Now we can see the claims of the Perron-Frobenius theorem holds for the irreducible matrix $ A $:\n",
    "\n",
    "1. The dominant eigenvalue is real-valued and non-negative.  \n",
    "1. All other eigenvalues have absolute values less than or equal to the dominant eigenvalue.  \n",
    "1. A non-negative and nonzero eigenvector is associated with the dominant eigenvalue.  \n",
    "1. As the matrix is irreducible, the eigenvector associated with the dominant eigenvalue is strictly positive.  \n",
    "1. There exists no other positive eigenvector associated with the dominant eigenvalue.  \n",
    "\n",
    "\n",
    "\n",
    "<a id='prim-matrices'></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5761a4fd",
   "metadata": {},
   "source": [
    "### Primitive matrices\n",
    "\n",
    "We know that in real world situations it’s hard for a matrix to be everywhere positive (although they have nice properties).\n",
    "\n",
    "The primitive matrices, however, can still give us helpful properties with looser definitions.\n",
    "\n",
    "Let $ A $ be a square nonnegative matrix and let $ A^k $ be the $ k^{th} $ power of $ A $.\n",
    "\n",
    "A matrix is called **primitive** if there exists a $ k \\in \\mathbb{N} $ such that $ A^k $ is everywhere positive."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "818667ed",
   "metadata": {},
   "source": [
    "### \n",
    "\n",
    "Recall the examples given in irreducible matrices:\n",
    "\n",
    "$$\n",
    "A = \\begin{bmatrix} 0.5 & 0.1 \\\\ \n",
    "                    0.2 & 0.2 \n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "$ A $ here is also a primitive matrix since $ A^k $ is everywhere nonnegative for $ k \\in \\mathbb{N} $.\n",
    "\n",
    "$$\n",
    "B = \\begin{bmatrix} 0 & 1 \\\\ \n",
    "                    1 & 0 \n",
    "\\end{bmatrix}\n",
    ", \\quad\n",
    "B^2 = \\begin{bmatrix} 1 & 0 \\\\ \n",
    "                      0 & 1\n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "$ B $ is irreducible but not primitive since there are always zeros in either principal diagonal or secondary diagonal.\n",
    "\n",
    "We can see that if a matrix is primitive, then it implies the matrix is irreducible but not vice versa.\n",
    "\n",
    "Now let’s step back to the primitive matrices part of the Perron-Frobenius theorem"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7f9b95a6",
   "metadata": {},
   "source": [
    "###  (Continous of Perron-Frobenius Theorem)\n",
    "\n",
    "If $ A $ is primitive then,\n",
    "\n",
    "1. the inequality $ |\\lambda| \\leq r(A) $ is **strict** for all eigenvalues $ \\lambda $ of $ A $ distinct from $ r(A) $, and  \n",
    "1. with $ v $ and $ w $ normalized so that the inner product of $ w $ and  $ v = 1 $, we have\n",
    "  $ r(A)^{-m} A^m $ converges to $ v w^{\\top} $ when $ m \\rightarrow \\infty $. The matrix $ v w^{\\top} $ is called the **Perron projection** of $ A $.  "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "55683d25",
   "metadata": {},
   "source": [
    "#### Example 1: primitive matrix\n",
    "\n",
    "Consider the following primitive matrix $ B $:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "553edf4f",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "B = np.array([[0, 1, 1],\n",
    "              [1, 0, 1],\n",
    "              [1, 1, 0]])\n",
    "\n",
    "np.linalg.matrix_power(B, 2)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "47d8b815",
   "metadata": {},
   "source": [
    "We compute the dominant eigenvalue and the corresponding eigenvector"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4bc8f8b3",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "eig(B)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "92c394f8",
   "metadata": {},
   "source": [
    "Now let’s give some examples to see if the claims of the Perron-Frobenius theorem hold for the primitive matrix $ B $:\n",
    "\n",
    "1. The dominant eigenvalue is real-valued and non-negative.  \n",
    "1. All other eigenvalues have absolute values strictly less than the dominant eigenvalue.  \n",
    "1. A non-negative and nonzero eigenvector is associated with the dominant eigenvalue.  \n",
    "1. The eigenvector associated with the dominant eigenvalue is strictly positive.  \n",
    "1. There exists no other positive eigenvector associated with the dominant eigenvalue.  \n",
    "1. The inequality $ |\\lambda| < r(B) $ holds for all eigenvalues $ \\lambda $ of $ B $ distinct from the dominant eigenvalue.  \n",
    "\n",
    "\n",
    "Furthermore, we can verify the convergence property (7) of the theorem on the following examples:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "01379dee",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "def compute_perron_projection(M):\n",
    "\n",
    "    eigval, v = eig(M)\n",
    "    eigval, w = eig(M.T)\n",
    "\n",
    "    r = np.max(eigval)\n",
    "\n",
    "    # Find the index of the dominant (Perron) eigenvalue\n",
    "    i = np.argmax(eigval)\n",
    "\n",
    "    # Get the Perron eigenvectors\n",
    "    v_P = v[:, i].reshape(-1, 1)\n",
    "    w_P = w[:, i].reshape(-1, 1)\n",
    "\n",
    "    # Normalize the left and right eigenvectors\n",
    "    norm_factor = w_P.T @ v_P\n",
    "    v_norm = v_P / norm_factor\n",
    "\n",
    "    # Compute the Perron projection matrix\n",
    "    P = v_norm @ w_P.T\n",
    "    return P, r\n",
    "\n",
    "def check_convergence(M):\n",
    "    P, r = compute_perron_projection(M)\n",
    "    print(\"Perron projection:\")\n",
    "    print(P)\n",
    "\n",
    "    # Define a list of values for n\n",
    "    n_list = [1, 10, 100, 1000, 10000]\n",
    "\n",
    "    for n in n_list:\n",
    "\n",
    "        # Compute (A/r)^n\n",
    "        M_n = np.linalg.matrix_power(M/r, n)\n",
    "\n",
    "        # Compute the difference between A^n / r^n and the Perron projection\n",
    "        diff = np.abs(M_n - P)\n",
    "\n",
    "        # Calculate the norm of the difference matrix\n",
    "        diff_norm = np.linalg.norm(diff, 'fro')\n",
    "        print(f\"n = {n}, error = {diff_norm:.10f}\")\n",
    "\n",
    "\n",
    "A1 = np.array([[1, 2],\n",
    "               [1, 4]])\n",
    "\n",
    "A2 = np.array([[0, 1, 1],\n",
    "               [1, 0, 1],\n",
    "               [1, 1, 0]])\n",
    "\n",
    "A3 = np.array([[0.971, 0.029, 0.1, 1],\n",
    "               [0.145, 0.778, 0.077, 0.59],\n",
    "               [0.1, 0.508, 0.492, 1.12],\n",
    "               [0.2, 0.8, 0.71, 0.95]])\n",
    "\n",
    "for M in A1, A2, A3:\n",
    "    print(\"Matrix:\")\n",
    "    print(M)\n",
    "    check_convergence(M)\n",
    "    print()\n",
    "    print(\"-\"*36)\n",
    "    print()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "22c8c22c",
   "metadata": {},
   "source": [
    "The convergence is not observed in cases of non-primitive matrices.\n",
    "\n",
    "Let’s go through an example"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "96893e6d",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "B = np.array([[0, 1, 1],\n",
    "              [1, 0, 0],\n",
    "              [1, 0, 0]])\n",
    "\n",
    "# This shows that the matrix is not primitive\n",
    "print(\"Matrix:\")\n",
    "print(B)\n",
    "print(\"100th power of matrix B:\")\n",
    "print(np.linalg.matrix_power(B, 100))\n",
    "\n",
    "check_convergence(B)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9457a529",
   "metadata": {},
   "source": [
    "The result shows that the matrix is not primitive as it is not everywhere positive.\n",
    "\n",
    "These examples show how the Perron-Frobenius theorem relates to the eigenvalues and eigenvectors of positive matrices and the convergence of the power of matrices.\n",
    "\n",
    "In fact we have already seen the theorem in action before in [the Markov chain lecture](https://intro.quantecon.org/markov_chains_I.html#mc1_ex_1).\n",
    "\n",
    "\n",
    "<a id='spec-markov'></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4b9c4318",
   "metadata": {},
   "source": [
    "#### Example 2: connection to Markov chains\n",
    "\n",
    "We are now prepared to bridge the languages spoken in the two lectures.\n",
    "\n",
    "A primitive matrix is both irreducible and aperiodic.\n",
    "\n",
    "So Perron-Frobenius theorem explains why both [Imam and Temple matrix](https://intro.quantecon.org/markov_chains_I.html#mc-eg3) and [Hamilton matrix](https://en.wikipedia.org/wiki/Hamiltonian_matrix) converge to a stationary distribution, which is the Perron projection of the two matrices"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "28b08bcb",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "P = np.array([[0.68, 0.12, 0.20],\n",
    "              [0.50, 0.24, 0.26],\n",
    "              [0.36, 0.18, 0.46]])\n",
    "\n",
    "print(compute_perron_projection(P)[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e45f4f74",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "mc = qe.MarkovChain(P)\n",
    "ψ_star = mc.stationary_distributions[0]\n",
    "ψ_star"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e8950819",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "P_hamilton = np.array([[0.971, 0.029, 0.000],\n",
    "                       [0.145, 0.778, 0.077],\n",
    "                       [0.000, 0.508, 0.492]])\n",
    "\n",
    "print(compute_perron_projection(P_hamilton)[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c1aaecd1",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "mc = qe.MarkovChain(P_hamilton)\n",
    "ψ_star = mc.stationary_distributions[0]\n",
    "ψ_star"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "015ebd31",
   "metadata": {},
   "source": [
    "We can also verify other properties hinted by Perron-Frobenius in these stochastic matrices.\n",
    "\n",
    "Another example is the relationship between convergence gap and convergence rate.\n",
    "\n",
    "In the [exercise](https://intro.quantecon.org/markov_chains_I.html#mc1_ex_1), we stated that the convergence rate is determined by the spectral gap, the difference between the largest and the second largest eigenvalue.\n",
    "\n",
    "This can be proven using what we have learned here.\n",
    "\n",
    "Please note that we use $ \\mathbb{1} $ for a vector of ones in this lecture.\n",
    "\n",
    "With Markov model $ M $ with state space $ S $ and transition matrix $ P $, we can write $ P^t $ as\n",
    "\n",
    "$$\n",
    "P^t=\\sum_{i=1}^{n-1} \\lambda_i^t v_i w_i^{\\top}+\\mathbb{1} \\psi^*,\n",
    "$$\n",
    "\n",
    "This is proven in [[Sargent and Stachurski, 2023](https://intro.quantecon.org/zreferences.html#id24)] and a nice discussion can be found [here](https://math.stackexchange.com/questions/2433997/can-all-matrices-be-decomposed-as-product-of-right-and-left-eigenvector).\n",
    "\n",
    "In this formula $ \\lambda_i $ is an eigenvalue of $ P $ with corresponding right and left eigenvectors $ v_i $ and $ w_i $ .\n",
    "\n",
    "Premultiplying $ P^t $ by arbitrary $ \\psi \\in \\mathscr{D}(S) $ and rearranging now gives\n",
    "\n",
    "$$\n",
    "\\psi P^t-\\psi^*=\\sum_{i=1}^{n-1} \\lambda_i^t \\psi v_i w_i^{\\top}\n",
    "$$\n",
    "\n",
    "Recall that eigenvalues are ordered from smallest to largest from $ i = 1 ... n $.\n",
    "\n",
    "As we have seen, the largest eigenvalue for a primitive stochastic matrix is one.\n",
    "\n",
    "This can be proven using [Gershgorin Circle Theorem](https://en.wikipedia.org/wiki/Gershgorin_circle_theorem),\n",
    "but it is out of the scope of this lecture.\n",
    "\n",
    "So by the statement (6) of Perron-Frobenius theorem, $ \\lambda_i<1 $ for all $ i<n $, and $ \\lambda_n=1 $ when $ P $ is primitive.\n",
    "\n",
    "Hence, after taking the Euclidean norm deviation, we obtain\n",
    "\n",
    "$$\n",
    "\\left\\|\\psi P^t-\\psi^*\\right\\|=O\\left(\\eta^t\\right) \\quad \\text { where } \\quad \\eta:=\\left|\\lambda_{n-1}\\right|<1\n",
    "$$\n",
    "\n",
    "Thus, the rate of convergence is governed by the modulus of the second largest eigenvalue."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "eb28aa23",
   "metadata": {},
   "source": [
    "## Exercises"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "163c3c5c",
   "metadata": {},
   "source": [
    "##  (Leontief’s Input-Output Model)Exercise 38.1\n",
    "\n",
    "[Wassily Leontief](https://en.wikipedia.org/wiki/Wassily_Leontief) developed a model of an economy with $ n $ sectors producing $ n $ different commodities representing the interdependencies of different sectors of an economy.\n",
    "\n",
    "Under this model some of the output is consumed internally by the industries and the rest is consumed by external consumers.\n",
    "\n",
    "We define a simple model with 3 sectors - agriculture, industry, and service.\n",
    "\n",
    "The following table describes how output is distributed within the economy:\n",
    "\n",
    "||Total output|Agriculture|Industry|Service|Consumer|\n",
    "|:--------------:|:--------------:|:--------------:|:--------------:|:--------------:|:--------------:|\n",
    "|Agriculture|$ x_1 $|0.3$ x_1 $|0.2$ x_2 $|0.3$ x_3 $|4|\n",
    "|Industry|$ x_2 $|0.2$ x_1 $|0.4$ x_2 $|0.3$ x_3 $|5|\n",
    "|Service|$ x_3 $|0.2$ x_1 $|0.5$ x_2 $|0.1$ x_3 $|12|\n",
    "The first row depicts how agriculture’s total output $ x_1 $ is distributed\n",
    "\n",
    "- $ 0.3x_1 $ is used as inputs within agriculture itself,  \n",
    "- $ 0.2x_2 $ is used as inputs by the industry sector to produce $ x_2 $ units,  \n",
    "- $ 0.3x_3 $ is used as inputs by the service sector to produce $ x_3 $ units and  \n",
    "- 4 units is the external demand by consumers.  \n",
    "\n",
    "\n",
    "We can transform this into a system of linear equations for the 3 sectors as\n",
    "given below:\n",
    "\n",
    "$$\n",
    "x_1 = 0.3x_1 + 0.2x_2 + 0.3x_3 + 4 \\\\\n",
    "    x_2 = 0.2x_1 + 0.4x_2 + 0.3x_3 + 5 \\\\\n",
    "    x_3 = 0.2x_1 + 0.5x_2 + 0.1x_3 + 12\n",
    "$$\n",
    "\n",
    "This can be transformed into the matrix equation $ x = Ax + d $ where\n",
    "\n",
    "$$\n",
    "x =\n",
    "\\begin{bmatrix}\n",
    "    x_1 \\\\\n",
    "    x_2 \\\\\n",
    "    x_3\n",
    "\\end{bmatrix}\n",
    ", \\; A =\n",
    "\\begin{bmatrix}\n",
    "    0.3 & 0.2 & 0.3 \\\\\n",
    "    0.2 & 0.4 & 0.3 \\\\\n",
    "    0.2 & 0.5 & 0.1\n",
    "\\end{bmatrix}\n",
    "\\; \\text{and} \\;\n",
    "d =\n",
    "\\begin{bmatrix}\n",
    "    4 \\\\\n",
    "    5 \\\\\n",
    "    12\n",
    "\\end{bmatrix}\n",
    "$$\n",
    "\n",
    "The solution $ x^{*} $ is given by the equation $ x^{*} = (I-A)^{-1} d $\n",
    "\n",
    "1. Since $ A $ is a nonnegative irreducible matrix, find the Perron-Frobenius eigenvalue of $ A $.  \n",
    "1. Use the [Neumann Series Lemma](https://intro.quantecon.org/eigen_I.html#la-neumann) to find the solution $ x^{*} $ if it exists.  "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b5502512",
   "metadata": {},
   "source": [
    "## Solution to[ Exercise 38.1 (Leontief’s Input-Output Model)](https://intro.quantecon.org/#eig_ex1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "35e95bd2",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "A = np.array([[0.3, 0.2, 0.3],\n",
    "              [0.2, 0.4, 0.3],\n",
    "              [0.2, 0.5, 0.1]])\n",
    "\n",
    "evals, evecs = eig(A)\n",
    "\n",
    "r = max(abs(λ) for λ in evals)   #dominant eigenvalue/spectral radius\n",
    "print(r)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "090879bd",
   "metadata": {},
   "source": [
    "Since we have $ r(A) < 1 $ we can thus find the solution using the Neumann Series Lemma."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6e19191d",
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "I = np.identity(3)\n",
    "B = I - A\n",
    "\n",
    "d = np.array([4, 5, 12])\n",
    "d.shape = (3,1)\n",
    "\n",
    "B_inv = np.linalg.inv(B)\n",
    "x_star = B_inv @ d\n",
    "print(x_star)"
   ]
  }
 ],
 "metadata": {
  "date": 1733269015.4057481,
  "filename": "eigen_II.md",
  "kernelspec": {
   "display_name": "Python",
   "language": "python3",
   "name": "python3"
  },
  "title": "The Perron-Frobenius Theorem"
 },
 "nbformat": 4,
 "nbformat_minor": 5
}