Abstract
Investigations into the acute effects of psychedelics on brain imaging have emphasised increased ‘brain entropy’ as a potential neural correlate. To date, 12 previous studies have reported brain entropy effects, each reporting a single and unique metric, none of which have been examined in an independent cohort. Here we evaluated acute psilocybin effects on these 12 brain entropy metrics in an independent cohort of 28 healthy participants. Following a single psilocybin dose, participants completed pre/post resting-state BOLD fMRI scans (28 pre-drug, 93 post-drug scans during the acute drug effects). Each scan was accompanied by a plasma sample to quantify plasma drug levels and estimate brain serotonin 2A receptor occupancy, as well as a rating of subjective drug intensity. We assessed relations between brain entropy and these measures with linear mixed-effects models. There was a significantly positive association for Shannon entropy of path-length and instantaneous correlation distributions and divergent associations of network-wise sample entropy at varying time-scales. We did not observe significant psilocybin effects for seven of 12 brain-entropy metrics. Whole-brain entropy metrics showed limited correlations between each other. Our observations suggest a nuanced acute effect of psychedelics on brain entropy, underscoring the need for reproducing effects. The variable effects and limited inter-metric correlation undermines the generalisability of ‘brain entropy’ as a singular construct. Future studies in clinical cohorts are crucial to elucidate the link between these metrics and therapeutic effects of psychedelics.
Introduction
Psychedelics
Psychedelic drugs induce profound altered states of consciousness including affective, sensory and cognitive effects mediated by biased agonism at the brain serotonin 2A receptor (Becker et al., 2022; Madsen et al., 2019; Nichols, 2016). In combination with psychological support, clinical studies up to phase 2b indicate promising clinical efficacy of psychedelics in the treatment of affective and behavioural neuropsychiatric disorders that may be associated with their acute effects (Anderson et al., 2020; Bogenschutz et al., 2022; Carhart-Harris et al., 2021, 2018; Goodwin et al., 2022; Griffiths et al., 2016, 2016; Holze et al., 2022; Johnson et al., 2017). Similarly, psychedelics induce acute and lasting effects on behaviour and personality in healthy participants (Griffiths et al., 2011; McCulloch et al., 2022a; Schmid and Liechti, 2018). In parallel to evaluating clinical treatment effects, human brain functional magnetic resonance imaging (fMRI) has begun to shed light on neural pathways affected by psychedelics (McCulloch et al., 2022b).
The Entropic Brain Hypothesis
Several prominent theories of critically relevant psychedelic effects on brain function have been advanced (reviewed in Doss et al., 2021). A prominent theory is the “Entropic Brain Hypothesis” (EBH), which posits that the ‘richness’ of the phenomenology of the acute psychedelic state reflects brain-wide increases in entropy of functional brain signals (Carhart-Harris, 2018; Carhart-Harris et al., 2014). The concept of “entropy” serves as a quantification of the degree of information or complexity contained in a system; it is typically expressed in bits and lacks physical dimensions. Entropy has been conceptualised as a metric of information content and is largely defined in terms of a probability distribution; the metric is commonly referred to by the eponymous name “Shannon entropy” 
where H(X) refers to the Shannon entropy of probability mass function X containing bins (x) with height (p) (Shannon, 1948). Other entropy metrics have been defined, e.g., Lempel-Ziv complexity and sample entropy (Delgado-Bonal and Marshak, 2019; Lempel and Ziv, 1976; Richman and Moorman, 2000; Ziv and Lempel, 1978), and adapted to characterise information contained in various objects, e.g., complex networks or graph structures (Rubinov and Sporns, 2010).
Previous Studies Evaluating the Entropic Brain Hypothesis
To date, 12 studies have evaluated either acute or lasting psychedelic effects on a brain entropy metric using blood oxygen level dependent (BOLD) fMRI data (Figure 1). Four papers analysed data from a study evaluating intravenous psilocybin administration (“PsiIV”) in up to 15 healthy participants (Carhart-Harris et al., 2014; Kringelbach et al., 2020; Lebedev et al., 2015; Tagliazucchi et al., 2014). Two papers reported effects from a study evaluating intravenous LSD administration (“LSDIV”) in 15 healthy participants (Lebedev et al., 2016; Luppi et al., 2021). Two papers evaluated data from both datasets (Singleton et al., 2022; Varley et al., 2020). Three papers reported effects from a study evaluating oral ayahuasca administration (“Ay”) in 9 healthy participants (Felippe et al., 2021; Viol et al., 2019, 2017). Finally, one paper reported effects from a study evaluating lasting effects of orally administered psilocybin (“PsiL”) in 11 healthy participants (Barrett et al., 2020). Notably, and highlighted in a recent review (McCulloch et al., 2022b), each of these reports quantified a distinctly different metric of brain entropy. Here we group these metrics into three categories: 1) “static connectivity” of a matrix using a graph theory framework, 2) “dynamic connectivity”, i.e., the time-varying relation between two or more time-series, 3) “dynamic activity”, i.e., the entropy of regional time-series (Figure 1). Nine of these metrics are based on the Shannon entropy of distributions, two are Lempel-Ziv complexity metrics of a time-series and one is sample entropy of a time-series. Taken together, although there is a clear interest in evaluating psychedelic effects on brain entropy, prominent limitations include that none of these measures have been evaluated in an independent cohort, the set of effects have been evaluated in only a few datasets, some with atypical modes of drug administration, and the inter-correlation between these metrics has not been considered. Previous studies compared pre-drug or placebo and a single post-drug scan. For the most part, these previous studies report increased brain entropy with only a few exceptions, e.g., two studies report non-significant psilocybin effects and one study reports no lasting changes (Figure 1, Table 1). See Supplementary Table S1 for details on the previous studies. Previous fMRI studies have reported reduced brain entropy in other states of “reduced consciousness”, including NREM sleep (Kung et al., 2022), minimal-consciousness (Maki-Marttunen et al., 2013), and anaesthesia (Pappas et al., 2019), whereas other studies have shown increased brain entropy following caffeine (Chang et al., 2018) and Salvinorin A intake (Doss et al., 2020).
Overview of previous psychedelic fMRI entropy quantifications. Methods are grouped according to the type of brain entropy metric. Symbols in the top left of each box represent the dataset that the original publication analysed. LSDIV (blue) refers to data collected in relation to a 75 µg intravenous LSD administration, PsiIV (green) refers to data collected in relation to a 2 mg intravenous psilocybin administration, Ay (orange) refers to data collected in relation to an oral ayahuasca administration, and PsiL (green oval) refers to data collected one-week before, one-week after, and one-month after a 25 mg/70 kg oral psilocybin administration. Brain images show an axial slice illustrating the parcellation/atlas used. Images in the top right describe the input into the entropy function as one of “region-wise connectivity matrix”, “ROI-wise BOLD signal” and “voxel-wise BOLD signal”. Illustrations in the lower part of each box graphically represent simplified analysis steps for each method. The original publication for each metric is denoted in the bottom left corner of each box. See Supplementary Table 1 for more details.
Summary of entropy quantification methods, previous findings, and findings reported within this manuscript. For the “Original findings” and “Our findings” columns, grey cells describe no association of effect, light grey refers to no reporting of acute effects on brain entropy, yellow describes marginal effects, and green statistically significant effects.
Aim of this study
In the current study, we sought to evaluate the acute effects of psilocybin administration on these 12 brain entropy metrics in an independent dataset of 28 healthy individuals. Based on the entropic brain hypothesis, we hypothesised that brain entropy metrics would be increased following psilocybin administration. Participants in the current study were scanned multiple times following psilocybin administration and each scan was accompanied by a self-report measure of subjective drug intensity (SDI) and blood sample to quantify plasma psilocin level (PPL). We characterised psilocybin effects by evaluating the relation between brain entropy and three metrics: SDI, PPL, and brain serotonin 2A (5-HT2A) receptor occupancy (Occ2A), estimated based on its relation to PPL established in a previous study from our lab (Madsen 2019). Additionally, we explored the inter-correlation between brain entropy metrics to characterise their associations.
Methods
Twenty-eight healthy volunteers participated in the study (10 female, mean age ± SD: 33 ± 8) and were recruited from a database of individuals interested in participating in a study involving psychedelics. A detailed description of the study design can be found in the Supplementary Text and has been reported previously (Madsen et al., 2021). The study protocol was approved by the ethics committee of the capital region of Copenhagen (H-16026898) and the Danish Medicines Agency (EudraCT no.: 2016-004000-61). The study was registered at clincaltrials.gov (NCT03289949). Data presented here were collected between 2018 and 2021. A subset of the functional brain imaging data presented here has been included in different studies reported previously (Madsen et al., 2021; Olsen et al., 2022). Details of recruitment, procedures during the psilocybin session, ethical approvals, MRI acquisition and quality control, are described in the Supplementary Text. Analyses were pre-registered on the 3rd of August 2022 (https://aspredicted.org//bw8y7.pdf).
Data Collection
After obtaining written informed consent and screening for neurological, somatic and psychiatric illness, participants completed a single-blind, cross-over study design wherein participants received a single 0.2-0.3 mg/kg dose of psilocybin (mean ± SD dose: 19.7 ± 3.6 mg, administered in units of 3 mg capsules) or 20 mg of ketanserin. Data from ketanserin scans are outside the scope of the current evaluation and not presented here. After drug administration, participants completed MRI scan sessions including resting-state fMRI (see Supplementary Text for details) approximately 40, 80, 130, and 300 minutes after administration. Following each scan, participants were asked, “On a scale from 0 to 10 how intense is your experience right now” to measure SDI and a venous blood-draw used to quantify PPL (see Supplementary Text for details). After each resting-state fMRI scan, participants were asked if they had fallen asleep (no participants reported doing so). Occ2A, i.e., occupancy of psilocybin at the 5-HT2A receptor is closely related to PPL and SDI (Madsen et al., 2019). Here we applied the previously reported parameter estimates relating PPL to occupancy based on the Hill-Langmuir equation: 
where Occmax refers to the maximum measurable occupancy, Cp refers to the measured concentration of the ligand in plasma (i.e., PPL), and EC50 refers to the concentration in plasma at which occupancy is equal to 50% of Occmax (fixed parameters used to compute Occ2A: EC50 = 1.95 µg/L and Occmax = 76.6%).
Preprocessing
Preprocessing and denoising was uniform across all entropy metrics despite differences in the pipelines of the original publications. Our pipeline included slice-timing correction (where applicable), unwarping, realignment, co-registration of structural scans to functional data, segmentation, normalisation, and smoothing. Two MR-scanners were used to acquire the data, and some functional data were temporally downsampled so that the sampling frequency was consistent across scan sessions. Denoising in CONN (Whitfield-Gabrieli and Nieto-Castanon, 2012) included linear detrending, aCompCor (Behzadi et al., 2007), 12-motion (three translations, three rotations and their first derivatives) and artefact-flagged volume regression, band-pass filtering and parcellation. Cerebellar ROIs were removed from included atlases as they were not consistently within the field of view. See Supplementary Text for more details.
Entropy of Static Connectivity
Four studies evaluated the entropy of static connectivity given by the matrix of Pearson correlation coefficients, R, computed from N-regional time-series data (Felippe et al., 2021; Lebedev et al., 2015; Viol et al., 2019, 2017), N being the number of ROIs in the atlas used by the study.
Out-network Connectivity Distribution
Following a graph-theory framework, ROIs from the 200-region Craddock-atlas (Craddock et al., 2012) were partitioned into “networks” using the Louvain modularity algorithm applied to the average connectivity matrix across scan sessions (Blondel et al., 2008). The “Out-network Connectivity”, referred to as “diversity coefficient” in the original publication and “Brain Connectivity Toolbox” (Bullmore and Sporns, 2009; Lebedev et al., 2015), of an ROI was calculated for each scan session as the Shannon entropy of the distribution of connectivity estimates between a given ROI and the set of ROIs assigned to a different network.
Degree Distribution
Degree refers to the number of non-zero elements in any given row of a thresholded matrix. ROI-specific degrees are computed based on R, the Pearson correlation matrix between ROIs, with N=105 using the Harvard-Oxford-105 atlas (Desikan et al., 2006). Both this analysis and Path-length distribution use the absolute correlation values. The thresholding for this analysis occurred in two steps. In the first step, any correlation for which the corresponding p-value was above 0.05 was set to 0. In the second step the goal is to reach a pre-specified mean degree across rows. In order to achieve this, a threshold below which all absolute values are set to 0 is gradually increased until the mean number of non-zero elements is at the desired level. Here we applied a scan-specific threshold that produced a mean degree of 27 because this was the threshold that produced the largest effect in the original publication (Viol et al., 2017). This means that each scan may have a different absolute threshold value, but identical mean degree. The final entropy quantification is simply the Shannon entropy of the distribution of degrees across ROIs. We also calculated entropy for mean degrees of 1 up to the point at which for any given scan session an increase in absolute threshold did not produce an increase in mean degree, i.e., 48. This also applies to the path-length entropy described below.
Path-length Distribution
Again using absolute correlation values, the matrix was thresholded using only the mean-degree criteria and not the p-value threshold. The matrix was then binarised, setting all non-zero elements to 1. The “shortest path length” was then computed as the fewest edges one must traverse to go from one node to another. The Shannon entropy of the distribution of path lengths for each node to all other nodes was then calculated (Viol et al., 2019). Path-length distribution was evaluated for correlation coefficient thresholds up to a mean degree of 53.
Von Neumann Entropy
Entropy of the Pearson correlation matrix, R, derived for the Harvard-Oxford-105 atlas, was calculated through the von Neumann entropy: 
, where λ are the eigenvalues of the scaled correlation matrix 
. The von Neumann entropy may also be defined as 
, where log represents the matrix logarithm (Felippe et al., 2021).
Entropy of Dynamic Connectivity
Intra-network Synchrony Distribution
Nine brain networks were defined according to a previous study (Smith et al., 2009): auditory, dorsal attention, default mode, left and right frontoparietal, motor, salience, visual 1 and visual 2. For a given network, for a given time point, the variance across voxels within the network was evaluated. The Shannon entropy was then calculated on the histogram of the variance estimates over time (Carhart-Harris et al., 2014).
Motif-connectivity Distribution
Dynamic functional brain connectivity was evaluated in four regions (10mm diameter spheres) located at bilateral hippocampi, MNI coordinates: right: (26, –21, –16), left: (–34, –22, –16), and anterior cingulate cortices, right: (4, 35, 18), left: (–2, 23, 28) using a non-overlapping sliding window approach with varying window lengths (15-150s). In each window, the partial correlation coefficient and corresponding p-value was calculated for every region pair, controlling for the remaining regions and the motion framewise displacement time-series. These time-series were standardised before windowing. The 4 x 4 partial correlation matrix was binarised for every window, according to a corrected significance threshold p=0.0083 (i.e., 0.05/6, where 6 is the number of region pairs). A probability distribution of the frequency of each of the 64 possible graph structures was established and the Shannon entropy was calculated (Tagliazucchi et al., 2014).
LEiDA-state Markov-rate
Notably, this entropy metric was not applied to evaluate psychedelic effects in the original paper. Rather, the authors provided a computational framework wherein parameters were learned by optimising this entropy measure. For each scan session, Leading Eigenvector Dynamics Analysis (LEiDA) (Cabral et al., 2017) was applied to the time-series of 90 AAL atlas regions (Tzourio-Mazoyer et al., 2002). The phase series was computed using the Hilbert transform and, for each time point, a phase coherence matrix was estimated based on the cosine of the difference between pairwise instantaneous phases. The phase coherence matrices were decomposed using the eigenvalue decomposition and the first eigenvector was retained for every time point. The set of eigenvectors was clustered using K-means with K = 3 states. Subsequently, the transition probability matrix was computed for each scan session. The entropy rate of the transition matrix, P(i, j), for each state, i, was calculated as 
. where, p is the eigenvector of P. The final entropy measure is given as 
(Kringelbach et al., 2020).
Dynamic Conditional Correlation Distribution
Regional time-series were evaluated for each of the regions described in the Shen 268 region atlas (Shen et al., 2013). Windowless framewise correlation coefficients were calculated for all edges using the Dynamic Conditional Correlation (DCC) toolbox (Lindquist et al., 2014). Subsequently, the probability distribution over each ROI-to-ROI DCC time-series was established, and the Shannon entropy was calculated. Each ROI was assigned to one of eight networks: default mode, fronto-parietal, medial-frontal, motor, subcortical-cerebellar, visual association, visual 1, and visual 2. Each ROI-to-ROI pair was assigned to its respective network-to-network association (e.g., motor-to-motor, default mode-to-motor) and the mean entropy of each network-to-network association was calculated. Although the original publication applies bin-width correction, they do not report an effect of bin width and we report findings using MATLAB’s histcounts function, which automatically calculates bin-width (Barrett et al., 2020). Thus, we did not implement bin-width correction.
Meta-state Complexity
Regional time-series were evaluated for each of the regions described in the Lausanne 463 region atlas (Daducci et al., 2012). BOLD time-series across all scan sessions were clustered using K-means into K = 4 states using the Pearson correlation distance metric. The clustering procedure was repeated 200 times with random initialisations and the best repeat in terms of K-means loss was extracted. The four states were grouped into two meta-states because the clustering procedure typically produces sign-symmetric states. Each volume was assigned to meta-state 0 or 1 and the Lempel-Ziv complexity (LZ76 exhaustive algorithm) of this binary sequence was calculated (Singleton et al., 2022).
Integration/Segregation-state Distribution
Regional time-series were evaluated for each region described in the Schaefer 200 region atlas (Schaefer et al., 2018), augmented with 32 subcortical regions from the Tian atlas (Tian et al., 2020). A sliding-window correlation analysis was performed using a window defined by convolving a rectangular window of size 44 seconds with a temporal Gaussian kernel (FWHM = 3s). The correlation matrix was established for each window (stride of 1), and the Louvain modularity algorithm (Blondel et al., 2008) was applied to estimate the module degree z-score and participation coefficient for each region. The Louvain modularity algorithm was repeated 100 times to ensure an optimal assignment. K-means clustering with K = 2 states was applied to a cartographic profile, i.e., a two-dimensional unnormalised histogram of these measures, using the correlation distance and 500 replications. The Shannon entropy was computed on the probability distribution of state occurrences (Luppi et al., 2021).
Entropy of Regional Dynamics
Multi-scale Sample Entropy
Networks were defined using the Yeo 17-network atlas (Yeo et al., 2011). Sample entropy is defined as the negative logarithm of the conditional probability that if two vectors with length m (set to 2) are dissimilar below a threshold distance r (set as 0.3), then vector pairs with length m + 1 will also have distance below the threshold (Delgado-Bonal and Marshak, 2019). Scales 1-5 were evaluated for each network, meaning that each time-series was split into non-overlapping windows of length (scale) s volumes and the means of each window were concatenated to form a condensed time-series upon which sample entropy was calculated (Lebedev et al., 2016).
BOLD Complexity
Regional time-series were evaluated for each of the regions described in the Schaefer 1000 region atlas (Schaefer et al., 2018). BOLD time-series for each ROI were first Hilbert-transformed. The amplitude of the Hilbert series was then binarised around the mean amplitude for that region, i.e., assigned as “1” if greater than the mean and “0” if less. These binarised time-series were combined into an T × N matrix, where N = 1000 is the number of regions and T and is the number of time points.This matrix was collapsed into a single vector to compute 1) the Lempel-Ziv complexity over time (LZct, LZ78 algorithm) wherein regional time-series were concatenated or 2) Lempel-Ziv complexity over space (LZcs) wherein time-adjacent “region series” were concatenated. LZct represents a calculation of the temporal entropy of each ROI, whereas LZcs represents a calculation of the spatial entropy at each timepoint. The original publication (Varley et al., 2020) reported only LZcs, but LZct is also described in (Schartner et al., 2017) whom Varley and colleagues reference as the source of their methods.
Statistical Model
Effects of psilocybin on brain entropy metrics were estimated using a linear mixed effects model with relevant R packages, i.e., predictmeans (v1.0.6), lme4 (v1.1.30), nlme (v3.1.157), lmerTest (v3.1.3) and LMMstar (v0.7.6). We regressed each metric against each of the three measures (PPL, SDI, or Occ2A) separately with a subject-specific random intercept and adjusting for motion, age, sex, and scanner. A test statistic for the association between metric and measure was obtained using the Wald statistic. To ensure adequate control of the family-wise error rate (FWER) across regions within each of the 12 metrics, (e.g., 17 networks for one time scale of multi-scale sample entropy), we calculate pFWER adjusted using the max T test method (Westfall and Troendle, 2008) in a permutation framework similar to (Lee and Braun, 2012), employing 10000 permutations. As such, if observed data superseded all permutations, the p-value is reported as p < 0.0001. “Motion” reflects the framewise displacement computed using the Artifact Detection Toolbox (ART) (see Supplementary Text) and “scanner” controls for MR scanner, of which there were two. We do not adjust p-values across metrics, nor across SDI, PPL and Occ2A; unadjusted p-values are reported for non-regional metrics as pperm. We defined findings as statistically significant if they were associated with all three psilocybin effects, SDI, PPL and Occ2A (collectively summarised “PsiFx”) at pperm < 0.05 for non-regional metrics or pFWER < 0.05 for regional metrics. Effect sizes are reported as Pearson’s correlation coefficient between the partial residuals of the entropy metrics (adjusted for covariates using the mixed-model described above) and each of PsiFx. The strength of Pearson’s correlation coefficients for significant associations are described as “weak” (≤0.3), “moderate” (>0.3 and ≤0.6), or “strong” (>0.6) as previously defined (Akoglu, 2018).
Correlation between metrics
Simple Pearson’s correlations were calculated between each whole-brain entropy metric pair as well as with motion. Of the two graph theory metrics requiring thresholding, the threshold producing a mean degree of 27 was used. For the motif-connectivity distribution, 15 and 100 second windows were selected to represent fast and slow dynamics, respectively. All scans remaining after preprocessing were used in these analyses.
Code
We shared relevant analysis scripts with original authors, hoping to ensure as much as possible that our computations aligned with original reports; we are thankful for the feedback we received. All functions used to derive entropy estimates from preprocessed data have been compiled into the “Copenhagen Brain Entropy Toolbox” (CopBET), a Matlab-based toolbox that can be found here: https://github.com/anders-s-olsen/CopBET. The permutation testing code is also available here. Code for other statistical analyses and figures can be made available upon request.
Results
Participants showed substantial SDI and PPL following drug administration as anticipated (Supplementary Figure S1). See Supplementary Table S1 and Table 1 for a summary of entropy metrics, previous findings, and our findings.
Entropy of Static Connectivity
Out-network Connectivity Distribution
Shannon entropy of regional out-network connectivity was not significantly associated with PsiFx in any of the 181 non-cerebellar brain regions after controlling for multiple comparisons (pFWER > 0.07 for all regions for at least one of PsiFx, Supplementary Table S2).
Degree Distribution
The entropy of degree distribution at a correlation coefficient threshold corresponding to a mean degree of 27 was not associated with any of PsiFx (pperm > 0.18, Figure 2A). We also did not observe significant effects for thresholds producing a mean degree between 1 and 48 (Supplementary Table S3).
Scatter plots describing the relation between static connectivity entropy and PPL. Y-axis values are partial residuals i.e., entropy values adjusted for age, sex, MR scanner and motion.
Path-length Distribution
The entropy of path-length distribution was significantly positively associated with PsiFx at the a priori described threshold producing mean degree 27 (pperm < 0.04, Figure 2B). The associations were weak to moderate (Pearson’s rho = 0.39, 0.27, and 0.23 for PPL, Occ2A and SDI, respectively). Significant weak to moderate positive associations with PsiFx were also observed across thresholds producing mean degrees from 22 to 38 (Supplementary Table S3).
Von Neumann Entropy
The Von Neumann entropy of correlation matrices was not significantly associated with PsiFx (pperm > 0.35, Figure 2C and Supplementary Table S3).
Entropy of Dynamic BOLD Connectivity
Intra-network Synchrony Distribution
Intra-network synchrony distribution was not significantly associated withPsiFx in any of the nine networks after controlling for multiple comparisons (pFWER > 0.98, Supplementary Table S2, Supplementary Figure S2).
Motif-connectivity Distribution
The four-ROI motif-connectivity state distribution was not significantly associated with PsiFx at any window length from 15 to 150s except a single weak association at window-length 100s (pperm < 0.05, Pearson’s rho 0.30, 0.25, and 0.24 for PPL, Occ2A, and SDI, respectively) surrounded by non-significant findings (Supplementary Figure S3 and Supplementary Table S3).
LEiDA-state Markov-rate
LEiDA-state Markov-rate was not significantly associated with PsiFx (pperm > 0.7 for all PsiFx, Figure 3A, Supplementary Table S3).
Scatter plots and linear models describing the relation between whole-brain dynamic connectivity entropy measures and plasma psilocin levels. Y-axis values are partial residuals i.e., entropy values adjusted for age, sex, MR scanner and motion. Panel B (Meta-state complexity) shows a non-significant association with PPL but this entropy metric does show a significant linear relation with Occ2A and SDI.
Meta-state Complexity
Meta-state entropy was positively associated with Occ2A (pperm = 0.03) and SDI (pperm = 0.003), but not PPL (pperm = 0.076, Figure 3B). Associations were weak (Pearson’s rho = 0.22, 0.33, and 0.20 for Occ2A, SDI, and PPL, respectively; Supplementary Table S3).
Integration/Segregation-state Distribution
Integration sub-state entropy was not significantly associated with PsiFx (pperm > 0.06 for all PsiFx, Figure 3C, Supplementary Table S3).
Dynamic Conditional Correlation Distribution
Dynamic conditional correlation entropy was significantly positively associated with PsiFx in 35 of 36 network-network connections (18/36 pFWER < 0.0001, i.e., observed data superseded all permutations, 29/36 pFWER < 0.001, 35/36 pFWER < 0.05; Figure 4A; Supplementary Table S2). Associations were moderate to strong (Pearson’s rho range: 0.35 to 0.78, Supplementary Table S2). The one association with at least one non-significant relation was for edges within the motor cortex.
Heatmap of the Pearson’s correlation and p-values for the association between Dynamic Conditional Correlation entropy and plasma psilocin level for each within and between network entropy estimate. *** represents pFWER < 0.0001, ** pFWER < 0.001, and * pFWER < 0.05 for associations with PPL. B: A scatter plot of the network edge with the strongest association between DCC entropy and PPL (Pearson’s Rho = 0.74). Y-axis values are partial residuals i.e., entropy values adjusted for age, sex, MR scanner and motion.
Entropy of Dynamic Activity
Multi-Scale Sample Entropy
At scale 1, (i.e., no time-series compression), sample entropy was significantly positively associated with PsiFx (pFWER < 0.05) in 7 of 17 networks (i.e., Central Visual, Dorsal Attention A, Control A, B and C, Default-Mode A and C). At scales 2, 3, and 4, no associations were significantly associated with PsiFx (pFWER > 0.05). At scale 5, sample entropy was significantly negatively associated with PsiFx in 14 of 17 networks; Control A and C and Default-Mode A (pFWER < 0.001), Somatomotor A, Dorsal Attention A and B, Salience-Ventral-Attention B, Limbic B, Control B, Temporal-Parietal, Default-Mode B and C (pFWER < 0.05, Figure 5, Supplementary Table S2). Scale 1 associations were weak to moderate (Pearson’s rho range: 0.26 to 0.47), as were Scale 5 associations (Pearson’s rho range: –0.27 to –0.49).
Forest plots representing the estimate of the association between sample entropy and PPL in each network of the Yeo-17 network parcellation and at each of scales 1 to 5. Colours represent networks and error bars represent the 95% confidence interval. *** pFWER < 0.0001, ** pFWER < 0.001 and * pFWER < 0.05 for associations with PPL.
Spatial and Temporal Dynamic BOLD Complexity
Temporal BOLD complexity (LZct) was not associated with PPL (pperm = 0.14, Figure 6B) but was associated with Occ2A (pperm = 0.03) and SDI (pperm = 0.009). Associations were weak (Pearson’s rho = 0.17, 0.23, and 0.30 for PPL, Occ2A, and SDI, respectively). Spatial BOLD complexity (LZcs) was not significantly associated with PsiFx (pperm > 0.6, Figure 6A, Supplementary Table S3).
Scatter plots and linear models describing the relation between spatial (A) and temporal (B) Lempel-Ziv entropy of dynamic BOLD activity and PPL. Y-axis values are partial residuals i.e., entropy values adjusted for age, sex, MR scanner and motion. Temporal BOLD complexity (LZct) was significantly, but weakly associated with both Occ2A and SDI despite not being significantly associated with PPL.
Correlation Between Whole-brain Entropy Quantifications
We estimated the correlation between whole-brain entropy metrics to explore their association with one another across all included scans. Some pairs of metrics were both positively correlated and negatively correlated and some pairs were effectively not correlated with one another (Figure 7). After correction for multiple comparisons, five entropy quantification pairs were positively related (LZcs & LZct, LZcs & Von-Neumann, Path-distribution & Degree-distribution, LEiDA state & Von-Neumann, Slow-4ROI-state & Von-Neumann) and four were negatively related (Path-distribution & LEiDA state, Path-distribution & Von-Neumann, Degree-distribution & LEiDA state, Degree-distribution & Von-Neumann). LZcs was significantly negatively associated with motion. See Supplementary Table S4 for pairwise correlation coefficients and p-values.
Heatmap showing the correlations between whole-brain entropy metrics and motion. Colours represent the Pearson’s correlation coefficient. * represents punc < 0.05 and ** pcorr < 0.05 where pcorr represents a Bonferroni correction for 55 tests.
Discussion
Summary
Recent studies have reported acute psychedelic effects on functional brain entropy, but to date none of these metrics have been evaluated in an independent sample. In this study we evaluated 12 previously reported entropy metrics in a novel sample of 28 healthy participants scanned with BOLD fMRI before and several times after psilocybin with concomitant measurements of subjective drug intensity and plasma psilocin level. We observed significant psilocybin effects that echoed previous reports for only two brain entropy metrics: path-length distribution, wherein we replicate increased entropy at previously reported thresholds; and sample entropy, wherein we replicate a previously observed increase in entropy at short scales and decrease at long scales. We observed a strong positive relation between psilocybin effects and brain entropy measured by Dynamic Conditional Correlation analyses that has not been previously reported. Two Lempel-Ziv complexity metrics showed some evidence for associations with psilocybin effects. For 7 of 12 brain entropy metrics previously reported, we did not observe a significant association with psilocybin measures and we see limited correlation between entropy metrics. These mixed findings underscore the importance of corroborating outcomes in independent datasets. Although we observe some evidence supporting the entropic brain hypothesis, these variable findings underscore the broadness of this theory and the need to more clearly establish which brain entropy metrics of functional brain imaging signals are acutely affected by psychedelics.
Brain Entropy Metrics Associated with Psilocybin Effects
We report a significant positive association between the Shannon entropy of the distribution of path lengths across the whole brain as previously reported by (Viol et al., 2019) and all three psilocybin metrics evaluated: PPL, Occ2A, and SDI. We observed statistically significant associations at a range of thresholds that produce graphs with a mean degree from 22 to 38, Viol and colleagues reported significant differences between conditions at thresholds producing mean degrees from 24 to 35. Characteristic path length is a description of the number of edges that must be traversed to get from any one brain region to another, a putative measure of capacity for information flow. Our results suggest that one of the effects of psilocybin on the brain can be described as a broadening of the histogram (i.e., distribution) of path lengths across region-to-region connections in the brain. Notably, this does not imply that the average path length is shorter or longer, only that there is a wider distribution of these across the whole brain i.e., it is more equally likely that the path-length between any two nodes is 1, 2 or 3 instead of being more likely to be one of these. Our convergent results are encouraging considering that the previously reported dataset used a different drug (ayahuasca) and different imaging settings, suggesting robustness of the metric. Thus, path-length entropy may be a useful candidate biomarker of neural psychedelic effects, though associations were weak-moderate and interpretation is not straightforward. We are not aware of other studies evaluating the entropy of path-length distribution so comparison to other drugs or psychiatric conditions are not yet possible and should be evaluated in future studies.
We observed a statistically significant positive relation between psilocybin effects and DCC entropy for all within or between network relations except within the motor network. Dynamic Conditional Correlation Distribution entropy is a measure of the width of the distribution of instantaneous connectivity values for any region-region edge across each scan. The previous study found no change in DCC at one-week and one-month post administration (Barrett et al., 2020) and importantly did not evaluate acute effects during psilocybin. We observed moderate to strong correlations with psilocybin effects (i.e., Pearson’s rho with PPL > 0.7 for three edges, all including DMN (DMN-DMN, DMN-Frontoparietal, DMN-Medial Frontal), and Pearson’s rho > 0.5 for 28/36 network edges). The strength of these associations is remarkable, perhaps as large as any previously reported fMRI effect of psychedelic action, suggesting that DCC entropy may be a strong candidate neural correlate of acute psychedelic effects and among the strongest correlations observed in pharmaco-fMRI. Our results suggest that psilocybin increases the variability of connectivity between regions across time across almost all region-region pairs, which are summarised into networks. Similar to entropy of path-length distribution, we are not aware of other pharmaco-fMRI studies evaluating DCC entropy. Notwithstanding, the sheer magnitude of the observed associations suggests DCC entropy may be a sensitive marker for acute psychedelic effects on the brain and so we encourage independent replication.
We observed a significant positive relation between psilocybin effects and scale-1 sample entropy (i.e., temporal resolution = 2-seconds) in seven out of 17 brain networks. Conversely, we observed a significant negative relation between scale 5 sample entropy (temporal resolution = 10-seconds) in 14 out of 17 networks. Multi-scale sample entropy measures the irregularity of a signal over its entire length. Increased sample entropy in most networks at scale-1 and decreased sample entropy at scale-5 align in both cases with the original observation (Lebedev et al., 2016). Unfortunately, we were unable to align the previously reported network labels with available versions of the atlas; thus, it is hard to resolve the spatial overlap between studies. fMRI-measured multi-scale sample entropy has been shown to be increased in the default-mode, visual, motor and lateral-prefrontal networks following caffeine (Chang et al., 2018) and decreased at scale 1 during sleep (Kung et al., 2022), although certain parameters used in their calculations were different to those employed here. As such, it is possible that the effects that we, and Lebedev and colleagues, observed may reflect differences in wakefulness and may thus be non-specific to psychedelic effects. Positive symptoms of schizophrenia have been positively associated with sample entropy at scales 1 and 2, and negatively associated in certain brain regions at scales 3, 4, and 5 (Yang et al., 2015). This is consistent with our observations and is also phenomenologically consistent, as the high-dose psychedelic state has some overlap with some positive symptoms of schizophrenia e.g., sense of veridicality, alterations in visual perception (though psychedelics do not produce ‘true’ hallucinations, i.e., sensory appearances indistinguishable from reality, as are present in schizophrenia), and sense of self. We are aware of the problematic history of psychedelic ‘psychotomimetic’ research and urge caution in overinterpretation of this apparent convergence (Nichols and Walter, 2021). Our convergent results with Lebedev and colleagues are intriguing considering that the original paper reported effects following intravenous LSD administration whereas we administer psilocybin orally. Taken together, divergent effects on sample entropy at short and long temporal resolutions may be a candidate biomarker for psychedelic effects, if they cannot be explained by, e.g., wakefulness state.
Brain Entropy Metrics Associated with Downstream Psilocybin Effects Only
Intriguingly, Lempel-Ziv complexity of two measures (meta-state complexity and temporal BOLD complexity) were significantly positively associated with Occ2A and SDI, but not PPL. The original study of meta-state time-series does not report a statistical analysis of intravenous LSD nor intravenous psilocybin effects on LZc of the meta-state time-series, but the data are publicly available and do not support a significant effect of either drug (Singleton et al., 2022). The original study of BOLD complexity reports an increase in spatial BOLD complexity following LSD, but not psilocybin and does not report any findings pertaining to temporal BOLD complexity. We speculate that the less significant relations with PPL reflect that the latter is non-linearly related to brain 5-HT2A receptor occupancy (Madsen et al., 2019) and therefore may be a less precise metric of acute psychedelic effects. MEG and EEG studies have reported increased LZc following psychedelic administration (Schartner et al., 2017; Singleton et al., 2023; Timmermann et al., 2023, 2019; Toker et al., 2022), providing convergence for its utility as a marker of psychedelic effects. Across previous studies analysing regional timeseries, there is inconsistency in the quantification of LZc in the temporal (LZct) and spatial (LZcs) domain. Our borderline statistically significant associations were observed for LZct only. It is our perspective that LZct is more sensible and should be used in future studies as it preserves region-specific temporal information, whereas LZcs is seemingly sensitive to arbitrary region order.
Brain Entropy Metrics not Associated with Psilocybin Effects
We did not observe a significant association between seven of the 12 brain entropy metrics considered here (Table 1). Of these seven metrics, the original studies reported either increased entropy following psychedelic drug administration (Carhart-Harris et al., 2014; Lebedev et al., 2015; Tagliazucchi et al., 2014; Viol et al., 2017), no effect (Luppi et al., 2021), or did not formally evaluate the effect of psychedelic drug administration (Felippe et al., 2021; Kringelbach et al., 2020) (Table 1, Supplementary Table S1). Our observed entropy estimates for hippocampal-ACC motif entropy are markedly different from those previously reported (Tagliazucchi et al., 2014). We are concerned that the originally reported values are not mathematically possible, see the Supplementary Text and Supplementary Figure S3 for a detailed consideration. Although our null findings with respect to these metrics does not establish that they have no relation to acute psychedelic drug effects, they imply a smaller relation that limits their utility as biomarkers of acute psychedelic effects. The discrepancy between our observations and those reported previously underscores the need to replicate or corroborate findings in independent cohorts to validate initial reports.
Our inability to replicate previous findings may be due to greater statistical power and different statistical models i.e., linear regression with PsiFx. Incongruence may also be attributed to differences in data collection or preprocessing parameters (e.g., orally administered psilocybin instead of intravenously administered LSD). However, if entropic brain effects are not consistent across drugs this would indicate that these metrics are not useful neural correlates of the psychedelic experience. We encourage all future fMRI studies evaluating psychedelic effects on brain function to measure subjective drug intensity at time of scanning and to collect plasma samples for quantification of plasma drug levels as described in the psychedelic fMRI consensus paper (McCulloch et al., 2022b).
Correlation Between Measures
Despite the large set of brain entropy metrics that have been reported previously, no studies have considered whether these measures are inter-correlated. We observed positive associations between path-length and degree distribution, which are based on the same graph-theory representation of connectivity, and between LZcs and LZct which are conceptually very much related. Notably, we observed four pairs of brain entropy metrics that were significantly negatively correlated. This highlights the importance of specificity in describing “brain entropy”. Many of these metrics represent distinctly different constructs, their individual meaning and collective representation of psychedelic effects is muddied by superficially considering them all metrics of “brain entropy”. Future studies should be cognisant of this variable relation in considering whether findings are consistent or convergent across studies.
Brain entropy outwith fMRI
Although we focus here on fMRI quantifications of entropy, it is worth noting that psychedelic effects on brain entropy have been applied to MEG and EEG (Schartner et al., 2017; Singleton et al., 2023; Timmermann et al., 2023, 2019; Toker et al., 2022). The entropy measures applied in these studies leverage the high temporal sampling rate that is not clearly applicable to temporally slower fMRI and were not evaluated here. Further, the methods capture different aspects of physiological response to psychedelics. Future work evaluating psychedelic effects on brain entropy using multimodal neuroimaging and evaluating relations between alternative quantifications of brain-entropy will contribute meaningfully to the field.
Limitations
Our study is not without its limitations. Brain imaging data were acquired on two different MRI scanners with different sequences (e.g., different TRs) requiring temporal downsampling of some data to match the other. However, each participant was scanned on only one scanner, enabling us to map within-subject changes onto psilocybin effects independent of scanner differences. Our study did not include a placebo condition, but we did acquire a pre-drug scan with which we estimated brain entropy metrics in the absence of psilocybin effects. Due to the large within-subject variability in fMRI outcomes in participants scanned several days apart, pre-drug vs post-drug scans performed on the same day may be superior to placebo scans performed weeks apart for evaluating drug effects because it limits this within-subject variance component (Noble et al., 2019). PPL and SDI were associated with increased motion in the scanner, see Supplementary Figure S4; although we included an estimate of motion as a covariate in our models, and employed scrubbing, motion correction and denoising strategies, we cannot rule out that motion confounds our reported effects. It has been reported in many groups that head motion is increased following psychedelic drug administration so this is not a limitation unique to our data (Timmermann et al., 2023). Most fMRI scans were 10 minutes long, though some were only five. This may not be long enough to derive stable estimates of brain entropy metrics, e.g., previous studies have recommended >13 minutes for single-echo fMRI (Birn et al., 2013). Future studies, e.g., using openly available data with longer scan durations, could inform recommended scan durations to establish stable brain entropy estimates or other methods that improve signal quality, e.g., multi-echo fMRI (Kundu et al., 2017). Approximately half of the scans analysed herein utilised a multi-band acceleration protocol that may negatively affect signal-to-noise (Risk et al., 2021), though these effects may be less pronounced for task-free imaging as performed here (Demetriou et al., 2018). There is a limited understanding of the impact of multi-band acceleration on entropy measures and future work exploring this effect is warranted. We corrected for estimated physiological noise using aCompCor but did not statistically model physiological effects such as changes in respiration, heart rate, or vasoconstriction which are affected by psilocybin and may have confounded our findings (Dyer and Gant, 1973; Holze et al., 2023). Our statistical models assume a close temporal relation between brain entropy and PsiFx (measured adjacent to scans), thus, if changes in brain entropy occur after PsiFx, they would not be well captured.
Conclusion
In conclusion, we observed acute effects of psilocybin on three of 12 previously reported brain entropy metrics. We report novel evidence for a strong association effect on DCC distribution entropy, implicating it as a potential biomarker of acute psychedelic effects as well as presenting convergent evidence for weak associations with increases in path-length distribution entropy. Lempel Ziv complexity showed marginal associations with psilocybin measures. We did not observe significant associations for seven of 12 metrics evaluated, suggesting nuanced support for the popular theory that psychedelics acutely increase brain entropy. Our observations implicate potential brain biomarkers of acute psychedelic effects and emphasise the need for both transparency in reporting brain entropy metrics and corroborating previously reported findings in independent datasets.
Data Availability
All data produced in the present study are available upon reasonable request to the authors. All functions used to derive entropy estimates from preprocessed data have been compiled into the "Copenhagen Brain Entropy Toolbox" (CopBET), a Matlab-based toolbox that can be found here: https://github.com/anders-s-olsen/CopBET. The permutation testing code is also available here. Code for other statistical analyses and figures can be made available upon request.
Supplementary Material
Supplementary Text
Study Design
Plasma Psilocin Level Quantification
BOLD scanning parameters
BOLD data preprocessing
fMRI Quality Control and missing data
Partial Residual Plots
In-depth Comparison with Tagliazucchi 2015
Pre-registration
Funding
Supplementary Figures
Time courses of subjective drug intensity and plasma psilocin level
Intra-network Synchrony Distribution regression plots
Tagliazucchi peak vs baseline plot replication
Relation between motion and PPL / SDI
Supplementary Tables
Supplementary Text
Study Design
All participants were screened to exclude individuals with significant somatic disease, neurological disorders and/or psychiatric disease using the Mini-International Neuropsychiatric Interview (Danish translation version 6.0.0). Prior to inclusion, participants received written and verbal descriptions of the study protocol, including the expected effects of psilocybin and potential side effects of administration; all participants provided written consent prior to study participation.
Exclusion criteria were 1) present or previous personal or family (first degree relatives) history of psychiatric disease (DSM axis 1 or WHO ICD-10 diagnostic classification) including substance misuse disorder; 2) Present or previous neurological condition/disease, significant somatic disease or any intake of drugs suspected to influence test results; 3) nonfluent Danish language skills; 4) vision or hearing impairment; 5) learning disability; 6) pregnancy; 7) breastfeeding; 8) MRI contraindications; 9) alcohol or drug abuse 10) allergy to test drugs; 11) significant exposure to radiation within the past year (e.g., medical imaging investigations); 12) intake of QT-prolonging medication or electrocardiogram (ECG) results indicative of heart disease, 13) history of significant adverse response to a hallucinogenic drug. 14) use of hallucinogenic drugs less than 6 months prior to inclusion; 15) blood donation less than 3 months before project participation; 16) body weight less than 50 kg; 17) low plasma ferritin levels (< 12 μg/L).
The order of psilocybin vs. ketanserin administration was balanced across participants as well as possible, drug sessions were separated by at least three weeks. The data presented here is part of a study approved by the ethics committee for the Capital Region of Copenhagen (journal identifier: H-16028698, amendments: 56023, 56967, 57974, 59673, 60437, 62255) and Danish Medicines Agency (EudraCT identifier: 2016-004000-61, amendments: 2017014166, 2017082837, 2018023295).
Participants were asked to refrain from alcohol the night before intervention days and to eat a light breakfast and refrain from caffeine intake on the morning of intervention. Prior to intervention days, a urine sample was drawn and screened for amphetamines, opioids, benzodiazepines, barbiturates, tetrahydrocannabinol, cocaine, ketamine, phencyclidine, and gamma hydroxybutyrate (Rapid Response, BTNX Inc., Markham, Canada), a blood sample and an electrocardiogram were taken to ensure that individuals were in good overall physiological health (e.g., absent evidence for a prolonged QT interval).
On the day upon which participants received psilocybin, they arrived in the morning, completed state questionnaires and a pre-drug MRI scan session. Following a brief conversation with the psychological support personnel (two people) present for the psilocybin administration day, participants received psilocybin in a room adjacent to the MRI scanner. Participants could lay on a bed or sit in a chair in this room and music was played after psilocybin administration but not during scanning.
Components of preparation, support and integration sessions
Preparation
The purpose of the session
The participant’s motivation for entering the study and expectancy management
The participant’s personal background
Previous experience with psychedelics or altered states of consciousness
Information about possible effects of psilocybin
Tools for navigating the psychedelic experience, with a particular focus on the MR environment
The participant’s social support-network and potential contact person
Practical road map of the intervention day, including rules for the study
Agreements regarding therapeutic touch and physical support
Intervention (before administration)
The purpose of the session, including any questions arising since preparation
Anchoring – and repetition of tools for navigating the experience
Support the participant in letting go of everyday tasks and preoccupations, including hand over of their mobile phone and helping with arrangements of various matters with the participant’s family/contact person or social network
Repeating the practicalities of the intervention, including rules and agreements made during preparation
Administration of psilocybin
Intervention (end of day)
Drawing a mandala
Information about afterglow effects and making a plan for self-care upon arrival at home
Confirmation of time for integration session
Integration
The purpose of the session
Integration wheel
Eliciting a full narrative about psychedelic experience
Witnessing and debriefing of difficult material
Relating the psychedelic experience to the participant’s life going forward
Plan for self-care
Thoughts about ending the study and information about the possibility of contacting the psychologists, if needed
Any questions that the participant may have
Plasma Psilocin Level Quantification
Plasma psilocin level was measured using ultra performance liquid chromatography and tandem mass spectrometry, reflecting free, unconjugated psilocin and not glucuronidated psilocin (Kolaczynska et al., 2021).
7 mL blood samples were drawn from an intravenous access in the antecubital vein and collected in EDTA vials, placed on ice, centrifuged, and plasma was aliquoted and stored at –20 or –80°C. Psilocin was obtained from Lipomed (Arlesheim, Switzerland) and Cerilliant (Round Rock, TX, USA) for calibrator and control batches, respectively, while the deuterated internal standard (IS) psilocin-d10 was from Cerilliant. Acetonitrile, methanol HPLC grade and water were obtained from Fishers Scientific (Loughborough, UK). Ascorbic acid was obtained from VWR (Hassrode, Belgium).
Stock solution (1000 mg/l) of psilocin and IS were prepared in acetonitrile and stored in amber ampoules at –20°C until use. Working standard solutions from 0.5 µg/l to 1000 µg/l were freshly prepared in 50% methanol in water for each analysis, and the IS-solution was 100 µg/l in 50% methanol. For preparation of calibration standards and quality controls blank plasma was preserved with 1% fluoride and stored at –20°C. Two quality controls (QC) were prepared at low (5.0 µg/kg) and high (50 µg/kg) levels and stored at –80°C. These two along with a freshly spiked blank plasma sample at 2.5 µg/kg were analysed in each run.
Protein precipitation was performed on a fully automated Tecan Freedom EVO 200 robotic platform (Tecan group Ltd, Männedorf, Switzerland) that included all pipetting, centrifugation, and evaporation steps. Each plasma sample (100 mg) was transferred to a 96-well 2.0 ml deep-well plate and 20 µl IS-solution was added to each well, followed by precipitation with 700 µl acetonitrile and shaking. The samples were centrifuged at 1000 g for 10 min, and the supernatant was evaporated to dryness under a stream of nitrogen at 35°C. Afterwards, the samples were reconstituted in 100 µl mixture of 12.5% methanol:12.5% acetonitrile:75% 0.05% formic acid in water, shaken and centrifuged again. Finally, the supernatant was transferred to a 96-well plate and 1 µl was injected into the chromatographic system.
Chromatographic separation was performed on a HSS T3 column (100 x 2.1 mm, 1.8um, Waters, Milford, MA, USA) using an ACQUITY Ultra Performance Liquid Chromatography system (UPLC) from Waters. The mobile phase was composed of solvent A: 1 mM ammonium formate in 0.1% formic acid in water and B: 0.1% formic acid in 1:1 mixture of acetonitrile:methanol. The column was maintained at 45°C with a flow 0.4 ml/min, and a gradient elution was applied from 2% to 100% B within 3.2 min with a total analysis time of 4.5 min. Detection was done by tandem mass spectrometry using an ACQUITY TQS from Waters. Ionisation was achieved by electrospray in positive mode, and the source temperature was set at 150°C and desolvation temperature at 600°C. Two transitions were used for psilocin, m/z 205 –> 58 and 205 –> 160, with a cone voltage of 20 V and collision energy at 14 and 18 eV, respectively. For the IS the transition was m/z 215 –> 164 with cone 20 V and collision energy of 18 eV. Argon was used as collision gas at 0.45 Pa, and desolvation and cone gasflow were fixed at 1000 L/hr and 150 L/hr, respectively. Data were acquired and processed with MassLynx 4.2 software (Waters).
Quantification was performed by an eight-point linear calibration curve (0.1, 0.5, 1.0, 5.0, 25, 50, 100, 200 µg/kg) with weighting 1/x. Limits of detection (LOD) and quantification (LOQ) were 0.1 and 0.5 µg/kg, respectively, while the upper limit of quantification was 200 µg/kg. The overall process efficiency was found to be 63% based on an obtained extraction efficiency of 81% and matrix effect of 25% that the stable isotope labelled IS adjusted for. QC plasma samples were measured in each series with a RSD of 5% and an accuracy of 78% for the low level and less than 5% and an accuracy 89% for the high level. Similar performance was obtained with QCs preserved with 1 mM ascorbic acid demonstrating that stored plasma QCs at –80°C without ascorbic acid were stable for at least 6 months. The freshly spiked plasma QC at 2.5 µg/kg analysed in each series had a RSD of 13% and an accuracy of 88% (n = 11). For the final presentation of results, units of concentration (µg/kg) were converted into µg/L using a conversion factor or 1.02 kg per litre of plasma.
BOLD Scanning Parameters
Due to logistical constraints, MRI scans were acquired on one of two Siemens 3T MAGNETOM Prisma scanners (“Scanner A”, “Scanner B”). The first 15 participants were scanned on Scanner A, the final 13 participants were scanned on Scanner B, all participants completed all scans on a single scanner.
A high-resolution, T1-weighted 3D MPRAGE structural image was acquired (Scanner A: inversion time = 900 ms, TE = 2.58 ms, TR = 1900 ms, flip angle = 9°, in-plane matrix = 256×256, resolution = 0.9×0.9 mm, 224 slices; slice thickness = 0.9 mm, no gap, Scanner B: inversion time = 920 ms, TE = 2.41 ms, TR = 1810 ms, flip angle = 9°, in-plane matrix = 288 × 288, in-plane resolution = 0.8×0.8 mm, number of slices = 224, slice thickness = 0.8 mm).
BOLD fMRI data was obtained using a T2*-weighted gradient echo-planar imaging (EPI) sequence.
Scanner A: 64-channel head coil, TR = 2000 ms, TE = 30 ms, flip angle = 90°, in-plane matrix = 64×64 voxels, in-plane resolution=3.6×3.6 mm, number of slices = 32, slice thickness = 3.0 mm, gap = 0.75 mm, phase-encoding direction = AP, scan time = 10 minutes, volumes = 300.
Scanner B: 32-channel head coil, TR = 800 ms, TE = 37 ms, flip angle = 52°, in-plane matrix = 104×104 voxels, in-plane resolution = 2×2 mm, slices = 72, slice thickness = 2.0 mm, gap = 0, multi-band acceleration factor = 8, scan time = 10 minutes, volumes = 750. Due to an acquisition error, 8 participants from Scanner B were scanned for only five minutes (i.e., 375 volumes).
BOLD Data Preprocessing
BOLD data were preprocessed in SPM12 (http://www.fil.ion.ucl.ac.uk/spm) and denoised using the CONN toolbox (Whitfield-Gabrieli and Nieto-Castanon, 2012) using MATLAB (The Mathworks, inc., version 2019b 9.7.0). Preprocessing pipelines were identical for both scanners unless specified. We applied slice-timing correction (scanner A only because scanner B used multi-band), data were unwarped and realigned, structural scans were co-registered to functional data. High-resolution structural images were segmented into grey-matter, white-matter and CSF-maps; functional data were normalised to MNI152 space based on parameters estimated during the segmentation procedure and smoothed using a 6mm FWHM gaussian kernel. Denoising included aCompCor regression of white-matter and CSF time-series (first five components for each and their first derivatives) (Behzadi et al., 2007) as well as regression of 6-motion parameters and their first derivatives, and scrubbing of volumes flagged using the artifact detection tool (ART, https://www.nitrc.org/projects/artifact_detect). Some analyses utilised unsmoothed data to align with previous methods (Intra-network synchrony, Motif-connectivity distribution, Multi-scale Sample Entropy.). A bandpass filter of 0.008 – 0.09 Hz was applied. Regional time-series were averaged in CONN using analysis-specific atlases. Following preprocessing and denoising and to align the temporal frequency between the two scanners, scanner B time-series were temporally downsampled to match the TR of scanner A (i.e., 2s) using the MATLAB command “downsample”.
fMRI Quality Control and Missing Data
All BOLD data were visually inspected to ensure accurate co-registration to MNI space. Connectivity matrices and regional time-series were visually inspected for artefacts or missing data.
Of the 130 resting-state fMRI scan sessions acquired during psilocybin intervention days for the 28 participants (76 Scanner A and 54 Scanner B), Nine scans and no participants were excluded from analyses: two scans from one participant because of excessive head motion through the scan session (>50% of volumes flagged for censoring), one scan from one participant due to technical problems with data acquisition, two scans from one participant due to 1) interrupted scan and 2) technical problems with data acquisition, one scan each from three participants due to interrupted scans. Thus, 121 scan sessions were included in analyses. One scan session did not have an accompanying measure of PPL, another was without a measure of SDI, these individual scans were excluded from these respective analyses.
Partial Residual Plots
To display the relation between PPL and each whole-brain entropy metric, plots were constructed that display PPL values (in μg/L) on the x-axis and adjusted residual entropy values on the y-axis. These are adjusted by removing the covariate effects. This was performed by taking the raw entropy score and deducting the beta of each covariate multiplied by the relevant value of that covariate for each scan.
In-depth Comparison with Tagliazucchi 2015
We do not observe a significant relation between psilocybin effects and hippocampal-ACC motif distribution entropy, calculated as described in (Tagliazucchi et al., 2014). The original publication describes a significant difference comparing both pre-psilocybin vs post-psilocybin and post-psilocybin vs post-placebo at longer (c.75-150 second) windows. Unfortunately, the entropy values presented by Tagliazucchi and colleagues leave significant cause for concern. They use non-overlapping windows, the shortest window-lengths they used were 10-seconds, and their scans were 300-seconds long. Thus, at the shortest window lengths they would observe 30 brain states whereas at the longest window length (150-seconds) only two events would be observed. Thus at the longest window-lengths, there are two possible Shannon entropy values, one (i.e., two different states –2 x 0.5 x log2(0.5) = 1) and zero (i.e., two of the same state –1 x 1 x log2(1) = 0). Thus, it is mathematically impossible to have Shannon entropy values in excess of 1, which are reported in this original publication. Additionally, one would expect the Shannon entropy values at short window-lengths to be higher than those at long window-lengths, as there are more events in the histogram of motifs (i.e., 5 of state 1, 3 of state 2, 10 of state 3, 6 of state 4, 6 of state 5, Shannon entropy = 2.2 bits) but they report window-lengths approaching zero as window length approaches 10-seconds. Although unlike the high entropy values at long window length, this is mathematically possible, as the same state may reoccur with increasing frequency at shorter window-lengths (i.e., 29 of state 1, 1 of state 2, Shannon entropy = 0.21), but this seems highly unlikely and is contrary to what we observe (Supplementary Figure S3).
Pre-registration
This study was pre-registered on the 3rd of August 2022 and the pre-registration document can be viewed here https://aspredicted.org/bw8y7.pdf. Notably, three of 12 papers were identified after pre-registration and therefore not described therein.
Funding
The work was supported by Innovation Fund Denmark (grant number 4108–00004B), Independent Research Fund Denmark (grant number 6110–00518B), and Ester M. og Konrad Kristian Sigurdssons Dyreværnsfond (grant number 850–22–55,166–17-LNG). MKM was supported through a stipend from Rigshospitalet’s Research Council (grant number R130-A5324). BO has received funding from the European Union’s Horizon 2020 research and innovation program under the Marie Sklodowska-Curie grant agreement No 746,850. DEM salary is supported by COMPASS Pathways Ltd. Funding agencies did not impact the study and played no role in manuscript preparation and submission.
Supplementary Figures
Time-course of (A) plasma psilocin levels and (B) subjective drug intensity per participant. Each colour represents a participant and the black line represents a locally estimated scatter-plot smoothing across the sample. Pre-drug PPL measurements are set at 0 mins. fMRI scans were performed before drug administration and then at approximately +40, +80, +120 and +300 minutes.
Scatter plots describing the relation between Intra-network synchrony of dynamic BOLD activity and plasma psilocin level, for each network. Y-axis values are adjusted for age, sex, MR scanner and motion. All pFWER < 0.98.
Scatter plot of motif-connectivity distribution entropy with overlaid locally estimated scatter-plot smoothing showing the relation between window-size and entropy derived from the linear model. Blue points and line represent the intercept value from the linear model comparing PPL and entropy (i.e., model estimated entropy when PPL = 0). The red points and line represent model estimated entropy when PPL = 15 μg/L, approximately PPL at peak drug effects. Note the discrepancy in y-axis values and that we observed decreasing entropy with increasing window size vis-a-vis Figure 5A from Tagliazucchi et al., 2015.
Scatter plots displaying the adjusted linear relation between motion (i.e., framewise displacement, FD) and (A) PPL or (B) SDI. Y-axis values are partial residuals adjusted for covariates, i.e., age, sex, scanner and subject. P-values are uncorrected. Pearson’s Rho are 0.33 for PPL and 0.43 for SDI.
Acknowledgements
A sincere thank you to Andrea Luppi, Enzo Tagliazucchi, Alexander Lebedev, Manoj Doss, Thomas Varley, and Parker Singleton for their advice and feedback regarding their entropy metrics, including advice from Drs Doss and Varley to not use their entropy metrics, as they no longer believed them to be valid, which we appreciate, but disregarded.
Footnotes
Conflict Statement DEWM salary is supported by an unrestricted grant from COMPASS Pathways who have no involvement in the preparation or conception of this manuscript or related data collection. MKM has received an honorarium as a speaker for H. Lundbeck. GMK has served as a consultant for Sanos, Gilgamesh, Onsero, Pangea, Abbvie, PureHealthTech, and has received honoraria as speaker for H. Lundbeck and Sage Therapeutics.
References
