Overview

We re-analysed MEG recordings from healthy subjects with open eyes, after taking a placebo and after taking a psychedelic drug. The data come from three different experiments; in each, a different drug was administered intravenously to different participants: lysergic acid diethylamide (LSD)27, ketamine (KET)10 and psilocybin (PSIL)28. After artefact removal and source modelling (see following sections for details), we analysed 2 sec segments of 90 source channels at 600 Hz: 5–7 min data for 15 participants for LSD, 6–10 min data for 19 participants for KET and 2–5 min data for 14 participants for PSIL, each time comparing resting state MEG for the drug condition with a placebo condition.

Ethics statement

All studies were approved by a UK National Health Service research ethics committee and participants gave informed consent. Experiments were performed in accordance with relevant guidelines and regulations.

Participants and drug dose

For all three datasets, participant exclusion criteria have been previously described in detail (PSIL28, KET10, LSD27). Briefly, participants were excluded if they were younger than 21, pregnant, had personal or immediate family history of psychiatric disorder, suffered from substance dependence, had cardiovascular disease, suffered from claustrophobia, blood or needle phobia, had ever had a significant adverse response to a hallucinogenic drug, or if they had a medically significant condition rendering them unsuitable for the study. All participants had previous experience with a hallucinogenic drug but not within 6 weeks of the study (for LSD and PSIL only). For KET, participants were additionally excluded if they smoked, were female, or had a body mass index outside the range of 18–30 kg/m2.

LSD and PSIL were each administered intravenously at a fixed single dose of 75 μg and 2 mg, respectively, over the course of less than one minute. By contrast KET was administered with an initial bolus of 0.25 mg/kg delivered over one minute followed by maintenance infusion at a rate of 0.375 mg/h for forty minutes. PSIL and KET data were obtained immediately after drug administration, whereas for LSD the data were obtained four hours after drug administration due to LSD’s slow pharmacodynamics.

Data acquisition and preprocessing

Participants lay in a supine position for KET and LSD but were seated for PSIL. Pulse rates and blood oxygenation levels were continually monitored throughout the experiment via a probe over the left hand index finger. Whole-head MEG recordings were made using a CTF 275-channel radial gradiometer system sampled at 1200 Hz (0–300 Hz band-pass). An additional 29 reference channels were recorded for noise cancellation purposes and the primary sensors were analysed as synthetic third-order gradiometers. For LSD and KET, in addition to the MEG channels, ECG, horizontal and vertical participant electro-oculograms, and electromyograms from the bilateral frontalis and temporal muscles were obtained and participant compliance was monitored via an eyetracking camera.

All MEG recordings were band-pass filtered (1–150 Hz), downsampled to 600 Hz and segmented into epochs of 2 sec length. Each epoch was then visually inspected, and those with gross artifacts (e.g. head movements, jaw clenches) were removed from the analysis. An automated algorithm was used to remove further epochs contaminated with muscle artefacts. In this algorithm, a set of 30 gradiometer sensors were predefined at the edge of the MEG dewar (vacuum flask), as these are most likely to be contaminated by muscle artefacts. Using Hanning windowed fourier transformations, the mean spectral power for these sensors in the 105–145 Hz frequency band for each epoch was calculated. If the resulting power averaged across these sensors exceeded 10 fT2, then that epoch was eliminated from subsequent analysis. On the remaining epochs, independent component analysis (ICA) was performed, as implemented in Fieldtrip/EEGLAB, to identify and remove ocular, muscle and cardiac artifacts from the data. For LSD and KET, any components that showed a correlation (r > 0.1) in the time domain with the EOG/EMG electrodes were automatically removed, whereas these were identified manually for the PSIL data. Likewise, any components that showed correlations (r > 0.1) with similarly filtered EOG/EMG channels after being bandpass filtered in the range 105–145 Hz were removed. Visual inspection was also used to remove artifact components.

Source modelling of the data was performed using the fieldtrip toolbox29. For each participant, individual forward models were generated from their individual structural MRI scan30. In order to reduce the data, an atlas-based beamformer approach was used31. Broadband virtual sensor time-series were constructed using a linearly constrained minimum variance beamformer32 at 90 cortical and subcortical seed locations as specified in the automated anatomical labelling atlas33.

Prior to computing signal diversity measures, the data were further low-pass filtered with cut-off at 30 Hz to assure that possible muscle artefacts were excluded. Some residual muscle artefacts were seen in the LSD dataset27 and so the conservative approach was to apply the 30 Hz low-pass filter to all datasets.

Measures

Lempel-Ziv complexity

(LZc, LZs, LZc N , LZs N ); We computed Lempel-Ziv complexity following our previous studies19,21. As schematically shown in Fig. 1, the instantaneous amplitude (obtained via Hilbert transform) of each source channel is binarised using its mean over observations as a threshold. I.e. the continuous signal of each source channel is transformed into a string of 1200 binary digits (for our case of 2 sec segments at 600 Hz sampling rate), resulting in a matrix with binary entries with a row for each channel and a column for each observation. To assess the signal diversity across all channels and observations, this binarised data matrix is concatenated observation-by-observation into one binary string. Then the encoding step of the Lempel-Ziv 1978 (LZ78) compression algorithm (implemented by adapting open source code) is applied to this binary string. The LZ78 algorithm divides the string into non-overlapping and unique binary substrings. The more diverse the binary string, the more substrings will be listed (a sequence containing only zeros or only ones would lead to the minimal number of substrings being obtained). The total number of these substrings is what we call Lempel-Ziv complexity LZc.

Figure 1: Schematic of the computation of LZc and SCE. LZc: (a) x i is the activity of the ith channel and a i is the (Hilbert) amplitude of x i . (b) b i is a i binarised, using the mean activity of a i as the binarisation threshold. (c) After binarisation of all n signals, (d) the multidimensional time series are concatenated observation-by-observation into one binary sequence k and then (e) repeated patterns are searched and listed into a dictionary of binary words via a Lempel-Ziv algorithm. Lempel-Ziv complexity LZc is proportional to the size of this dictionary. SCE: (a) Two time series. (b) The analytic signals of these two, which are complex signals with the real part being the original signal and the imaginary part being the Hilbert transform of the original signal. (c) A binary synchrony time series is created for this pair of signals; a 1 indicates that the phases of the complex values of the analytic signals are similar (difference of less than 0.8 modulo 2π). (d) Such time series are obtained to represent each channel’s synchrony with seed channel i. e) SCE(i) is the entropy over observations in the resulting data matrix. The overall SCE is then the mean value of SCE(i) across choices of seed channel i. Full size image

Given the observation-by-observation concatenation of the binarised data matrix, LZc captures temporal signal diversity of single channels as well as spatial signal diversity across channels. In order to assess temporal signal diversity only, we further applied this procedure to single source channels independently, and we denote the resulting single channel Lempel-Ziv complexity by LZs. I.e. LZs quantifies the temporal signal diversity of single channels, with a high score for a uniformly random binary string and a low score for a string of zeros only.

We normalize LZc (also LZs) by dividing the raw value by the value obtained for the same binary input sequence randomly shuffled. Since the value of LZc (also LZs) for a binary sequence of fixed length is maximal if the sequence is entirely random, the normalized values indicate the level of signal diversity on a scale from 0 to 1.

In order to test if changes in the diversity measures in the drug vs placebo contrast can be explained away by previously characterised changes in the overall power spectrum, we employed the following control procedure (as in ref. 19). We generate surrogate time series through phase-shuffling of the data, maximising signal diversity as measured by non-normalised LZs (analogously non-normalised LZc) under the constraint of preserving the spectral power profile of the data. If the observed difference in non-normalised LZs between drug and placebo was completely preserved in the equivalent contrast for the surrogate data, we would conclude that the observed difference is entirely due to changes of the power spectrum and not due to changes in signal diversity beyond spectral changes. If, on the other hand, the difference was completely absent for the surrogate data, we would conclude that the observed difference in non-normalised LZs was entirely due to changes in signal diversity that were not expected from the spectral power profiles.

We apply this rationale by comparing the ratio of diversity scores for data and surrogate data as follows. Assume that D 1 is the LZs (LZc) score for the placebo condition (state 1) and D 2 that for the drug condition (state 2). Assume D 1 < D 2 and let N 1 be the average LZs score of the phase-shuffled data of state 1 (analogously N 2 for state 2). If D 1 /N 1 > D 2 /N 2 , then it must be that N 1 < N 2 . Thus, this outcome implies that the difference in signal diversity between the states, D 1 < D 2 , can be entirely explained by the difference in power spectra alone, as the maximal signal diversity given the spectral power profile of state 1 is much smaller than that of state 2. Conversely, without such reversal of results we can exclude that the observed signal diversity changes are entirely due to changes in the spectral power profile, i.e. some of the difference in diversity between the states must be due to the specific properties of the data captured by LZs (analogously LZc).

For a given participant and condition we thus computed LZc (LZs) for all 2 sec data segments. Then we phase-randomised each channel’s time series for each 2 sec segment and then recalculated LZc (LZs). Phase randomisation was performed by first applying a discrete Fourier transform, then randomising phases of the frequencies while keeping their amplitudes fixed, and then applying the inverse Fourier transform (i.e. only the relative temporal positions [phases] of the Fourier sinusoids, whose superposition equals the signal, are randomly changed). For a given participant and condition, measuring the average diversity across many phase-shuffled data segments gives the highest signal diversity possible for the spectral power profile for that participant in that condition, i.e. N 1 . We denote Lempel Ziv complexity measures normalised by average scores for phase-randomised data (D 1 /N 1 ) as LZc N (LZs N ).

Coalition entropy

(ACE, SCE, ACE N , SCE N ); Synchrony coalition entropy (SCE) is a measure of the entropy over time of the constitution of the set of channels that are in synchrony, see ref. 21 for details and Fig. 1 for a schematic. SCE normalized by its score for phase-randomised data is denoted as SCE N .

Amplitude coalition entropy (ACE) is defined as in ref. 21 as the entropy (over time) of the constitution of the set (coalition) of channels that are ‘active’, given the binarization scheme described above for LZc for defining ‘active’/‘inactive’ channels. As for LZc, we normalise ACE by its value for a random shuffling of the data. We further consider renormalisation by the mean value obtained from phase-randomised surrogate data (ACE N ), as described above for LZc.

Note that ACE and SCE are only defined for multidimensional time series, capturing signal diversity over channels and time.

Normalized spectral power and phase coherence (PC)

The behaviour of the signal diversity measures is compared to that of normalized spectral power and phase coherence. We defined power bands as: δ = 1–4 Hz, θ = 4–8 Hz, α = 8–15 Hz, β = 15–30 Hz and γ = 30–70 Hz. For the computation of spectral power, the data was not low-pass filtered at 30 Hz, as was the case prior to the computation of the signal diversity measures. The power of a spectral band was computed using Welch’s method34 for each 2 sec segment of each of the 90 sources, normalised by the sum of the power of all 5 bands, then averaged across sources and trials per subject.

As a measure of overall synchrony we use the mean phase coherence (PC) across all pairs of channels. Let and describe the analytic signals of two channels at time t. Then

with T being the length of the segment (1200 observations = 2 sec) and N = 90 the number of sources. Reported PC scores are averages across trials. Scores of PC lie between 0 and 1, with 1 indicating perfect synchrony and 0 indicating no synchrony at all (phase differences uniformly randomly distributed).

Questionnaire scores

Across the three experiments, different questionnaires were employed for participants to retrospectively evaluate their psychedelic experience. Here, we analyse a subset of the questionnaire items that were common across all three experiments, i.e. for PSIL, KET and LSD:

strange: “Things looked strange”.

geom: “I saw geometric patterns”.

vivid: “My imagination was extremely vivid”.

time: “My perception of time was distorted”.

space: “My sense of size and space was distorted”.

ego: “I experienced a disintegration of my ‘self’ or ‘ego’”.

muddle: “My thinking was muddled”.

merge: “I experienced a sense of merging with my surroundings”.

control: “I feared losing control of my mind”.

spirit: “The experience had a spiritual or mystical quality”.

peace: “I felt a profound inner peace”.

float: “I felt like I was floating”.

past: “I saw events from my past”.

sounds: “Sounds influenced things I saw”.

For all three experiments subjective questionnaires were completed retrospectively on the day of the experiments after most drug effects had subsided, asking participants to recall experiences they had at the time of peak drug effects. For KET and PSIL this was typically an hour after drug delivery had ceased but was approximately 8 hours after LSD due to its relatively prolonged pharmacodynamic profile. MEG scans were taken for KET and PSIL at peak drug effects, whereas for LSD, MEG scans were typically taken one to two hours post peak effects, i.e. immediately post drug administration for KET and PSIL, and around four hours post drug administration for LSD. Participants answered each question using a visual analogue scale format with a bottom anchor of “no, not more than usually” and a top anchor of “yes, much more than usually”. In addition, we consider the mean score over all these questions as an index of overall intensity of the psychedelic state. We call this index “total”.

Further, we obtained a single subjective rating of the overall drug-effect intensity, acquired while participants were inside the MEG scanner. This rating is for each drug denoted as “InScanner” and was obtained for LSD27 and PSIL28 by asking participants “please rate the intensity of the drug effects during the last scan” and for KET10 by asking “please rate your subjective high on a scale between 0 and 40”, using a two-digit button box.

Statistics

Analyses were performed using non-overlapping segments of length 2 sec for a total length between 2 min and 10 min of MEG recording per participant and state. For each segment, the signal diversity measures ACE, SCE and LZc were computed for 30 random picks of 10 channels, and the mean across these 30 scores was considered the score for the segment. We chose 10 channels since in our previous study19 this was the smallest channel number for which we still found a robust decrease of ACE, SCE and LZc for EEG signals in propofol-anaesthesia. To verify that the results obtained were not dependent on the particular random channel selection, we performed a re-run of the complete analysis and indeed obtained almost identical results. For LZs, the mean across all 90 channels was set as the score for the segment. The mean and standard error of the diversity measures’ scores were computed across segments. At the single participant level, the effect size of differences between states was measured using Cohen’s d35. We call an effect size high if d > 0.7. For group level comparisons, a two-sided t-test was applied, with Bonferroni correction (by the number of measures) where indicated.

Computation of correlation between measures