Working memory load modulates oscillatory activity and the distribution of fast frequencies across frontal theta phase during working memory maintenance

Neurobiology of Learning and Memory 183 (2021) 107476 Available online 2 June 2021 1074-7427/© 2021 The Authors. Published by Elsevier Inc. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/). Working memory load modulates oscillatory activity and the distribution of fast frequencies across frontal theta phase during working memory maintenance Alba Fern´ andez a , * , 1 , Diego Pinal b , 1 , Fernando Díaz a , Montserrat Zurr´ on a a Cognitive Neuroscience Laboratory, Department of Clinical Psychology and Psychobiology, Universidade de Santiago de Compostela, Santiago de Compostela, Galiza, Spain b Psychological Neuroscience Lab, Escola de psicologia, Universidade do Minho, Portugal ARTICLE INFO Keywords: Neural oscillations Working memory Cross-frequency synchronization ABSTRACT Working memory (WM) is a keystone of our cognitive abilities. Increasing load has been shown to dampen its performance and affect oscillatory neural activity in different frequency bands. Nevertheless, mixed results regarding fast frequencies activity and a lack of research on WM load modulations of cross-frequency phaseamplitude coupling mechanisms preclude a better understanding of the impact of increased WM load levels on brain activity as well as inter-regional communication and coordination supporting WM processes. Hence, we analyzed the EEG activity of 25 participants while performing a delayed-matching-to-sample (DMS) WM task with three WM load levels. Current density power and distribution at the source level for theta, beta, and gamma frequencies during the task’s delay period were compared for each pair of WM load conditions. Results showed maximal increases of theta activity in frontal areas and of fast frequencies’ activity in posterior regions with WM load, showing the involvement of frontal theta activity in WM maintenance and the control of attentional resources and visual processing by beta and gamma activity. To study whether WM load modulates communication between cortical areas, posterior beta and gamma amplitudes distribution across frontal theta phase was also analysed for those areas showing the largest significant WM load modulations. Higher beta activity amplitude at bilateral cuneus and right middle occipital gyrus, and higher gamma activity amplitude at bilateral posterior cingulate were observed during frontal theta phase peak in low than high memory load conditions. Moreover, greater fast beta amplitude at the right postcentral gyrus was observed during theta phase trough at right middle frontal gyrus in high than low memory load conditions. These results show that WM load modulates whether interregional communication occurs during theoretically optimal or non-optimal time windows, depending on the demands of frontal control of posterior areas required to perform the task successfully. 1. Introduction Working memory (WM) is the capacity that allows us to retain and manipulate for brief periods of time small amounts of information no longer available in the environment (Baddeley, 1998, 2003). WM involves different subprocesses, namely, the initial encoding of information in temporary stores, the maintenance during a few seconds of its representation in those stores and, finally, the retrieval of this information for its use. Different research has studied WM using delayed–matching–to–sample (DMS), n-back, and Sternberg tasks. Among these, one of the most used tasks to study WM subprocesses is the DMS task, in which participants must memorize non-verbalizable patterns of items and, after a brief delay, correctly match them to a presented probe pattern. This task can also be used to study WM load’s effect by increasing the number of items to be memorized. Such increases in WM load are accompanied by decreased task performance, reflected in longer reaction times and reduced response accuracy rates (Alvarez & Cavanagh, 2004; Vogel, Woodman, & Luck, 2001). Electroencephalography (EEG) is a technique that allows recording neuroelectric oscillatory activity. Thus, it enables the exploration of the neural bases of WM subprocesses, such as maintenance, and the effects * Corresponding author. E-mail address: [email protected] (A. Fern´ andez). 1 Contributed equally to this work. Contents lists available at ScienceDirect Neurobiology of Learning and Memory journal homepage: www.elsevier.com/locate/ynlme https://doi.org/10.1016/j.nlm.2021.107476 Received 2 November 2020; Received in revised form 13 May 2021; Accepted 30 May 2021 Neurobiology of Learning and Memory 183 (2021) 107476 2 of WM load on them. Studies that employ EEG oscillatory activity analyses have demonstrated that oscillatory neural activity in slow frequency bands such as theta (4 – 7/8 Hz) may play a key role in temporarily maintaining information in WM stores. Probably theta activity recorded in frontal electrodes reflects brain activity related to attentional control over the stored information or another kind of central executive functions (for a review, see Sauseng, Griesmayr, Freunberger, & Klimesch, 2010). Regarding fast oscillations (i.e., frequencies over 13 Hz), many studies have observed their presence throughout various posterior cortical regions during the maintenance period in WM tasks. However, the results of these studies are inconsistent as to how beta and gamma activity are implicated in WM processes (for a review see Pavlov & Kotchoubey, 2020). Hence, some studies have reported increases in their amplitude relative to a baseline period in visual and occipitotemporal areas (Honkanen, Rouhinen, Wang, Palva, & Palva, 2015; Tallon- Baudry, Bertrand, Peronnet, & Pernier, 1998), while other work reports suppression of these frequencies in those same areas during the execution of similar tasks (Brookes et al., 2012; Proskovec, Wiesman, Heinrichs-Graham, & Wilson, 2018). These results, however, are obtained in paradigms with a single WM load condition. WM load increases the demands made on neural resources to achieve a successful task performance. Thus, variations in its level are expected to impact oscillatory activity related to WM processing. Indeed, theta oscillatory activity has been shown to be modulated by this factor. Therefore, theta amplitude in a spatial version of an n-back task (in which participants must retain a sequence of items to judge whether the current item matches an item n places back in the sequence) was larger at frontal midline sites as WM load demands (i.e., the number of items back) increased (Gevins, 1997). Theta activity has also been found to increase with WM load over frontal midline sites during the delay period of a DMS task (Eschmann, Bader, & Mecklinger, 2018) as well as of a Sternberg task (participants are asked to memorize a list of items to decide upon presentation of a probe if it was included on the previous list) with nonverbal symbols (Maurer et al., 2015). Neural sources of this theta activity have been identified at the medial prefrontal cortex (Kaplan et al., 2016; van Ede, Jensen, & Maris, 2017) However, with beta and gamma oscillations, overall the picture is not as clear (for a review see Pavlov & Kotchoubey, 2020). When looking at the influence of WM load in beta and gamma activity during the maintenance period of WM tasks, the results reported in previous studies are inconsistent, showing both increases and suppression of this activity depending on the specific work. Some studies have observed increased beta and gamma amplitude with WM load. For example, in a study with intracranial recordings, increased beta amplitude during the maintenance period of a Sternberg task was found over prefrontal, parietal, and temporal areas for longer than shorter lists of letters presented sequentially (Howard et al., 2003). Similarly, using a DMS task where participants had to memorize a sample stimulus with 1 to 6 squares, gamma amplitude during the maintenance period was also observed to increase with greater WM loads over frontal, parietal, and temporal regions in MEG data analyses (Palva, Kulashekhar, Hamalainen, & Palva, 2011). On the other hand, a similar number of studies have found either no significant change in fast oscillatory activity with WM load changes or even a decrease of beta and/or gamma activity with WM load modulations. For example, Proskovec, Wiesman, Heinrichs-Graham, and Wilson (2019) found a reduction of beta activity with increased WM load in posterior regions during the maintenance period of a visual DMS task where participants had to memorize the location of 2 or 4 black squares located in a 7 x 9 grid. Also, regarding gamma oscillations, Poch, Campo, and Barnes (2014) found no effect of WM load on its amplitude during the maintenance period of a visual retro-cueing WM task where participants had to memorize the location and orientation of 1, 2, or 4 rectangles to later match them with a probe. Likewise, Pahor and Jauˇ sovec (2017) found no effect of WM load on gamma activity during the maintenance period of a Sternberg task where participants had to memorize a visual array of 4, 6, or 8 colored squares. Fast EEG frequencies are less well studied than theta activity, especially with non-invasive techniques such as scalp recorded EEG, which may be one contributing factor to these mixed results. Further, even when the same or similar tasks are used in different studies, there are still differences in the particular parameters selected (i.e., the number of and the features of the items to retain in memory for each WM load condition, the type of stimuli used, the duration of the delay period, etc.). Given that both frequencies are involved in many cognitive processes beyond WM maintenance, it is difficult to dissect what is precisely driving the increase or decrease of this activity. Additionally, only a few of the previously mentioned studies have attempted to identify the cortical origins of such scalp-recorded EEG oscillatory activity. In summary, several studies report that WM load modulates these fast frequencies. However, the direction of this modulation remains unclear, and it seems to be more widely distributed in the scalp than in the case of theta oscillations (Honkanen et al., 2015; Howard et al., 2003; Palva et al., 2011; Proskovec et al., 2019). Besides the independent response of theta, beta, and gamma bands to WM processes and their WM load modulation, there has been a growing interest in the interaction between slow (theta) and fast oscillatory (beta and gamma) activity during WM processing. This interest stems from recent theoretical models highlighting the relevance of the interaction between oscillatory brain activity in different frequencies to neural communication and activity coordination (e.g., Canolty & Knight, 2010; Roux & Uhlhaas, 2014). Experimental evidence in support of these models includes reports of slow frequencies phase and fast oscillatory activity phase or amplitude synchronization in n-back (Yang & Huang, 2018) and DMS tasks as recorded with EEG (Sauseng et al., 2009) as well as with MEG (Siebenhühner, Wang, Palva, & Palva, 2016). For example, it has been observed that with higher memory loads, there is an increase of theta – gamma phase synchronization during WM maintenance that is unaffected by irrelevant items and is positively correlated with the individuals’ WM capacity (Sauseng et al., 2009). Regarding the synchronization produced between the phase of a slow frequency (e.g., theta) and the amplitude of a faster frequency (e.g., gamma or beta) in the same or a different cortical area; theoretical models of brain function suggest that communication between neuronal populations depends on these groups being phase synchronized and thus having their windows for input and output open simultaneously (Fries, 2005). Therefore, it is hypothesized that the slow frequency cycle will determine time windows of enhanced or decreased receptivity to establish functional links (Fries, 2005). These functional links will be established whenever the periods of higher amplitude of the faster oscillation co-occur with a determined phase angle of the slower frequency in the same or different brain regions. In this regard, it is believed that neuronal populations will be more or less likely to be excited, depending on the phase angle of the slow oscillatory activity. Specifically, they are supposed to be more likely to be excited during the trough and less likely during the peak (Fries, 2005). When this phenomenon is produced between distant brain regions, it is known as interregional phase – amplitude coupling (PAC). This phenomenon has been proposed to play an essential role in memory processes in general and in WM in particular (Roux & Uhlhaas, 2014). Furthermore, regarding WM load’s influence in this mechanism, intracranial recordings in the human hippocampus have shown that the specific theta frequency at which gamma amplitude was coupled decreased as a function of WM load (Axmacher et al., 2010). Evidence for this theta-gamma mechanism in memory processes, and in WM specifically, was provided by a study in which transcranial alternating current stimulation (tACS) was used to alter the individual ratio of theta to gamma frequencies. To that end, the authors induced slower theta rhythms than the individuals’ natural theta frequency showing maximum coupling to gamma amplitudes in a subset of the experimental participants. Those participants receiving tACS would then had more A. Fern´ andez et al. Neurobiology of Learning and Memory 183 (2021) 107476 3 gamma cycles produced within a theta cycle, which was hypothesized to increase their WM span. Indeed, those participants demonstrated an increased digit span during stimulation than before it, while the performance of those participants receiving sham stimulation did not change (Vosskuhl, Huster, & Herrmann, 2015). Nevertheless, the observed modulation of frequency within theta range in PAC mechanisms does not provide any information or evidence regarding the theta activity phase angle at which fast oscillatory activity is coupled. Indeed, to the best of our knowledge, it has not been previously studied whether WM load could modulate the theta frequency phase window at which fast oscillatory activity (beta and gamma) is nested. This would shed new insights on how cortical areas communicate depending on WM load. In this work, our main goal is studying whether WM load influences interregional cortical communication between the brain regions showing the largest modulation of oscillatory activity by WM load. To do so, we will analyse how WM load modulates the relationship between amplitudes of fast frequencies (beta and gamma) with the phase of a slower frequency (i.e., theta) into source space. With this aim in mind, in the first part of the manuscript we will use source estimation to investigate the neural sources that show the largest effects of WM load manipulations in the theta, beta and gamma frequency bands during the maintenance of information in working memory. In the second part of the manuscript, we will study the cortical communication between these regions. In order to do this, oscillatory brain activity was recorded during the maintenance period of a visual DMS task with three different memory load levels, and source estimation was used to transform scalp recorded activitiy into source space activity. A DMS task was selected because it allowed us to accurately delimit the maintenance period in WM with no interfering processes. 1.1. The working hypotheses were: (1) WM load should modulate theta and fast frequencies (beta and gamma) amplitudes during the maintenance period of a DMS task; and for theta, it was expected that amplitude would increase as a function of WM load. (2) The distribution at the source space level of the largest differences between WM load conditions should differ among frequencies. For theta, it was expected that the maximum increase in power with WM load would be located in prefrontal areas. (3) Communication between areas with the largest difference in power between WM load conditions should also be modulated by WM load. It was expected that WM load modulates the distribution of fast frequencies amplitude across theta phase, in areas where the modulation of WM load was maximal. Specifically: for each pair of communicating areas, we expected that beta and gamma amplitudes during theta trough were greater in higher WM load conditions compared to their amplitudes in lower WM load conditions. Consequentially, it was also expected that beta and gamma amplitudes during theta peak, were greater in lower WM load conditions compared with higher WM load conditions. 2. Methods 2.1. Participants A total of thirty-two healthy volunteers (20 women, mean age 21 ± 3.3 years old) participated in the study. All of them were right-handed as assessed by the Edinburgh Handedness Inventory (Oldfield, 1971) and had normal or corrected-to-normal vision and no history of psychiatric or neurological disorders. Seven women had to be excluded from data analysis due to excessive artefacts in the EEG data (see section 2.4 for details). The remaining 25 participants had a mean age of 21.1 ±2.70 years old. Participants reported no recent consumption of alcohol or recreational drugs and were instructed to abstain from caffeine and other stimulants at least 30 min before the EEG recording. All participants signed an informed consent form before the experimental session. The study was approved by the ethics on investigation autonomic committee of Galicia (CAEIG, Galicia, Spain; code: 82017/498). 2.2. Task Participants performed a delayed-matching-to-sample visuospatial task, preceded by a 6-trials training block in the same task. The task contained 120 trials divided into three blocks with a 90-second interval between blocks. Each trial began with a warning sound (1000 Hz, 150 ms duration). After a jittered pre-stimulus interval of 450–550 ms, a sample stimulus (encoding stage) appeared, and remained 1750 ms on screen. The sample stimulus consisted of three domino tiles, and participants were instructed to memorize the number of dots and their position within each tile. After a maintenance period (3250–3750 ms), a probe stimulus (information retrieval stage) appeared on screen for 1750 ms. The probe stimulus also comprised three domino tiles that could match the sample stimulus or have one of the dots in a different location within one of the domino tiles. Participants had to decide by button-press whether the probe stimulus matched the previously encoded sample stimulus or not. Between trials, there was a 900–1100 ms jittered interval (Fig. 1). Each tile was an 8 ×4 cm rectangle formed by two white squares, presented on the centre of a flat-screen located at a distance of 1 m from the participant, covering a visual angle of 4.58◦x 2.29◦. Tiles could have between 3 and 4 dots randomly distributed across the four corners of each square (i.e., eight different possible locations). The dots were black and were located 0.5 cm from the tile edges and 1 cm from each other. Furthermore, the number of dotted tiles in the sample and probe stimuli was manipulated to create three different working memory load conditions. Given that only EEG epochs with correct answers would be analysed and that a larger number of errors was expected in the higher load conditions, the number of trials in each condition was adjusted. Hence, more trials were included for conditions with higher loads to compensate for these facts and still retain a similar number of valid epochs in each condition for posterior analysis. Thus, on the low working memory load (LL) condition, both tiles at the sides were blank -i.e., empty of dots- (30 trials); on the medium WM load (ML) condition, the centre tile was blank (40 trials), and on the high WM load (HL) condition all tiles had dots inside (50 trials). Stimulus’ presentation was semi-randomized so that a maximum of three “match” or “non-match” consecutive trials could appear. 50% of trials were “match”, and 50% were “non-match” trials. 2.3. Behavioural data Percentage of correct responses (hits) and reaction times measured from probe stimulus onset to button press only from trials with correct responses were calculated for each participant and WM load condition. 2.4. EEG recording and signal processing Participants sat in a comfortable armchair inside a noise and light attenuated Faraday chamber. EEG activity was recorded between 0.001 and 100 Hz with a 50 Hz notch filter and digitized at 500 Hz from 60 Ag- AgCl active scalp electrodes positioned according to the 10–10 system, with nose tip reference and the ground electrode at Fp1 location. Besides, vertical and horizontal EOG were recorded with two electrodes positioned on the outer canthi of both eyes (HEOG); and two electrodes placed above and below the right eye (VEOG). Electrode impedances were kept below 10 kΩ. All 32 continuous data recordings were inspected visually and with the Raw Data Inspection utility in Brain Vision Analyzer to determine A. Fern´ andez et al. Neurobiology of Learning and Memory 183 (2021) 107476 4 their suitability for further analysis. The following parameters were use in Raw Data Inspection: for all 60 scalp electrodes, it would mark as bad sections of data from 200 ms before to 200 ms after a difference larger than ±200 μ V in any 200 ms window. It would also mark voltage steps larger than 50 μ V and segments with 0.5 μ V or lower amplitudes for 100 ms or longer intervals. Recordings marked with one or more of these artefacts in more than 20% of the entire data for more than 6 channels were excluded from further analysis. As a result, recordings of 7 female participants were discardaded from further analyses leaving a sample size of 25 participants. Data from the remaining 25 participants were then offline resampled to 512 Hz to optimize them for posterior time–frequency analysis. A phase shift-free Butterworth filter between 0.5 and 80 Hz (12 dB/octave roll-off) was applied. Ocular artefacts were corrected using independent component analysis (ICA) in Brain Vision Analyzer, and only components compatible with ocular movements were rejected. During this step, all components and topographic maps for each participant were also revised in order to control and verify that components picking up low amplitude muscular noise (i.e., extracephalic components with frequency peaks above 30 Hz) or other noise sources (e.g. noisy channels) were not entered in signal reconstruction for further analyses. EEG data were segmented in 5500 ms epochs from sample stimulus onset (i.e., comprising the maintenance period) for those trials with correct answers. Semi-automatic artefact rejection was also applied to exclude EEG segments with data points exceeding ±125 µV. Epochs were further grouped into low, medium, and high working memory load (LL, ML, and HL) trials for each participant (mean epochs for LL =27.64 ±2.14; mean epochs for ML =32.44 ±3.20; mean epochs for HL = 32.84 ±4.38). 3. Data analyses and results Due to part of the analyses being dependent on the results of previous steps, it was decided to group both analysis and results in the same section in order to facilitate their reading. 3.1. Behavioural analyses and results To test for potential effects of WM load on task performance, percentage of hits and reaction times were subjected to repeated-measures ANOVAs with the within-subjects factor WM Load (LL, ML, HL). Greenhouse–Geisser correction was used when the sphericity assumption was violated, and Bonferroni correction was employed for posthoc pairwise comparisons. Visual working memory load significantly dampens task performance, for both percentage of hits F(2,48) =134; p <0.001, and reaction times F(1.5, 36.3) =139.8; p =0.01. Participants’ reaction times were significantly longer for HL (mean =1115.6 ms, SD =234.5 ms) than both, ML (mean =980.7 ms, SD =210 ms) and LL (mean =709.4 ms, SD =130.2 ms) conditions (p <0.001 and p <0.001, respectively); as well as for ML than LL (p <0.001). In addition, percentage of hits significantly decreased with higher WM load (LL: mean =95.1%, SD = 5.5; ML: mean =84.1%, SD =8.8; HL: mean =67%, SD =9). Thus, percentage of hits was significantly lower for HL than ML and LL (p < 0.001 and p <0.001, respectively) and lower for ML than LL (p <0.001). 3.2. Current density power analyses and results The initial epoch was further limited to a 2000 ms epoch following the offset of the sample stimulus, so that only the maintenance period was included. In order to assess potential WM load-related differences on oscillatory activity power at the cortical sources of scalp-recorded EEG signals, standardized low-resolution electromagnetic tomography analysis (sLORETA) software (Pascual-Marqui, 2002) was used. First, for each single frequency and electrode pair in the epoched data, sLORETA uses the mean voltage values to calculate the crossspectra, which is then averaged across epochs. Thus, a full crossspectra matrix is obtained for each participant. These cross-spectra matrices are feed to sLORETA algorithms to calculate the current source density power estimates for 6430 voxels of 5 ×5 mm representing cortical grey matter and hippocampus, based on a 3-shell spherical head model registered to the Talairach brain atlas (Pascual- Marqui, 2002). Finally, differences in mean current density power between each pair of WM load conditions (i.e. HL vs LL, HL vs ML, and ML vs LL) were evaluated for each of the following frequency bands: theta (4–7 Hz), beta 1 (13 – 20 Hz), beta 2 (20–30 Hz), gamma 1 (30–50 Hz) and gamma 2 (50–70 Hz) using statistical non-parametric mapping (SnPM) as implemented in sLORETA. Hence, independent analyses were run to compare current density power maps for each pair of WM conditions and frequency band (i.e., 15 independent SnPM analyses) using 5000 randomizations of paired samples t-tests on log transformed data. Each one of these randomizations calculated 6430 paired samples t-tests (one for each cortical gray matter voxel) and stores the two most extreme values (i.e., one per tail of the distribution) of the t-statistics to build a distribution of extreme t-statistics. The t-statistic obtained with the observed Fig. 1. DMS task used in the experiment. A. Fern´ andez et al. Neurobiology of Learning and Memory 183 (2021) 107476 5 data was then compared against this distribution of extreme statistics to assess its statistical significance using an alpha threshold set at p <0.05. This procedure allows us to correct for multiple comparisons when establishing the significance of the observed t score and do not require the assumption of Gaussianity (for details see Nichols & Holmes, 2002). All the gyri containing significant voxels in the Brodmann area (BA) showing the most extreme observed t-statistic for each SnPM test (see Table 1 and Fig. 2) were selected as regions of interest (ROIs) for the interregional cross-frequency coupling analyses. A comprehensive list of the results of these current density power maps comparisons for each WM load conditions pair and frequency can be found in the supplementary information file. 3.3. Analysis and results of the distribution of fast frequencies amplitude across theta cycle angles sLORETA current density power analysis, showed that WM load modulations of oscillatory activity were largest in posterior areas for beta and gamma frequency bands and in frontal regions for the theta frequency. Thus, we assessed whether the distribution of posterior beta and gamma amplitudes across the frontal theta cycle angles was modulated by WM load. To do this, we used the initial 5500 ms epoch starting at sample stimulus onset. Based on the original voltage values, the LORETA algorithm, as implemented in Brain Vision Analyzer 2.1 (BVA), was used to estimate for each data point in the epoch the current density at the regions that showed significant WM load-related activity modulations in the previous analysis (see Table 1). In this step, gyri were used instead of voxels due to the source space in BVA LORETA transformation being restricted to 2394 voxels at 7 ×7 mm spatial resolution, making it difficult to have a one–one correspondence with the voxels form sLORETA analyses. A continuous complex Morlet wavelet transformation with Gabor normalization and a Morlet parameter c of 5 was applied to these current density estimates in order to calculate their time–frequency decomposition between 1 and 70 Hz in 30 frequency steps. After this step, to avoid edge artefacts, the data was segmented in 3000 ms epochs starting at sample stimulus offset. These epochs were averaged for each participant and experimental condition. For each experimental condition and participant average, instantaneous theta phase values were extracted from the bilateral anterior cingulate cortex, and right inferior and middle frontal gyri ROIs. In addition, instantaneous current density amplitude values for beta 1, beta 2, gamma 1, and gamma 2 were extracted for the following posterior ROIs: left fusiform gyrus, bilateral cuneus, precuneus and posterior cingulate, right middle occipital gyrus, superior parietal lobe, and postcentral gyrus. Note that for medial regions showing current density power modulations by WM load, such as cuneus, precuneus, and cingulate gyrus, we decided to analyse both hemispheres due to the sLORETA estimation of significant activity being too close to the midline. A detailed list of the voxels comprising each ROI is available at the Supplementary Information file. To calculate the distribution of fast frequencies amplitude along frontal theta cycle, the same procedure as used by Sauseng et al. (2009) and Pinal, Zurr´ on, Díaz, and Sauseng (2015) was used. Beta 1, beta 2 gamma 1, and gamma 2 instantaneous current density amplitude values from each of the aforementioned posterior ROIs were z-transformed for each single epoch, frequency band, and ROI independently. The posterior ROIs’ z-transformed current density amplitude values were then grouped for each frontal ROIs’ theta cycle angle and sorted according to these instantaneous theta phase values. This step was done independently for each of the studied frontal ROIs. These theta-phase-sorted ztransformed beta and gamma instantaneous current density amplitude values were independently averaged for each frontal ROI into 16 frontal theta phase bins, each bin spanning 22.5◦of a theta cycle of 360◦. The equivalence between bins and the phase angles in degrees can be seen in Fig. 3. As a result, for each of the three frontal ROIs, there are 16 theta phase bins with their correspondent current density beta 1, beta 2, gamma 1, and gamma 2 z-transformed amplitude values from each of the seven posterior ROIs. This procedure was done separately for each WM load condition using custom made MATLAB R2016a scripts. The number of trials included in these computations was equated between groups and working memory load conditions. Finally, for each pair of frontal and posterior ROIs, theta – phase sorted z-transformed current density amplitude values for beta 1, beta 2, gamma 1, and gamma 2 in each WM load condition were entered into repeated – measures ANOVAs. The focus of this study was exclusively the interaction between the two within – subjects factors: WM load (LL, ML, and HL) and Frontal Theta Phase (16 phase bins). Therefore, independently for each pair of frontal and posterior regions, a repeated – measure ANOVA was made with Frontal Theta Phase (16 phase bins) and WM Load (LL, ML, HL) as within-subject factors for beta 1, beta 2, gamma1 and gamma 2. The main interaction effects for the four frequency bands were controlled for multiple comparisons with the Holms- Bonferroni correction; and Greenhouse – Geisser correction were applied in all cases the condition of sphericity was not met. For all the significant interaction effects, post hoc comparisons were computed with the Bonferroni correction. The post hoc comparisons focused on contrasting the amplitude of a given fast frequency at each single frontal theta-phase bin between WM load conditions. Furthermore, an associated p-value of 0.01 was used as statistical threshold Table 1 Brodmann areas, gyri and voxel coordinates with extreme observed t-statistics from the current density power analyses. For each SnPM test with statistically significant results, it is listed the Brodmann Area (BA) with the most extreme t statistic, the label of the gyri showing significant voxels in that BA, the coordinates in the Montreal Neurologic Institute space of the peak voxels, and the t and associated p -values for those peak voxels. Frequency band Compared conditions Brodmann Area X(MNI) Y(MNI) Z(MNI) Gyrus t Extreme - p Theta (4–7 Hz) ML >LL 47 40 40 −5 Right middle frontal gyrus 4.23 0.0055 40 30 −5 Right inferior frontal gyrus 4.08 HL >LL 25 0 10 −5 Right anterior cingulate 5.20 0.0001 Beta 1 (13–20 Hz) ML >LL 18 & 30 5 −80 5 Right cuneus 5.23 0.0008 18 5 −85 0 Right middle occipital gyrus 5.09 30 20 −65 10 Right posterior cingulate 4.95 HL >LL 30 5 −70 10 Right posterior cingulate 4.75 0.0007 Beta 2 (20–30 Hz) HL >LL 18 −25 −70 15 Left posterior cingulate 3.77 0.0156 −20 −70 15 Left cuneus 3.77 ML >HL 20 −35 −40 −25 Left fusiform gyrus 3.57 0.0482 Gamma 1 (30–50 Hz) No significant results were found for any comparison Gamma 2 (50–70 Hz) ML >LL 7 25 −60 65 Right superior parietal lobe 3.49 0.0413 HL >LL 7 30 −50 65 Right superior parietal lobe 3.69 0.0279 25 −55 70 Right postcentral gyrus 3.50 ML >HL 7 & 31 20 −65 30 Right precuneus 3.45 0.0188 A. Fern´ andez et al. Neurobiology of Learning and Memory 183 (2021) 107476 6 Due to the number of analyses, only statistically significant interactions between the factors will be reported. The ANOVAs showed significant interactions between WM Load and Frontal Theta Phase for right inferior frontal gyrus theta phase and gamma 1 amplitude at bilateral posterior cingulate cortex (Table 2). Additionally, significant interactions were found between WM Load and Frontal Theta Phase for right middle frontal gyrus theta phase and beta 1 amplitude at bilateral cuneus; and beta 2 amplitude at right postcentral gyrus and right middle occipital gyrus (Table 2). Finally, a significant interaction between WM load and frontal theta phase for bilateral anterior cingulate theta phase and beta 2 amplitude values at right middle occipital gyrus was found (Table 2). Post-hoc analysis showed the following results: 3.3.1. Theta phase at Right Inferior Frontal gyrus Gamma 1 (30 – 50 Hz) Gamma 1 amplitude at bilateral posterior cingulate was significantly greater at theta pre-peak (bin 14) in the ML condition compared to the HL condition (Fig. 4A). 3.3.2. Theta phase at Right Middle Frontal gyrus Beta 1 (13 – 20 Hz) Beta 1 amplitude at bilateral cuneus was significantly greater at theta post-peak (bin 4) in the LL condition than the ML condition. (Fig. 4C). Beta 2 (20 – 30 Hz) Beta 2 amplitude at the right postcentral gyrus was significantly Fig. 2. Brain areas showing statistically significant extreme t-statistics in the current density power analyses. In each image, it is plotted the brain sites showing the extreme significant t-sta- tistic in the SnPM contrast showing statistically significant between WM load condition differences for the studied frequency bands. Colour bars are based on t-values (see Table 1), with warm colours indicating increases of current source density power with memory load and cold colours decreases. Note that comparisons within the Gamma 1 frequency band (30–50 Hz) are not included since no statistically significant results were observed between WM load conditions. Fig. 3. Equivalence between bins and phase angles in degrees. A. Fern´ andez et al. Neurobiology of Learning and Memory 183 (2021) 107476 7 greater in the HL condition compared to the LL condition at theta trough (bin 7) (Fig. 4B). Additionally, beta 2 amplitude at the right middle occipital gyrus was significantly higher at theta peak (bin 15) in the LL condition than the ML condition (Fig. 4D). 3.3.3. Theta phase at Anterior Cingulate Beta 2 (20 – 30 Hz) Beta 2 amplitude at the right middle occipital gyrus was significantly greater at theta post-peak (bin 3) for the ML condition compared to the HL condition (Fig. 4E). 4. Discussion In this study, the effect of WM load on the spectral activity and location of theta and fast (i.e., beta and gamma) frequencies was assessed during the maintenance period of a DMS task with three different levels of WM load to select the regions of interest for an interregional cross-frequency phase-amplitude coupling analysis. Thus, it was, for the first time, explored whether WM load modulates the distribution of fast frequencies amplitude across theta cycle angles. 4.1. Behavioural results Participants’ performance was modulated by WM load regarding both, percentage of hits and reaction time, with lower correct responses and slower reaction time with higher WM loads. These results were expected as it has been previously observed in studies with similar tasks that increases in WM load worsen participants’ performance (Alvarez & Cavanagh, 2004; Vogel et al., 2001). 4.2. Current density power results Besides the behavioural effects of an increasing WM load, the results showed that during the maintenance period of the DMS task, theta and beta/gamma activity increased as a function of WM load, and the largest modulation in WM load was located at different sources for fast and slow frequencies. Besides, we found that the beta/gamma amplitude synchronization with theta cycle phase angles appears to be modulated by WM load. Firstly, regarding the current density power of theta and fast frequencies, it was found that theta activity (4 – 7 Hz) had the largest increases with WM load in prefrontal areas of the right hemisphere (anterior cingulate cortex, right inferior frontal, and right middle frontal gyri). These findings are in line with those from previous studies that found frontal theta activity increases with WM load (Gevins, 1997; Jensen & Tesche, 2002; Meltzer, Zaveri, Goncharova, Distasio, Papademetris, Spencer, Spencer, & Constable, 2008). Furthermore, greater prefrontal theta activity has been found in association with an increase of executive control demands (Griesmayr, Gruber, Klimesch, & Sauseng, 2010; Cavanagh, Eisenberg, Guitart-Masip, Huys, & Frank, 2013). Taking these findings together, it has been proposed that theta activity has a key role in WM maintenance and executive control (Maurer et al., 2015; Sauseng et al., 2010). Regarding fast frequencies, the present results showed that fast beta and gamma were modulated by WM load in posterior brain regions. Specifically, beta activity (13 – 30 Hz) had the maximum increase with WM load in occipital areas, and fast gamma (50 – 70 Hz) activity had the maximum increase with WM load in right parietal regions. Beta activity has been found to increase in occipital and occipitotemporal regions with WM load during the delay period of DMS tasks when memorizing different features of geometrical shapes, such as color, shape or location (Palva et al., 2011; Honkanen et al., 2015). This lends support to the idea that beta activity is related to the maintenance of object representations in WM (Palva et al., 2011). Gamma activity has also be found to increase with WM load in the parietal cortex during Sternberg (Howard et al., 2003) and DMS tasks (Palva et al., 2011). In this regard, gamma involvement in WM has been frequently associated with the maintenance of visual representations in WM (Jokisch & Jensen, 2007) and WM capacity (Tallon-Baudry et al., 1998). Furthermore, findings of fast gamma activity in the parietal cortex during WM maintenance link this activity with increasing attentional demands (Fries, 2009; Palva et al., 2011). These results, however, contradict the studies that have observed decreases in beta and gamma activity on spatial WM tasks when the focus is on memorizing the stimulus location (Poch et al., 2014; Proskovec et al., 2018; 2019). Interestingly, when analyzing the activity related to the memorization of the different stimulus features, Honkanen et al. (2015) found that fast activity increased with WM load when participants were memorizing shape or colour but failed to find a WM load effect when participants memorized location alone. Increases in fast activity when memorizing distinct symbols such as letters (Howard et al., 2003; Michels et al., 2010) or irregular shapes (Honkanen et al., 2015; Morgan et al., 2011; Tallon-Baudry et al., 1998) seem to be common. Even though the participants in the current study had to detect whether a dot in the domino tile had changed its location, it is possible that, at least for higher WM loads, they relied on the overall shape pattern formed by the dots in the three tiles, rather than the individual location of each dot. This fact would explain why posterior increases in beta and gamma activity with higher WM loads were found, since they would be processing irregular shapes rather than the dots’ location per se. Additionally, it was found that in two specific brain regions, this modulation of fast frequency oscillatory activity by WM load was reversed: fast beta and fast gamma decreased with higher WM load, with the largest decreases in the left fusiform gyrus and in the right precuneus respectively. Although increasing gamma activity has been previously observed in the left precuneus during the manipulation of features of two stimuli in WM (Morgan et al., 2011), the present task only required participants to maintain without mentally manipulating any feature of the stimuli (i.e., stimulus’ configuration). It could be the case that gamma activity in the precuneus is specifically associated with more complex WM manipulations than those demanded in the present study. Thus, in the HL condition, it could be preferred to allocate processing resources to the superior parietal lobe and right postcentral gyrus, where increases of gamma activity with WM load were found, while activity in the right precuneus might be no longer required, allowing a reduction of gamma Table 2 Significant interactions between WM Load Factor (HL, ML, and LL conditions) and Frontal Theta Phase Factor (16 phase bins) on beta and gamma normalized current density amplitude values. F (with degrees of freedom in parenthesis) and p-values are indicated alongside the frontal and posterior regions for theta and beta/gamma frequencies, respectively. Frontal region Frontal region frequency band Posterior region Posterior region frequency band F (30/720) p Right inferior frontal gyrus Theta (4–7 Hz) Posterior cingulate (bilateral) Gamma 1 (30–50 Hz) 1.859 0.004 Right middle frontal gyrus Theta (4–7 Hz) Cuneus (bilateral) Beta 1 (13–20 Hz) 1.890 0.003 Right postcentral gyrus Beta 2 (20–30 Hz) 1.797 0.006 Right middle occipital gyrus Beta 2 (20–30 Hz) 1.816 0.005 Anterior cingulate (bilateral) Theta (4–7 Hz) Right middle occipital gyrus Beta 2 (20–30 Hz) 2.235 <0.001 A. Fern´ andez et al. Neurobiology of Learning and Memory 183 (2021) 107476 8 activity to be observed. Summarizing the present study findings, increasing theta activity with WM load in prefrontal areas supports the involvement of theta activity in WM maintenance and executive control. Regarding fast frequencies, higher beta activity with increasing WM load in areas related to perceptual visual processing, combined with gamma activity increases in attention–related regions may indicate that, for this task, beta activity was enough for visual representations, while parietal fast gamma activity reflects the attentional demands of the task. This increase of fast activity in higher WM load conditions is in line with that described in previous studies when processing and memorizing global shapes, even though some spatial information was required to solve the present task. These results highlight the need to take into account the methodological aspects of the task, as well as the memorizing strategies used by the participants in future studies of WM load effects on brain activity. Fig. 4. Pairwise comparisons for fast frequencies normalized amplitude in the three memory load conditions across the theta cycle angles at the right inferior frontal gyrus (rIFG), right middle frontal gyrus (rMFG), and anterior cingulate cortex (ACC). For each fast frequency and region, it is indicated in the table on the right the task conditions showing significant differences, the bin of theta frequency band where it is localized, and the associated p-values. In the line graph, the fast frequencies’ mean normalized amplitude values are shown (Standard Error -SE- with bars) for each region. The frontal theta cycle bins showing significant differences are grey shaded. A) theta phase at rIFG and gamma 1 normalized amplitudes at bilateral posterior cingulate. B) theta at rMFG and beta 2 amplitudes at right postcentral gyrus; C) theta at rMFG and beta 1 amplitudes at bilateral cuneus; D) theta at rMFG and beta 2 amplitudes at right middle occipital gyrus; E) theta at ACC and beta 2 amplitudes at right middle occipital gyrus. A. Fern´ andez et al. Neurobiology of Learning and Memory 183 (2021) 107476 9 4.3. Distribution of fast frequencies amplitude across theta cycle angles results Lastly, regarding the distribution of posterior fast frequencies amplitude across frontal theta phase angles, it was found that the distribution of posterior beta and gamma amplitudes in the course of the phase cycle of frontal theta is influenced by WM load, providing novel evidence about WM load effects on interregional cross-frequency phase: amplitude mechanisms. Specifically, greater amplitudes of beta activity at the bilateral cuneus and right middle occipital gyrus as well as higher gamma activity at the bilateral posterior cingulate were observed for lower memory load conditions when compared with higher memory load conditions during or near frontal theta peak. Conversely, greater fast beta amplitude at the right postcentral gyrus was observed during theta trough at rMFG in the high memory load condition compared to the low memory load condition. These results evidence that WM load modulates communication between cortical areas. Following the model by Fries (Fries, 2005), it is expected that optimal neural communication is reflected in a concentration of fast frequencies activity around specific angles of theta phase. Therefore, for demanding situations, the preferred phase angle will be theta trough and, therefore, ideally: i) greater amplitudes in the higher load conditions (compared to lower loads) should be observed during theta trough; ii) as well as greater amplitudes in the lower load conditions (compared to higher loads) during theta peak. Supporting the first prediction based on the model, fast beta amplitude at the postcentral gyrus was greater during theta trough at the rMFG for high than low WM load. According to Fries (2005) model, this coupling occurs at the considered optimal window for communication. The postcentral gyrus is a somatosensory region. However, it has also been found to play a role in attention (Shulman, Astafiev, McAvoy, d’Avossa, & Corbetta, 2007), as well as to be part of a posterior attentional network (Tomasi, Ernst, Caparelli, & Chang, 2006). Further, its lesion has been shown to impair attentional skills (Bajaj, Dailey, Rosso, Rauch, & Killgore, 2018). Given that the rMFG is thought to be related to the maintenance of information in WM, especially in non-verbal tasks (Daniel, Katz, & Robinson, 2016; Linden, 2007), it could be argued that for higher memory load, rMFG and postcentral gyrus established communication during the model’s optimal window for communication as a reflection of the increased demands on attentional control over WM contents. Regarding the second prediction, in the current study, the posterior cingulate has been found to communicate during a non-optimal window for communication (i.e., frontal theta peak) with the rIFG, an area known to be implicated in executive control and inhibition in WM (Kasahara et al., 2011). The posterior cingulate cortex is considered part of the default mode network (Andrews-Hanna, Reidler, Sepulcre, Poulin, & Buckner, 2010) and has also been related to spatial processing and spatial navigation (Rolls, 2018; 2019). Specifically, it has been observed to play a role in the visual processing dorsal stream, related to spatial information (Kravitz, Saleem, Baker, & Mishkin, 2011). Since the task used in the present study requires the maintenance of the stimuli’ global configuration, it may be the case that there was a minimal demand for communication to allocate cognitive control resources on spatial processing-related areas during low load trials. Likewise, the right middle occipital gyrus is a visual processing area. Still, it has also been found to be related to visual memorization strategies. For instance, it has been observed that simulating a lesion on this area with transcranial magnetic stimulation impairs digit span task performance in participants who use visualization strategies (Hilbert et al., 2019). This area appears to communicate during non-optimal windows with both the rMFG and the ACC, which are related to nonverbal maintenance of information and to cognitive control, respectively (Daniel et al., 2016; MacDonald, Cohen, Stenger, & Carter, 2000). This non–optimal communication has been found for lower WM loads in both cases, following the same pattern as indicated above, where lower WM loads demand less cognitive resources, and thus there is a lesser need for frontal control of occipital areas. It is noticeable that most areas that yielded a significant result are lateralized to the right hemisphere. It is believed that coordinate spatial relation processing is lateralized to the right hemisphere and that this pattern is also generalized to visual working memory (van der Ham, van Wezel, Oleksiak, & Postma, 2007; van der Ham, Postma, and Laeng, 2014). In conclusion, the present study results offered novel evidence showing that, during the maintenance period of a DMS task, WM load modulates interregional communication between frontal and posterior areas. This shows that interregional communication depends on the cognitive demands required by the task for a successful performance, with WM load acting as a switch that enables or hinders the optimal communication between frontal and posterior areas depending on the need for frontal control in visual and attention processing regions. CRediT authorship contribution statement Alba Fern´ andez: performed the data analysis and interpretation under the supervision of Diego Pinal and wrote the first draft of the manuscript. Diego Pinal: designed the study and the DMS task alongside Fernando Díaz and Montserrat Zurr´ on, supervised the data analysis and interpretation, and provided critical revisions to the manuscript. Fernando Díaz: designed the study and the DMS task alongside Diego Pinal and Montserrat Zurr´ on, and provided critical revisions to the manuscript. Montserrat Zurr´ on: designed the study and the DMS task alongside Fernando Díaz and Diego Pinal, and provided critical revisions to the manuscript. Acknowledgements This study was supported by grants from the Spanish Government, Ministerio de Economía y Competitividad (PSI2014-55316-C3-3-R; PSI2017-89389-C2-2-R), with FEDER Funds; the Galician Government, Consellería de Cultura, Educaci´ on e Ordenaci´ on Universitaria, Axudas para a Consolidaci´ on e Estruturaci´ on de Unidades de Investigaci´ on Competitivas do Sistema Universitario de Galicia: GRC (GI-1807-USC); Ref: ED431-2017/27, with FEDER funds. D. P. was also supported by the Portuguese Foundation for Science and Technology (FCT) through the fellowship SFRH/BPD/120111/2016. Declaration of Competing Interest Authors declare no conflict of interest. Appendix A. Supplementary material Supplementary data to this article can be found online at https://doi. org/10.1016/j.nlm.2021.107476. References Alvarez, G. A., & Cavanagh, P. (2004). The capacity of visual short-term memory is set both by visual information load and by number of objects. Psychological Science, 15 (2), 106–111. https://doi.org/10.1111/j.0963-7214.2004.01502006.x. Andrews-Hanna, J. R., Reidler, J. S., Sepulcre, J., Poulin, R., & Buckner, R. L. (2010). Functional-anatomic fractionation of the brain’s default network. Neuron, 65(4), 550–562. https://doi.org/10.1016/j.neuron.2010.02.005. Axmacher, N., Henseler, M. M., Jensen, O., Weinreich, I., Elger, C. E., & Fell, J. (2010). Cross-frequency coupling supports multi-item working memory in the human hippocampus. Proceedings of the National Academy of Sciences, 107(7), 3228–3233. https://doi.org/10.1073/pnas.0911531107. Baddeley, A. (1998). Recent developments in working memory. Current Opinion in Neurobiology, 8(2), 234–238. Baddeley, Alan (2003). Working memory: Looking back and looking forward. Nature Reviews Neuroscience, 4(10), 829–839. https://doi.org/10.1038/nrn1201. A. Fern´ andez et al.