\section{Introduction}

Together with the Global Neuronal Workspace, Integrated Information Theory (IIT) \cite{Balduzzi2008,Tononi2008,Oizumi2014,tononi2016integrated,Tononi2017} is currently the most prominent approach to obtaining a scientific theory of consciousness, i.e, a theory that can identify the neural correlates of consciousness through some direct or indirect measure. Difficulties in agreeing on experimental setups that may discriminate between different theories have recently led to a sort of adversarial collaboration between them. \cite{Melloni2021} Within said collaborative framework, each theory provides the conditions for it to be experimentally falsified. However, despite such partnership and the obtaining of promising results in clinical tests, \cite{Casali2013,Bodart2017,Haun2017} reaching definitive answers may prove more difficult than initially expected because of the multifactorial complexity of the biological substrate of consciousness, namely the living brain. \cite{Arsiwalla2018a}
\\

The appeal of IIT for computation science and, in general, the AI community relies on its straightforward definition of consciousness as a formal property of any physical system, composed of nodes and interactions, that instantiates it. Consciousness {\it{is}} maximally integrated information $\Phi$; \cite{Tononi2008} the latter quantity being amenable to a strict mathematical definition in a network of interconnected elements that can be thought of as a probabilistic causal graph. Were said statement true, IIT would be providing not only a quantifiable prediction for a system to be conscious but also a heuristic approach to generate conscious systems ---the holy grail of the AI program. Yet it could also happen that IIT offers just a necessary, but not sufficient, condition for the emergence of consciousness in a system. As a logical consequence, a system implementing an artificial intelligence algorithm need not be conscious whatsoever \cite{garrido2022independence} and, hence, need not {\it{understand}} the meaning of the problems that it solves. \cite{garrido2022artificial}
\\

Interestingly, IIT derives its mathematical definition of consciousness from five axioms that fully account for what, phenomenologically, being conscious implies. To wit: The present experience of consciousness exists, is informative, composed, integrated, and excludes other conscious experiences \cite{Oizumi2014}. Said axioms should naturally determine the emergence of a physical scale in which the system, from its own perspective, presents a maximally integrated structure of intrinsic information, in which the information of the whole is maximally irreducible to the sum of information in subsets of any partition of the system. Additionally, the physical state of the substrate instantiating maximal $\Phi$, for the same token, presents a maximal cause-and-effect composition (for the details of the link between measures of integrated information and causality in an interconnected network, one may see \cite{Albantakis2019a,Albantakis2019}).
\\

Whereas IIT continues developing as a research program that strives for finding the best definition of distance between informational contents \cite{Barbosa2021} and including in its framework many other features of consciousness, \cite{Haun2019}
it has been criticized on diverse grounds \cite{Searle,Pautz2019a,Merker2022}. More specifically for our interests in this paper, some authors \cite{Morch2019,Sanchez-Canizares2022} have raised the issue of the so-called 'intrinsicality problem' (IP), namely, the fact that consciousness cannot be an intrinsic property of the system because maximal $\Phi$ crucially depends on the possible existence of bigger values of $\Phi$ if the initial system is appropriately linked to or embedded in larger systems. IP's threat, if true, could manifest that the causality involved in the constitution of the system goes beyond the conceptualization of causality provided by IIT in its current form. \cite{Sanchez-Canizares2022}
\\

Even though the literature on the computation of $\Phi$ in different networks is rapidly increasing and gaining insights on the logical architectures favoring maximal $\Phi$ (basically those consisting of modular, homogeneous, and specialized networks with feedback connections \cite{Oizumi2014}),
the prohibitively large computation times for big networks, scaling with the number of nodes ($n$) as
$\mathcal{O}(n53^n)$,
\cite{mayner2018pyphi} make it necessary the search for alternative heuristics and/or new methods for measuring $\Phi$ in larger systems. \cite{Nilsen2019,Toker2019}
Not much is known about the actual distribution of $\Phi$ over different network types and topologies.
This paper thus explores the benefit of random search algorithms in order to compute $\Phi$.
If, according to IIT's creators "the quantity and quality of consciousness are what they are," and consequently "the cause–effect structure cannot vary arbitrarily with the chosen measure of intrinsic information" \cite{Barbosa2021}, how  maximal $\Phi$ varies with the number of systems nodes must be a relevant matter to get a glimpse of how consciousness could eventually scale up in networks. Hence, our paper studies the expected increase of maximal $\Phi$ with the number of nodes, paving the way for gaining insights on the network's architectures that favor consciousness and assessing the seriousness of the IP, if IIT holds, in future works.
In such an endeavor, our paper will also show the limits of some optimization techniques for IIT as the number of network nodes increases.
\\

This paper is organized as follows. First, we provide a cursory definition of integrated information, a problem definition of $\Phi$ optimization, and a technical description of the methodology that we have implemented to optimize $\Phi$. Most critically, we provide arguments to justify why our proposed algorithm does not use smarter black-box global optimization techniques such as Bayesian optimization. Then, we show in an experiments section empirical evidence to support the claim that our method is able to optimize $\Phi$. Finally, we present some concluding remarks and refer to possible research lines to optimize $\Phi$ in further work.

\section{Proposed Methodology}
In this section, we begin formally defining integrated information and the problem of optimizing integrated information in the space of transition probability matrices. Then, we show how this problem cannot be solved with common metaheuristics or techniques such as vanilla Bayesian optimization \cite{garrido2021advanced} and introduce a random search algorithm to solve the problem that can be used as a baseline for further research.

\subsection{Integrated information}
Throughout this paper, we use the definitions of integrated information employed in \cite{mayner2018pyphi}. Details regarding the differences between integrated information ($\phi$) for a mechanism ---a set of nodes and interactions in a graph forming a concept--- and conceptual integrated information ($\Phi$) for a constellation of concepts forming a maximally irreducible conceptual structure, can be consulted in \cite{Oizumi2014}. Computations in our paper always refer to  $\Phi$.
\\

We can define integrated information for a mechanism in state $s$ as

\begin{align}
\phi (s) = D [ p(s_{0} \to s) \| p^{MIP} (s_0 \to s)
]
\,,
\end{align}

where $D$ is a measure of the distance between the probability distribution $p$ of obtaining state $s$ as the cause (effect) of future (past) states $s_{0}$ and the correspondent probability distribution $p^{IMP}$ when the system has been optimally partitioned in its weakest links, i.e. decomposed into those parts that are most independent (least integrated), providing a minimum of effectively integrated information for the mechanism as a whole.
\\

When one finds maximal $\phi$ for a mechanism in state $s$ in a system, one has a concept. After determining the structure of concepts in a candidate system, one is allowed to define the conceptual integrated information ($\Phi$) as
$\Phi(C) = D(C \| C^{MIP})$, where $C$ refers to the structure of concepts,  $C^{MIP}$ refers to the optimal partition of the structure of concepts that makes $\Phi(C)$ minimal, and $D$ is an appropriately adapted measure of the distance between the corresponding probability distributions. Different versions of IIT have employed different definitions of the distance between probability distributions ---the Kullback-Leibler divergence for $\phi$ \cite{Balduzzi2008}, the Wasserstein Distance, also known as Earth Mover’s Distance, for both $\phi$ and $\Phi$ \cite{Oizumi2014}, or the Intrinsic Difference \cite{Barbosa2021}. As the PyPhi algorithm favors the Earth Mover’s Distance as the measure of integrated information \cite{mayner2018pyphi}, which quantifies how much two distributions differ by taking into account the distance between system states, we will use it here too.

\subsection{Problem definition}
We now provide a formal definition of the global black-box optimization problem that we need to solve. Concretely, we need to tackle the optimization of a function $\Phi : \mathbb{R}^{2^D,D} \to \mathbb{R}$, from a D-$2^D$-dimensional transition probability matrix, with conditional independencies, space to a real number.

\begin{align}
\omega^\star = \arg \max_{\omega \in \Omega} \Phi(\omega) \quad s.t. \quad \omega \in \Upsilon \,,
\end{align}

where $\omega^\star$ is the TPM that maximizes $\Phi(\cdot)$, $\omega$ is a TPM of the D-dimensional TPM space, $\Omega \in \mathbb{R}^{2^D,D}$ is the D-dimensional TPM space, $\Upsilon$ is the space of TPM with conditional independencies and $D \in \mathbb{N}^+$ is the number of dimensions of the TPM, whose range in this problem is $\mathbb{N}$. Recall that the objective function $\Phi$ is expensive to evaluate when the number of nodes of the TPM is high. If $\omega \notin \Upsilon$, then $\Phi(\omega) = \emptyset$. Consequently, $\mathbb{R}^{2^D,D}$ is constrained to $\Upsilon$. This is especially problematic, because of the curse of dimensionality. In particular, as $D$ grows, the $\Omega_D$ space will contain higher $\Phi$ TPMs, which is desirable. However, the number of valid TPMs will be lower with respect to all the possible TPMs than in lower dimensional TPM spaces. Consequently, as the evaluation of $\Phi$ is expensive, we have a tradeoff between risking to find high $\Phi$ solutions in a higher dimensional TPM space where the volume of valid TPMs is low or searching in a lower dimensional TPM space where we have a high probability of obtaining a valid TPM but its $\Phi$ is potentially low. More formally, if $\delta(\omega)$ is a delta function whose value is $1$ if $\omega \in \Upsilon$, $D_{min}$ is the lowest dimensionality considered for the number of nodes of the TPM, $D_{max}$ is the highest dimensionality considered for the number of nodes of the TPM and $p(\omega|\Omega_D)$ is the probability of sampling each possible $\omega \in \Omega_D$, we have that:

\begin{align}
 \int \delta(\omega) p(\omega|\Omega_{D_{max}}) d\omega \leq ... \leq \int \delta(\omega) p(\omega|\Omega_{D_{min}}) d\omega \,.
\end{align}

Moreover, as $D$ grows, the cost of evaluating $\Phi$ also grows. Hence, if the budget is constrained by time, and not by iterations, we must be sure about searching in a high dimensional $\Omega$ space. We tackle this issue and define a new algorithm to solve this problem, among others that we further explain, in the following subsection.

\subsection{Proposed algorithm}
As we have seen in the previous section, the global optimization problem of the $\Phi$ function $\omega^\star = \arg \max_{\omega \in \Omega} \Phi(\omega)$ is not smooth with respect to the smallest variations $\nabla \omega$ that could be performed in the transition probability matrix space $\Omega$. For example, in the case of binary nodes, switching just one node from zero to one can make the $\Phi$ from $\mathbb{R}$ to $\emptyset$. In particular, this is troublesome for some advanced black-box global optimization techniques, as we describe in this section.
\\

Bayesian optimization \cite{garrido2021advanced,shahriari2015taking} is the state-of-the-art class of methods to solve black-box global optimization problems that iteratively recommend a suggestion $\mathbf{x}$ making a tradeoff between exploration of promising unknown regions of the search space $\mathcal{X}$ and potential good explored solutions based on the predictions made by a probabilistic surrogate model as a Gaussian process \cite{schulz2018tutorial}, which formally is defined as:

\begin{align}
\mathbf{x}^\star = \arg\min_{\mathbf{x} \in \mathcal{X}} f(\mathbf{x}) \,,
\end{align}

where $\mathbf{x} \in \mathbb{R}^D$ is a candidate point of the hypercube $\mathcal{X} \subseteq \mathbb{R}^D$ and $f(\cdot)$ is the multi-modal function to be optimized. In order to be considered a black-box global optimization problem, $f$ must have no known analytical expression, and hence we are not able to extract gradients; it can be noisy and it is expensive to evaluate. Most critically, if the function $f(\cdot)$ is cheap to evaluate, then it is recommended to use metaheuristics as genetic algorithms \cite{kumar2010genetic}, that also perform an exploration-exploitation tradeoff. The reason is that their execution time to suggest recommendations $\mathbf{x}$ is faster than in the case of Bayesian optimization, whose complexity is a cubic on the number of previous recommendations $\mathcal{O}(N^3)$, where $N$ is the number of previous recommendations.
\\

\SetKwComment{Comment}{/* }{ */}

\begin{algorithm}[htb!]
\caption{Prior guided random search algorithm to optimize $\Phi$ on the D-dimensional TPM space $\Omega$}
\label{alg:phi}
\KwData{$p(\omega), D_{min}, D_{max}, \mu, \epsilon, T$}
\KwResult{$\Phi^\star, \omega^\star$}
$i \gets 0$\;
$\Phi^\star \gets 0$\;
\While{$i \leq T$}{
$j \gets 0$\;
$\Phi_{list} = []$\;
\While{$j \leq \epsilon$}{
    $D_i \gets multinomial\_sample(p(\omega), D_{min}, D_{max})$\;
    $\omega_i \gets uniform\_TPM\_sample(D_i)$\;
    $\Phi^i = \Phi(\omega_i)$\;
  \If{$\Phi^\star \geq \Phi_i$}{
    $\Phi^\star = \Phi_i$\;
    $\omega^\star = \omega_i$\;
  }
  $\Phi_{list} \gets append\_result(\Phi_i, D_i, \Phi_{list})$\;
  }
  $p(\omega) \gets update\_prior(p(\omega), \mu, \Phi_{list})$\;
  }
\end{algorithm}

The problem definition presented in the previous section involves the computation of $\Phi$, which it is very expensive if the number of nodes of the TPM is high. Consequently, this scenario does not present a good fit for metaheuristic algorithms, as they require the evaluation of the objective function to be cheap. Moreover, vanilla Bayesian optimization requires that the function does not vary with respect to its parameters up to a certain $K$-Lipschitz value. \cite{hager1979lipschitz} Critically, in this problem the smallest variation on any elements of the TPM matrix can make $\Phi$ jump from $\mathbb{R}$ to $\emptyset$. Hence this assumption, required by Bayesian optimization, is violated and would need to be specifically modified to tackle this particular problem, whose complexity is even bigger as the space being optimized is $\mathbb{R}^{2^D,D}$, making Bayesian optimization also hard to use. Having reviewed these alternatives, there also exists another reason why those techniques may be inefficient. Concretely, if the $K$-Lipschitz value assumption of Bayesian optimization is violated, the algorithm cannot perform exploitation. Hence, any exploitative behavior on the D-dimensional TPM space $\Omega$ would be inefficient. Consequently, the only alternative that remains is to perform a prior guided random search algorithm, adapted to this particular problem, and encode into the prior our knowledge about the $\Phi$ optimization problem.
In particular, we enhance a random search procedure instead of a grid search because random search significantly outperforms grid search when the percentage of explained variance on the variable to optimize is not uniform over the explanatory variables \cite{bergstra2012random}. As assuming that the variance of $\Phi$ is uniformly explained by all the nodes for all the possible TPMs is a heavy assumption, we can not rely on uniform grid search, as random search would be a smarter option in this case. Recall that random search outperforms grid search when some variables of the input space explain the objective function better than others. Consequently, we decided to enhance the random search procedure. The prior guided random search illustrated on Algorithm \ref{alg:phi} will obtain $\Phi^\star$, which is an approximation to $\Phi^{best}$, that is, $\Phi^\star \approx \Phi^{best}$, which could only be found using an exhaustive search throughout all the TPM D-dimensional space $\Omega$. Where, in Algorithm $\ref{alg:phi}$, $\mu$ is the adaptation rate of $p(\omega)$, the multinomial a priori distribution of the probability of obtaining $\Phi^\star$ in a certain $D$-dimensional space, $\epsilon$ is the number of iterations that are needed to update $p(\omega)$, $T$ is the budget of total iterations that the user can afford, $D_{min}$ is the lower dimension of the TPM D-dimensional space and $D_{max}$ is the higher dimension of the TPM D-dimensional space. The sampling functions define a probabilistic graphical model where the TPM is modeled as a random variable that depends on the multinomial random variable that encodes the probability distribution assigned to the dimensions of the D-dimensional TPM space. More formally, we obtain $\omega_i$ sampling from the following joint probability distribution:

\begin{align}
\omega = uniform\_matrix(\omega \mid D) \quad multinomial(D \mid \theta, D_{min}, D_{max}) \,.
\end{align}

Sampling from the previous joint distribution of both random variables we can sample TPMs conditioned to the probability of obtaining a TPM belonging to the space of valid matrices $\Upsilon$. The initial state of the TPM necessary to compute $\Phi$ is tested in an iterative fashion until a valid state is found.
\\

At the end of a batch of $\epsilon$ iterations, we update our beliefs about the $\theta$ parameter values of our multinomial distribution $p(\omega|\theta)$ up to a learning rate value $\mu$. We do so in a Bayesian learning fashion, basically defining a likelihood function $f(\Phi_{list}, \mu)$ that penalizes the dimensions whose results are worse and rewards those whose values are better up to $\mu$. This likelihood function is multiplied by the prior beliefs $p(\omega|\theta)$ and the posterior is obtained by normalizing the previous operation with a normalization constant $\sum_{D_{min}}^{D_{max}} p(\omega|\theta) f(\Phi_{list})$ to ensure that the resulting distribution is a valid multinomial distribution and the prior is conjugate. Finally, we smooth the maximum and the minimum of the vector of probabilities by a quantity $\kappa=0.02$ that could be also modified by the user. The specific parametric likelihood function $f(\Phi_{list}, \mu)$ can be set by the user. Concretely, we just provide as default to linearly penalize or reward each dimension in $0.2 \mu r$ points where $r$ is the ranking of each dimension $D$ that consists of a $[-1,1]$ linear space vector. More formally, if $p(\omega|\theta)_{new}$ is the posterior distribution obtained after the batch of $\epsilon$ iterations, then, to obtain this distribution we perform the following update:

\begin{align}
p(\omega|\theta)_{new} = \frac{f(\Phi_{list}, \mu) p(\omega|\theta)}{\sum_{D_{min}}^{D_{max}} f(\Phi_{list}, \mu) p(\omega|\theta)}
\end{align}

Lastly, we remark that the values of the set of prior beliefs $\theta$, the learning rate $\mu$, and the batch size of iterations $\epsilon$ are hyper-parameters that the user can set and will modify the behavior of the algorithm. Additionally, they also allow the user to set the value of $D_{min}, D_{max}$ to define the TPM search space $\Omega$. In the next section, we illustrate the algorithm performance, being compared with a grid search method and a random search method as baselines.

\section{Experiments}

In order to estimate the integrated information of a transition probability matrix, we use the PyPhi tool implemented on Python3 \cite{mayner2018pyphi}. In the computation of $\Phi$, to measure distances between probability distributions, we use the Earth Mover's Distance \cite{andoni2008earth}, which is equivalent to the Hamming Distance if the nodes of the graph are binary. We split the experiment section into two parts. First, we show the results obtained by the prior guided random search algorithm with respect to other baselines, just to illustrate how we can obtain graphs whose $\Phi$ score is bigger algorithmically. Driven by the results of these experiments and because of the size of the TPM space that makes it computationally unfeasible to explore, we append another section where we try to statistically make inference about $\Phi$ as the number of nodes grows. We exhaustively search for an initial state for the TPM and retrieve the first feasible state for a randomly sampled TPM, making a point estimation about the mean $\Phi$ of the TPM with the first feasible state found, for computational reasons. All the code, written on Python3, and the experiments are uploaded into the Github repository (https://github.com/EduardoGarrido90/iit\_opt).
\\

\subsection{Prior Guided Random Search experiments}

The following experiments model the obtained $\Phi^\star$ by the various methods as a Gaussian random variable with mean $\mu$ and standard deviation $\sigma$. All the experiments are repeated $25$ times to infer their respective values of $\mu$ and standard deviation $\sigma$. Namely, we compare our proposed prior guided random search method with respect to a grid search method and a random search method. The random search method samples uniformly both the biggest dimension of the TPM and the values of the matrix. However, the grid search method uniformly splits the dimensions of the TPM across all the space range, and the binary TPMs are selected by obtaining the binary representations of a linear space of integers from $0$ to $2^n$ where $n$ is the number of nodes which are the columns of the matrices.
\\

First, we begin with the following experiment, to verify whether all the methods work and provide valid $\Phi$ results. We perform 25 repetitions of the experiment with each method to ensure statistical significance. We set the lower dimension, or the number of nodes, to 3 and the higher dimension to 4. We execute a total of 50 iterations with $\epsilon=5$ for the prior guided search algorithm. We set $\mu=0.1$ and the prior distribution values on the dimensions to $p(\omega) = [0.2, 0.8]$. With these values, we obtain the $\Phi$ results across the iterations shown in Figure \ref{fig:exp1}. The best $\Phi$ found in the search was $3.034082$.
\\

\begin{figure}[h]
    \centering
    \includegraphics[width = 0.99\textwidth]{exp1_final.png}
\caption{Average results of $\Phi$ with the toy experiment configuration of the Prior Guided Random Search method (blue) compared to the Grid (green) and Random search (orange) baselines. Better seen in color.}
    \label{fig:exp1}
\end{figure}

As we can see in Figure \ref{fig:exp1}, the prior guided random search method, plotted in blue, outperform the baselines. The results are statistically significant as the plotted results are the means of $25$ different experiments. As the prior guided random search method performs almost twice as better as the grid search method, a statistical hypothesis testing method will provide a p-value lower than $.01$ to reject the null hypotheses that the differences are random, accepting the alternative hypothesis that both methods have different distributions. Concretely, random search also outperforms grid search. And lastly, the more interesting conclusion is that the prior guided random search method outperforms random search. Regarding grid search, we have found that the restriction of the matrices of the linear grid seems to provide bad results, resulting mostly on invalid TPMs and not being a method that is worth using to optimize $\Phi$ in the $\Omega$ space. Hence, for the next experiments, and also for computational reasons due to the complexity of $\Phi$ with respect to the number of nodes, we will only compare the prior guided random search method with the random search method.
\\
We have performed a little experiment, which we will call Graphical Experiment, of one repetition just to show graphically the obtained transition probability matrix, the connectivity matrix, and the initial state. This experiment involves 100 repetitions with the minimum dimension of nodes being $3$ and the maximum number of nodes $4$. The initial prior is $[0.3, 0.7]$. We obtained, with the prior guided random search algorithm, a graph with $\Phi=2.6823$, with a final prior with value $[0.12, 0.88]$, the connectivity matrix shown on Figure \ref{fig:cm}, the transition probability matrix shown on Figure \ref{fig:tpm} and the state shown on Figure \ref{fig:state}.

\begin{figure}[h]
    \centering
    \includegraphics[width = 0.99\textwidth]{cm.png}
\caption{Connectivity matrix associated to the best result obtained in the Graphical Experiment.}
    \label{fig:cm}
\end{figure}

\begin{figure}[h]
    \centering
    \includegraphics[height = 0.6\textheight]{tpm.png}
\caption{Transition probability matrix associated to the best result obtained in the Graphical Experiment. 0 is painted in purple and 1 in yellow.}
    \label{fig:tpm}
\end{figure}

\begin{figure}[h]
    \centering
    \includegraphics[width = 0.99\textwidth]{state.png}
\caption{Initial state associated with the best result obtained in the Graphical Experiment. 0 is painted in purple and 1 in yellow.}
    \label{fig:state}
\end{figure}

In the following experiment, we want to verify whether the behavior differences in the method generalize to more than 2 dimensions in the space of TPMs, this time, from 2 nodes to 5. With a prior distribution with weights $p(\omega) = [0.1,0.1,0.3,0.5]$ and only $5$ repetitions to save computational resources and $50$ iterations. We obtain the results plotted in Figure \ref{fig:exp2} where we see how the Prior Guided Random Search algorithm is also outperforming the random search method.
\\

\begin{figure}[h]
    \centering
    \includegraphics[width = 0.99\textwidth]{exp2.png}
\caption{Average results of $\Phi$ with the small experiment configuration of the Prior Guided Random Search method (blue) compared to the Random search (orange) baseline. Better seen in color.}
    \label{fig:exp2}
\end{figure}

We can see, once again, that both searches are able to find better TPMs according to the $\Phi$ measure once the iterations are executed and that prior guided random search once again outperforms random search, not only on average but in absolute value. In particular, the maximum $\Phi$ found by the prior guided random search method is $3.1681$ and the best result obtained by random search is $2.0832$.
\\

As we have stated previously, the complexity of evaluating $\Phi$ with respect to the nodes is high, in particular, $\mathcal{O}(n53^n)$, which must be multiplied by $2^n$ in these searches as we look for a feasible state for a randomly sampled TPM. Consequently, and having shown empirical evidence of the usefulness of the random searches method and how our proposed method even outperforms random search for this specific task, our last experiment only involves $1$ repetition of each method, for computational reasons. Being unable to test $\Phi$ in TPMs with $6$ or more nodes due to the complexity of the algorithm, in this experiment, we will check whether better results of $\Phi$ are obtained if the number of nodes and iterations is bigger. The minimum dimension is now $4$ and the maximum dimension is $5$. The prior distribution is $[0.2,0.8]$ and the number of iterations is $500$, with $\mu=1$ and $\epsilon=5$, obtained the results plotted in Figure \ref{fig:exp3}.
\\

\begin{figure}[h]
    \centering
    \includegraphics[width = 0.99\textwidth]{exp3.png}
\caption{Average results of $\Phi$ with the small experiment configuration of the Prior Guided Random Search method (blue) compared to the Random search (orange) baseline. Better seen in color.}
    \label{fig:exp3}
\end{figure}

We can see how better results are obtained if we have a bigger budget of iterations or time. Moreover, this experiment also shows that the Prior Guided Random Search algorithm outperforms Random search once the prior distribution fits the optimal number of nodes. As the $\Omega$ space is very sparse, a much higher number of iterations would be required to obtain a better approximation to the optimal $\Phi$. More critically, for further work, with more computation power that would enable searches for graphs of more than $6$ nodes, the space would even be more sparse, so much more computation would be needed to find graphs with bigger $\Phi$. However, we have shown with these experiments how our algorithm would adapt to the sparsity of every dimension and keep searching for solutions in the number of nodes that represents the better trade-off between the probability of finding better $\Phi$ and the sparsity of the dimension that makes difficult to find a feasible $\Phi$.

\subsection{Experiments to explore the graph space and make statistical inference of $\Phi$}
The complexity of the space of transition probability matrices makes it computationally intractable to be explored exhaustively. Moreover, in order to measure $\Phi$, we must not only take into account the matrix that encodes the transitions but also the connectivity matrix, that represents the connections, or edges, of the graph. Consequently, the complexity, or size measure, of the space is even higher.

However, we can frame this problem into the classical statistical inference framework. The population is the total of possible combinations for a given number of graphs and we can extract a random sample from it. The goal is to make inferences about the mean $\Phi$ parameter $\mu$ associated with each number of nodes. Then, we can perform statistical hypothesis testing to determine whether the difference of $\Phi$ is statistically significant when $N$ grows. By doing this procedure, we can validate, with the previous set of experiments, that our algorithm is going to work even better when the number of nodes grows, assuming that the same inference will hold when the number of nodes grows, which we can not empirically test since the complexity of computing $\Phi$ when the number of nodes grows is exponential.

The first toy experiment involves 100 graphs of 3 nodes and 100 graphs of 4 nodes. We obtain an empirical $\Phi$ mean of $0.0517$ in the case of 3 nodes and, in the case of 4 nodes, an empirical mean of $0.1711$. The confidence intervals at 95\% are $[0.0212,0.0822]$ in the case of 3 nodes and $[0.1093,0.2328]$ in the case of 4 nodes. We also perform a 2 sample t-test whose null hypothesis is that both population means are equal and the alternative hypothesis is that they are different. We obtain a t-statistic of $-3.38$ whose associated p-value is $0.0008$. At a confidence level of $\alpha=0.01$ we can reject the null hypothesis that both means are equal, hence accepting the alternative hypothesis that they are different, being the $\Phi$ population mean of graphs of 4 nodes bigger than the $\Phi$ population mean of graphs of 3 nodes. Also, the percentage of infeasible solutions found in the case of $3$ nodes has been $72.6775\%$ and in the case of $4$ nodes has been $84.7094\%$, providing empirical evidence of Eq. 2. To test consistency, we performed the same experiment with a population of $200$ in both groups. In this case, we obtained an even bigger, in absolute terms, a t-test statistic of value $-4.6420$ with associated p-value $4.69e-6$, rejecting the null hypothesis. The obtained empirical means are $0.0569$ and $0.1847$, pretty similar to the ones obtained in the previous experiment with similar confidence intervals.

The second experiment involves a population of 100 graphs of 4 nodes and 100 graphs of 5 nodes. We obtain an empirical $\Phi$ mean of $0.1432$ in the case of 4 nodes and, in the case of 5 nodes, an empirical mean of $0.3267$. The confidence intervals at 95\% are $[0.0822,0.2043]$ in the case of 4 nodes and $[0.2217,0.4317]$ in the case of 5 nodes. We also perform a 2 sample t-test whose null hypothesis is that both population means are equal and the alternative hypothesis is that they are different. We obtain a t-statistic of $-2.94$ whose associated p-value is $0.0036$. At a confidence level of $\alpha=0.01$ we can reject the null hypothesis that both means are equal, hence accepting the alternative hypothesis that they are different, being the $\Phi$ population mean of graphs of 5 nodes bigger than the $\Phi$ population mean of graphs of 3 nodes. Also, the percentage of infeasible solutions found in the case of $4$ nodes has been $84.98\%$ and in the case of $5$ nodes has been a $91.84\%$

\section{Philosophical discussion}

This work presents an improvement on how to compute $\Phi$ for arbitrary networks and TPMs. Whereas the present calculations still work with toy models with a maximum of 6 nodes, the obtained results with a guided random method strongly suggest that the latter outperforms other methods to explore optimal TPMs as the number of nodes increases. In this sense, it reveals a promising procedure for further progress in the time computation of $\Phi$ and, eventually, for the ascertainment of the validity of IIT as the scientific theory of consciousness in realistic networks. We aim at answering questions like these: Can we quantitatively provide a method to generate consciousness? May this methodology be adapted to deep neural networks if the assumptions of IIT were true, providing deep learning machine consciousness architectures \cite{merchan2020machine}?
\\

On the other hand, one should keep in mind the intrinsic limitations of the PyPhi algorithm, which obviously affect our results. According to its creators, ``[T]he analysis can only be meaningfully applied to a system that is Markovian and satisfies the conditional independence property. These are reasonable assumptions for the intended use case of the software: analyzing a causal TPM derived using the calculus of perturbations. However, there is no guarantee that these assumptions will be valid in other circumstances, such as TPMs derived from observed time series (e.g., EEG recordings)." \cite{mayner2018pyphi}
\\

The ideas of the last paragraph allow us to dive into additional philosophical problems related to the overall IIT procedures. One of those problems is the origin of networks. In the case of artificial networks (1), the unity of the network is artificially created by network designers, avoiding external influences that can disrupt the workings of the network. In the case of natural networks (2), however, one may wonder who or what causes the unity of the system. Moreover, whereas in case (1) possible consciousness as defined by IIT remains extrinsic, as a bigger maximal $\Phi$ can be reached by adding new nodes to the initial network, it is no longer clear ---indeed very controversial--- that the same procedure works in case (2). Case (2) could be the only scenario in which, in keeping with IIT tenets, consciousness is intrinsic. But if  that is the case, IP still holds for case (2).
\\

Last but not least, our procedure to obtain bigger $\Phi$ as the number of nodes is increased and the system's architecture becomes more complex might not reflect the ontology of case (2). The bottom-line problem is that an optimization procedure might just be an epistemic tool to obtain the magnitude which maximizes an extant, meaningful, quantity. However, the procedure itself need not unravel the ontological processes underlying the reality of said quantity.

\section{Conclusions and further work}

Throughout this paper, we have dealt with integrated information theory as a relatively novel approach to quantifying truly integrated information in a network of interacting nodes. Beyond its technical interest in the field, research in IIT is exploding as it aims to become a scientific theory of consciousness by validating its axioms in artificial and natural networks. Nevertheless, at present, procedures to compute $\Phi$ remain too lengthy. Even if the creation of the PyPhi software has provided IIT researchers with a magnificent tool, computation times increase too quickly with the number of nodes in a network. Yet, studying how maximal $\Phi$ varies when new nodes and interactions are added to the system may prove crucial to see how, allegedly, consciousness scales up and whether the intrinsicality problem presents a true threat to IIT.
\\

Our paper aims to provide a new method to evaluate maximal $\Phi$ with the help of state-of-the-art optimization procedures of network quantities. However, as stated in the technical sections, the requirements of conditional independence for the TPMs charge our procedure with a burden that may not be eased out while keeping all its computational strengths. For the time being, ours remain a toy model that seems promising to pursue in future works.
\\

This work thus opens up several lines of research. To enhance the behavior of the method, we hypothesize that several actions could be taken. First, it may be possible to optimize on a smooth, continuous space of conditional TPMs to use Bayesian optimization or genetic algorithms. This approach would need to formulate an analytical transformation from the space of TPMs into a smooth and continuous space that can be modeled with a Gaussian process, that is the probabilistic surrogate model of Bayesian optimization. By doing so, we can exploit the exploitation-exploration behavior of Bayesian optimization, which may outperform the current prior-guided random search algorithm. If this transformation cannot be eventually found, we can extract the best-observed result by the prior-guided random search algorithm and execute a combinatorial optimization algorithm with the best-observed result as the initial point. This may also provide good results if $\Phi$ varies smoothly as a result of varying the TPM or the CM in a single element. Lastly, as we have seen, $\Phi$ is computationally costly. An interesting line of research would be to find a cheaper approximation of $\Phi$ with respect to a particular graph with a regression graph machine learning algorithm. In order to fit this algorithm it would be necessary to create a huge dataset of graphs annotated by its $\Phi$ values and study whether a regression graph machine learning algorithm is able to accurately approximate $\Phi$.

\bibliographystyle{acm}


\section{Introduction}
\label{intro}

Back-propagation (BP) has been the mainstay of training deep neural networks over the past decade. However, it is not considered biologically plausible \cite{10.1145/3354265.3354275, Lillicrap2020BackpropagationAT}. Hence, research interest has always remained in finding learning algorithms that the human brain might be using. Over the last few years a myriad of Bio-inspired algorithms have been proposed. These algorithms have now achieved comparable accuracy as Backprop (BP) and also have been known to be applied to a vast array of deep learning tasks such as computer vision, reinforcement learning, language modelling, amongst others.

Bio-algorithms differ from BP in a few crucial aspects. Firstly, weight symmetry is required by BP between the forward and backward connections, that is, the same weight pathways must be used for both the forward and backward passes. This is not biologically plausible in the brain due to the uni-directional nature of brain synapses and is called the weight transport problem \cite{GROSSBERG198723}. Secondly, the problem of update locking occurs in BP. That is all layers must complete their processing before the backward pass can begin and the weights can be updated. This does not appear to be biological in nature as many plasticity algorithms in the brain learn in a local fashion, that is, learn between two pairs of neurons itself, without depending upon the activity of other layers \cite{BALDI2017110}. Thus, the brain supports multiple sets of neurons to learn in parallel. This also brings us to the third crucial difference, which is that BP learns in a global fashion, by calculating an error at the end of the last layer and then feeding it back to each layer. On the contrary, the brain learns in a local fashion, dependent only upon the activity of two connected neurons without a global loss to be learnt and backpropagated \cite{MAASS19971659}.

The motivation of using and developing Bio-algorithms, though, still mostly remains to make these algorithms more biologically plausible than BP. Tangible benefits over BP have not yet been demonstrated and therefore, adoption of Bio-learning algorithms outside the Bio-learning community remains small. We take a step in this direction and conduct a holistic benchmarking of Bio algorithms vs. BP in various experimental settings, to understand the practical benefits of Bio-algorithms. We look at popular Bio-learning algorithms like Feedback Alignment (FA) \cite{Lillicrap2016RandomSF}, Direct Feedback Alignment (DFA) \cite{DFA}, Hebbian Learning (HB) \cite{Krotov7723, Amato2019HebbianLM, 9207242, hebbnet, miconi}, Predictive Coding (PC) \cite{PC} and Difference Target Propagation (DTP) \cite{TargetProp}. Many other Bio-algorithms exist but we choose these as they have gained traction over the years and serve as representative samples from the different categories of Bio-algorithms.

For our benchmarking scenarios, we include a range of aspects important for practical usage of a learning algorithm. We measure performance of the algorithms under i) Limited availability of training data, ii) Limited training budget (in terms of the number of training epochs), iii) Application of noise to input samples, and iv) Application of sparsity to the parameters of the base neural network models (i.e. pruning the parameters of the model). In this paper, we have the following key contributions-

\begin{enumerate}
    \item We show that Bio-algorithms perform surprisingly well in the scenario of limited training data. Four of the five Bio-algorithms tested, namely, FA, DFA, HB and PC outperform BP when only 20\% of data is available for testing on multiple datasets.
    \item Another interesting finding is that some Bio-inspired algorithms perform exceedingly well in scenarios of limited training budgets. Specifically, HB and DTP learn much faster than BP and arrive at high stable accuracies in far lesser training epochs than BP. HB especially, learns within 5 epochs on multiple datasets where BP takes around 100 epochs to learn.
    \item To further understand the superior performance of Bio-algorithms, especially Hebbian learning, we plot the features learnt by BP vs. HB. We find that HB learns much more fine-grained and explainable features whereas BP seems to learn features that are much more random looking. We postulate that this may be the underlying reason for the better learning capability of HB.
    \item We also find that Bio-learning algorithms are robust against different kinds of noises and also adept under the application of parameter sparsity to the neural networks, achieving similar accuracy as BP.
\end{enumerate}

The above findings, especially points 1 and 2, are important advantages of Bio-learning algorithms that make them very useful in learning with less data applications or resource-constrained training environments like real-time devices. They also motivate the need to work on Bio-learning algorithms for practical advantages and beyond the traditional motivation of solving the biological implausibility of Backprop.

\begin{figure*}
\includegraphics[right]{Figures/Biolearning.pdf}
\caption{Illustration of how different learning algorithms. Feedback Alignment and Direct Feedback Alignment back-propagate gradients over a fixed random matrix instead of the original weights, thereby relaxing the condition of having the same pathway (weights) for both the forward and backward pass. Hebbian learning does not need any gradients to be back-propagated as it learns in an unsupervised manner. Only the last layer, which is the classification layer, learns based on gradients.}
\label{fig:algos}
\end{figure*}
\section{Related Work}

Hebbian learning has been around since 1949 when it was first proposed by Donald Hebb \cite{Brain_Function}. It has since been applied in many different areas including continual learning \cite{pmlr-v80-miconi18a, thangarasa2020enabling}, principal components analysis \cite{Oja1982SimplifiedNM, Sanger, Gorrell2006GeneralizedHA} and biology \cite{Oja1989NeuralNP, Keysers, Thiagarajan} amongst others. Hebbian learning has recently been adapted for modern neural networks and applied to train deep learning models \cite{Krotov7723, Amato2019HebbianLM, 9207242, hebbnet, miconi}. However, the focus for most of these studies has been to get Hebbian learning successfully running on deep learning architectures. A gap remains in benchmarking Bio-learning vs. Backprop on multiple aspects like training with fewer epochs, utilizing less data, evaluating the performance under network sparsity, etc. We focus on these aspects and present novel findings on the performance of Bio-learning vs. Backprop in this work.

Another Bio-inspired algorithm that has become popular is Feedback Alignment \cite{Lillicrap2016RandomSF}. It works by having a different set of weights for the forward and backward pass. In particular, the weights for the backward pass are a set of fixed random matrices over which the gradients are back-propagated. So, it is similar to BP in that an error is still calculated and back-propagated but it is different in that the weights for the backward pass do not have to be the same as the forward pass. This is biologically desirable because it breaks the symmetry condition requiring the forward synapses to also act as backward synapses.

A variant of FA, Direct Feedback Alignment (DFA) \cite{DFA}, relaxes this condition further whereby the error does not need to be back-propagated across every layer but can directly flow from the last layer to the layer in question through a shortcut connection. This makes the model even more biologically plausible and has also seen great success with many works derived from it \cite{10.3389/fnins.2019.00525, scaling_dfa}.

The fourth broad category of Bio-algorithms corresponds to Target Propagation algorithms that include Target Propagation \cite{TargetProp}, Difference Target Propagation \cite{theoretical_framework}, Equilibrium Propagation \cite{EquiProp}, amongst others. Target Propagation (TP) tries to solve the problem of credit assignment and computes targets rather than gradients at each layer. A variant of TP, Difference Target Propagation (DTP) actually works well in practice and differs from TP in that it linearly corrects for the imperfectness of TP. Equilibrium Propagation (EP) is similar in idea but is used for energy based networks like Hopfield networks.

Lastly, other Bio-algorithms include Predictive Coding \cite{PC}, Direct Random Target Projection \cite{10.3389/fnins.2021.629892}, Weight Mirroring \cite{WM}, Local Representation Alignment \cite{LRA}, amongst others. Predictive coding is a theoretical neuroscience inspired algorithm approximating BP using local learning rules. It attempts to arrive at similar gradient updates as BP but using only biologically plausible local updates. It has gained traction in recent years and presents a promising new direction for bridging Bio-learning and BP.

We find that while many interesting Bio-algorithms have been proposed in the literature, a thorough benchmarking between BP and various Bio-algorithms has not been undertaken on aspects outside the standard measurement metrics like test accuracy. Some papers propose theoretical frameworks for comparison \cite{cricket} or provide brief comments on other aspects, for instance, number of training episodes \cite{BackpropFreeRL}, but do not conduct in-depth experiments and ablations tantamount to a usable research benchmark. In particular, areas of i) Limited training data ii) Limited training epochs and iii) Fewer network parameters have not been benchmarked before in an in-depth manner. We bridge this gap and benchmark five popular and varied Bio-algorithms, namely, Hebbian Learning, FA, DFA, DTP and Predictive coding to give a wide-ranging coverage of the various groups of Bio-algorithms. We conduct experiments on multiple datasets in carefully setup experiments, to demonstrate empirically the areas where Bio-algorithms surpass BP and can be used by practitioners outside the Bio-learning community for their specific deep learning tasks.

\begin{figure*}[h]
\centering
\begin{minipage}{7cm}
\centering
\includegraphics[width=7cm]{Figures/cifar10-20-0.2.pdf}
\centerline{(a) Using 20\% of training set and 20 epochs}\medskip
\end{minipage}%
\begin{minipage}{7cm}
\centering
\includegraphics[width=7cm]{Figures/cifar10-20-1.pdf}
\centerline{(b) Using 100\% of training set and 20 epochs}\medskip
\end{minipage}
\begin{minipage}{1.8cm}
\centering
\centerline{}\medskip
\end{minipage}
\begin{minipage}{7cm}
\centering
\includegraphics[width=7cm]{Figures/cifar10-100-0.2.pdf}
\centerline{(c) Using 20\% of training set and 100 epochs}\medskip
\end{minipage}%
\begin{minipage}{7cm}
\centering
\includegraphics[width=7cm]{Figures/cifar10-100-1.pdf}
\centerline{(d) Using 100\% of training set and 100 epochs}\medskip
\end{minipage}
\begin{minipage}{1.8cm}
\centering
\includegraphics[width=1.8cm]{Figures/cifar10-legend.png}
\centerline{}\medskip
\end{minipage}
\caption{Experiments using limited data or limited training budget or both. Bio-algorithms perform much better than BP with Hebbian in particular surpassing BP by 20\% for the case of lesser data and fewer epochs (Fig. (a)). Hebbian also converges in around 5 epochs. Using CIFAR-10 dataset. Mean values plotted and standard deviation denoted by shaded regions over 10 runs.}
\label{fig:less_data_and_fewer_epochs}
\end{figure*}

\begin{figure*}[h]
\centering
\begin{minipage}{7cm}
\centering
\includegraphics[width=7cm]{Figures/cifar100-20-0.2.pdf}
\centerline{(a) Using 20\% of training set and 20 epochs}\medskip
\end{minipage}%
\begin{minipage}{7cm}
\centering
\includegraphics[width=7cm]{Figures/cifar100-20-1.pdf}
\centerline{(b) Using 100\% of training set and 20 epochs}\medskip
\end{minipage}
\begin{minipage}{1.8cm}
\centering
\centerline{}\medskip
\end{minipage}
\begin{minipage}{7cm}
\centering
\includegraphics[width=7cm]{Figures/cifar100-100-0.2.pdf}\textbf{s}
\centerline{(c) Using 20\% of training set and 100 epochs}\medskip
\end{minipage}%
\begin{minipage}{7cm}
\centering
\includegraphics[width=7cm]{Figures/cifar100-100-1.pdf}
\centerline{(d) Using 100\% of training set and 100 epochs}\medskip
\end{minipage}
\begin{minipage}{1.8cm}
\centering
\includegraphics[width=1.8cm]{Figures/cifar100-legend.png}
\centerline{}\medskip
\end{minipage}
\caption{Experiments using limited data or limited training budget or both. Bio-algorithms perform much better than BP with Hebbian in particular surpassing BP by 16\% for the case of lesser data and fewer epochs (Fig. (a)). Hebbian also converges in around 5 epochs. Using CIFAR-100 dataset. Mean values plotted and standard deviation denoted by shaded regions over 10 runs.}
\label{fig:less_data_and_fewer_epochs_cifar100}
\end{figure*}

\begin{figure*}[h]
\centering
\begin{minipage}{8cm}
\centering
\includegraphics[width=8cm]{Figures/Hebbian_filters.png}
\centerline{(a) Filters learnt by Hebbian learning}\medskip
\end{minipage}%
\hspace{1cm}
\begin{minipage}{8cm}
\centering
\includegraphics[width=8cm]{Figures/BP_filters.png}
\centerline{(b) Filters learnt by Backprop}\medskip
\end{minipage}
\caption{Visualizing the filters learnt by Hebbian learning vs. Backprop. Hebbian learning learns much more fine-grained representations compared to Backprop. Hebbian learns filters which are orientation-sensitive (vertical, horizontal and diagonal grayscale edges), color-sensitive (blue, green, red filters) or both (filters with combination of orientation and two or more colors). Plotting filters from the first CNN layer on the CIFAR-10 dataset after 100 epochs of training.}
\label{fig:representations}
\end{figure*}

\begin{figure*}[!h]
\centering
\begin{minipage}{.48\textwidth}
    \centering
    \includegraphics[width=\textwidth]{Figures/cifar10-random_noise.pdf}
    \caption{Experiment on adding Random noise to the input images. Bio-learning (especially DFA) is quite robust to Random noise. Using CIFAR-10 dataset. Mean values plotted and standard deviation denoted by shaded regions over 10 runs.}
    \label{fig:random_noise}
\end{minipage}
\hspace{0.5cm}
\begin{minipage}{.48\textwidth}
    \centering
    \includegraphics[width=\textwidth]{Figures/cifar10-pepper_noise.pdf}
    \caption{Experiment on adding Pepper noise to the input images. Bio-learning (especially DFA) is quite robust to Pepper noise as well. Using CIFAR-10 dataset. Mean values plotted and standard deviation denoted by shaded regions over 10 runs.}
    \label{fig:pepper_noise}
\end{minipage}
\end{figure*}

\begin{figure*}[!h]
\centering
\begin{minipage}{.48\textwidth}
    \centering
    \includegraphics[width=\textwidth]{Figures/cifar100-random_noise.pdf}
    \caption{Experiment on adding Random noise to the input images. Bio-learning (especially DFA) is robust, surpassing BP. Using CIFAR-100 dataset. Mean values plotted and standard deviation denoted by shaded regions over 10 runs.}
    \label{fig:random_noise_cifar100}
\end{minipage}
\hspace{0.5cm}
\begin{minipage}{.48\textwidth}
    \centering
    \includegraphics[width=\textwidth]{Figures/cifar100-pepper_noise.pdf}
    \caption{Experiment on adding Pepper noise to the input images. Bio-learning (especially DFA) is robust. Using CIFAR-100 dataset. Mean values plotted and standard deviation denoted by shaded regions over 10 runs.}
    \label{fig:pepper_noise_cifar100}
\end{minipage}
\end{figure*}

\section{Approach}

In this section, we introduce the three Bio-inspired learning algorithms used in our experiments, that is, Hebbian learning (Grossberg’s Instar Rule), Feedback Alignment (FA) and Direct Feedback Alignment (DFA).

\subsection{Hebbian Learning}

\quad \textbf{Vanilla Hebbian learning:} We first describe the vanilla Hebbian rule which states that when a neuron is fired by the stimulation of another neuron connected to it, the strength of the weighted connection between the two is enhanced, and vice versa \cite{10.1007/978-3-642-70911-1_15, lowel1992selection}. Hebbian learning is a completely local, unsupervised learning algorithm and  satisfies biological plausibility. The weight update is given by

\begin{equation}
\Delta{w_{ij}} = \eta z_jx_i
\label{eq1}
\end{equation}

where $w_{ij}$ is the weight connecting neurons \textit{i} and \textit{j}, $x_{i}$ and $z_{j}$ are the input and output (after activation) of neuron \textit{j} and $\eta$ is the learning rate.

However, this vanilla rule has a significant drawback that can cause the weight to grow infinitely. Therefore, variants of the vanilla Hebbian rule have been proposed to address this problem.

\textbf{Grossberg’s Instar Rule:} Instar rule \cite{grossberg1976adaptive} introduces the mechanism of weight decay for solving the problem of unbounded growth of weights. From a mathematical point of view, the introduction of the weight decay mechanism leads the weight vector $w_{ij}$ to cluster towards the input vector $x_i$ \cite{lagani2021hebbian}. Giving the same input to the neuron repeatedly will eventually make the angle between the weight vector and the input vector to be zero. In other words, the maximum inner product between $x_i$ and $w_{ij}$ will occur in the above case.

\begin{equation}
\Delta{w_{ij}} = \eta z_j(x_i-w_{ij}) \label{eq2}
\end{equation}

where the symbol interpretation is same as Equation \ref{eq1}.

\subsection{Feedback Alignment}

\quad FA was proposed to eliminate the weight transport mechanism relied on by Backprop (BP), as it is not considered biologically plausible in the brain. Weight transport implies that the weight matrix between neurons needs to be involved in operations for both forward and backward paths \cite{akrout2019deep}. Unlike the BP algorithm, FA replaces the transpose of the weight matrix that is required in the backward pass, with a fixed random matrix instead \cite{Lillicrap2016RandomSF}.

In the forward pass, the feed-forward calculation follows:

\begin{equation}
a_{i+1} = w_{i+1}z_i+b_{i+1} \label{eq3}
\end{equation}

\begin{equation}
z_{i+1} = \sigma (a_{i+1}) \label{eq4}
\end{equation}

where $z_i$ is the output of the $i$-th layer, $w_{i+1}$ is the feed-forward weight, $b_{i+1}$ is the bias vector, and $\sigma$ is the non-linear activation function.

In the backward pass of FA, the error $e$ follows:

\begin{equation}
e_{i} = (B_{i+1}e_{i+1}) \odot \sigma'(a_i) \label{eq5}
\end{equation}

where $B_{i+1}$ is the fixed random matrix to pass the error through, $\odot$ is an element-wise multiplication operator and $\sigma'$ is the derivative of the non-linear activation function.

Thus, the weight update becomes

\begin{equation}
\Delta{w_{i+1}} = -\eta e_{i+1} z_{i}^{T}
\label{eq6}
\end{equation}

\subsection{Direct Feedback Alignment}

\quad DFA is a modified version of FA that breaks the mechanism whereby errors need to be propagated layer by layer. Instead, it passes the errors directly from the final layer to the layer in question, circumventing any intermediate layers \cite{DFA}. This allows the algorithm to have greater biological plausibility.

The backward pass of DFA is slightly different from FA and follows:

\begin{equation}
e_{i} = (B_{i+1}e_{f}) \odot \sigma'(a_i) \label{eq7}
\end{equation}

where $e_{f}$ is the error after the final layer.

The weight update for DFA remains the same as that for for FA, and follows Equation \ref{eq6}.
\section{Experiments}
\label{results}

We conduct experiments on Hebbian Learning (HB), Feedback Alignment (FA), Direct Feedback Alignment (DFA) and Backprop (BP). Multiple datasets including CIFAR-10 \cite{Krizhevsky09learningmultiple} and CIFAR-100 \cite{Krizhevsky09learningmultiple} are used for testing. For CIFAR-100, two other Bio-algorithms Difference Target Propagation (DTP) and Predictive Coding (PC) are also benchmarked for a more in-depth comparison. We use Convolutional Neural Networks (CNNs) for the experiments owing to their higher accuracies than Fully Connected Networks (FCNs). For all experiments, 4 layered networks are used with the first 3 layers being convolutional layers having 100, 196 and 400 neurons respectively and the last layer being a linear layer for the case of BP, FA and DFA and ridge classifier for the case of Hebbian. For FA and DFA the implementation from \cite{10.3389/fnins.2021.629892} is used and for Hebbian learning the implementation from \cite{miconi} is used. For DTP the implementation from \cite{theoretical_framework} is used, called DDTP, while for PC the implementation from \cite{PC} is used. For fair comparison between Bio-learning and BP, vanilla networks are used for all experiments without drop-out, batch-norm or any other training tricks. Hyper-parameters used are provided in the Appendix. All experiments are conducted for 100 epochs unless stated otherwise.

We first conduct experiments on the limited data and limited training budget scenarios where only a subset of the original dataset or a smaller training budget is available for training (Fig. \ref{fig:less_data_and_fewer_epochs} (a)). As can be seen, all the three Bio-learning algorithms perform much better than BP. Hebbian learning, in particular, performs really well achieving 20\% higher accuracy than BP for the case of less data and fewer training epochs. Also, notably, Hebbian learning converges much faster than BP, usually achieving convergence in just around 5 epochs. A similar trend is seen when the epochs is increased to 100 (Fig. \ref{fig:less_data_and_fewer_epochs} (c)), where FA achieves the top accuracy amongst all the algorithms. When the condition of less data availability is relaxed (Fig. \ref{fig:less_data_and_fewer_epochs} (b)), even then Bio-algorithms surpass BP by a sizable margin. Only in the case of availability of the full dataset and the full training budget (Fig. \ref{fig:less_data_and_fewer_epochs} (d)), is BP able to get comparable accuracy as the top performing Bio-algorithm. Even then BP requires the full 100 epochs to learn whereas all the three Bio-algorithms finish learning within 40 epochs and within 5 epochs for Hebbian learning.

Similar results are seen for the case of CIFAR-100 (Fig. \ref{fig:less_data_and_fewer_epochs_cifar100}). Hebbian learning (HB) and Difference Target Propagation (DTP) perform the best and much better than BP (around 16\% better in terms of accuracy) for the case of limited data and limited training budget (Fig. \ref{fig:less_data_and_fewer_epochs_cifar100} (a)). Even with increasing the training budget, Bio-algorithms still surpass BP (Fig. \ref{fig:less_data_and_fewer_epochs_cifar100} (c)). Only when the full dataset is shown to BP is it able to achieve comparable accuracy to the top performing Bio-algorithm DTP (Fig. \ref{fig:less_data_and_fewer_epochs_cifar100} (b)). DTP also outperforms BP for the case of the full dataset and full training budget (Fig. \ref{fig:less_data_and_fewer_epochs_cifar100} (d)). Notably, Hebbian learning continues to learn extremely fast and converges to stable accuracy within 5 epochs where other algorithms take around 60 epochs to converge.

To understand this superior performance of Hebbian learning, we visualize the features learnt by Hebbian learning vs. BP. Referring to Fig. \ref{fig:representations}, we note that the filters learnt by Hebbian learning are much more fine-grained compared to BP. Hebbian learns filters which are orientation-sensitive (vertical, horizontal and diagonal grayscale edges), color-sensitive (blue, green, red filters) or both (filters with combination of orientation and two or more colors). BP on the other hand seems to learn more random-like features with no clear orientation or color-sensitivity. We believe this fine-grained representation ability of Hebbian might be the underlying reason for its better performance.

We also note that Bio-learning is quite robust to different kinds of noises. Referring to Fig. \ref{fig:random_noise} and Fig. \ref{fig:pepper_noise}, we see that Bio-learning (especially DFA) performs better than Backprop under the addition of random noise and pepper noise, respectively for CIFAR-10. For CIFAR-100, DFA performs the best and even better than BP (Fig. \ref{fig:random_noise_cifar100} and Fig. \ref{fig:pepper_noise_cifar100}). Bio-learning is also robust to network sparsification achieving similar results as Backprop even at very high sparsity rates of 95\% of all weights pruned, see Table \ref{table:sparsity_cifar10} and Table \ref{table:sparsity_cifar100}. We don't include Hebbian learning here since as per the implementation of Hebbian learning \cite{miconi}, it is already pruned by default and hence, we are unable to perform a fair comparison with the other networks.

\begin{table}
\centering
{\begin{tabular}{|l|c|c|}
\hline
\multicolumn{1}{|c|}{Model}& \multicolumn{2}{c|}{Accuracy} \\
\hline
\multicolumn{1}{|c|}{FA}& \multicolumn{2}{c|}{62.17\% $\pm$ 0.24\%} \\
\hline
\multicolumn{1}{|c|}{DFA}&  \multicolumn{2}{c|}{57.75\% $\pm$ 0.16\%} \\
\hline
\multicolumn{1}{|c|}{BP}&  \multicolumn{2}{c|}{65.37\% $\pm$ 0.65\%} \\
\hline
\end{tabular}}
\caption{Experiment on the application of 95\% weight sparsity to the networks. We find that Bio-learning is robust to network sparsification and achieves similar accuracy as Backprop. Using CIFAR-10 dataset. Mean and standard deviations reported over 10 runs.}
\label{table:sparsity_cifar10}
\end{table}

\begin{table}
\centering
{\begin{tabular}{|l|c|c|}
\hline
\multicolumn{1}{|c|}{Model} & \multicolumn{2}{c|}{Accuracy}\\
\hline
\multicolumn{1}{|c|}{FA}& \multicolumn{2}{c|}{{30.93 $\pm$ 0.16}}\\
\hline
\multicolumn{1}{|c|}{DFA}& \multicolumn{2}{c|}{30.29 $\pm$ 0.25}\\
\hline
\multicolumn{1}{|c|}{BP}& \multicolumn{2}{c|}{{37.68 $\pm$ 0.48}}\\
\hline
\end{tabular}}
\caption{Experiment on the application of 95\% weight sparsity to the networks. We find that Bio-learning is quite robust to network sparsification. Using CIFAR-100 dataset. Mean and standard deviations reported over 10 runs.}
\label{table:sparsity_cifar100}
\end{table}

\section{Discussion, Limitations and Future Work}
\label{Disc}

We have noted the superior performance of Bio-inspired learning compared to Backprop on several aspects. Under limited data or limited training budget Bio-inspired learning does much better and also converges much faster. We believe higher quality representations learnt by the Bio-algorithms are the underlying cause for this. This is an exciting direction and the research landscape for Bio-learning is still nascent. One of the biggest bottlenecks for Bio-learning is that the training tricks and techniques to maximize performance of Bio-learning are still under-developed. Backprop benefits from many years of combined research on optimizers, hyper-parameter tuning techniques and heuristics, custom layers like Batchnorm \cite{pmlr-v37-ioffe15} and other tricks like Dropout \cite{JMLR:v15:srivastava14a}. Such progress is severely lacking for Bio-learning making the discovery of the full potential of Bio-learning difficult. More research in these areas will help Bio-learning substantially.

Another limitation of Bio-learning is the lack of in-depth understanding on why it works so well. We have shown that the representations learnt by Hebbian learning are of higher quality than Backprop. However, the next question that arises is why does Hebbian learning learn higher quality representations? Is it because the representations learnt by Hebbian are more invariant to input orientations and configurations \cite{hebb_representation}?  Having theoretical and more in-depth understanding of these questions will help design better use-cases for Bio-learning. Currently, empirical results are the only guiding factor in deciding which areas might Bio-learning be most suitable for. However, a richer understanding of the learning dynamics will automatically expose a wider set of research areas and applications that Bio-learning might be useful for.
\section{Conclusion}

We have investigated Bio-inspired learning algorithms like Hebbian learning, Feedback Alignment, Direct Feedback Alignment, Difference Target Propagation and Predictive Coding in this work, and benchmarked them against Backprop on multiple aspects. We measure the performance of these algorithms under conditions of i) Limited availability of training data ii) Limited training budget iii) Application of different kinds of noises on the inputs and iv) Application of parameter sparsity on the base neural network models. Hebbian learning in particular, surpasses the accuracy performance of Backprop by upto 20\% under limited data and limited training budget scenarios. Bio-algorithms also achieve much faster convergence (upto 20$\times$ for the case of Hebbian learning) compared to Backprop. These insights help motivate the usefulness of Bio-algorithms in important areas like learning with less data or learning in resource-constrained training environments like real-time devices. They highlight the need to work on Bio-learning algorithms for tangible practical benefits and beyond the purely theoretical motivations of solving the biological implausibility of Backprop. Furthermore, we also plot the representations learnt by the Bio-algorithms to understand their superior performance, and find that they learn higher quality representations compared to Backprop. Lastly, we find that Bio-learning is robust to noise and parameter sparsity (even at extreme sparsities of 95\%), achieving similar or higher accuracy than Backprop.
\section{Appendix}

\subsection{Experiment design and hyper-parameters}

We describe below in detail the hyper-parameters and experimental setup for the various experiments. For Hebbian, we use the same experimental setup as the original and utilize ZCA, pruning and triangle method for computing Hebbian updates. We also list below tables capturing the hyper-parameters used for each experiment.

\vspace{1cm}
\centering
\small
\begin{tabular}{p{1.5cm}p{1.5cm}p{1.5cm}p{2.5cm}p{2.8cm}p{2cm}}
\hline
\multicolumn{6}{c}{\centering{Hyper-parameters for Fig. 2}} \\ \hline
Algorithm & Epochs & Batch-size & Initial LR & LR Schedule & Weight Decay\\
\hline
HB&20 or 100&100&1e-5&Not used&0.95\\
FA&20 or 100&100&5e-5&Not used&Not used\\
DFA&20 or 100&100&5e-5&Not used&Not used\\
BP&20 or 100&100&1e-5&Step LR, Step Size 1&Not used\\
\hline
\end{tabular}

\vspace{0.5cm}
\centering
\small
\begin{tabular}{p{1.5cm}p{1.5cm}p{1.5cm}p{2.5cm}p{2.8cm}p{2cm}}
\hline
\multicolumn{6}{c}{\centering{Hyper-parameters for Fig. 3}} \\ \hline
Algorithm & Epochs & Batch-size & Initial LR & LR Schedule & Weight Decay\\
\hline
HB&20 or 100&100&1e-5&Not used&0.95\\
FA&20 or 100&100&5e-5&Not used&Not used\\
DFA&20 or 100&100&5e-5&Not used&Not used\\
DDTP&20 or 100&100&3e-4, 9e-4, 1e-4, resp. for each layer&Not used&Not used\\
PC&20 or 100&100&1e-5&Not used&Not used\\
BP&20 or 100&100&1e-5&Step LR, Step Size 1&Not used\\
\hline
\end{tabular}

\vspace{0.5cm}
\centering
\small
\centering
\small
\begin{tabular}{p{1.5cm}p{1.5cm}p{1.5cm}p{2.5cm}p{2.8cm}p{2cm}}
\hline
\multicolumn{6}{c}{\centering{Hyper-parameters for Fig. 5}} \\ \hline
Algorithm & Epochs & Batch-size & Initial LR & LR Schedule & Weight Decay\\
\hline
HB&20 or 100&100&1e-5&Not used&0.95\\
FA&20 or 100&100&5e-5&Not used&Not used\\
DFA&20 or 100&100&5e-5&Not used&Not used\\
BP&20 or 100&100&1e-5&Step LR, Step Size 1&Not used\\
\hline
\end{tabular}

\vspace{0.5cm}
\centering
\small
\begin{tabular}{p{1.5cm}p{1.5cm}p{1.5cm}p{2.5cm}p{2.8cm}p{2cm}}
\hline
\multicolumn{6}{c}{\centering{Hyper-parameters for Fig. 6}} \\ \hline
Algorithm & Epochs & Batch-size & Initial LR & LR Schedule & Weight Decay\\
\hline
HB&20 or 100&100&1e-5&Not used&0.95\\
FA&20 or 100&100&5e-5&Not used&Not used\\
DFA&20 or 100&100&5e-5&Not used&Not used\\
BP&20 or 100&100&1e-5&Step LR, Step Size 1&Not used\\
\hline
\end{tabular}

\vspace{0.5cm}
\centering
\small
\begin{tabular}{p{1.5cm}p{1.5cm}p{1.5cm}p{2.5cm}p{2.8cm}p{2cm}}
\hline
\multicolumn{6}{c}{\centering{Hyper-parameters for Fig. 7}} \\ \hline
Algorithm & Epochs & Batch-size & Initial LR & LR Schedule & Weight Decay\\
\hline
HB&20 or 100&100&1e-5&Not used&0.95\\
FA&20 or 100&100&5e-5&Not used&Not used\\
DFA&20 or 100&100&5e-5&Not used&Not used\\
BP&20 or 100&100&1e-5&Step LR, Step Size 1&Not used\\
\hline
\end{tabular}

\vspace{0.5cm}
\centering
\small
\begin{tabular}{p{1.5cm}p{1.5cm}p{1.5cm}p{2.5cm}p{2.8cm}p{2cm}}
\hline
\multicolumn{6}{c}{\centering{Hyper-parameters for Fig. 8}} \\ \hline
Algorithm & Epochs & Batch-size & Initial LR & LR Schedule & Weight Decay\\
\hline
HB&20 or 100&100&1e-5&Not used&0.95\\
FA&20 or 100&100&5e-5&Not used&Not used\\
DFA&20 or 100&100&5e-5&Not used&Not used\\
BP&20 or 100&100&1e-5&Step LR, Step Size 1&Not used\\
\hline
\end{tabular}

\end{minipage}