## 4 Discussion  
We examined the reproducibility of theta-gamma coupling between repeated measurements using a set of EEG data for healthy subjects performing the 3-back working memory task as well as a publicly available benchmark LFP recording with previously reported strong PAC. The analysis was concerned with the overall structure of the comodulograms and the frequencies of maximal TGC, in view of the possibility of using these frequencies in customizing noninvasive neuromodulation techniques to optimally impact the neural dynamics of every individual subject. For studying repeatability, the results of processing the two halves of each recording were treated as repeated measurements. Thus, the findings pertain to TGC reproducibility at short time intervals comparable to the session length (a few minutes). Several PAC calculation pipelines (Tables 1 and 2) were compared in terms of their influence on the obtained coupling patterns, peak frequencies, and their temporal stability. The algorithms were constructed on the basis of an analysis of the approaches to TGC estimation used in the literature and the available options at each step of the computation. A number of these options have previously been compared on synthetic data, while there is insufficient information about their influence on PAC estimates for actual neural signals. The benchmark LFP data from the rat hippocampus showed a large theta peak in the power spectrum that was repeatable between the half-sessions (Figure 2). The comodulograms were largely similar between the considered algorithms, although there were noticeable peak differences in certain methods (Figure 3). A considerable degree of reproducibility between the half-sessions was found for both the general shape of the comodulograms (measured by the Spearman correlations) and the positions of PAC maxima in them (quantified by the peak frequency changes, Figures 4 and 5, with absolute values 0-2.5 Hz for phase and 0-5 Hz for amplitude). Each considered algorithm of PAC computation produced one bright peak, all of them having roughly similar locations but varying shapes, with relatively broad and 'patchy' peaks for the two methods employing the normalization of the PAC index by surrogate data. It may be possible to smooth these comodulograms using a larger number of surrogates, but, in addition to increased computation time, there is no guarantee that the results will be more repeatable or accurate. In the EEG data, the presence of theta peaks in electrode Fz was not universal, and many of them were only minimally noticeable against the background (Figure 6). This is largely consistent with previous reports providing variable evidence for theta oscillations during human episodic memory, with some studies reporting memory-related theta increases and others reporting decreases [2]. Thus, a test-retest study [63] of spectra during a modified Sternberg task showed theta peaks in the AFz signal of some but not all subjects. This raises concerns about the range of applicability of TGC analysis because the existence of a clear peak at the frequency for phase is recognized as a prerequisite for a meaningful interpretation of any CFC pattern [25]. The TGC comodulograms generally showed a considerable variation both between the algorithms and between half-sessions (Figures 7 and 8). The Spearman correlations of the comodulograms between half-sessions did not significantly differ on average from zero for the two methods that include normalization, while a significant positive average comodulogram similarity was found for the other methods (Figure 9). However, to a large degree, this was not specific to the particular individual. Indeed, the similarity of comodulograms from different subjects was comparable to that of comodulograms from the same participant (Figure 10). The fact that the methods using unnormalized KL-MI, unlike the algorithms with normalization, produced comodulograms with significant 'unspecific' similarity (both within and between individuals) presents interest for further research.