Enhanced Repertoire of Brain Dynamical States

During the Psychedelic Experience

Enzo Tagliazucchi,

1

Robin Carhart-Harris,

2

Robert Leech,

3

David Nutt,

2

and Dante R. Chialvo

4

1

Neurology Department and Brain Imaging Center, Goethe University, Frankfurt am Main,

Germany

2

Imperial College London, Centre for Neuropsychopharmacology, Division of Experimental

Medicine, London, United Kingdom

3

Computational, Cognitive and Clinical Neuroimaging Laboratory (C3NL), Division of Brain

Sciences, Imperial College London, United Kingdom

4

Consejo Nacional de Investigaciones Cientificas y Tecnologicas (CONICET), Buenos Aires,

Argentina

r

r

Abstract:

The study of rapid changes in brain dynamics and functional connectivity (FC) is of increasing

interest in neuroimaging. Brain states departing from normal waking consciousness are expected to be
accompanied by alterations in the aforementioned dynamics. In particular, the psychedelic experience
produced by psilocybin (a substance found in “magic mushrooms”) is characterized by unconstrained
cognition and profound alterations in the perception of time, space and selfhood. Considering the spon-
taneous and subjective manifestation of these effects, we hypothesize that neural correlates of the psyche-
delic experience can be found in the dynamics and variability of spontaneous brain activity fluctuations
and connectivity, measurable with functional Magnetic Resonance Imaging (fMRI). Fifteen healthy sub-
jects were scanned before, during and after intravenous infusion of psilocybin and an inert placebo.
Blood-Oxygen Level Dependent (BOLD) temporal variability was assessed computing the variance and
total spectral power, resulting in increased signal variability bilaterally in the hippocampi and anterior
cingulate cortex. Changes in BOLD signal spectral behavior (including spectral scaling exponents)
affected exclusively higher brain systems such as the default mode, executive control, and dorsal atten-
tion networks. A novel framework enabled us to track different connectivity states explored by the brain
during rest. This approach revealed a wider repertoire of connectivity states post-psilocybin than during
control conditions. Together, the present results provide a comprehensive account of the effects of psilo-
cybin on dynamical behavior in the human brain at a macroscopic level and may have implications for
our understanding of the unconstrained, hyper-associative quality of consciousness in the psychedelic
state. Hum Brain Mapp 00:000–000, 2014.

V

C

2014

Wiley Periodicals, Inc.

Key words:

Psilocybin; fMRI; functional connectivity; resting state; psychedelic state

r

r

Contract grant sponsors: CONICET (Argentina) and LOEWE
Neuronale

Koordination

Forschungsschwerpunkt

Frankfurt—

NeFF (Germany).
*Correspondence to: Enzo Tagliazucchi, Neurology Department
and Brain Imaging Center, Goethe University, Frankfurt am
Main, Germany. E-mail: tagliazucchi.enzo@googlemail.com

Received for publication 11 January 2014; Revised 21 May 2014;
Accepted 23 May 2014.
DOI: 10.1002/hbm.22562
Published online 00 Month 2014 in Wiley Online Library
(wileyonlinelibrary.com).

r

Human Brain Mapping 00:00–00 (2014)

r

V

C

2014

Wiley Periodicals, Inc.

INTRODUCTION

Psilocybin (phophoryl-4-hydroxy-dimethyltryptamine) is

the phosphorylated ester of the main psychoactive com-
pound found in magic mushrooms. Pharmacologically
related to the prototypical psychedelic LSD, psilocybin has
a long history of ceremonial use via mushroom ingestion
and, in modern times, psychedelics have been assessed as
tools to enhance the psychotherapeutic process [Grob
et al., 2011; Krebs et al., 2012; Moreno et al., 2006]. The
subjective effects of psychedelics include (but are not lim-
ited to) unconstrained, hyperassociative cognition, dis-
torted sensory perception (including synesthesia and
visions of dynamic geometric patterns) and alterations in
one’s sense of self, time and space. There is recent prelimi-
nary evidence that psychedelics may be effective in the
treatment of anxiety related to dying [Grob et al., 2011]
and obsessive compulsive disorder [Moreno et al., 2006]
and there are neurobiological reasons to consider their
potential as antidepressants [Carhart-Harris et al., 2012,
2013]. Similar to ketamine (another novel candidate antide-
pressant) psychedelics may also mimic certain psychotic
states such as the altered quality of consciousness that is
sometimes seen in the onset-phase of a first psychotic epi-
sode [Carhart-Harris et al., 2014]. There is also evidence to
consider similarities between the psychology and neurobi-
ology of the psychedelic state and Rapid Eye Movement
(REM) sleep [Carhart-Harris, 2007; Carhart-Harris and
Nutt, 2014], the sleep stage associated with vivid dreaming
[Aserinsky and Kleitman, 1953].

The potential therapeutic use of psychedelics, as well as

their capacity to modulate the quality of conscious experi-
ence in a relatively unique and profound manner, empha-
sizes the importance of studying these drugs and how
they act on the brain to produce their novel effects. One
potentially powerful way to approach this problem is to
exploit human neuroimaging to measure changes in brain
activity during the induction of the psychedelic state. The
neural correlates of the psychedelic experience induced by
psilocybin have been recently assessed using Arterial Spin
Labeling (ASL) and BOLD fMRI [Carhart-Harris et al.,
2012]. This work found that psilocybin results in a reduc-
tion of both CBF and BOLD signal in major subcortical
and cortical hub structures such as the thalamus, posterior
cingulate (PCC) and medial prefrontal cortex (mPFC) and
in decreased resting state functional connectivity (RSFC)
between the normally highly coupled mPFC and PCC.
Furthermore, our most recent study used magnetoence-
phalography (MEG) to more directly measure altered neu-
ral activity post-psilocybin and here we found decreased
oscillatory power in the same cortical hub structures
[Muthukumaraswamy et al., 2013, see also Carhart-Harris
et al., 2014 for a review on this work).

These results establish that psilocybin markedly affects

BOLD, CBF, RSFC, and oscillatory electrophysiological
measures in strategically important brain structures, pre-
sumably involved in information integration and routing

[Carhart-Harris et al., 2014; de Pasquale et al., 2012; Hag-
mann et al., 2008; Leech et al., 2012]. However, the effects
of psilocybin on the variance of brain activity parameters
across time has been relatively understudied and this line
of enquiry may be particularly informative in terms of
shedding new light on the mechanisms by which psyche-
delics elicit their characteristic psychological effects. Thus,
the main objective of this article is to examine how psilo-
cybin modulates the dynamics and temporal variability of
resting state BOLD activity. Once regarded as physiologi-
cal noise, a large body of research has now established
that resting state fluctuations in brain activity have enor-
mous neurophysiological and functional relevance [Fox
and Raichle, 2007]. Spontaneous fluctuations self-organize
into internally coherent spatiotemporal patterns of activity
that reflect neural systems engaged during distinct cogni-
tive

states

(termed

“intrinsic”

or

“resting

state

networks”—RSNs) [Fox and Raichle, 2005; Raichle, 2011;
Smith et al., 2009]. It has been suggested that the variety
of spontaneous activity patterns that the brain enters dur-
ing task-free conditions reflects the naturally itinerant and
variegated quality of normal consciousness [Raichle, 2011].
However, spatio-temporal patterns of resting state activity
are globally well preserved in states such as sleep [Boly
et al., 2009, 2012; Brodbeck et al., 2012; Larson-Prior et al.,
2009; Tagliazucchi et al., 2013a,b,c] in which there is a
reduced

level

of

awareness—although

very

specific

changes in connectivity occur across NREM sleep, allow-
ing the decoding of the sleep stage from fMRI data
[Tagliazucchi et al., 2012c; Tagliazucchi and Laufs, 2014].
Thus, if the subjective quality of consciousness is markedly
different in deep sleep relative to the normal wakeful state
(for example) yet FC measures remain largely preserved,
this would suggest that these measures provide limited
information about the biological mechanisms underlying
different conscious states. Similarly, intra-RSN FC is
decreased under psilocybin [Carhart-Harris et al., 2013]
yet

subjective

reports

of

unconstrained

or

even

“expanded” consciousness are common among users (see
Carhart-Harris et al. [2014] for a discussion). Thus, the
present analyses are motivated by the view that more sen-
sitive and specific indices are required to develop our
understanding of the neurobiology of conscious states, and
that measures which factor in variance over time may be
particularly informative.

A key feature of spontaneous brain activity is its

dynamical nature. In analogy to other self-organized sys-
tems in nature, the brain has been described as a system
residing in (or at least near to) a critical point or transition
zone between states of order and disorder [Chialvo, 2010;
Haimovici et al., 2013; Tagliazucchi and Chialvo, 2011;
Tagliazucchi et al., 2012a]. In this critical zone, it is
hypothesized that the brain can explore a maximal reper-
toire of its possible dynamical states, a feature which
could confer obvious evolutionary advantages in terms of
cognitive and behavioral flexibility. It has even been pro-
posed that this cognitive flexibility and range may be a

r

T

agliazucchi et al.

r

r

2

r

key property of adult human consciousness itself [Tononi,
2012]. An interesting research question therefore is
whether changes in spontaneous brain activity produced
by psilocybin are consistent with a displacement from this
critical point—perhaps towards a more entropic or super-
critical state (i.e. one closer to the extreme of disorder than
normal waking consciousness) [Carhart-Harris et al., 2014].
Further motivating this hypothesis are subjective reports
of hyper-associative cognition under psychedelics, indica-
tive of unconstrained brain dynamics. Thus, in order to
test this hypothesis, it makes conceptual sense to focus on
variability in activity and FC parameters over time, instead
of the default procedure of averaging these over a pro-
longed period. In what follows, we present empirical data
that tests the hypothesis that brain activity becomes less
ordered in the psychedelic state and that the repertoire of
possible states is enhanced. After the relevant findings
have been presented, we engage in a discussion to suggest
possible strategies that may further characterize quantita-
tively where the “psychedelic brain” resides in state space
relative to the dynamical position occupied by normal
waking consciousness.

MATERIALS AND METHODS

Study Overview and Design

This was a within-subjects placebo-controlled study. The

study was approved by a local NHS Research Ethics Com-
mittee and University of Bristol Research and Develop-
ment department, and conducted in accordance with
Good Clinical Practice guidelines. A Home Office License
was obtained for storage and handling of a Schedule 1
drug and the University of Bristol sponsored the research.

Participants

This is a new analysis on previously published data

[Carhart-Harris et al., 2012, 2013]. Fifteen healthy subjects
took part: 13 males and 2 females (mean age 5 32,
SD 5 8.9). Recruitment was via word of mouth. All sub-
jects were required to give informed consent and undergo
health screens prior to enrollment. Study inclusion criteria
were: at least 21 years of age, no personal or immediate
family history of a major psychiatric disorder, substance
dependence, cardiovascular disease, and no history of a
significant adverse response to a hallucinogenic drug. All
of the subjects had used psilocybin at least once before
(mean number of uses per subject 5 16.4, SD 5 27.2) but
not within 6 weeks of the study.

Anatomical Scans

Image acquisition was performed on a 3T GE HDx sys-

tem. Anatomical scans were performed before each
functional scan and thus during sobriety. These were

three-dimensional fast spoiled gradient echo scans in an
axial orientation, with field of view 5 256 3 256 3 192 and
matrix 5 256 3 256 3 192 to yield 1 mm isotropic voxel
resolution (repetition time/echo time [TR/TE] 5 7.9/3.0
ms; inversion time 5 450 ms; flip angle 5 20



).

Drug Infusion

All subjects underwent two 12 min eyes closed resting

state BOLD fMRI scans on two separate occasions at least
7 days apart: placebo (10 ml saline, 60 s intravenous injec-
tion) was given on one occasion and psilocybin (2 mg dis-
solved in 10 ml saline, 60 s infusion) on the other. Seven
of the subjects received psilocybin in scan 1, and 8
received it in scan 2. Injections were given manually by a
doctor within the scanning suite. The 60 s infusions began
exactly 6 min after the start of the 12 min scans. The sub-
jective effects of psilocybin were felt almost immediately
after injection and were sustained for the remainder of the
RS scan. For more details on the subjective effects of intra-
venous psilocybin see [Carhart-Harris et al., 2012, 2013].

fMRI Data Acquisition and Preprocessing

BOLD-weighted fMRI data were acquired using a gradi-

ent echo planar imaging sequence, TR/TE 3,000/35 ms,
field-of-view 5 192 mm, 64 3 64 acquisition matrix, paral-
lel acceleration factor 5 2, 90



flip angle. Fifty-three oblique

axial slices were acquired in an interleaved fashion, each
3 mm thick with zero slice gap (3 3 3 3 3 mm voxels). A
total of 240 volumes were acquired, with infusion taking
place in the middle of the session. Data were motion cor-
rected using FSL MCFLIRT and a high-pass filter of 100 s
was applied. Data were smoothed using a Gaussian kernel
of 5 mm FWHM. Motion time courses were regressed
from the data during the preprocessing step, together with
average CSF and white matter time courses. As an addi-
tional control to exclude epochs of large movement ampli-
tude,

all

volumes

associated

with

a

mean

head

displacement larger than 1 SD were erased from the
analysis.

Analysis of Spatiotemporal Dynamics

The methods here employed are based on a statistical

physics framework useful for characterizing fluctuations in
systems composed of large number of coupled degrees of
freedom. From this perspective, fluctuating activity in a
particular region and correlations between different
regions are interdependent or related [Ross, 1966]. There-
fore, any change in the cortical dynamics due to a given
intervention is expected to be reflected both at the level of
the variance of the activity at one region and in the
strength of the interactions between different regions. In
what follows, the methodology accounting for these two
aspects (i.e. the temporal variability and the dynamical

r

E

nhanced Repertoire of Brain Dynamical States

r

r

3

r

changes in correlations) are explained in detail. All the
numerical analyses for these calculations were performed
using in-house MATLAB scripts.

Analysis of Temporal Variability

The analysis referred to here as “temporal variability” is

concerned with the variance in the amplitude of the BOLD
signal and it can be expressed both in a temporal (i.e., the
standard variance) and frequency domain, as will be
explained below. As a time domain measure of variability
in the amplitude of the BOLD signal, straightforward
voxel-wise computation of BOLD variance was performed,
resulting in whole-brain BOLD variance maps for each
condition (psilocybin/placebo before and after infusion)
and participant. An additional evaluation of variance can
be performed in the frequency domain, noting that by Par-
seval’s theorem:

ð

x

2

t

ð Þdt5

ð

jAðf Þj

2

df

(1)

Then, the following series of equalities holds (for a sig-

nal with zero mean):

r

2

5

lim

T!1

1

T

ð

T

2

2

T

2

x

2

t

ð Þdt5 lim

T!1

1

T

ð

jAðf Þj

2

df 5

ð

U

f

ð Þdf

(2)

In Eqs. (1) and (2), A(f) represents the Fourier transform

of x(t) and U(f) is the power spectral density. Therefore,
the variance can also be obtained by integrating U(f) across
the whole range of frequencies.

Further Spectral Analysis of BOLD Fluctuations

Equations (1) and (2) show that the integral of the

power spectral density equals the variance of the signal. If
this integral is not performed across all frequencies but in
a certain range, the result will correspond to the contribu-
tion of the frequencies to the variance in the specified
ranges. It is generally considered that slow (0.01–0.1 Hz)
BOLD frequencies carry neural significance in the resting
state fMRI signal [Cordes et al., 2001], therefore, we eval-
uated separately the low frequency power (LFP) in this
range and compared it before and after the psilocybin
infusion.

Given the fact that the spectral content of spontaneous

BOLD fluctuations, like many other complex systems in
nature, follows a power law of the form

1

f

a

[Expert et al.,

2011], its power spectrum density can also be character-
ized by a single parameter (a). This parameter condenses
the scaling behavior and is demonstrative of the long-
range temporal correlations (or memory) of any given sig-
nal. Thus, for uncorrelated noise, a 5 0 (i.e. a flat spectrum,
or so called “white noise”) is obtained, whereas for signals
presenting some degree of long-term correlations, values

of a > 0 (i.e., the so called “colored noise”) are observed.
To obtain the scaling exponent a, the first derivative of
BOLD time series was first computed (in order to mini-
mize the influence of nonstationarities) and the fit was
performed on the spectrum of the derivative. Note that for
a power spectrum of the form

1

f

a

, differentiation decreases

the value of a by 2, so the exponent of the original time
series can be recovered from that of the derivative.

Point-Process Analysis

A recent series of studies demonstrated that the continu-

ous BOLD signal can be transformed into a discrete point-
process encoding the timings of the most functionally rele-
vant events [Petridou et al., 2012, Davis et al., 2013; Taglia-
zucchi et al., 2010a,b, 2012a]. In this approach, relevant
events are defined by a threshold crossing (e.g. whenever
the signal departs 11 SD from its mean value). Besides
allowing a dramatic compression of the information pres-
ent in a BOLD timeseries, this approach has the particular
benefit of allowing a simple criterion to eliminate motion
artifacts: i.e. “points” (threshold crossings) occurring dur-
ing high movement epochs (i.e. larger than a certain value)
are simply ignored in the analysis.

Two interdependent observables are defined once the

point-process is obtained: the rate (number of crossings
divided by the series length) and inter-event intervals
(average temporal separation between two consecutive
points). On average, there is a clear inverse relationship
between these two variables. Furthermore, the rate is
expected to increase (or the interval to decrease) for a sig-
nal with a high contribution of fast frequencies.

As an efficient alternative to the spectral analysis

described in the previous section, we computed the voxel-
wise distribution of rates and intervals and compared it
between the different conditions.

Computation of Dynamical Functional

Connectivity States and its Associated Entropy

As already discussed above, for spatio-temporal fluctua-

tions arising from the dynamics of large scale systems, the
relation between the temporal fluctuations of the average
signal and its spatial correlation function is well defined
[Ross, 1966]. It is known that under general assumptions,
the mean field peak-to-peak amplitude of a signal (or its
variance) is directly proportional to its mean correlation
value. Intuitively, this is analogous to the principle that
synchronized/desynchronized clapping produces stron-
ger/weaker collective effects. Therefore, to complement
the investigation of regional changes in the variance of the
BOLD signal in various regions of the cortex, it is logical
to also look at measures of the mean functional correla-
tions between participating regions. Given the highly tran-
sient nature of these correlations, we term these indices
“dynamical functional connectivity states”. They are

r

T

agliazucchi et al.

r

r

4

r

computed as follows (see Fig. 4 for an illustration of the
general procedure).

The BOLD time series is first divided into M non-

overlapping windows of length L. Then, the partial corre-
lation values between a set of N brain regions of interest is
computed inside each temporal window as follows,

R

C

X; Y

ð

Þ5min RðX; YjZÞ

(3)

where R(X,Y]Z) is defined as follows,

R X; YjZ

ð

Þ5

R X; Y

ð

Þ2RðY; ZÞ

ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
12RðX; ZÞ

2

q

(4)

In Eq. (3), the minimum is taken across all signals Z dif-

ferent from X,Y. In Eq. (4), R(X, Y) represents the linear cor-
relation between variables X and Y. Thus, the partial
correlation RC(X,Y) measures the correlation of both varia-
bles removing the effect of a set of controlling variables. In
what follows, the set of variables will include all BOLD sig-
nals from the N brain regions under study plus the time
series of absolute head displacement (to remove spurious
correlations), as estimated during the motion correction
step. Finally, for each pair of regions a link is established
between them if the correlation P value is significant at the
level of P < 0.05, Bonferroni corrected (N (N 2 1)/2). Per-
forming this computation for all temporal windows gives
the temporal evolution of the connectivity graph.

For N brain regions, the total number of dynamical

functional connectivity graphs (i.e. all possible connectivity
motifs) equals 2

(N(N 2 1)/2)

. For example, for N 5 3, 4, 5, the

number of states is 8, 64, 1,024. Once the dynamical func-
tional connectivity states are obtained by the procedure
above described, a symbolization procedure can be used
to map the limited amount of states into discrete symbols
(a bijection or one-to-one mapping of each temporal

sequence of states into “words”—whose letters represent
the different connectivity graphs). Then, the entropy of
this sequence can be computed using the straightforward
definition of Shannon’s entropy:

H5

X

NðN21Þ=2

i51

p

i

log ð1=p

i

Þ

(5)

where p

i

is the probability of the ith state occurring in the

sequence. We evaluate the entropy of dynamical functional
connectivity states in the network of regions of peak statisti-
cal significance presented in Table I and Figure 1 (i.e., the
left and right hippocampus and ACC). This choice is moti-
vated by the well-known observation that, for a system of
coupled elements, the temporal variability of the activity
time series is driven by variability in the collective interac-
tions [Ross, 1966]. This process resulted in four regions and
a total of 64 connectivity states (i.e. the number of possible
dynamical interaction motifs between the four nodes).

Statistical Testing

To perform statistical significance testing, a paired t-test was

used, as implemented in SPM8. Only clusters passing a thresh-
old of P < 0.05 Family Wise Error (FWE) corrected for multiple
comparisons were considered significant. For display pur-
poses, significance maps were thresholded at P < 0.005, show-
ing only clusters passing the above mentioned criterion.
Significance testing for the entropy differences reported in
Figure 5 was also evaluated using a paired t-test.

RESULTS

The statistical significance maps for the BOLD signal

variance and total spectral power can be found in

TABLE I. Regions corresponding to local maxima of statistical significance (provided they are more than 8 mm

apart) for increased cortical BOLD variance and total spectral power after psilocybin infusion (compared to before

psilocybin infusion)

MNI coordinates

AAL region

Hemisphere

t-value

BOLD variance (after psilocybin > before psilocybin)

(234, 222, 216)

Hippocampus

Left

4.51

(26, 222, 216)

Hippocampus

Right

3.54

(22, 22, 28)

Anterior cingulate

Left

3.72

(4, 34, 18)

Anterior cingulate

Right

3.74

BOLD total spectral power (after psilocybin > before psilocybin)

(32, 214, 213)

Hippocampus

Right

4.03

(226, 218, 215)

Hippocampus

Left

4.39

(226, 1, 38)

Parahippocampal Gyrus

Left

4.77

(37, 22, 238)

Parahippocampal gyrus

Right

4.55

(7, 27, 17)

Anterior cingulate

Right

4.02

(27, 43, 22)

Anterior cingulate

Left

4.12

Montreal Neurological Institute (MNI) coordinates and Automatic Anatomical Labeling (AAL) regions are provided.

r

E

nhanced Repertoire of Brain Dynamical States

r

r

5

r

Figure 1.

Psilocybin infusion modifies the temporal variability of BOLD signal in
a network comprising anterior cingulate cortex (ACC) and bilateral
hippocampus. A) Maps of statistical significance for variance (r

2

) and

total spectral power (TP) increases after psilocybin infusion. Results
are shown overlaid separately into an anatomical MNI152 template,
overlaid together and also rendered together into a three-
dimensional anatomical image. In all cases only clusters surviving a
threshold of < 0.05, Family Wise Error (FWE) cluster corrected
(after passing an uncorrected threshold of < 0.005) are shown.

B) BOLD Variance time courses (obtained over a 1 min. sliding win-
dow) into the four regions of peak statistical significance defined in
Table I for the psilocybin and the placebo infusion. C) Variance of the
time course of intra-hippocampal correlations (computed over a
range of non-overlapping window lengths) before and after psilocybin
infusion. Results for placebo infusion are shown as insets (before and
after in blue and red, respectively). [Color figure can be viewed in the
online issue, which is available at wileyonlinelibrary.com.]

Figure 1A and the significance peaks are summarized in
Table I (in this and the following sections, unless explic-
itly stated no differences were found by examining the
opposite contrasts). Both measures show increased vari-
ability following psilocybin administration—both in the
temporal and spectral domain—with peaks in the ante-
rior cingulate cortex and bilateral hippocampus (total
spectral power is also increased in the bilateral parahip-
pocampal gyri). As shown in Figure 1 there is a large
overlap in the regions affected by both measures, as can
be expected from Eq. (2). The differences found for the
power spectral estimation in the parahippocampal gyrus
are not apparent for the variance. Most likely this dis-
crepancy is due to the numerical evaluation of these two
equivalent quantities.

To find the temporal evolution of the variance

increases, we selected the four significance peaks for
BOLD variance as ROIs (Table I) and used a sliding win-
dow analysis (1 min. window length) to compute the evo-
lution of the variance. Since these ROI are selected as
regions of high variance, it is already known that fluctua-
tions in their activity will have an increased variance rel-
ative to baseline. However, this analysis provides
additional

information

about

when

those

increases

occurred (e.g. just after the psilocybin infusion or later).
In Figure 1B the results of this analysis are presented. A
homogeneous increase in signal variance is oberved dur-
ing the first 3 min postpsilocybin infusion, with a tend-
ency to plateau or slightly decrease afterwards. A longer
experiment would be needed to confirm the persistence
of this effect, but it does seem to correspond well with
the known pharmaco-dynamics of intravenous psilocybin,
i.e. subjects reported the most intense subjective effects
within the first 3 min of the infusion and the effects were
persistent for the duration of the scan [Carhart-Harris
et al., 2012, 2013].

We also studied the variance of intra-hippocampal con-

nectivity time courses, computed as the linear correlation
of all hippocampal voxels over nonoverlapping windows
of different lengths. Results are shown in Figure 1C. It can
be seen that the variance is higher in both hippocampi
after the psilocybin infusion.

Next, we studied the spectral content of spontaneous BOLD

fluctuations. This was done by computing the low frequency
(0.01–0.1 Hz) power (LFP) and the power spectrum scaling
exponent (a). Statistical significance maps are presented in
Figure 2A and statistical significance peaks in Table II.

After psilocybin infusion, diffuse widespread decreases in

LFP and the scaling exponent a were observed in frontal and
parietal regions. Changes in LFP and the scaling exponent (a)
were found consistently in the same spatial locations, which is
expected since a scaling closer to that of uncorrelated noise (i.e.
an a value closer to 0) will result in weaker low frequency
spectral power. These effects were confined to the postpsilocy-
bin period. No differences were found when performing the
same comparison in the placebo condition.

To further characterize BOLD fluctuations, we trans-

formed whole brain signals into a spatio-temporal point-
process and extracted two statistics: the rate (average num-
ber of events) and the interval (average separation
between two events). By the very definition of the point-
process, a signal with high power in fast frequencies will
give a high rate and therefore a low separation between
events i.e., the average interval will be small. Results
shown in Figure 2B confirm this observation. Statistical
significance peaks are presented in Table III. A rate
increase and interval decrease after psilocybin was
observed in parietal and frontal regions largely overlap-
ping with those of Figure 2A. No differences were found
when performing the same comparison in the placebo
condition.

To identify the brain networks associated with the

changes shown in Figure 2A,B, we computed the overlap
of the statistical significance maps with a set of well-
established Resting State Networks (RSNs) [Beckmann and
Smith, 2004; Beckmann et al., 2005]. These included: two
visual (medial and lateral), an auditory, sensori-motor,
default mode, executive control and two lateralized dorsal
attention networks.

To compute the overlap for a given map, we counted

the number of voxels included in each RSN and normal-
ized by the total number of voxels in the RSN mask. Using
this approach, larger maps present a higher chance of hav-
ing large overlaps, therefore, we constructed a null-
hypothesis by randomizing the phases of the maps (after
transforming to Fourier space), resulting in images with
the same second order statistics [see also Tagliazucchi
et al., 2013c]. A total of 100 randomizations were per-
formed for each comparison and an empirical p-value was
constructed, counting the ratio of instances in which the
real overlap exceeded the overlap computed with the
randomized map. The results of these analyses are pre-
sented in Figure 3. For LFP and a, significant overlaps
were detected in the default mode, control and attention
networks, whereas all sensory networks remained unaf-
fected. For the point-process rate (PPR) and point-process
interval (PPI) the changes were locally confined to the
default mode network.

Next, the dynamic functional connectivity states were

obtained from the set of regions showing increased tempo-
ral variability (Fig. 1 and Table I) by applying the analysis
outlined in Figure 4. Using different window lengths (rang-
ing from 15 to 150 s.) we evaluated the entropy of the dis-
tribution of connectivity states in the network comprised by
two ACC ROIs and the bilateral hippocampi. An entropy
increase was found when comparing the results between
the periods before and after the psilocybin infusion and
after the psilocybin infusion vs. after the placebo infusion,
but no changes were observed when comparing before and
after placebo (see Fig. 5A). The entropy increase was not
seen for very short window lengths but was manifest for all
lengths larger than approximately 1 min.

r

E

nhanced Repertoire of Brain Dynamical States

r

r

7

r

Figure 2.

Psilocybin infusion modifies BOLD spectral content in a distrib-
uted fronto-parietal network. A) Maps of statistical significance
for decreased low frequency power (LFP) and power spectrum
scaling exponent a after psilocybin infusion. Results are shown
overlaid separately into an anatomical MNI152 template, over-
laid together and also rendered together into a three-
dimensional anatomical image. Notice that in all cases only

clusters surviving a threshold of < 0.05, FWE cluster corrected
(after passing an uncorrected threshold of < 0.005) are shown.
B) Maps of statistical significance of increased power point rate
(PPR) and decreased point process interval (PPI) after psilocybin
infusion (same renderings and statistical thresholds as in A).
[Color figure can be viewed in the online issue, which is
available at wileyonlinelibrary.com.]

The most frequent states in each condition are shown in

Figure 5B. These have sparse connectivity and represent
functional connections between homologous ROIs and
between ROIs in the same hemisphere, but cross-
hemisphere connections between hippocampal and ACC

ROIs also appear in the psilocybin condition. The psilocybin
state is also characterized by a larger repertoire of states:
i.e. novel motifs that are exclusive to the psychedelic state
and which are shown in the last row of Figure 5B. These
motifs are among the most interconnected states possible.

TABLE II. Regions corresponding to local maxima of statistical significance (provided they are more than 8 mm
apart) for decreased cortical BOLD low frequency power and power spectrum scaling exponent after psilocybin

infusion (compared to before psilocybin infusion)

MNI coordinates

AAL region

Hemisphere

t-value

BOLD low frequency power (after psilocybin < before psilocybin)

(11, 260, 40)

Precuneus

Right

4.92

(38, 267, 44)

Angular gyrus

Right

4.60

(234, 269, 44)

Angular gyrus

Left

5.05

(239, 242, 42)

Inferior parietal cortex

Left

4.66

(239, 45, 23)

Middle frontal gyrus

Left

5.10

(47, 44, 23)

Inferior frontal gyrus

Right

5.08

(5, 56, 29)

Middle frontal gyrus (orbital)

Right

4.28

BOLD power spectrum scaling exponent a (after psilocybin < before psilocybin)

(4, 259, 33)

Precuneus

Right

4.86

(53, 261, 33)

Angular gyrus

Right

4.31

(248, 255, 34)

Angular gyrus

Left

4.20

(240, 234, 60)

Postcentral gyrus

Left

4.80

(6, 222, 60)

Paracentral lobule

Right

4.79

(50, 230, 60)

Superior parietal cortex

Right

5.73

(241, 14, 33)

Middle frontal gyrus

Left

4.32

(52, 7,28)

Precentral gyrus

Right

5.12

(38, 41, 28)

Middle frontal gyrus

Right

6.22

(227, 43, 28)

Middle frontal gyrus

Left

5.10

(10, 55, 26)

Middle frontal gyrus (orbital)

Right

4.48

TABLE III. Regions corresponding to local maxima of statistical significance (provided they are more than 8 mm
apart) for increased cortical point-process rate (PPR) and decreased point-process interval (PPI) after psilocybin

infusion (compared to before psilocybin infusion)

MNI coordinates

AAL region

Hemisphere

t-value

Point-process rate (after psilocybin > before psilocybin)

(5, 252, 25)

Precuneus

Right

4.47

(47, 266, 29)

Angular gyrus

Right

3.99

(241, 266, 24)

Angular gyrus

Left

2.84

(241, 47, 25)

Middle frontal gyrus (orbital)

Left

2.82

(241, 41, 24)

Middle frontal gyrus

Left

2.97

(35, 40, 24)

Middle frontal gyrus

Right

3.14

Point-process interval (after psilocybin < before psilocybin)

(5, 257, 25)

Precuneus

Right

4.06

(50, 258, 36)

Angular gyrus

Right

3.98

(234, 266, 36)

Angular gyrus

Left

3.00

(237, 232, 62)

Postcentral gyrus

Left

3.35

(41, 239, 62)

Paracentral lobule

Right

3.18

(226, 29, 60)

Postcentral gyrus

Left

3.59

(28, 17, 54)

Middle frontal gyrus

Right

3.72

(220, 15, 54)

Superior frontal gyrus

Left

2.87

(222, 59, 12)

Superior frontal gyrus

Left

3.54

(4, 54, 6)

Anterior cingulated

Right

3.29

r

E

nhanced Repertoire of Brain Dynamical States

r

r

9

r

Figure 3.

BOLD spectral changes after psilocybin infusion are located in
higher order brain networks and leave primary sensory areas
unaffected. A) Overlap between statistical significance maps pre-
sented in Figure 2 and a group of well-established cortical RSN
(from Beckmann et al. [2005]) together with the average overlap
obtained after 100 spatial randomizations preserving first order
statistics (image phase shuffling). (*) indicates an empirical
value smaller than 0.05, Bonferroni corrected. This value is

defined as the ratio of instances in which the real maps yield a
higher overlap than the randomized versions. B) Whole brain
grey matter average probability distributions for a, LFP, PPR, and
PPI, before and after psilocybin infusion. In the inset, the same
distributions are shown before and after the placebo infusion.
[Color figure can be viewed in the online issue, which is avail-
able at wileyonlinelibrary.com.]

r

T

agliazucchi et al.

r

r

10

r