Massive relaxation

We now consider methods to construct a positive definite or at least invertible learning matrix. For example, far from a minimum the Hessian ${\bf H}$ may not be positive definite and like a differential operator ${{\bf K}}$ with zero modes, not even invertible. Massive relaxation can transform a non-invertible or not positive definite operator ${\bf A}_0$, e.g. ${\bf A}_0 = {{\bf K}}$ or ${\bf A}_0 = - {\bf H}$, into an invertible or positive definite operators:

$${\bf A} = {\bf A}_0 + m^2 {\bf I} : {\rm Massive \,\, relaxation}$$ (638)

A generalization would be to allow $m = m(x,y)$. This is, for example, used in some realizations of Newton`s method for minimization in regions where ${\bf H}$ is not positive definite and a diagonal operator is added to $-{\bf H}$, using for example a modified Cholesky factorization [19]. The mass term removes the zero modes of ${{\bf K}}$ if $-m^2$ is not in the spectrum of ${\bf A}_0$ and makes it positive definite if $m^2$ is larger than the smallest eigenvalue of ${\bf A}_0$. Matrix elements $(\, \phi,\, ({\bf A}_0 - z {\bf I})^{-1} \,\phi\,)$ of the resolvent ${\bf A}^{-1} (z)$, $z=-m^2$ representing in this case a complex variable, have poles at discrete eigenvalues of ${\bf A}_0$ and a cut at the continuous spectrum as long as $\phi$ has a non-zero overlap with the corresponding eigenfunctions. Instead of multiples of the identity, also other operators may be added to remove zero modes. The Hessian ${\bf H}_L$ in (156), for example, adds a $x$-dependent mass-like, but not necessarily positive definite term to ${{\bf K}}$. Similarly, for example ${\bf H}_P$ in (182) has $(x,y)$-dependent mass ${\bf P}^{-2} {\bf N}$ restricted to data points.

While full relaxation is the massless limit $m^2\rightarrow 0$ of massive relaxation, a gradient algorithm with $\eta^\prime$ can be obtained as infinite mass limit $m^2\rightarrow \infty$ with $\eta \rightarrow \infty$ and $m^2/\eta = 1/\eta^\prime$.

Constant functions are typical zero modes, i.e., eigenfunctions with zero eigenvalue, for differential operators with periodic boundary conditions. For instance for a common smoothness term $-\Delta$ (kinetic energy operator) as regularizing operator ${{\bf K}}$ the inverse of ${\bf A} = {{\bf K}} +m^2 {\bf I}$ has the form

$${\bf A}^{-1} (x^\prime , y^\prime ; x,y ) = \, \, \frac{1}{-\Delta + m^2}.$$ (639)

$$=\int_{-\infty}^\infty \! \frac{d^{d_X}\!k_x \, d^{d_Y}\!k_y... ...i k_x (x-x^\prime) + ik_y (y-y^\prime )} }{k_x^2 + k_y^2+m^2},$$ (640)

with $d=d_X+d_Y$, $d_X$ = dim($X$), $d_Y$ = dim($Y$). This Green`s function or matrix element of the resolvent kernel for ${\bf A}_0$ is analogous to the (Euclidean) propagator of a free scalar field with mass $m$, which is its two-point correlation function or matrix element of the covariance operator. According to $1/x = \int_0^\infty \! dt\, e^{-xt}$ the denominator can be brought into the exponent by introducing an additional integral. Performing the resulting Gaussian integration over $k = (k_x,k_y)$ the inverse can be expressed as

$${\bf A}^{-1} (x^\prime , y^\prime ; x,y ; m) = m^{d-2} {\bf A}^{-1} (m (x-x^\prime) , m (y-y^\prime) ; 1)$$

$$=(2 \pi)^{-d/2} \left( \frac{m}{\vert x-x^\prime\vert+\vert ... ...{(d-2)/2}( m \vert x-x^\prime\vert+ m \vert y-y^\prime\vert ),$$ (641)

in terms of the modified Bessel functions $K_\nu (x)$ which have the following integral representation

$$K_\nu (2 \sqrt{\beta \gamma}) = \frac{1}{2} \left(\frac{\gam... ... \int_0^\infty \! dt\, t^{\nu-1} e^{\frac{\beta}{t}-\gamma t}.$$ (642)

Alternatively, the same result can be obtained by switching to $d$-dimensional spherical coordinates, expanding the exponential in ultra-spheric harmonic functions and performing the integration over the angle-variables [117]. For the example $d=1$ this corresponds to Parzen´s kernel used in density estimation or for $d=3$

$${\bf A}^{-1} (x^\prime , y^\prime ; x,y ) = \frac{1}{4 \pi ... ...ert} e^{- m \vert x-x^\prime\vert - m \vert y-y^\prime\vert }.$$ (643)

The Green's function for periodic, Neumann, or Dirichlet boundary conditions can be expressed by sums over ${\bf A}^{-1} (x^\prime , y^\prime ; x,y )$ [77].

The lattice version of the Laplacian with lattice spacing $a$ reads

$$\hat \Delta f(n) = \frac{1}{a^2} \sum_j^d [ f(n-a_j)-2 f(n) + f(n+a_j) ],$$ (644)

writing $a_j$ for a vector in direction $j$ and length $a$. Including a mass term we get as lattice approximation for ${\bf A}$

$$\hat{\bf A} (n_x,n_y ; m_x,m_y ) = -\frac{1}{a^2}\sum_{i=1}^... ...{n_x+a_i^x,m_x} -2 \delta_{n_x,m_x} + \delta_{n_x-a^x_i,m_x})$$

$$-\frac{1}{a^2}\sum_{j=1}^{d_Y} \delta_{n_x,m_x} (\delta_{n_... ...\delta_{n_y-a^y_j,m_y}) +m^2 \delta_{n_x,m_x} \delta_{n_y,m_y}$$ (645)

Inserting the Fourier representation (101) of $\delta (x)$ gives

$$\hat{\bf A} (n_x,n_y ; m_x,m_y ) = \frac{2d}{a^2} \int_{-\p... ...\, d^{d_Y}\!k_y}{(2\pi)^d} e^{ik_x (n_x-m_x) + ik_y (n_y-m_y)}$$

$$\times \left( 1 + \frac{m^2 a^2}{2d} -\frac{1}{d} \sum_{i=1}... ...s k_{x,i} -\frac{1}{d} \sum_{j=1}^{d_Y} \cos k_{y,j} \right),$$ (646)

with $k_{x,i}$ = $k_x a^x_i$, $\cos k_{y,j}$ = $\cos k_y a^y_j$ and inverse

$$\hat {\bf A}^{-1} (n_x , n_y ; m_x,m_y) = \int_{-\pi}^\pi \!... ...t {\bf A}^{-1} (k_x , k_y) e^{ik_x (n_x-m_x) + ik_y (n_y-m_y)}$$

$$=\frac{a^2}{2d} \int_{-\pi}^\pi \!\! \frac{d^{d_X}\!k_x \, ... ...os k_{x,i} \!-\! \frac{1}{d} \sum_{j=1}^{d_Y} \cos k_{y,j} }.$$ (647)

(For $m=0$ and $d\le 2$ the integrand diverges for $k\rightarrow 0$ (infrared divergence). Subtracting formally the also infinite $\hat {\bf A}^{-1} (0,0 ; 0,0)$ results in finite difference. For example in $d=1$ one finds $\hat {\bf A}^{-1} (n_y ; m_y) - \hat {\bf A}^{-1} (0;0)$ = $-(1/2) \vert n_y-m_y \vert$ [103]. Using $1/x = \int_0^\infty \! dt\, e^{-xt}$ one obtains [195]

$$\hat {\bf A}^{-1} (k_x , k_y) = \frac{1}{2} \int_0^\infty \!... ...um_i^{d_X} \cos k_{x,i} + \sum_j^{d_Y} \cos k_{y,j} \right) },$$ (648)

with $\mu = d/a^2 + m^2/2$. This allows to express the inverse $\hat {\bf A}^{-1}$ in terms of the modified Bessel functions $I_\nu (n)$ which have for integer argument $n$ the integral representation

$$I_\nu (n) = \frac{1}{\pi} \int_0^\pi \!d\Theta\, e^{n \cos \Theta}\cos (\nu \Theta).$$ (649)

One finds

$$\hat {\bf A}^{-1} (n_x , n_y ; m_x,m_y) = \frac{1}{2} \int_0... ..._{j =1}^{d_Y} K_{\vert m_{y,j} - m^\prime_{y,j}\vert} (t/a^2).$$ (650)

It might be interesting to remark that the matrix elements of the inverse learning matrix or free massive propagator on the lattice $\hat {\bf A}^{-1}(x^\prime,y^\prime;x,y)$ can be given an interpretation in terms of (random) walks connecting the two points $(x^\prime, y^\prime)$ and $(x,y)$ [56,195]. For that purpose the lattice Laplacian is splitted into a diagonal and a nearest neighbor part

$$- \hat \Delta = \frac{1}{a^2} \left( 2d {\bf I} - {\bf W} \right),$$ (651)

where the nearest neighbor matrix ${\bf W}$ has matrix elements equal one for nearest neighbors and equal to zero otherwise. Thus,

$$\left( - \hat \Delta + m^2\right)^{-1} = \frac{1}{2 \mu} \le... ...um_{n=0}^\infty \left(\frac{1}{2 \mu a^2}\right)^n {\bf W}^n ,$$ (652)

can be written as geometric series. The matrix elements ${\bf W}^n (x^\prime,y^\prime;x,y)$ give the number of walks $w[(x^\prime,y^\prime)\rightarrow (x,y)]$ of length $\vert w\vert$ = $n$ connecting the two points $(x^\prime, y^\prime)$ and $(x,y)$. Thus, one can write

$$\left( - \hat \Delta + m^2\right)^{-1}(x^\prime,y^\prime;x,y... ...rrow (x,y)]} \left(\frac{1}{2 \mu a^2}\right)^{\vert w\vert} .$$ (653)