\chapter{Supplementary: Recurrent Predictive Coding for Associative Memory}\label{app:chap3}

\section{Expression for retrieved patterns in the rPCNs (Eq~\ref{eq:retrieval})}
\label{S1_Appendix}
Here we provide the details of derivation of Eq.~\ref{eq:retrieval}, showing that all three models perform the same retrieval at convergence. At convergence, the retrieval dynamics of the explicit PCN satisfy:
\begin{eqnarray}
    \bfS^{-1}(\mathbf{x} - \bar{\mathbf{x}}) = 0
\end{eqnarray}

This equation can be written into its block matrix form:
\begin{eqnarray}
    \begin{bmatrix}\bfS_{kk} & \bfS_{km} \\ \bfS_{mk} & \bfS_{mm}\end{bmatrix}^{-1}\begin{bmatrix}\mathbf{x}_{k} - \bar{\mathbf{x}}_k\\
\mathbf{x}_{m} - \bar{\mathbf{x}}_m
\end{bmatrix} = 0
\end{eqnarray}

The inverse of $\bfS$ can thus be written as:
\begin{eqnarray}
    \begin{bmatrix}\bfS_{kk} & \bfS_{km} \\ \bfS_{mk} & \bfS_{mm}\end{bmatrix}^{-1} = \begin{bmatrix}\bfQ^{-1} & -\bfQ^{-1}\bfS_{km}\bfS_{mm}^{-1}\\
-\bfS_{mm}^{-1}\bfS_{km}^{\top}\bfQ^{-1} & \bfS_{mm}^{-1}+\bfS_{mm}^{-1}\bfS_{km}^{\top}\bfQ^{-1}\bfS_{km}\bfS_{mm}^{-1}
\end{bmatrix}
\end{eqnarray}
\noindent where $\bfQ$ is called the Schur complement of $\bfS_{mm}$ in $\bfS$ \citep{haynsworth1968schur, boyd2004convex} and $\bfQ=\bfS_{kk}-\bfS_{km}\bfS_{mm}^{-1}\bfS_{km}^{\top}$. Since only $\mathbf{x}_k$ is relaxed during retrieval, we get:
\begin{eqnarray}
    \bfQ^{-1}(\mathbf{x}_{k} - \bar{\mathbf{x}}_k)-\bfQ^{-1}\bfS_{km}\bfS_{mm}^{-1}(\mathbf{x}_{m} - \bar{\mathbf{x}}_m) = 0
\end{eqnarray}
\noindent which immediately gives us the retrieval dynamics in Eq.~\ref{eq:retrieval}.

Now we show that the retrieval of the implicit/dendritic models also follows Eq.~\ref{eq:retrieval} at the convergence of inference. Notice that the sufficient and necessary condition for the convergence is $\bm{\epsilon} = \mathbf{x}-\bfW\mathbf{x}-\bm{\nu} = 0$. Splitting it into blocks corresponding to $\mathbf{x}_k$ and $\mathbf{x}_m$ this becomes:
\begin{eqnarray}
    \begin{bmatrix}\bfW_{kk} & \bfW_{km}\\
    \bfW_{mk} & \bfW_{mm}\end{bmatrix}\begin{bmatrix}\mathbf{x}_{k}\\
    \mathbf{x}_{m}\end{bmatrix} + \begin{bmatrix}
    \bm{\nu}_k \\ \bm{\nu}_m
    \end{bmatrix}
    =\begin{bmatrix}\mathbf{x}_{k}\\
    \mathbf{x}_{m}
    \end{bmatrix}
\end{eqnarray}

Since we only relax the top $k$ corrupted entries in $\mathbf{x}$, we have:
\begin{eqnarray}\label{eq:implicit-retrieval-temp}
    \bfW_{kk}\mathbf{x}_{k}+\bfW_{km}\mathbf{x}_{m}+\bm{\nu}_k=\mathbf{x}_{k}\Rightarrow(\bfI_{kk}-\bfW_{kk})\mathbf{x}_{k}=\bfW_{km}\mathbf{x}_{m}+\bm{\nu}_k
\end{eqnarray}

We now investigate the parameter values in the above equation. Notice that we assumed the convergence of learning at the time of retrieval, which gives us Eqs.~\ref{eq:nu} and~\ref{eq:W}. By splitting Eq.~\ref{eq:nu} into blocks we have:
\begin{eqnarray}
    \bm{\nu}_k = (\bfI_{kk} - \bfW_{kk})\bar{\mathbf{x}}_k - \bfW_{km}\bar{\mathbf{x}}_m
\end{eqnarray}

Substituting the $\bm{\nu}_k$ in Eq.~\ref{eq:implicit-retrieval-temp} with this relationship we get:
\begin{eqnarray}\label{eq:implicit-retrieval-final}
    \bfW_{km}(\mathbf{x}_m - \bar{\mathbf{x}}_m) = (\bfI_{kk} - \bfW_{kk})(\mathbf{x}_k - \bar{\mathbf{x}}_k)
\end{eqnarray}

Notice that Eq.~\ref{eq:W} connects the value of $W$ and $S$, which can also be written into the block-matrix form:
\begin{eqnarray}
    \left(\begin{bmatrix}\bfW_{kk} & \bfW_{km}\\ \bfW_{mk} & \bfW_{mm}\end{bmatrix}\begin{bmatrix}\bfS_{kk} & \bfS_{km}\\ \bfS_{mk} & \bfS_{mm}\end{bmatrix}-\begin{bmatrix}\bfS_{kk} & \bfS_{km}\\ \bfS_{mk} & \bfS_{mm}\end{bmatrix}\right)\Bigr\rvert_{diag = 0} = 0
\end{eqnarray}

\noindent which gives us two useful relationships:
\begin{eqnarray}
    \bfW_{kk}\bfS_{km}+\bfW_{km}\bfS_{mm}=\bfS_{km}\Rightarrow \bfW_{km}=(\bfI_{kk}-\bfW_{kk})\bfS_{km}\bfS_{mm}^{-1}
\end{eqnarray}
\begin{eqnarray}\label{eq:digonal-WS}
    \left(\bfW_{kk}\bfS_{kk} + \bfW_{km}\bfS_{mk} - \bfS_{kk}\right)\rvert_{diag = 0} = 0
\end{eqnarray}

Substituting the expression of $\bfW_{km}$ in terms of $\bfW_{kk}$ into Eq.~\ref{eq:implicit-retrieval-final} we have:
\begin{eqnarray}\label{eq:implicit-equi-final}
    (\bfI_{kk}-\bfW_{kk})\bfS_{km}\bfS_{mm}^{-1}(\mathbf{x}_m - \bar{\mathbf{x}}_m) = (\bfI_{kk} - \bfW_{kk})(\mathbf{x}_k - \bar{\mathbf{x}}_k)
\end{eqnarray}

The above expression is already close to Eq.~\ref{eq:retrieval} which we seek to prove, i.e., we could obtain Eq.~\ref{eq:retrieval} by cancelling $\bfI_{kk}-\bfW_{kk}$ on both sides of Eq.~\ref{eq:implicit-equi-final}, hence we will now show that $\bfI_{kk}-\bfW_{kk}$ is invertible. To do it we first observe that according to Eq.~\ref{eq:digonal-WS}, $\bfW_{kk}\bfS_{kk} + \bfW_{km}\bfS_{mk} - \bfS_{kk}$ is a diagonal matrix, which we call $D$. Therefore we have:
\begin{eqnarray}
    \begin{aligned}
        \mathbf{D} &= (\bfI_{kk}-\bfW_{kk})\bfS_{kk} - \bfW_{km}\bfS_{mk} \\&= (\bfI_{kk}-\bfW_{kk})\bfS_{kk} - (\bfI_{kk}-\bfW_{kk})\bfS_{km}\bfS_{mm}^{-1}\bfS_{mk} \\&= (\bfI_{kk}-\bfW_{kk})(\bfS_{kk} - \bfS_{km}\bfS_{mm}^{-1}\bfS_{mk}) \\&= (\bfI_{kk}-\bfW_{kk})\bfQ
    \end{aligned}
\end{eqnarray}

Since $\bfQ$ is the Schur complement of $\bfS_{mm}$ in $\bfS$, a sample covariance matrix that we assume to be positive definite, it is also positive definite and thus invertible \cite{boyd2004convex}. Therefore $\bfI_{kk}-\bfW_{kk} = \mathbf{D}\bfQ^{-1}$. Notice that since $\bfW_{kk}$ has all its diagonal elements equal to 0, the matrix $\mathbf{D}\bfQ^{-1}$ must have $1$'s on its diagonal, which prevents the diagonal matrix $\mathbf{D}$ from having $0$'s on its diagonal (for $\left(\mathbf{D}\bfQ^{-1}\right)_{ii}=\mathbf{D}_{ii}\left(\bfQ^{-1}\right)_{ii}=1$, $\mathbf{D}_{ii}$ cannot be $0$). Therefore, $\mathbf{D}$ is also invertible, which makes $\bfI_{kk}-\bfW_{kk}$ an invertible matrix as well. This enables us to ``cancel out" the $\bfI_{kk}-\bfW_{kk}$ on both sides of Eq.~\ref{eq:implicit-equi-final} and establish the equivalence to the retrieval of explicit rPCN (Eq.~\ref{eq:retrieval}).

\section{Derivation of rPCNs and Hopfield Networks as metric learning}
\paragraph{Hopfield Networks}
Given a set of memorized patterns $\{\mathbf{x}(i)\}_{i=1}^N$ and a query $\mathbf{q}$, the energy function of a (vanilla) Hopfield Network can be written as:
\begin{eqnarray}\label{HN_energy}
    \mathcal{F}_{HN} = -\sum_{i = 1}^N \left(\mathbf{q}^{\top}\mathbf{x}{(i)}\right)^2 
\end{eqnarray}

If we denote the collection of all patterns as a matrix $X\in \mathbb{R}^{N\times d}$ where $d$ denotes the dimension of each pattern, this energy can also be rewritten in terms of the sample covariance matrix of patterns $\bfS$:
\begin{eqnarray}\label{eq:EHN}
    \mathcal{F}_{HN} = -\frac{1}{2}\mathbf{q}^{\top}\bfX^{\top}\bfX\mathbf{q} \propto -\frac{1}{2}\mathbf{q}^{\top}\bfS \mathbf{q} = -\frac{1}{2}\mathbf{q}^{\top}(\bfS^{\frac{1}{2}})^{\top}\bfS^{\frac{1}{2}} \mathbf{q}= -\frac{1}{2}\| \bfS^{\frac{1}{2}} \mathbf{q} \|^2_2
\end{eqnarray}
which follows the general expression in Eq.~\ref{eq:general_metric} with $c=-\frac{1}{2}$ and $\bfL=\bfS^{\frac{1}{2}}$.

\paragraph{Explicit rPCN} The parameters $\bm{\mu}$ and $\bm{\Sigma}$ of the energy of the explicit rPCN will converge to the sample statistic $\bar{\mathbf{x}}$ and $\bfS$ after learning. Thus, the energy given a query $\mathbf{q}$ would become:
\begin{eqnarray}
    \mathcal{F}_{\tn{exp}}=\frac{1}{2}\log|\bfS| + \frac{1}{2}
    (\mathbf{q}-\bar{\mathbf{x}})^{\top}\bfS^{-1}(\mathbf{q}-\bar{\mathbf{x}})
\end{eqnarray}

Assuming zero-mean memories ($\bar{\mathbf{x}}$=0) and dropping the constant term $\frac{1}{2}\log|\bfS|$ unrelated to $\mathbf{q}$ after training, this energy can be simplified as:
\begin{eqnarray}
    \mathcal{F}_{\tn{exp}} = \frac{1}{2}
    \mathbf{q}^{\top}\bfS^{-1}\mathbf{q} = \frac{1}{2}
    \mathbf{q}^{\top}(\bfS^{-\frac{1}{2}})^{\top}\bfS^{-\frac{1}{2}}\mathbf{q} = \frac{1}{2}\| \bfS^{-\frac{1}{2}} \mathbf{q} \|^2_2
\end{eqnarray}
which follows the general expression in Eq.~\ref{eq:general_metric} with $c=\frac{1}{2}$ and $\bfL=\bfS^{-\frac{1}{2}}$.

\paragraph{Implicit rPCN} To start, note that the training phase of implicit rPCN can be seen as a constrained optimization problem. By Eq.~\ref{eq:implicit-obj}, assuming the entire set of patterns $\bfX$ and zero bias ($\bm{\nu} = 0$), we have
\begin{eqnarray} \label{constrained-optimation-problem}
    \min_{\bfW} \frac{1}{2}\|\bfX-\bfX\bfW\|_F^2 \quad \text{  s.t. } \quad \textnormal{diag}(\bfW) = \mathbf{0}
\end{eqnarray}
where for simplicity of notation, $\bfW$ in the appendix is the transpose of $W$ used in the main text. Then, we can equivalently write the constraints into the Lagrangian
\begin{eqnarray}
\mathcal{F}_{\tn{imp}} = \frac{1}{2} \|\bfX-\bfX \bfW\|_F^2 + \bm{\lambda}^{\top} \textnormal{diag}(\bfW)
\end{eqnarray}
where $\bm{\lambda} = (\lambda_1,...,\lambda_d)$ is a vector of Lagrangian multipliers. Taking gradient with respect to $\bfW$ yields that
\begin{eqnarray}
\begin{aligned}
    \frac{\partial \mathcal{F}_{\tn{imp}}}{\partial \bfW} &= \frac{\partial}{\partial \bfW} \left(\frac{1}{2}\Tr(\bfX^{\top}\bfX-\bfX^{\top}\bfX \bfW-\bfW^{\top}\bfX^{\top}\bfX+\bfW^{\top}\bfX^{\top}\bfX \bfW)+\bm{\lambda}^{\top}\tdiag(\bfW) \right) \\ &= - \frac{\partial}{\partial \bfW}\Tr(\bfW \bfX^{\top}X) + \frac{\partial}{\partial \bfW} \Tr(\bfW^{\top}\bfX^{\top}\bfX \bfW) + \tdiagMat(\bm{\lambda})\\&= \bfX^{\top}\bfX(\bfW-\bfI) + \tdiagMat(\bm{\lambda})
\end{aligned}
\end{eqnarray}

Similarly, taking gradient with respect to $\bm{\lambda}$ yields:
\begin{eqnarray}\label{langrange-lambda}
    \frac{\partial \mathcal{F}_{\tn{imp}}}{\partial \bm{\lambda}} = \tdiag(\bfW)
\end{eqnarray}
Setting the gradient $\frac{\partial \mathcal{F}_{\tn{imp}}}{\partial \bfW}$ to $\mathbf{0}$ yields:
\begin{eqnarray}\label{lagrange-W}
    \hat{\bfW} = I - \bfS^{-1} \tdiagMat(\bm{\lambda})
\end{eqnarray}
By substituting $\hat{\bfW}$ into Eq.~\ref{langrange-lambda} and setting it to $\mathbf{0}$ we get:
\begin{eqnarray}
    \hat{\bm{\lambda}} = \mathbf{1} \oslash \tdiag(\bfS^{-1})
\end{eqnarray}
where $\oslash$ is the element-wise division. Finally, by substituting $\hat{\bm{\lambda}}$ back into Eq.~\ref{lagrange-W} we get the expression of the optimal $\bfW$:
\begin{eqnarray}
    \hat{\bfW} = \bfI - \bfS^{-1}\tdiagMat(\mathbf{1} \oslash \tdiag(\bfS^{-1}))
\end{eqnarray}

It can also be verified that $(\hat{\bfW}, \hat{\bm{\lambda}})$ is indeed the global minimum by substituting it in Eq.~\ref{constrained-optimation-problem}.\par
Now, to express rPCN as performing metric learning in the form of Eq.~\ref{eq:general_metric} with a query $\mathbf{q}$, note that by replacing $W$ with the optimal $\hat{\bfW}$:
\begin{eqnarray}
\begin{aligned}
    \mathcal{F}_{\tn{imp}} &= \frac{1}{2} \|\mathbf{q} - \bfW^{\top}\mathbf{q}\| \\&= \frac{1}{2}\| (\bfI - {\bfW})^{\top}\mathbf{q}\|^2_2 \\&= \frac{1}{2}\|\tdiagMat(\mathbf{1} \oslash \tdiag(\bfS^{-1}))\bfS^{-1}\mathbf{q}\|_2^2
\end{aligned}
\end{eqnarray}
which follows the general expression of Eq.~\ref{eq:general_metric}.

\section{Experimental results with MNIST}
\label{S1_Figure}


\begin{figure}[ht!]

\includegraphics[width=\textwidth]{figs/chap3/Appendix-MNIST.pdf}
\caption{A: Retrieval MSEs of the single-layer models across multiple $N$s. B: Retrieval MSEs of the multi-layer models across multiple $N$s. }

\label{fig:performance-mnist}
\end{figure}

In this section of the Appendix we show the quantitative results obtained with the MNIST dataset. Fig~\ref{fig:performance-mnist}A shows the memory retrieval MSEs when the explicit, implicit and dendritic models were trained to memorize various numbers of MNIST images, and were given half-covered memories to retrieve from. The performances of these models are similar to those with grayscale CIFAR10 (Fig~\ref{fig:chap3-mnist}), where the implicit model has the lowest retrieval MSE, while the dendritic model's retrieval performance will start to deteriorate when $N$ becomes large, due to its unstable neural dynamics. The performance of the explicit model is the worst due to the need to compute the inverse term. In fact, we noticed that the learning of the explicit models with MNIST can hardly converge with some subsets of data, which resulted in the non-monotonically increasing retrieval MSEs as we increase $N$. We suspect that this is due to the innate dark backgrounds of the MNIST images: the black pixels in the MNIST images have perfect collinearity (e.g., the top-right and top-left black pixels), which makes the covariance matrix rank-deficient and thus non-invertible. This will prevent the learning of $\Sigma$ from converging.

Fig~\ref{fig:performance-mnist}B shows the performance of the multi-layer models in the same retrieval task from half-covered MNIST images. For comparison the same line from Fig~\ref{fig:performance-mnist}A for the implicit model is also added. The number of parameters in these models were also controlled to be the same. The results are mostly similar to those with the grayscale CIFAR10 images, where the performance different between the implicit and dendritic models are remedied by the hierarchical pre-processing, and the implicit model performs the best as it represents the linear solution with the lowest MSE. Interestingly though, the purely hierarchical model from \cite{salvatori2021associative} performs slightly better than the hybrid models. We leave the investigation into this phenomenon for future works.

\section{Examples of retrieved color CIFAR10 images with different levels of corruptions.}
\label{S2_Figure}
We also performed preliminary experiments with the dendritic model on colored CIFAR10 images, and show some examples here. Particularly, we found that when presenting the network a cue consisting of $1/2$ (panel A), $1/4$ (panel B) and $1/8$ (panel C) of the original pixels, the model has successfully recovered, respectively and on average, 89, 79, and 49 of the original memories in the experiments using colored CIFAR10.

\begin{figure}[ht!]
\centering
\includegraphics[width=0.9\textwidth]{figs/chap3/Appendix-CIFAR.png}
\label{fig:appendix-cifar}
\caption{Experiments with colored CIFAR10 imags. See main text above for details.}
\end{figure}


\newpage
\section{Implementation details of experiments with grayscale CIFAR10}
\label{S1_Table}
For all sample size of memories $N$, we use a batch size of $N/8$. For the inference iterations with the multi-layer models, the first number $400$ is the number of inference iterations during training and within each training iteration. The second number $100000$ is the number of inference iterations used to retrieve the original patterns.

\begin{table}[!ht]
\scriptsize
\centering
\begin{tabular}{|c|c|c|c|c|c|}
\hline
model & corruption & $\alpha$ & learning iters &  $\beta$ & inference iters \\
\hline
explicit & cover half & 1.00E-05 & 800 &  1.00E-01 & 5000 \\
\hline
implicit/dendritic & cover half & 1.00E-04 & 800 &  1.00E-01 & 50000 \\
\hline
implicit/nonlinear implicit & noise (var=0.1) & 1.00E-04 & 800 &  1.00E-03 & 20000 \\
\hline
nonlinear implicit & cover half & 1.00E-04 & 800 &  1.00E-01 & 50000 \\
\hline
hybrid/hierarchical & cover half & 1.00E-03 & 200 &  1.00E-02 & 400/100000 \\   
\hline
\end{tabular}
\label{tab:cifar10-hparams}
\end{table}

\section{Implementation details of experiments with MNIST}
\label{S2_Table}
For all sample sizes of memories $N$, we use a batch size of $N/8$. For the inference iterations with the multi-layer models, the first number $400$ is the number of inference iterations during training and within each training iteration. The second number $100000$ is the number of inference iterations used to retrieve the original patterns.

\begin{table}[!ht]
\scriptsize
\centering
\begin{tabular}{|c|c|c|c|c|c|}
\hline
model & corruption & $\alpha$ & learning iters & $\beta$ & inference iters \\
\hline
explicit & cover half & 1.00E-05 & 1000 & 1.00E-01 & 5000 \\
\hline
implicit/dendritic & cover half & 2.00E-04 & 1000 & 2.00E-02 & 50000 \\
\hline
hybrid/hierarchical & cover half &  1.00E-03 & 200 & 1.00E-02 & 400/100000 \\
\hline
\end{tabular}
\label{tab:mnist-hparams}
\end{table}