General anesthesia reduces complexity and temporal asymmetry of the informational structures derived from neural recordings in Drosophila

TL;DR

This study uses computational mechanics to show that general anesthesia in fruit flies decreases neural complexity and temporal asymmetry, indicating a less structured and more reversible brain state.

Contribution

It introduces a novel application of computational mechanics to neural data, revealing how anesthesia reduces complexity and temporal asymmetry in Drosophila.

Findings

Neural complexity decreases under anesthesia across various brain regions.

Anesthesia reduces temporal asymmetry in neural signals.

Complexity differences are observable at short timescales.

Abstract

We apply techniques from the field of computational mechanics to evaluate the statistical complexity of neural recording data from fruit flies. First, we connect statistical complexity to the flies' level of conscious arousal, which is manipulated by general anesthesia (isoflurane). We show that the complexity of even single channel time series data decreases under anesthesia. The observed difference in complexity between the two states of conscious arousal increases as higher orders of temporal correlations are taken into account. We then go on to show that, in addition to reducing complexity, anesthesia also modulates the informational structure between the forward- and reverse-time neural signals. Specifically, using three distinct notions of temporal asymmetry we show that anesthesia reduces temporal asymmetry on information-theoretic and information-geometric grounds. In contrast to prior work, our results show that: (1) Complexity differences can emerge at very short timescales and across broad regions of the fly brain, thus heralding the macroscopic state of anesthesia in a previously unforeseen manner, and (2) that general anesthesia also modulates the temporal asymmetry of neural signals. Together, our results demonstrate that anesthetized brains become both less structured and more reversible.

Figures

Tables

Table 1. Table 1: Significant χ 2 superscript 𝜒 2 \chi^{2} and p 𝑝 p values of effects of channel c 𝑐 c , memory length λ 𝜆 \lambda , and condition ψ 𝜓 \psi , for informational quantities Q 𝑄 Q obtained via LME analysis. First order effects correspond to significant channel, memory, or condition responses on an informational quantity, while second and third-order effects correspond to interactions between these effects. χ 2 superscript 𝜒 2 \chi^{2} values are reported with n − 1 𝑛 1 n-1 degrees of freedom in the parentheses, corresponding to the number of effects removed under the null model, described in Sec.  III.1

$Q$	$1^{st}$ Order		$2^{nd}$ Order		$3^{rd}$ Order
$C_{\mu }^{+}$	$\lambda :{\chi }^2\left(1\right)=443.64$		$\lambda c:{\chi }^2\left(14\right)=42.876$	$p=8.97\times {10}^{-5}$	$\lambda c\psi :{\chi }^2\left(14\right)=24.00$	$p=0.0458$
			$\lambda \psi :{\chi }^2\left(1\right)=14.63$	$p=1.31\times {10}^{-4}$
$\Delta C_{\mu }^{+}$	$\lambda :{\chi }^2\left(1\right)=20.97$	$p=4.65\times {10}^{-6}$	$\lambda c:{\chi }^2\left(14\right)=37.19$	$p=6.90\times {10}^{-4}$
$\Xi$	$\psi :{\chi }^2\left(1\right)=4.870$	$p=0.0273$	$\lambda \psi :{\chi }^2\left(1\right)=5.565$	$p=0.0183$	$\lambda c\psi :{\chi }^2\left(14\right)=31.79$	$p=4.29\times {10}^{-3}$
	$\lambda :{\chi }^2\left(1\right)=6.725$	$p=9.51\times {10}^{-3}$
$d$	$\psi :{\chi }^2\left(1\right)=5.896$	$p=0.0152$	$\lambda \psi :{\chi }^2\left(1\right)=6.119$	$p=0.0134$
	$\lambda :{\chi }^2\left(1\right)=460.5$		$\lambda c:{\chi }^2\left(14\right)=35.86$	$p=1.09\times {10}^{-3}$
${𝒟}_{KLS}$	$\lambda :{\chi }^2\left(1\right)=127.4$		$\lambda c:{\chi }^2\left(14\right)=85.81$
			$\lambda \psi :{\chi }^2\left(1\right)=127.4$

Equations

$$\displaystyle\reflectbox{$\vec{\reflectbox{$r$}}$}_{\ell}\sim_{\epsilon}\reflectbox{$\vec{\reflectbox{$r$}}$}_{\ell}^{\prime}\quad\text{if}\quad P(\vec{r}\;|\;\reflectbox{$\vec{\reflectbox{$r$}}$}_{\ell})=P(\vec{r}\;|\;\reflectbox{$\vec{\reflectbox{$r$}}$}_{\ell}^{\prime}),$$

$$C_{\mu }:=H\left[\mathcal{S}\right]=-\sum _{i}P(S_i)\log P(S_i).$$

$$\Xi :=C_{\mu }^{+}-C_{\mu }^{-},$$

$$d:=2C_{\mu }^{\pm }-C_{\mu }^{+}-C_{\mu }^{-}.$$

$${\mathcal{D}}_{KLS}={\mathcal{D}}_{KL}(M^{+}\Vert M^{-})+{\mathcal{D}}_{KL}(M^{-}\Vert M^{+}),$$

$$\Delta Q:=Q^w-Q^a,$$

$$\displaystyle\mathcal{Q}=\mathbf{F}\mbox{\boldmath$\beta$}+\mathbf{R}\mathbf{b}+\mathcal{E}.$$

$$t=\frac{⟨\Delta Q{⟩}_f}{s_f/\sqrt{|f|}},$$

$$[\mathcal{C}{]}_{(i,j,k)}=C_{\mu }^{(i,j,k)}.$$

$$\displaystyle\mathbf{F}_{f}=\begin{pmatrix}\vec{\mbox{\boldmath$\lambda$}}&\mathbf{D}&\vec{\mbox{\boldmath$\Psi$}}_{W}&\lambda\mathbf{D}&\lambda\vec{\mbox{\boldmath$\Psi$}}_{W}&\mathbf{D}_{\Psi_{W}}&\lambda\mathbf{D}_{\Psi_{W}}\\ \vec{\mbox{\boldmath$\lambda$}}&\mathbf{D}&\vec{\mbox{\boldmath$\Psi$}}_{A}&\lambda\mathbf{D}&\lambda\vec{\mbox{\boldmath$\Psi$}}_{A}&\mathbf{D}_{\Psi_{A}}&\lambda\mathbf{D}_{\Psi_{A}}\end{pmatrix},$$

$${\Psi }_{W(A)}=\{\begin{array}{ll}1& \text{if ψ= wakeful (anaesthetised)}\\ 0& \text{otherwise}.\end{array}$$

Full text

I Introduction

Complex phenomena are everywhere in the physical world. Typically, these emerge from simple interactions among elements in a network, such as atoms making up molecules or organisms in a society. Despite their diversity, it is possible to approach these subjects with a common set of tools, using numerical and statistical techniques to relate microscopic details to emergent macroscopic properties Thurner et al. (2018). There has long been a trend of applying these tools to the brain, the archetypal complex system, and much of neuroscience is concerned with relating electrical activity in networks of neurons to psychological and cognitive phenomena Gazzaniga (2009). In particular, there is a growing body of experimental evidence Boly et al. (2013), that neural firing patterns can be strongly related to the level of conscious arousal in animals.

In humans, level of consciousness varies from very low in coma and under deep general anaesthesia, to very high in fully wakeful states of conscious arousal Laureys and Schiff (2012). With the current technology, precise discrimination between unconscious vegetative states and minimally conscious states are particularly challenging and remains a clinical challenge Laureys et al. (2015). Therefore, substantial improvement in accuracy of determining such conscious states using neural recording data will have significant societal impacts. Towards such a goal, neural data has been analysed using various techniques and notions of complexity to try to find the most reliable measure of consciousness Engemann et al. (2018); Sitt et al. (2014).

One of the most successful techniques to date in distinguishing levels of conscious arousal is the perturbational complexity index Massimini et al. (2005); Casali et al. (2013); Casarotto et al. (2016), which measures the neural activity patterns that follows a perturbation of the brain through magnetic stimulation. The evoked patterns are processed through a pipeline then finally summarised using Lempel-Ziv complexity Casali et al. (2013). This method is inspired by a theory of consciousness, called integrated information theory (IIT) Tononi (2004); Tononi et al. (2016), which proposes that a high level of conscious arousal should be correlated with the amount of so-called integrated information, or the degree of differentiated integration in a neural system (see Ref. Oizumi et al. (2014) for details). While there are various ways to capture this essential concept Mediano et al. (2019); Barrett and Seth (2011), one way to interpret integrated information is as the amount of loss of information a system has on its own future or past states based on its current state, when the system is minimally disconnected Tegmark (2016); Oizumi et al. (2016a, b).

These complexity measures, inspired by IIT, are motivated by the fundamental properties of conscious phenomenology, such as informativeness and integratedness of any experience Tononi (2004). While there are ongoing efforts to accurately translate these phenomenological properties into mathematical postulates Oizumi et al. (2014), such translation often contains assumptions about the underlying process which are not necessarily borne out in reality. For example, the derived mathematical postulates in IIT assume Markovian dynamics, i.e., that the future evolution of a neural system is determined statistically by its present state Barrett and Seth (2011). Moreover, IIT requires computing the correlations across all possible partitions between subsystems, which is computationally heavy Tegmark (2016) in relation to methods which do not require such partitioning to work. Assuming that the hierarchical causal influences in the brain would manifest as oscillations across a range of frequencies and spatial regions Buzsaki (2006), non-Markovian temporal correlations likely play a significant role in explaining any experimentally measurable behaviours, including the level of conscious arousal. There is therefore, scope for applying more general notions of complexity to meaningfully distinguish macroscopic brain states that support consciousness.

A conceptually simple approach to quantifying the complexity of time series data, such as the fluctuating potential in a neuron, is to construct the minimal model which statistically reproduces it. Remarkably, this minimal model, known as an epsilon machine ($\epsilon\text{-machine}$), can be found via a systematic procedure which has been developed within the field of computational mechanics Crutchfield and Young (1989); Travers and Crutchfield (2011); Crutchfield (2017). Crucially, $\epsilon\text{-machines}$ account for multiple temporal correlations contained in the data and can be used to quantify the statistical complexity of a process – the minimal amount of information required to specify its state. As such they have been applied over various fields, ranging from neuroscience Heinz Haslinger et al. (2009); Klinkner et al. (2006) and psychology Shalizi and Klinkner (2004) to crystallography Varn and Crutchfield (2004) and ecology Boschetti (2008), to the stock market Park et al. (2007). Lastly, unlike IIT the $\epsilon\text{-machine}$ analysis can be performed for data coming from a single channel.

In this paper, we use the statistical complexity derived from an $\epsilon\text{-machine}$ analysis of neural activity to distinguish states of conscious arousal in fruit flies (D. melanogaster). We analyse neural data collected from flies under different concentrations of isoflurane Cohen et al. (2016, 2018). By analysing signals from individual electrodes and disregarding spatial correlations, we find that statistical complexity distinguishes between the two states of conscious arousal through temporal correlations alone. In particular, as the degree of temporal correlations increases, the difference in complexity between the wakeful and anaesthetised states becomes larger. In addition to measuring complexity, the $\epsilon\text{-machine}$ framework also allows us to assess the temporal irreversibility of a process- the difference in the statistical structure of the process when read forwards vs. backwards in time. This may be particularly important for wakeful brains which are thought to be sensitive to the statistical structure of the environment which runs forward in time Cohen et al. (2018); Hohwy (2013); Tononi (2010). Using the nuanced characterisation of temporal information flow offered by the $\epsilon\text{-machine}$ framework Crutchfield et al. (2009), we then analyse the time irreversibility and crypticity of the neural signals to further distinguish the conscious states. We find that the asymmetry in information structure between forward and reverse-time neural signals is reduced under anaesthesia.

The present approach singularly differentiates between highly random and highly complex information structure; accounts for temporal correlations beyond the Markov assumption; and quantifies temporal asymmetry of the process. None of the standard methods possesses all of these features within a single unified framework. Before presenting these results in detail in Sec. III and discussing their implications in Sec. IV, we begin with a brief overview of the $\epsilon\text{-machine}$ framework we will use for our analysis.

II Theory: $\epsilon$-Machines and statistical complexity

To uncover the underlying statistical structure of neural activity that characterises a given conscious state, we treat the measured neural data, given by voltage fluctuations in time, as discrete time series. To analyse these time series, we use the mathematical tools of computational mechanics, which we outline in this section. We start with a general discussion on the ways to use time series data to infer a model of a system while placing $\epsilon\text{-machines}$ in this context. Next, we explain how we construct $\epsilon\text{-machines}$ in practice. Finally, we show how this can be used to extract a meaningful notion of statistical complexity of a process.

II.1 From time series to $\epsilon$-Machines

In abstract terms, a discrete-time series is a sequence of symbols $\mathbf{r}=(r_{0},\ldots,r_{k},\ldots)$ that appear over time, one after the other Rabiner (1989). Each element of $\mathbf{r}$ corresponds to a symbol from a finite alphabet $\mathcal{A}$ observed at the discrete time step labelled by the subscript $k$. The occurrence of a symbol, at a given time step, is random in general and thus the process, which produces the time series, is stochastic Doob (1953). However, the symbols may not appear in a completely independent manner, i.e., the probability of seeing a particular symbol may strongly depend on symbols observed in the past. These temporal correlations are often referred to as memory, and they play an important role in constructing models that are able to predict the future behaviour of a given stochastic process Gu et al. (2012).

Relative to an arbitrary time $k$, let us denote the future and the past partitions of the complete sequence as \mathbf{r}=(\reflectbox{\vec{\reflectbox{$r$}}},\vec{r}), where the past and the future are \reflectbox{\vec{\reflectbox{$r$}}}=(\ldots,r_{k-2},r_{k-1}) and $\vec{r}=(r_{k},r_{k+1},\ldots)$ respectively. In general, for the prediction of the immediate future symbol $r_{k}$, knowledge of the past $\ell$ symbols \reflectbox{\vec{\reflectbox{$r$}}}_{\ell}:=(r_{k-\ell},\ldots,r_{k-2},r_{k-1}), may be necessary. The number of past symbols we need to account for in order to optimally predict the future sequence is called the Markov order Gagniuc (2017).

In general, the difficulty of modelling a time series increases exponentially with its Markov order. However, not all distinct pasts lead to unique future probability distributions, leaving room for compression in the model. In a seminal work, Crutchfield and Young showed the existence of a class of models, which they called $\epsilon$-machines, that are provably the optimal predictive models for a non-Markovian process under the assumption of statistical stationarity Crutchfield and Young (1989); Travers and Crutchfield (2011). Constructing the $\epsilon$-machine is achieved by partitioning sets of partial past observations \reflectbox{\vec{\reflectbox{$r$}}}_{\ell} into causal states. That is, two distinct sequences of partial past observations \reflectbox{\vec{\reflectbox{$r$}}}_{\ell} and \reflectbox{\vec{\reflectbox{$r$}}}_{\ell}^{\prime} belong to the same causal state $S_{i}\in\mathcal{S}$, if the probability of observing a specific $\vec{r}$ given \reflectbox{\vec{\reflectbox{$r$}}}_{\ell} or \reflectbox{\vec{\reflectbox{$r$}}}_{\ell}^{\prime} is the same; that is

[TABLE]

where $\sim_{\epsilon}$ indicates that two histories correspond to the same causal state. The conditional probability distributions in Eq. (1) may always be estimated from a finite set of statistically stationary data via the naive maximum likelihood estimate, given by P(r_{k}|\reflectbox{\vec{\reflectbox{$r$}}}_{\ell})=\nu(r_{k},\reflectbox{\vec{\reflectbox{$r$}}}_{\ell})/\nu(\reflectbox{\vec{\reflectbox{$r$}}}_{\ell}), where $\nu(X)$ is the frequency of occurrence of sub-sequence $X$ in the data. For the case of non-stationary data, the probabilities obtained by this method will produce a non-minimal model that corresponds to a time-averaged representation of the time series. We now discuss how to practically construct an $\epsilon\text{-machine}$ for a given time series.

II.2 Constructing $\epsilon\text{-machines}$ with the CSSR algorithm

Several algorithms have been developed to construct $\epsilon\text{-machines}$ from time series data Tino and Dorffner (2001); Crutchfield and Young (1989, 1990). Here, we briefly explain the Causal State Splitting Reconstruction (CSSR) algorithm Shalizi and Klinkner (2004), which we use in this work to infer $\epsilon\text{-machines}$ predicting the statistics of neural data we provide as input.

The CSSR algorithm proceeds to iteratively construct sets of causal states accounting for longer and longer sub-sequences of symbols. In each iteration, the algorithm first estimates the probabilities P(r_{k}|\reflectbox{\vec{\reflectbox{$r$}}}_{\ell}) of observing a symbol conditional on each length $\ell$ prior sequence and compares them with the distribution $P(r_{k}|\mathcal{S}=S_{i})$ it would expect from the causal states it has so far reconstructed. If P(r_{k}|\reflectbox{\vec{\reflectbox{$r$}}}_{\ell})=P(r_{k}|\mathcal{S}=S_{i}) for some causal state, then \reflectbox{\vec{\reflectbox{$r$}}}_{\ell} is identified with it. If the probability is found to be different for all existing $S_{i}$, then a new causal state is created to accommodate the sub-sequence. By constructing new causal states only as necessary, the algorithm guarantees a minimal model that describes the non-Markovian behaviour of the data (up to a given memory length), and hence the corresponding $\epsilon\text{-machine}$ of the process.

The CSSR algorithm compares probability distributions via the Kolmogorov-Smirnov (KS) test Massey (1951); Hollander et al. (2013). The hypothesis that P(r_{k}|\reflectbox{\vec{\reflectbox{$r$}}}_{\ell}) and $P(r_{k}|\mathcal{S}=S_{i})$ are identical up to statistical fluctuations is rejected by the KS test at the significance level $\sigma$ when a distance $\mathcal{D}_{KS}$ 111The distance \mathcal{D}_{KS}=\max|F(r_{k}|\mathcal{S}=S_{i})-F(r_{k}|\reflectbox{\vec{\reflectbox{$r$}}}_{\ell})|, where $F(r_{k}|\mathcal{S}=S_{i})$ and F(r_{k}|\reflectbox{\vec{\reflectbox{$r$}}}_{\ell}) are cumulative distributions of $P(r_{k}|\mathcal{S}=S_{i})$ and P(r_{k}|\reflectbox{\vec{\reflectbox{$r$}}}_{\ell}) respectively. is greater than tabulated critical values of $\sigma$ Miller (1956). In other words, $\sigma$ sets a limit on the accuracy of the history grouping by parametrising the probability that an observed history \reflectbox{\vec{\reflectbox{$r$}}}_{\ell} belonging to a causal state $S_{i}$, is mistakenly split off and placed in a new causal state $S_{j}$. Our analysis, in agreement with Ref. Shalizi and Klinkner (2004), found that the choice of this value does not affect the outcome of CSSR within the tested range of $0.001<\sigma<0.01$. As a result, we set $\sigma=0.005$.

As it progresses, the CSSR algorithm compares future probabilities for longer sub-sequences, up to a maximum past history length of $\lambda$, which is the only important parameter that must be selected prior to running CSSR in addition to $\sigma$. If the considered time series is generated by a stochastic process of Markov order $\ell$, choosing $\lambda<\ell$ results in poor prediction because the inferred $\epsilon\text{-machine}{}$ cannot capture the long-memory structures present in the data. Despite this, the CSSR algorithm will still produce an $\epsilon\text{-machine}$ that is consistent with the approximate future statistics of the process up to order-$\lambda$ correlations Shalizi and Klinkner (2004). Given sufficient data, choosing $\lambda\geq\ell$ guarantees convergence on the true $\epsilon\text{-machine}{}$. One important caveat to note is that the time complexity of the algorithm scales asymptotically as $\mathcal{O}(|\mathcal{A}|^{2\lambda+1})$, putting an upper limit to the longest history length that is computationally feasible to use. Furthermore, the finite length of the time series data implies an upper limit on an ‘acceptable’ value of $\lambda$. Estimating P(r_{k}|\reflectbox{\vec{\reflectbox{$r$}}}_{\lambda}) requires sampling strings of length $\lambda$ from the finite data sequence. Since the number of such strings grows exponentially with $\lambda$, a value of $\lambda$ that is too long relative to the size $N$ of the data, will result in a severely under-sampled estimation of the distribution. A distribution P(r_{k}|\reflectbox{\vec{\reflectbox{$r$}}}_{\lambda}) that has been estimated from an under-sampled space is almost always never equal to $P(r_{k}|\mathcal{S}=S_{i})$, resulting in the algorithm creating a new causal state for every string of length $\lambda$ it encounters. A bound for the largest permissible history length is $L(N)\geq\log_{2}N/\log_{2}|\mathcal{A}|$, where $L(N)$ denotes maximum length for a given data size of $N$ Marton and Shields (1994); Cover and Thomas (1991). Once these considerations have been taken into account, the $\epsilon\text{-machine}$ produced by the algorithm provides us with a meaningful quantifier of the complexity of the process generating the time series, as we now discuss.

II.3 Measuring the complexity and asymmetry of a process

The output of the CSSR algorithm is the set of causal states and rules for transitioning from one state to another. That is, CSSR gives a Markov chain represented by a digraph Crutchfield and Young (1989); Gagniuc (2017) $G(V,E)$ consisting of a set of vertices $v_{i}\in V$ and directed edges $\{i,j\}\in E$, e.g. Figs. 1(c) and (d). Using these rules, one can find $P(S_{i})$, which represents the probability that the $\epsilon\text{-machine}$ is in the causal state $S_{i}$ at a any time. The Shannon entropy of this distribution quantifies the minimal number of bits of information required to optimally predict the future process; this measure, first introduced in Ref. Crutchfield and Young (1989), is called the statistical complexity:

[TABLE]

Formally, the causal states of a time series depend upon the direction in which the data is read Crutchfield et al. (2009). The main consequence of this result is that the set of causal states obtained by reading the time series in the forward direction $\mathcal{S}^{+}$, are not necessarily the same as those obtained by reading the time series in the reverse direction $\mathcal{S}^{-}$. Naturally, this corresponds to potential differences in forward and reverse-time processes and the associated complexities, which is known as causal irreversibility

[TABLE]

capturing the time-asymmetry of the process.

Another (stronger) measure of time-asymmetry is crypticity:

[TABLE]

This quantity measures the amount of information hidden in the forwards and reverse $\epsilon\text{-machines}$ that is not revealed in the future or past time series, respectively. Specifically, it combines the information that must be supplemented to determine the forwards $\epsilon\text{-machine}$ given the reverse $\epsilon\text{-machines}$ and the information to determine reverse $\epsilon\text{-machines}$ given the forwards $\epsilon\text{-machine}$. In each case, this is equivalent to the difference between the complexity of a bidirectional $\epsilon\text{-machine}$, denoted $C_{\mu}^{\pm}$ Crutchfield et al. (2009), and that of the corresponding unidirectional machine. Throughout this manuscript, we implicitly refer to the usual forward-time statistical complexity $C_{\mu}^{+}$ when writing $C_{\mu}$, unless otherwise stated.

Finally, an operational measure for time-asymmetry is defined by the microscopic irreversibility, which quantifies how statistically distinguishable the forwards and reverse $\epsilon\text{-machines}$ are, in terms of the sequences of symbols they produce. If the forward-time $\epsilon\text{-machine}$ produces the same sequences with similar probabilities to the reverse-time $\epsilon\text{-machine}$, then the process is reversible. Should a sequence available to $M^{+}$ be impossible for $M^{-}$ to produce, then the process is strictly irreversible. Here, we assess the distinguishability between two $\epsilon\text{-machines}$ by estimating the asymptotic rate of (symmetric) Kullback-Leibler (KL) divergence $\mathcal{D}_{KLS}$ between long output sequences; this measure is commonly applied to stochastic models Yang et al. (2019). Specifically

[TABLE]

where $\mathcal{D}_{KL}$ is the regular, non-symmetric estimated KL divergence rate Rached et al. (2004). The KL divergence can be proved to be a unique measure that satisfies all of the theoretical requirements of information-geometry Amari (2016); Oizumi et al. (2016a); Amari et al. (2018).

A few remarks are in order: in general, any one of the above measures vanishing does not imply that the other measures must also vanish. For instance, consider the case where the structures of the forward ($M^{+}$) and reverse-time ($M^{-}$) $\epsilon\text{-machines}$ are different but they happen to have the same complexities, i.e., $C_{\mu}^{+}=C_{\mu}^{-}$. Then, clearly we have $\Xi=0$ but $d\neq 0$ and $\mathcal{D}_{KLS}\neq 0$. On the other hand, consider the case when $M^{+}$ and $M^{-}$ are the same; here, we have $\Xi=\mathcal{D}_{KLS}=0$, yet may not have $d=0$. This means that vanishing $\mathcal{D}_{KLS}$ implies that $\Xi=0$ (but not the converse, and not $d=0$). This turns out be an interesting extremal case because, while the forward and reverse processes are identical, the non-vanishing crypticity accounts for the information required to synchronise the corresponding $\epsilon\text{-machines}$, i.e., producing the joint statistics of the paired $\epsilon\text{-machines}$. Moreover, we can conclude that microscopic irreversibility is a stronger measure than causal irreversibility; this comes at the expense of computational cost, i.e., the former is harder to compute than the latter. In essence, each measure above represents a different notion of temporal asymmetry, with its own operational significance. Causal irrversibility and crypticity are information-theoretic constructs, while microscopic irrversibility is a information-geometric construct.

In the next section, we describe the experimental and analytical methods, as well as the results: that the statistical complexity and temporal asymmetry of the neural time series, taken from fruit flies, significantly differ between states of conscious arousal.

III Experimental results and analysis

III.1 Methods

We analysed local field potential (LFP) data from the brains of awake and isoflurane-anaesthetised D. melanogaster (Canton S wild type) flies. Here, we briefly provide the essential experimental outline that is necessary to understand this paper. The full details of the experiment are presented in Refs. Cohen et al. (2018, 2016). LFPs were recorded by inserting a linear silicon probe (Neuronexus 3mm-25-177) with 16 electrodes separated by 25 $\mu$m. The probe covered approximately half of the fly brain and recorded neural activity as illustrated in Fig.  1(a). A tungsten wire inserted into the thorax acted as the reference. The LFPs at each electrode were recorded for 18s while the fly was awake and 18s more after the fly was anaesthetised (isoflurane, 0.6% by volume, through an evaporator). Flies’ unresponsiveness during anaesthesia was confirmed by the absence of behavioural responses to a series of air puffs, and recovery was also confirmed after isoflurane gas was turned off Cohen et al. (2016).

We used data sampled at 1kHz for the analysis Cohen et al. (2016), and to obtain an estimate of local neural activity, the 16 electrodes were re-referenced by subtracting adjacent signals giving 15 channels which we parametrise as $c\in[1,15]$. Line noise was removed from the recordings, followed by linear de-trending and removing the mean. The resulting data is a fluctuating voltage signal, which is time-binned (1ms bins) and binarised by splitting over the median, leading to a time series, see Fig. 1(b).

For each of the 13 flies in our data set, we considered 30 time series of length $N=18,000$. These correspond to the 15 channels, labelled numerically from the central to peripheral region as depicted in Fig. 1(a), and the two states of conscious arousal. Using the CSSR algorithm Shalizi and Klinkner (2004), we constructed $\epsilon\text{-machines}$ for each of these time series as a function of maximum memory length within the range $\lambda\in[2,11]$, measured in milliseconds. This is below the memory length $L(N)\sim 14$ beyond which we would be unable to reliably determine transition probabilities for a sequence of length $N$ (see Sec. II.2) 222$L(N)\sim 14$ only serves as a lower bound on $\lambda$, past which CSSR is guaranteed to return incorrect causal states for the neural data. In practice, this may occur at even lower memory lengths than this limit. We observed this effect marked by an exponential increase in the number of inferred causal states for $\lambda>11$, and thus excluded these memory lengths from the study.. For a given time direction $\xi\in\{+:\text{forward},-:\text{reverse}\}$, we recorded the resulting $3,900$ $\epsilon\text{-machine}$ structures and their corresponding statistical complexities $C_{\mu}^{(\xi,\psi)}$, and grouped them according to their respective level of conscious arousal, $\psi\in\{w,a\}$ for awake and anaesthesia, channel location, $c$, and maximum memory length, $\lambda$. Thus, the statistical complexity we computed in a given time direction is a function of the set of parameters $\{\psi,c,\lambda\}$ for each fly, $f$. We also determined the irreversibility $\Xi$, crypticity $d$, and symmetric KL divergence rate $\mathcal{D}_{KLS}$ for each fly and again grouped them over the same set of parameters $\{\psi,c,\lambda\}$. While we found that not all the data is strictly stationary, in that the moving means of the LFP signals were not normally distributed, the conclusions we draw from them are still broadly valid. As mentioned in Sec. II.1, $\epsilon\text{-machines}$ reconstructed from approximately stationary data are time-averaged models, and are likely to underestimate the true statistical complexity of the corresponding neural processes.

We are principally interested in the differences the informational quantities $Q^{\psi}\in\{C_{\mu}^{(\xi,\psi)},\;\Xi^{\psi},\;d^{\psi},\;\mathcal{D}^{\psi}_{KLS}\}$ have over states of conscious arousal and thus consider

[TABLE]

for fixed values of $\{f,c,\lambda\}$. Positive values of $\Delta Q$ indicate higher complexities observed in the wakeful state relative to the anaesthetised one. Finally, we use the notation $\langle Q^{\psi}\rangle_{x}$ to denote taking an average of any information quantity $Q^{\psi}$, over a specific parameter $x\in\{f,c,\lambda\}$. For example $\langle\Delta C_{\mu}^{+}\rangle_{f}$ means taking the fly-averaged difference in forward-time statistical complexity.

To assess the significance of each of the parameters $\psi$, $c$, $\lambda$, and $\xi$, or some combination of them, have on the response of the elements in the set $Q$ across flies, we conducted a statistical analysis using linear mixed effects modelling Harrison et al. (2018) (LME). The LME analysis describes the response of a given quantity $\mathcal{Q}$ by modelling it as a multidimensional linear regression of the form

[TABLE]

The resulting model in Eq. (7) consists of a family of equations where $\mathcal{Q}$ is the vector allowing for different responses of a quantity $Q$ for each specific fly, channel location, level of conscious arousal, and time direction where applicable. Memory length $\lambda$, channel location $c$, state of conscious arousal $\psi$, and time direction $\xi$ (again, where applicable) are the parameters that $\mathcal{Q}$ responds to. To account for variations in the response caused by interactions between parameters (e.g. between memory length and channel location), we included them in the model. Letting $X=\{\lambda,c,\psi,\xi\}$ be the set of the parameters which may induce responses, we can write all the non-empty $k$-combinations between them as $\mathcal{F}=\{\lambda,c,\psi,\xi,\lambda c,\lambda\psi,...,\lambda c\psi\xi\}=\binom{X}{k}\backslash\emptyset$. The elements in $\mathcal{F}$ are known as the fixed effects of Eq. (7), and are contained as elements within the matrix $\mathbf{F}$. The vector $\beta$, contains the regression coefficients describing the strength of each of the fixed effects $\mathcal{F}$.

In addition to fixed effects affecting the response of an element of $Q$ in our experiment, we also took into account any variation in response caused by known random effects. In particular, we expected stronger response variations to be caused by correlations occurring between the channels within a single fly, compared to between channels across flies. These random effects are contained as elements of the matrix $\mathbf{R}$, and the vector $\mathbf{b}$ encodes the regression coefficients describing their strengths. Finally, the vector $\mathcal{E}$ describes the normally-distributed unknown random effects in the model. The regression coefficients contained in the vectors $\beta$ and $\mathbf{b}$, were obtained via maximum likelihood estimation such that $\mathcal{E}$ are minimised. The explicit form of Eq. (7) used in this analysis is detailed in the Appendix IV.1.

With the full linear mixed effects model given by Eq. (7), we tested the statistical significance of a fixed effect in $\mathcal{F}$. This was accomplished by comparing the log-likelihood of the full model with all fixed effects, to the log-likelihood of a reduced model with the effect of interest removed Bates et al. (2015) (regression coefficients associated with the effect are removed). This comparison between the likelihood models is given by $\Lambda=2(h_{\text{full}}-h_{\text{reduced}})$, where $\Lambda$ is the likelihood ratio, $h_{\text{full}}$ is the log-likelihood of full model, and $h_{\text{reduced}}$ is the log-likelihood of the model with the effect of interest removed.

Under the null hypothesis, when a fixed effect does not have any influence on an informational quantity $Q$, i.e., the regression coefficients for the effect are vanishing, the likelihood ratio $\Lambda$ is $\chi^{2}$ distributed with degrees of freedom equal to the difference in the number of coefficients between the models. Therefore, we considered any fixed effect in the set $\mathcal{F}$ to have a statistically significant effect on a quantity if the probability of obtaining the likelihood ratio given the relevant $\chi^{2}$ distribution was less than 5% ($p<0.05$). Thus, for each significant effect we report the fixed effect being tested, i.e., an element of $\mathcal{F}$, the obtained likelihood ratio $\chi^{2}(n-1)$ with $n$ associated degrees of freedom, and the associated probability $p$ of obtaining the statistic under the null hypothesis.

In addition to all the quantities in the set $Q^{\psi}$, the LME and likelihood ratio test was also performed for $\Delta Q$, in order to find the significant interaction effects of the parameters. Here, we also modelled $\Delta Q$ as dependent on a fixed effect in $\mathcal{F}$ as in Eq. (7), but excluding the parameter $\psi$ as it was already implicitly considered. Once the significant effects of memory length, level of conscious arousal, and channel location were characterised with our statistical analysis, we followed with post-hoc, paired $t$-tests for elements in $Q^{\psi}$ given by

[TABLE]

where $s_{f}$ is the standard deviation of $\langle\Delta Q\rangle_{f}$, and $|f|=13$ is the sample size. The paired $t$-tests examine the nature of interactions between the parameters on a given quantity over the two states of conscious arousal. Positive $t$-scores indicate a quantity is larger for the wakeful state. We present the results of these analyses in the following sections, sorted categorically by whether time-direction is considered.

III.2 Results

III.2.1 Forward-time complexity results

In order to observe the effects of isoflurane on neural complexity, we began by visually inspecting the structure of the reconstructed $\epsilon\text{-machines}$ for the two levels of conscious arousal for the forward-time direction. We took special interest in observing the differences in the characteristics of the two groups of $\epsilon\text{-machines}$ heralding the two levels of conscious arousal. Here, memory length $\lambda$ plays an important role. At a given $\lambda$, the maximum number of causal states that may be generated scales according to $|\mathcal{A}|^{\lambda}$ Shalizi and Klinkner (2004). In our case, the alphabet is binary, $\mathcal{A}=\{0,1\}$. This greatly restricts the space of $\epsilon\text{-machine}$ configurations available for short history lengths Johnson et al. (2010). For $\lambda=2$ we can observe up to four distinct configurations, which is unlikely to reveal the difference based on conscious states. Given the previous findings Cohen et al. (2018), we generally expected that the data from the wakeful state would present more complexity than those from the anaesthetised state.

Visual inspection of the directed graphs indeed suggested higher $\epsilon\text{-machine}$ complexity during the wakeful state compared to the anaesthetised state, at a given set of parameters $\{f,c,\lambda\}$. In particular, the data from the anaesthetised state tended to result in fewer causal states and overall reduced graph connectivity. Panels (c) and (d) of Fig. 1 are examples of $\epsilon\text{-machines}$ (channel 1 data recorded from fly 1, at maximum memory length $\lambda=3$), where a simpler $\epsilon\text{-machine}$ is derived from the data under the anaesthetised condition. Differentiating between two conscious arousal states by visual inspection quickly becomes impractical because of the large number of $\epsilon\text{-machines}$. Moreover, for large values of $\lambda$ the number of causal states is exponentially large and it becomes difficult to see the difference in two graphs. To overcome these challenges, we looked at a simpler index, the statistical complexity $C_{\mu}$, to differentiate between conscious arousal states. To systematically determine the relationships between $C_{\mu}$ and the set of variables $\{c,f,\psi\}$ we employed the LME analysis outlined in Sec. III.1. We first tested whether $\lambda$ significantly affects $C_{\mu}$. We found $\lambda$ to indeed have a significant effect on $C_{\mu}$ ($\lambda$, ${\chi}^{2}(1)=443.64$, $p<10^{-16}$). Fig. 2 shows that independent of the conscious arousal condition or channel location, $C_{\mu}$ increases with larger $\lambda$. This indicates that the Markov order of the neural data is much larger than the largest memory length ($\lambda=11$) we consider. Nevertheless, we have enough information to work with.

We then sought to confirm if the complexity of $\epsilon\text{-machines}$ during anaesthesia are reduced, as suggested from visual inspection. Our statistical analysis indicates that $C_{\mu}$ is not invariably reduced during anaesthesia ($\psi$, ${\chi}^{2}(1)=0.212$, $p=0.645$) at all levels of $\lambda$ and all channel locations. This means that $C_{\mu}$ cannot simply indicate the causal arousal state without some additional information about time ($\lambda$) or space ($c$). In addition, we found that neither $c$ alone nor $c\psi$ strongly affects $C_{\mu}$. However, we found significant reductions in complexity when either the level of conscious arousal or the channel location, interacted with memory length ($\lambda\psi$, ${\chi}^{2}(1)=14.63$, $p=1.31\times 10^{-4}$) and (${c\lambda}$, $\chi^{2}(14)=42.876$, $p=8.97\times 10^{-5}$) respectively. Moreover, the three-way interaction also had a strong effect ($\lambda\psi c$, ${\chi}^{2}(14)=24.00$, $p=0.0458$).

As the three-way interaction between $\lambda$, $\psi$, and $c$ complicates interpretation of their effects, we performed a second LME analysis where we modelled $\Delta C_{\mu}$ instead of $C_{\mu}$, thus accounting for $\psi$ implicitly. In doing so, we investigated whether the change in statistical complexity due to anaesthesia is affected by memory length $\lambda$ or channel location $c$. Using this model, we found a non-significant effect of $c$ on $\Delta C_{\mu}$, while a significant effect of $\lambda$ on $\Delta C_{\mu}$ was seen ($\lambda$, ${\chi}^{2}(1)=20.97$, $p=4.65\times 10^{-6}$), indicating that $\Delta C_{\mu}$ overall changes with $\lambda$. Specifically, $\Delta C_{\mu}$ tended to increase with larger $\lambda$ when ignoring channel location, as is evident in Fig. 3 (top). Further, explaining our previous interaction between $\lambda$ and $\psi$, $\Delta C_{\mu}$ was not clearly larger than [math] for small memory length ($\lambda=2$; the top panel of Fig. 3). This suggests that the information to differentiate between states of conscious arousal is contained in higher order correlations. We also found that the interaction between ${\lambda}$ and channel location has a significant effect on $\Delta C_{\mu}$ ($\lambda c$, ${\chi}^{2}(14)=37.19$, $p=6.90\times 10^{-4}$), indicating that the effect of $\lambda$ is not constant across channels. Given that $\Delta C_{\mu}$ overall increases with $\lambda$, we considered that that the largest $\Delta C_{\mu}$ should occur at the largest $\lambda$. Fig. 3 (bottom) examines $\Delta C_{\mu}$ across channels at $\lambda=11$.

To further break down the interaction between $\lambda$ and $c$, we performed a one sample $t$-test at each value of memory length and channel location to find regions in the parameter space $(\lambda,c)$ where $C_{\mu}$ reliably differentiates wakefulness from anaesthesia across flies. We plot the $t$-statistic at each parameter combination in the top-left panel of Fig. 4, outlining regions in the parameter space where $\Delta C_{\mu}$ is significantly greater than [math] (with $p<0.05$, uncorrected, two-tailed), finding that the majority of the significance map is directed towards positive values of the $t$-statistic. However, only a subset of $(\lambda,c)$ cells contain values which are significantly different from [math]. Interestingly, we observed that for $\lambda=2$, $\Delta C_{\mu}$ is actually significantly negative, corresponding to greater complexity during anaesthesia, not during wakefulness. This marks $\lambda=2$ as anomalous relative to other levels of $\lambda$, and this reversal of the direction of the effect of anaesthesia likely contributed to the interaction between $\lambda$ and $\psi$.

Disregarding $\lambda=2$, we find $\Delta C_{\mu}$ to be significantly greater than [math] for channels 1, 3, 5-7, 9, 10, and 13, at varying levels of $\lambda$. As expected from our reported interaction between $\lambda$ and $c$, we observe $\Delta C_{\mu}$ to already be significantly greater than [math] at small $\lambda$ for channels 5-7, while $\Delta C_{\mu}$ only became significantly greater at larger $\lambda$ for channels 1, 3, 10 and 13. Further, other channels such as the most peripheral channel ($c=15$) did not have $\Delta C_{\mu}$ significantly greater than [math] at any $\lambda$. All significance results, due to LME tests, are reported in Table 1.

The above results suggest that the measured difference in complexity is present across various brain regions (top-left panel of Fig. 4), and that it grows as longer temporal correlations are taken into account (up to the largest value $\lambda=11$ tested). While Fig. 3 shows a continued increase in the difference of statistical complexity, $\Delta C_{\mu}$, as a function of $\lambda$, we did not analyse longer history lengths, due to limitations in the amount of the data and stability of the estimation of $C_{\mu}$. In addition to this general observation of increasing $\Delta C_{\mu}$ with $\lambda$, we observe that, remarkably, some brain regions discriminate the conscious arousal states with a history length of only 3. One trivial explanation for this effect is that under anaesthesia, the required memory length is indeed $\lambda=2$, while the optimal $\lambda$ for awake is much larger. However, a quick observation of Fig. 2 rules out this simple possibility; under both wakeful and anaesthetised states, $C_{\mu}$ continues to increase.

It is likely, however, that the tested range for $\lambda$ remains below the Markov order of the neural data; this is clearly indicated by the lack of a plateau in statistical complexity in Fig. 3. This suggests that we are far from saturating the Markov order of the process, and with more data we would be able to further distinguish between the two states. Future analyses with longer time series would also contribute to our understanding of the Markov order (maximum memory length) differences between the two states of conscious arousal. Nevertheless, our results, in Figs. 3 and 4, demonstrate that saturation of Markov order is not required for discrimination between conscious arousal states. This finding has a practical implication about the empirical utility of $\epsilon\text{-machines}$; even if the history length is too low, the inferred $\epsilon\text{-machine}$ and its statistical complexity can be useful. We now discuss the temporal asymmetry of neural processes.

III.2.2 Temporal asymmetry

Unlike other complexity measures, we obtain a distinct $\epsilon\text{-machine}$ from each given time series, and for each direction we read the time series, i.e., forward or backward in time. Based on the notion that wakeful brains should be better at predicting the next sensory input Cohen et al. (2018), we expect that anaesthesia should alter the information structures depending on the time direction. Our expectation translates to the following three hypotheses:

Causal irreversibility ($\Xi:=C_{\mu}^{+}-C_{\mu}^{-})$, which is purely based on the summary measure of statistical complexity, should be higher for awake but lower for anaesthetised brains; 2. 2.

Crypticity ($d:=2C_{\mu}^{\pm}-C_{\mu}^{+}-C_{\mu}^{-}$) should be higher for wakeful than anaesthetised brains; 3. 3.

Symmetric KL divergence rate ($\mathcal{D}_{KLS}$) should behave similarly.

On visual inspection of the variation in $\Xi$ for the wakeful and anaesthetised conditions, both appeared close to zero, suggesting that $\Xi$ would not have significant dependence on the condition. This impression was confirmed statistically with two-tailed $t$-tests against zero with corrections for multiple comparisons, as shown in the top-right panel of Fig. 4. Thus, Hypothesis 1 above, that irreversibility should be higher for wakeful over anaesthetised brains, is not supported by the data. However, as mentioned earlier, vanishing $\Xi$ does not imply that either $d=0$ or $\mathcal{D}_{KLS}=0$. To rule out the possibility that the information structure of $\epsilon\text{-machines}$ are different when read forwards, as opposed to backwards, depending on the condition, we also tested the latter two hypotheses.

With respect to crypticity, first, visual inspection of the two-tailed $t$-score map, which compares crypticity for the wakeful $d^{w}$ and anaesthetised $d^{a}$ conditions (bottom-left panel of Fig. 4) strongly implies that crypticity is larger in the former compared to the latter. This difference is largest over channels 5-7 and 9. To systematically evaluate this impression, we used LME statistical analysis (described in Sec. III.1) to determine the relationships between crypticity, $d$, and the set of variables $\{c,\lambda,\psi\}$ we employ. As expected, we found that both memory length ($\lambda$) and level of conscious arousal ($\psi$) significantly affected crypticity ($\lambda$, $\chi^{2}(1)=470.5$, $p<10^{-16}$) and ($\psi$, $\chi^{2}(1)=5.896$, $p=0.0152$) respectively. Crypticity also depended on a significant interaction between memory length and condition ($\lambda\psi$, $\chi^{2}(1)=6.119$, $p=0.0134$). Specific increases in crypticity around the middle brain region (bottom-left panel of Fig. 4) were also evident, with a strong interaction between channel location and memory length ($\lambda c$, $\chi^{2}(14)=35.86$, $p=1.09\times 10^{-3}$), which is similar to the result obtained for $\Delta C_{\mu}$. This LME analysis, together with the direction of effects in Fig. 4 (bottom-left) strongly confirms our Hypothesis 2.

Furthermore, as a more direct measure of microscopic structure, we also analysed the symmetric KL divergence rate, $\mathcal{D}_{KLS}$. Again, the two-tailed $t$-score map (Fig. 4, bottom-right panel) showed support for our hypothesis. Our formal statistical analysis with LME confirmed a critical interaction between memory length and condition ($\lambda\psi$, $\chi^{2}(1)=15.37$, $p<10^{-16}$), meaning that time-asymmetric information structure is lost due to anaesthesia, especially when a long memory length is taken into account. (We also note other significant effects: mainly the effect of memory length ($\lambda$, $\chi^{2}(1)=127.4$, $p<10^{-16}$) and interaction between memory length and channel location ($\lambda c$, $\chi^{2}(14)=85.81$, $p<10^{-16}$). Again, all significant results, due to LME tests, are reported in Table 1.

Taken together, these results show that the relative complexity of the forward versus reverse direction, as measured by causal irreversibility, does not distinguish between the wakeful and anaesthetised states. However, our crypticity results demonstrate that, under anaesthesia, the structures of the forward and reverse processes are relatively similar, whereas during wakefulness their structures differ. Fig. 5 demonstrates this effect with exemplar $\epsilon\text{-machines}$ reconstructed from a representative channel, from which we derived six distinct $\epsilon\text{-machines}$: three for wakeful (a, c, d), and three for anaesthetised (e-g) flies. Panel (b) shows how the time series and the transitions in the causal states of forward, reversed, and bidirectional $\epsilon\text{-machines}$ are related.

Our finding, that causal irreversibilities were not above zero for wakeful brains, corresponds to the fact that complexities of forward and reverse $\epsilon\text{-machines}$ were not significantly different. However, the bidirectional $\epsilon\text{-machines}$ for the wakeful condition were substantially more complex than those for the anaesthetised condition. The statistical complexity of bidirectional $\epsilon\text{-machines}$ should equal that of forward or reverse $\epsilon\text{-machines}$ if the process is completely time symmetric and deterministic Crutchfield et al. (2009), resulting in zero crypticity. However, for non-deterministic processes, additional information for synchronising the forward and reversed process may be needed, which would mean $d>0$. For instance, in Fig. 5(e-f), if we are told that the forward machine is in causal state $A$, we need extra information to determine whether the reversed machine is in causal state $X$ or $Y$. Yet, the detailed structure of the forward and reverse machines are the same in this example. Our analysis is supplemented with a study of the symmetric KL divergence rate between forwards and backwards processes, which measures the distance between the reconstructed $\epsilon\text{-machines}$. In other words, crypticity and symmetric KL divergence rate quantify two different notions of temporal asymmetry; the former is information theoretic, and the latter is information-geometric. Indeed, in general we find that the processes in the two directions are different in both ways and, further, their difference varies significantly between conditions, as shown in the bottom two panels of Fig. 4.

IV Discussion

Discovering a reliable measure of conscious arousal in animals and humans remains one of the major outstanding challenges of neuroscience. The present study addresses this challenge by connecting a complexity measure to the degree of conscious arousal, taking a step forward to strengthening the link between physics, complexity science, and neuroscience. Here, we have taken tools from the former and have applied them to a problem in the latter. Namely, we have studied the statistical complexity and time asymmetry of neural recordings in the brains of flies over two states of conscious arousal: awake and anaesthetised. We have demonstrated that differences between these macroscopic states can be revealed by both the statistical complexity of local electrical fluctuations in various brain regions, and various measures of temporal asymmetry of hidden models that explain their behaviour. Specifically, we have analysed the single-channel signals from electrodes embedded in the brain using the $\epsilon\text{-machine}$ formalism, and quantified the statistical complexity $C_{\mu}$, causal and microscopic reversibility $\Xi$ & $\mathcal{D}_{KLS}$, and crypticity $d$ of the recorded data for 15 channels in 13 flies over two states of conscious arousal. We find the statistical complexity is larger on average when a fly is awake than when the same fly is anaesthetised ($\Delta C_{\mu}>0$; Figs. 3 and 4), and that the structural complexity of information and its time reversibility captured by crypticity and KL rate are also reduced under anaesthesia ($\Delta d>0$ and $\Delta\mathcal{D}_{KLS}>0$; Fig. 4).

As we have demonstrated in this study, the local information contained within a single channel can contain information about global conscious states, which are believed to arise from interactions among many neurons. Theoretically, single channels can reflect the complexity of the multiple channels due to the concept of Sugihara causality Sugihara et al. (2012). This arises due to any one region of the brain causally interacting with the rest of the brain, making the temporal correlation in a single channel time series contain information about the spatial correlations, i.e., information that would be contained in multiple channels. With this logic, Ref. Tajima et al. (2015) infers the complexity of the multi-channel interactions from a single channel temporal structure of the time series. This is often known as the backflow of information in non-Markovian dynamics Breuer et al. (2009). The periodic structure of statistical complexity observed across channels in Fig. 2, demonstrates an unexpected example of spatial effects present in our study – one that was not observed with conventional LFP analyses. This observation may provide a motivation for multi-channel analyses.

While we already find differences between conscious states in the single channel based $\epsilon\text{-machine}$ analysis, it would be beneficial to extend the present analysis to the multi-channel scenario, in which $\epsilon\text{-machine}$ can be contrasted with the methods of IIT Casali et al. (2013); Casarotto et al. (2016); Tononi (2004); Tononi et al. (2016); Oizumi et al. (2014); Mediano et al. (2019); Barrett and Seth (2011); Tegmark (2016); Oizumi et al. (2016a, b). Formal comparison of the distinguishing power of conscious states among proposed methods (such as those in Ref. Engemann et al. (2018); Sitt et al. (2014)) will contribute to refining models and theories of consciousness.

Our results can be informally compared with a previous study, where the power spectra of the same data in the frequency domain Cohen et al. (2018) was analysed. There, a principal observation was the power in low-frequency signals in central and peripheral regions, which was more pronounced in the central region (corresponding to channel 1-6 in this study). Our $\epsilon\text{-machine}$ analysis here reveals that the region between periphery and centre (channels 5-7) shows most consistent difference in $C_{\mu}$ across history length $\lambda>2$. Ultimately, the reason for this difference is due to our distinct approach, insofar as $\epsilon\text{-machines}$ are provably the optimal predictive models of a large class of time series that take into account higher order correlations memory structure Crutchfield and Young (1989); Travers and Crutchfield (2011). Thus, our application of $\epsilon\text{-machines}$ contrasts with the power spectra analysis, by considering these higher order correlations for very high-frequency signals, instead of only two-point correlations in both high- and low-frequency signals. Finally, the top-left panel of Fig. 4 shows that in regions corresponding to channels 1 and 13, the differences in the conditions are only seen at high values of $\lambda$.

Our multi-time analysis further reveals an interesting effect when we look more closely at, e.g., the anaesthetised $\epsilon\text{-machine}$ example shown in Fig. 1(c). When we examine the binary strings belonging to each causal state, we find a clear split between active (consecutive strings of ones) and inactive (consecutive strings of zeros) neural behaviour corresponding to the left and right hand sides of Fig. 6 respectively. Previous studies have demonstrated an increase in low-frequency LFP and EEG power for mammals and birds during sleep and anaesthesia, mediated by similar neural states of activity and inactivity known as ‘up’ and ‘down’ states Sarasso et al. (2015); Lewis et al. (2012). A similar phenomenon has recently been observed in sleep deprived flies Raccuglia et al. (2019). Consistent with other studies, our study, using general anaesthesia, does not observe this slow oscillations. Future studies with more formal comparisons between up and down states and $\epsilon\text{-machines}$, in both theory and computer simulations, may be a fruitful avenue for further research in this regard.

An analysis in terms of $\epsilon\text{-machines}$ has also allowed us to discriminate between levels of conscious arousal by examining causal structures found in both forward and reverse time directions. Based on our previous finding Cohen et al. (2018) as well as related concepts in temporal predictive Hohwy (2013); Friston (2010); Tani and Nolfi (1999); Bastos et al. (2012) and causal matching Tononi (2010), we hypothesised that the wakeful brain may be tuned to causal structures of the world, which run forward in time, and thus $\epsilon\text{-machines}$ would be more complex for forward than reverse readings. Further, we hypothesised that such temporally tuned structural matching will be lost under anaesthesia. Our results (Sec. III.2.2) are highly intriguing in three ways. First, near-zero causal irreversibility indicates that reducing the structural complexity to a simple index is not enough to capture effects on the information structure that are sensitive to the direction of time. This is the case regardless of the level of consciousness (at least at the timescales of this study). Second, nonzero crypticity indicates that the underlying information structure is not symmetric in time. More precisely, the signals themselves encode different amounts of information when run forwards as opposed to backwards. Third, the KL divergence rate analysis definitively demonstrated the existence of greater temporal irreversibility in the wakeful as opposed to the anaesthetised state. Having said this, we are limited in drawing strong conclusions due, in part, to the relatively small observed effect size of $\Xi$, likely a consequence of our relatively small data set. Despite this, even at millisecond time scales, our study successfully identifies significant differences in the time direction of the neural recordings.

Identifying the decrease in temporal-reversibility due to anaesthesia in tandem with complexity is of broad interest in neuroscience. While some physicists and neuroscientists have conjectured links among physics, the brain, and even consciousness through the lens of the direction of the time, their accounts have remained rather speculative, and not built on any solid theoretical foundations (for related and alternative theoretical foundations, see the work by Cofré and colleagues Cofré et al. (2019); Cofré and Maldonado (2018)). For example, using reversely played movies, the sensitivity to the direction of time is shown to differ across brain regions in humans Hasson et al. (2008). In animal studies, some populations of neurons (in the hippocampus) are known to become activated in a particular sequential order while the animal experiences a particular event. For example, in anticipation of the event, the neurons activate in a forward direction, but in retrospection, the neurons activate in reverse order Diba and Buzsáki (2007). While direct links between these empirical findings and the $\epsilon\text{-machine}$ framework remains elusive, we foresee that our unified theoretical and analytical framework can potentially bridge this gap in the future.

Our study is not the first to apply complexity measure in consciousness research. Indeed, many definitions and measures of complexity have been proposed in the literature (see Ref. Edmonds (1997) for a list). Moreover, there is a flow of ideas going the other way as well Valdez et al. (2017); Sundar et al. (2018); Zanardi et al. (2018). However, many, if not most, of these measures cannot account for temporal correlations (memory), temporal asymmetry, or differentiate between random and structured processes. Our interdisciplinary study, based on $\epsilon\text{-machines}$, opens up new possibilities; physics can improve its theoretical constructs through the application of tools to empirical data, while neuroscience can benefit from rigorous quantitative tools that have proven their physical basis across different spatio-temporal scales. Among those complexity measures, $C_{\mu}$ can be easily interpreted in terms of temporal structure Feldman and Crutchfield (1998), as it has a direct relation to process predictability and memory requirements. We emphasise that statistical complexity $C_{\mu}$ derived from $\epsilon\text{-machines}$, drastically differentiates itself from other scalar complexity indices such as Lempel-Ziv complexity Schartner et al. (2015). For one Lempel-Ziv complexity is maximal for a random noise process whereas statistical complexity for the same process is zero (see Eq. (2)). In addition, the notion of temporal reversibility available in the $\epsilon\text{-machine}$ framework has no counterpart in Lempel-Ziv complexity. This is a critical difference since it is known that a low-complexity forward-time $\epsilon\text{-machine}$ consisting of only two causal states can have a very high-complexity reverse-time $\epsilon\text{-machine}$ with countably infinite states Ellison et al. (2011). Thus, explicitly considering the influence of time is critical for addressing questions about complexity. When coupled with our results, we can conclude that anaesthetised brains become less structured, more random, more reversible, and approach a stochastic process with a smaller memory capacity compared to the wakeful brains.

Overall, our results suggest that measures of complexity extracted from $\epsilon\text{-machines}$ might be able to identify further structures that are affected by anaesthesia at different spatial and temporal scales. It is also likely that applying a similar analysis to other data sets, in particular, human EEG data will lead to new discoveries regarding the relationship between consciousness and complexity that can be retrieved simply at the single channel level.

Acknowledgements.

RNM, FAP, NT, KM acknowledge support from Monash University’s Network of Excellence scheme and the Foundational Questions Institute (FQXi) grant on Agency in the Physical World. AZ was supported through Monash University’s Science-Medicine Interdisciplinary Research grant. DC was funded by an Overseas JSPS Postdoctoral Fellowship. NT was funded by Australian Research Council Discovery Project grants (DP180104128, DP180100396). NT and CD were supported by a grant (TWCF0199) from Templeton World Charity Foundation, Inc. We thank Felix Binder, Alec Boyd, Mile Gu, Rhiannon Jeans, and Jayne Thompson for valuable comments. KM is supported through Australian Research Council Future Fellowship FT160100073.

Appendix

IV.1 Linear mixed-effects model

In this section, we demonstrate an example of an LME analysis for the case of statistical complexity $C_{\mu}$ in the forward time direction. For the case when time direction $\xi$ is included as an effect, the only change this makes to the process is increasing the dimensions of the effects matrix $\mathbf{F}$. Performing an LME analysis on other quantities used in this study like crypticity or KL rate follow the same procedure outlined here.

The main goal of the LME analysis we perform in this study is to determine the degree of contributions each and combinations of memory length ($\lambda$), channel location ($c$), and level of conscious arousal ($\psi$) have on statistical complexity $C_{\mu}$. LME accomplishes this by modelling statistical complexity as a general linear regression equation (Eq. (7)), whose response is predicted by the aforementioned parameters $\lambda$, $c$, and $\psi$. In this Appendix, we show the exact form of the linear regression equation used in this analysis, while referring to the terminology introduced in the methods (Sec. III.1).

We begin by restating Eq. (7) for the case of statistical complexity, \mathcal{C}=\mathbf{F}\mbox{\boldmath\beta}+\mathbf{R}\mathbf{b}+\mathcal{E}, which has the form of a general multidimensional linear equation. We will set aside the right hand side of the equality for now. On the left hand side, statistical complexity takes the form of a column vector $\mathcal{C}$. Each row corresponds to the unique response of $C_{\mu}$, at a specific selection of parameters. There is a general freedom of choice associated with the number of parameters one would like to assign to the elements $\mathcal{C}$. We index the rows with fly number $f$, channel location $c$, and the conscious arousal state $\psi$. That is, the $(i,j,k)$th element is

[TABLE]

In other words, it is the $i$th fly’s $j$th channel in $k$th condition. Thus, $\mathcal{C}$ has length of $|f|\times|c|\times|\psi|=390$. Each $C_{\mu}$ in this vector is a function of $\lambda$.

The matrix $\mathbf{F}$ introducing the set of fixed effects $\mathcal{F}=\{\lambda,c,\psi,\lambda c,\lambda\psi,c\psi,\lambda c\psi\}$ into the model (known in the context of general linear models as the design matrix) can then be represented as $\mathbf{F}=(\mathbf{F}_{1},\dots,\mathbf{F}_{13})^{T}$, with each element corresponding to the design matrix of a specific fly. These individual fly response matrices can be explicitly expressed as

[TABLE]

where \vec{\mbox{\boldmath\lambda}}=(\lambda,\ldots,\lambda)^{T} and \vec{\mbox{\boldmath\Psi}}_{X}=(\Psi_{X},\ldots,\Psi_{X})^{T} are column vectors of length $15$ containing the predictor variables of memory length and level of conscious arousal respectively, $\mathbf{D}$ is the $15\times 15$ identity matrix which “selects out" the channel of interest, $\mathbf{D}_{\Psi_{X}}=\text{diag}(\Psi_{X},\ldots,\Psi_{X})$ is the $15\times 15$ matrix which “selects out" the condition of interest correlated with the level of conscious arousal, where

[TABLE]

In a similar fashion, the expression for the matrix containing the random effects $\mathbf{R}$ can be determined. For the case of our study, we only consider random effects arising due to correlations between channels within a specific fly. The result of this is an adjustment to the intercept of the linear model for each fly and channel combination. Therefore, the random effects matrix $\mathbf{R}$ is simply an identity matrix of dimension $390$. The accompanying elements of the random effects vector $\mathbf{b}$ consist of regression coefficients $b_{fc}$ describing the strength of each intercept adjustment.

The execution of the LME analysis which included coefficient fitting, and log-likelihood estimations was facilitated by running fitlme.m in MATLAB R2108b.