using Symbolics

, LinearAlgebra

, PlutoUI

, ForwardDiff

, Zygote

1.9 s

TableOfContents(title="Differentiating the Determinant and the Inverse", indent=true, depth=4, aside=true)

19.1 μs

Math: The gradient of the determinant

md"""

# Math: The gradient of the determinant

"""

152 μs

Theorem: $\nabla (det A) = cofactor (A) = (det A) A^{- T} = adj (A^{T})$
$d (det A) = tr (det (A) A^{- 1} d A) = tr(adj (A) d A) = tr(cofactor (A)^{T} d A)$
Note:the adjugate is defined as $(det A) A^{- 1}$ ; it is the transpose of the cofactor.

md"""

Theorem:

``\nabla(\det A) = \text{cofactor}(A) = (\det A)A^{-T} = \text{adj}(A^T)`` $(br)

``d(\det A) = \text{tr}( \det(A)A^{-1}dA) = \text{tr(adj}(A) dA) = \text{tr(cofactor}(A)^T dA)`` $(br)

Note:the adjugate is defined as ``(\det A)A^{-1}``; it is the transpose of the cofactor.

"""

228 μs

adj (generic function with 1 method)

adj(A) = det(A) * inv(A) # adjugate function

246 μs

cofactor (generic function with 1 method)

cofactor(A) = adj(A)'

231 μs

Formulas at a glance:
$A^{- 1}$ = adj( $A$ )/det( $A$ ) = cofactor( $A$ ) $^{T}$ /det( $A$ )
adj( $A$ ) = det( $A$ ) $A^{- 1}$ = cofactor( $A$ ) $^{T}$
cofactor( $A$ ) = $A^{- T}$ /det( $A$ ) = adj( $A$ ) $^{T}$

md"""

**Formulas at a glance:** $(br)

``A^{-1}`` = adj(``A``)/det(``A``) = cofactor(``A``)``^T``/det(``A``) $(br)

adj(``A``) = det(``A``)``A^{-1}`` = cofactor(``A``)``^T`` $(br)

cofactor(``A``) = ``A^{-T}``/det(``A``) = adj(``A``)``^{T}``

"""

352 μs

The adjugate is the transpose of the cofactor matrix. You may remember that the cofactor matrix has as (i,j) entry $(- 1)^{i + j}$ times the determinant obtained by deleting row i and column j.

md"""

The adjugate is the transpose of the cofactor matrix.

You may remember that the cofactor matrix has as (i,j) entry ``(-1)^{i+j}`` times the determinant obtained by deleting row i and column j.

"""

164 μs

Symbolics.Num

$a$

$b$

$c$

$d$

@variables a b c d

868 ns

$[\begin{array}{cc} a & c \\ b & d \end{array}]$

M = [ a c;b d]

1.1 μs

$[\begin{array}{cc} d & - b \\ - c & a \end{array}]$

simplify.(cofactor(M))

8.5 ms

$[\begin{array}{cc} d & - c \\ - b & a \end{array}]$

simplify.(adj( M))

7.8 ms

$[\begin{array}{cc} \frac{d}{a d - b c} & \frac{- c}{a d - b c} \\ \frac{- b}{a d - b c} & \frac{a}{a d - b c} \end{array}]$

simplify.(inv(M))

5.2 ms

Numerical Demonstration

md"""

## Numerical Demonstration

"""

5.3 ms

5×5 Matrix{Float64}:
 6.23738e-6  9.80317e-6  9.49748e-6  2.16775e-6  7.13224e-6
 4.29949e-7  1.04401e-6  9.43453e-6  3.49108e-6  3.30058e-6
 2.43275e-6  6.93099e-6  9.77647e-6  4.73139e-6  1.52485e-6
 4.08123e-6  5.79363e-6  5.63512e-6  9.05209e-6  5.3378e-6
 3.20542e-7  9.13855e-6  1.20064e-6  2.504e-6    3.50599e-6

A = rand(5,5); dA = rand(5,5) * .00001

1.3 μs

-1.88723e-6

-1.88722e-6

det(A+dA) - det(A), tr( adj(A) * dA )

13.6 μs

Forward Diff Autodiff

md"""

## Forward Diff Autodiff

"""

137 μs

5×5 Matrix{Float64}:
  0.227566    -0.231878    0.00127616   0.0537115   -0.113963
 -0.441303     0.453903    0.109316    -0.188827     0.132947
  0.0485198   -0.0217404  -0.0335765    0.068198    -0.0835962
 -0.00220219  -0.0833639   0.0166832    0.0266121   -0.0155698
  0.124369    -0.124389   -0.105723     0.00347312   0.0503901

ForwardDiff.gradient(A->det(A), A)

11.4 μs

5×5 Matrix{Float64}:
  0.227566    -0.231878    0.00127616   0.0537115   -0.113963
 -0.441303     0.453903    0.109316    -0.188827     0.132947
  0.0485198   -0.0217404  -0.0335765    0.068198    -0.0835962
 -0.00220219  -0.0833639   0.0166832    0.0266121   -0.0155698
  0.124369    -0.124389   -0.105723     0.00347312   0.0503901

adj(A')

6.9 μs

Reverse Mode Autodiff

md"""

## Reverse Mode Autodiff

"""

135 μs

Zygote does reverse ad

md"""

Zygote does reverse ad

"""

163 μs

5×5 Matrix{Float64}:
  0.227566    -0.231878    0.00127616   0.0537115   -0.113963
 -0.441303     0.453903    0.109316    -0.188827     0.132947
  0.0485198   -0.0217404  -0.0335765    0.068198    -0.0835962
 -0.00220219  -0.0833639   0.0166832    0.0266121   -0.0155698
  0.124369    -0.124389   -0.105723     0.00347312   0.0503901

Zygote.gradient(A->det(A), A )

21.9 μs

Symbolic Demonstration

md"""

## Symbolic Demonstration

"""

134 μs

$[\begin{array}{cc} d & - b \\ - c & a \end{array}]$

reshape(Symbolics.gradient(det(M),vec(M)), 2, 2)

204 μs

$[\begin{array}{cc} d & - b \\ - c & a \end{array}]$

simplify.(adj(M'))

371 ms

Direct Proof

A direct proof where you just differentiate the scalar with respect to every input can be obtained as a simple consequence from the cofactor expansion a.k.a. the Laplace expansion of the determinant based on the $i$ th row.

$det (A) = A_{i 1} C_{i 1} + A_{i 2} C_{i 2} + \dots A_{i n} C_{i n}$

Recall that the determinant is a linear (affine) function of any element with slope equal to the cofactor.

We then readily obtain $\frac{\partial det (A)}{\partial A_{i j}} = C_{i j}$ from which we conclude $\nabla (det A) = C$ .

Example:

md"""

# Direct Proof

A direct proof where you just differentiate the scalar with respect to every input

can be obtained as a simple consequence from the cofactor expansion a.k.a. the [Laplace expansion](https://en.wikipedia.org/wiki/Laplace_expansion) of the determinant based on the ``i``th row.

``\det(A) = A_{i1} C_{i1} + A_{i2}C_{i2} + \ldots A_{in}C_{in}``

Recall that the determinant is a linear (affine) function of any element with slope equal to the cofactor.

We then readily obtain

``\frac{\partial \det(A)}{\partial A_{ij}}=C_{ij}`` from which we conclude

``\nabla (\det A ) = C``.

$(br) **Example:**

"""

463 μs

$[\begin{array}{ccc} 4 & 1 & 5 \\ 3 & 4 & 5 \\ 6 & 7 & a \end{array}]$

begin

i,j = 3,3

Z= [Num(4) 1 5; 3 4 5 ; 6 7 8]

Z[i,j] = a

end

31.3 μs

$- 125 + 13 a$

det(Z)

3.7 ms

$13.0$

simplify( cofactor(Z)[i,j] )

1.8 s

Fancy Proof

Figure out linearization near the identity I

det(I+dA) = trace(dA) (think of the n! terms in the determinant and drop higher terms)

md"""

# Fancy Proof

Figure out linearization near the identity I $(br)

det(I+dA) = trace(dA) (think of the n! terms in the determinant and drop higher terms) $(br)

"""

257 μs

$det (A + A (A^{- 1} d A)) = det (A) det (I + A^{- 1} d A) = det (A) tr (A^{- 1} d A) = tr (det (A) A^{- 1} d A)$

md"""

``\det(A+ A(A^{-1}dA)) = \det(A)\det(I+A^{-1}dA) =

\det(A) \text{tr}( A^{-1}dA) =

\text{tr}(\det(A) A^{-1}dA)``

"""

159 μs

Application to derivative of the characteristic polynomial

md"""

## Application to derivative of the characteristic polynomial

"""

163 μs

Direct derivative (freshman calculus):

$\frac{d}{d x} \prod_{i} (x - λ_{i}) = \sum_{i} \prod_{j \neq i} (x - λ_{j}) = \prod_{i} (x - λ_{i}) {\sum_{i} (x - λ_{i})^{- 1}}$

Perfectly good simple proof, but if you want to show off...

md"""

### Direct derivative (freshman calculus):

``\frac{d}{dx} \prod_i (x-\lambda_i) \ = \ \sum_i \prod_{j\ne i} (x-\lambda _j) =

\prod_i (x-\lambda_i) \left\{ \sum_i (x-\lambda_i)^{-1} \right\}``

$(br)

Perfectly good simple proof, but if you want to show off...

"""

304 μs

With our new technology:

$d (det (x I - A)) = det (x I - A) tr ((x I - A)^{- 1} d (x I - A))$
$= det (x I - A) tr (x I - A)^{- 1} d x$

Note: $d (x I - A) = d x I$ when $M$ is constant and
$tr (M d x) = tr (M) d x$ since $d x$ is a scalar.

md"""

### With our new technology:

``d(\det (xI-A))= \det(xI-A) \text{tr}( (xI-A)^{-1} d(xI-A) ) ``$(br)

`` \hspace{1.2in}= \det(xI-A) \text{tr}(xI-A)^{-1} dx``

Note: ``d(xI-A)=dxI`` when ``M`` is constant and $(br) ``\text{tr}(Mdx) = \text{tr}(M)dx`` since ``dx`` is a scalar.

"""

357 μs

Check:

md"""

### Check:

"""

182 μs

-0.649343

-0.649343

begin

x = rand()

ForwardDiff.derivative( x-> det(x*I-A), x), # Automatic Differentiation

det(x*I-A)*tr(inv(x*I-A)) # Analytic Answer

end

13.6 μs

Application d(log(det(A)))

$d (\log (det (A))) = det (A^{- 1}) d (det (A)) = tr (A^{- 1}) d A$

The logarithmic derivative shows up a lot in applied mathematics. For example the key term in Newton's method: $f (x) / f^{'} (x)$ may be written $1 / {\log f (x)}^{'}$

md"""

## Application d(log(det(A)))

``d(\log(\det(A))) = \det(A^{-1}) d(\det(A)) = \text{tr}( A^{-1} ) dA`` $(br)

The logarithmic derivative shows up a lot in applied mathematics. For example the key term in Newton's method: ``f(x)/f'(x)`` may be written ``1/ \left\{ \log f(x) \right\}'``

"""

299 μs

Math: The Jacobian of the Inverse

$A^{- 1} A = I \to d (A^{- 1} A) = 0 = d (A^{- 1}) A + A^{- 1} d A$ from the product rule
$\to d (A^{- 1}) = - A^{- 1} (d A) A^{- 1} = - (A^{- T} \otimes A^{- 1}) d A$

md"""

# Math: The Jacobian of the Inverse

``A^{-1}A= I \rightarrow d(A^{-1}A)=0=d(A^{-1})A + A^{-1}dA`` from the product rule $(br)

``\rightarrow d(A^{-1}) = -A^{-1} (dA) A^{-1} = -(A^{-T} \otimes A^{-1}) dA``

"""

138 μs

true

let

A = rand(3,3)

ForwardDiff.jacobian(A->inv(A),A) ≈-kron(inv(A'),inv(A))

end

754 ms