Auditory and tactile frequency representations are co-embedded in modality-defined cortical sensory systems

Sensory information is represented and elaborated in hierarchical cortical systems that are thought to be dedicated to individual sensory modalities. This traditional view of sensory cortex organization has been challenged by recent evidence of multimodal responses in primary and association sensory areas. Although it is indisputable that sensory areas respond to multiple modalities, it remains unclear whether these multimodal responses reflect selective information processing for particular stimulus features. Here, we used fMRI adaptation to identify brain regions that are sensitive to the temporal frequency information contained in auditory, tactile, and audiotactile stimulus sequences. A number of brain regions distributed over the parietal and temporal lobes exhibited frequency-selective temporal response modulation for both auditory and tactile stimulus events, as indexed by repetition suppression effects. A smaller set of regions responded to crossmodal adaptation sequences in a frequency-dependent manner. Despite an extensive overlap of multimodal frequency-selective responses across the parietal and temporal lobes, representational similarity analysis revealed a cortical "regional landscape" that clearly reflected distinct somatosensory and auditory processing systems that converged on modality-invariant areas. These structured relationships between brain regions were also evident in spontaneous signal fluctuation patterns measured at rest. Our results reveal that multimodal processing in human cortex can be feature-specific and that multimodal frequency representations are embedded in the intrinsically hierarchical organization of cortical sensory systems.


Introduction
How the brain is organized to support perception of unisensory and multisensory signals is a fundamental question. Traditional organization models posit an initial processing of information in unisensory cortical systems before subsequent processing in multimodal systems residing in association cortices. Yet, it is now incontrovertible that brain regions traditionally considered to be dedicated to a single modality can be modulated or driven by sensory inputs delivered via multiple sensory modalities (Driver and Noesselt, 2008;Ghazanfar and Schroeder, 2006). These discrepancies motivate a debate between the traditional "modality-based scheme" for sensory cortex organization and a "function-based scheme" in which sensory regions perform particular operations regardless of input modality (Driver and Noesselt, 2008;Pascual-Leone and Hamilton, 2001). A function-based scheme is further supported by evidence from crossmodal neuroplastic changes in the brain associated with sensory deprivation (Heimler et al., 2015;Merabet and Pascual-Leone, 2010). Here we characterized BOLD fMRI responses to auditory and tactile stimulus sequences to investigate how sensory cortex is organized to support multimodal temporal frequency processing.
In humans, the relative sensitivities of the auditory (0.02-20 kHz) and somatosensory (2-1000 Hz) systems to temporal frequency information allow for the possibility that the two senses overlap in the processing of signals below 1000 Hz. Environmental cues in this frequency range may be particularly relevant for signaling our interactions with objects (Miller et al., 2018), fine texture information (Manfredi et al., 2014;Yau et al., 2009a), and speech (Lattner et al., 2005;Titze, 1994). Indeed, with both sensory modalities, the perception of sensory cues comprising repeating events can be organized into two domains that depend on the repetition rate of the events. At low repetition ratesthe flutter rangesignals are perceived as a stream of temporally discrete events. At higher repetition rates (>50 Hz), stimuli are perceived as a single continuous signal. There are clear correspondences between auditory and tactile frequency perception (Bensmaia et al., 2005) as well as extensive evidence for frequency-specific audiotactile interactions with flutter (Convento et al., 2019) and supra-flutter signals (Crommett et al., 2019Yau et al., 2010Yau et al., , 2009b. The intimate relationship between auditory and tactile frequency perception presumably reflects the existence of shared or interactive neural representations. The neural codes for low-frequency auditory (Bendor and Wang, 2007) and tactile stimulation (Romo and Salinas, 2003) in the flutter range (i.e., how neural spiking activity relates to the repetition rates of the sensory cues) are well established and neural populations that analogously signal flutter in both modalities have been identified (Vergara et al., 2016). In contrast, while neural populations are clearly tuned for supra-flutter auditory frequency (Bendor and Wang, 2005;Wang and Walker, 2012), an explicit coding of vibration frequencieswhich can be represented in spike timing (Harvey et al., 2013) has not been established. Accordingly, the relationship between auditory and somatosensory supra-flutter representations remains unclear (Saal et al., 2016). Consistent with a feature-based organization scheme in which auditory cortex is specialized for frequency processing irrespective of input modality, some auditory areas respond to tactile stimulation alone (Foxe et al., 2002;Fu et al., 2003;Kayser et al., 2005;Schurmann et al., 2006). Even with the demonstration of tactile responses in the auditory cortical system and the definition of extensive cortico-cortical and thalamo-cortical pathways connecting the somatosensory and auditory systems (Cappe et al., 2009b;Cappe and Barone, 2005;Hackett et al., 2007), questions remain regarding whether somatosensory inputs exclusively drive higher-order auditory regions and the nature of somatosensory influences on primary auditory cortex . Moreover, other studies (Liang et al., 2013;Perez-Bellido et al., 2017) have revealed that auditory information can also be decoded from putative somatosensory areas. Thus, questions remain regarding the specificity and distribution of auditory and tactile frequency representations in the human brain and the organization of the cortical systems supporting these representations.
We used a fMRI adaptation paradigm to probe auditory and tactile frequency responses in the human brain. Adaptation describes the tendency for repeated sensory stimulation to alter the responses of neurons that are tuned for features of the sensory stimuli (Solomon and Kohn, 2014). fMRI adaptation paradigms leverage this process to infer population-level (i.e., voxel-wise) tuning properties based on BOLD signal changes (Grill-Spector et al., 2006). We scanned healthy human participants performing a frequency monitoring task on brief auditory, tactile, and audiotactile stimulus sequences that contained repeating or changing frequencies. In univariate analyses, we quantified BOLD signal changes related to unimodal (auditory-only or tactile-only) and crossmodal adaptation effects to identify brain regions exhibiting activation patterns consistent with population-level frequency tuning. According to a strict modality-based organizing scheme, we predicted that unimodal auditory and tactile repetition effects would be largely confined to distinct brain regions. Alternatively, according to a strict function-based organizing scheme, we predicted that unimodal repetition effects would occur in overlapping regions which should also exhibit crossmodal adaptation effects. We additionally performed multivariate pattern analysis to establish the representational spaces occupied by the unimodal and crossmodal events in different regions spanning the parietal and temporal lobes. We used representational similarity analysis (RSA) to determine the relationships between the neural subspaces of different brain regions. We predicted that the representational spaces occupied by the unimodal and crossmodal events would be similar if auditory and tactile frequency information were represented by a common neural substrate. Lastly, we compared the cortical landscape defined by the event-related spatial activation patterns to the connectivity profile computed from temporal correlations between spontaneous signal fluctuations measured in the same participants in separate resting state scans. This comparison allowed us to test whether the cortical landscape defined by regional frequency-and modality-dependent responses emerged only with sensory stimulation or if the profile was reflected in intrinsic activity patterns. Collectively, our results indicate that frequency-selective auditory and tactile processing is broadly distributed over sensory systems which are organized according to the traditional modality-based scheme.

Participants
Twenty-four healthy adult participants were recruited from Baylor College of Medicine and the Houston metropolitan area. Participants indicated via self-report the absence of any current or past psychiatric or neurological conditions and were not taking any centrally acting medications. Participants reported normal tactile and auditory sensibilities. Testing procedures were performed in compliance with the policies and procedures of the Baylor College of Medicine Institutional Review Board. All participants provided informed written consent and were paid for their participation. A total for 4 participants were excluded from all analyses: 2 participants were excluded for below chance performance on the frequency monitoring task during the main adaptation scans (see below) and 2 participants were excluded for above chance performance on a control detection task (see below). All analyses thus included data from 20 participants (10m10f; mean age AE standard deviation: 23.2 AE 4.4 years). All participants were right handed according to the Edinburgh Handedness Inventory (mean score AE standard deviation: 79 AE 23.8) (Oldfield, 1971).

Overview: neuroimaging
Each participant underwent a functional localizer (Mapping) scan, 8 functional (Adaptation) scans involving continuous performance of a frequency monitoring task, 1 structural scan, and 2 resting state scans. Each participant also completed a control behavioral task in the scanner (Detection control).

Image acquisition
All scans were conducted in the Core for Advanced MRI (CAMRI) at Baylor College of Medicine. MRI data were acquired on a Siemens MAGNETOM Trio 3 T system using a Siemens 32-channel head coil. Structural images were acquired with a sagittal magnetization prepared rapid gradient echo (MP-RAGE) T1-weighted sequence (time echo (TE) ¼ 3.02 ms; time repetition (TR) ¼ 2600 ms; time to inversion ¼ 900 ms; flip angle ¼ 8 ; GRAPPA factor ¼ 2; 176 slices with 1 Â 1 Â 1 mm voxels). Functional images for task runs were obtained using an axial echo-planar imaging (EPI) sequence (TE ¼ 25 ms; TR ¼ 2750 ms; flip angle ¼ 80 ; GRAPPA factor ¼ 4; 56 slices; 2 mm 3 voxels; 0 mm gap). A total of 151 vol were acquired for the mapping scan, 121 vol were acquired for each of the 8 adaptation scans, and 121 vol were acquired for each of the 2 resting state scans.

Auditory and tactile stimulation
Auditory stimuli were generated in Matlab (v2011b; Mathworks, Inc.) running on a Macbook and sounds were delivered via MRI-compatible inear headphones (Model S14, Sensimetrics) after amplification (PCA1, Pyle). Vibrotactile stimulation was delivered simultaneously to the distal finger pads on digits 2-5 on the right hand using a piezoelectric stimulator (CM3, Cortical Metrics). Humans can perceive and compare the frequency of sinusoidal stimulation in within-and across-modality discrimination paradigms (Convento et al., 2018). All sounds and vibrations were matched for perceived intensity: In preliminary experiments with an independent rater group, intensity matching for the 100vs 300-Hz stimuli was performed for the auditory and tactile modalities separately. Subsequently, at the outset of each participant's session, the gain on the audio amplifier was manually adjusted to match the subjective intensity of the 100-Hz sound to that of the 100-Hz vibration and this setting was used with the 100-and 300-Hz sounds. To account for individual differences in perceptual sensitivity and to ensure further that stimulus amplitude did not systematically vary with stimulus frequency, we introduced a random jitter (AE5%) to the intensity-matched amplitude levels on every stimulus presentation.

Modality-mapping localizer scan
Each subject underwent a functional scan designed to identify brain regions responsive to tactile or auditory stimulation (Mapping scan). These data 1) provided an opportunity to replicate previous reports of overlapping auditory and tactile activations and 2) defined masks used to constrain the analyses of the frequency adaptation scans (see below). Using a block design paradigm, we presented unimodal blocks (16.5s) of tactile stimulation (25-400 Hz), auditory pure-tone stimulation (40-700 Hz), or auditory band-pass noise stimulation (center frequencies: 40-700 Hz; bandwidths: 0.45*center frequency). Each unimodal block consisted of 24 randomly ordered stimuli (stimulus duration: 437.5 ms; interstimulus interval (ISI): 250 ms). Each block type was repeated 5 times in pseudorandom order. Blocks were separated by 11s of no stimulation. Participants passively experienced the stimulation while maintaining fixation on a visual display.

Frequency adaptation scans
Each subject underwent 8 functional scans testing frequency adaptation. Each scan comprised 48 adaptation events. Each event (Fig. 1A) consisted of 3 adapting stimuli followed by 1 probe stimulus. Stimuli (500-ms duration; 250-ms ISI) were either 100 Hz or 300 Hz. We tested with these frequencies because of hardware constraints and because earlier work indicated that these frequencies should be reliably discriminable (Convento et al., 2018;Yau et al., 2010Yau et al., , 2009b. We tested 16 adaptation event types in a 2 Â 2 Â 2 Â2 factorial design (Fig. 1B) that varied the modality of the adapting stimuli (auditory or tactile), the modality of the probe stimulus (auditory or tactile), the frequency of the adaptors, and the frequency relationship between the adaptors and the probe (repeat or change). Thus, we tested within-modality adaptation sequences (auditory adaptorsauditory probes, AA; tactile adaptorstactile probes, TT) and across-modality adaptation sequences (auditory adaptorstactile probes, AT; tactile adaptorsauditory probes, TA). Each event type was presented 6 times per scan (yielding a total of 48 repeats per event type over all scans). Half of the REPEAT trials contained an adapting frequency of 100 Hz and the adapting frequency in the other half were 300 Hz. CHANGE trials also contained an equal number of 100and 300-Hz adaptors. This design ensured that the frequency and modality of the adapting stimuli did not predict the frequency or modality of the probe stimulus on any given event.
Subjects fixated a central visual crosshair during each trial. To maintain participants' attention on stimulus frequency throughout the scans, participants were asked to compare the frequencies of the adapting and probe stimuli on each trial and judge whether they were the same or different. However, participants only reported their judgments by button press when prompted on occasional catch events (3 per scan) with their left hand using an MRI-compatible button box (Current Designs). This design enabled us to monitor compliance on the frequency monitoring task without contaminating our estimates of sensory processing activity with motor response activity. Two participants failed to achieve abovechance performance over all of the conditions, so we excluded their data from additional analyses. The ordering and timing of events in each scan was optimized and pseudorandomized using Optseq2 (http://surfe r.nmr.mgh.harvard.edu/optseq).

Tactile detection control task
Because we were interested in comparing responses to auditory and tactile stimulation, it was imperative to verify that participants could not judge tactile stimulus frequency based on acoustic cues produced by the piezoelectric stimulator in the scanner. Accordingly, we tested each participant on a 2-interval forced-choice vibration detection task while the stimulator and scanner were operating but with the stimulator positioned next to the body without contacting the hand or any other body part. On each trial, a 500-ms tactile stimulus (25, 100, or 300 Hz) was generated by the piezoelectric stimulator in one of two intervals (ITI ¼ 1000 ms). Participants were required to report the interval in which they detected the stimulus. Trial and response intervals were cued visually, and participants were instructed to guess if uncertain. We reasoned that participants would perform at chance level if they were unable to perceive the tactile stimulus without physical contact. Conversely, we reasoned that participants would perform above chance level if they exploited alternative cues (e.g., stimulator-generated sounds) to perform the tactile detection task. To replicate the acoustic conditions of the adaptation and mapping runs during the tactile detection task, we operated the scanner using the same EPI sequence. Two participants achieved above-chance performance on the detection task, Fig. 1. fMRI adaptation paradigm and behavioral results. A, Event-related design comprising within-modality and across-modality adaptation events. Each 2.75-s event is a stimulus sequence comprising 3 adaptor stimuli that are matched in frequency and a probe stimulus whose frequency either repeats or changes with respect to the adapting frequency. The modality of the adaptors and probes can be matched (auditory-auditory, AA; tactile-tactile, TT) or different (auditory-tactile, AT; tactile-auditory, TA). Participants are tasked with monitoring the frequency of each stimulus and covertly determining if the frequency of the probe matches that of the adaptors. Participants report their decisions when prompted visually during catch events (CE). B, Adapting and probe frequencies could be either 100 Hz or 300 Hz. The experiment comprised a full parametric design crossing adapting frequency, probe frequency relation to adapting frequency (REPEAT or CHANGE), and stimulus modality. C, Average task performance for each event type. Circles indicate single participant performance levels (n ¼ 20). Error bars indicates S.E.M. so we excluded their data from additional analyses.

Data analysis: overview
We first defined brain regions that exhibited significant responses to auditory or tactile stimulation in analysis of the functional localizer scans. Regions identified from the localizer data then constrained the data space in the analyses of the frequency adaptation scans, which included univariate and multivariate tests. All analyses were performed in surface space. As a first step in the analysis of the frequency adaptation scans, we conducted group-level statistical tests of parameters estimated in general linear models (GLM) fitted to the time series of each surface node. By leveraging within-modality repetition effects, these univariate analyses were aimed at 1) identifying brain regions whose responses were consistent with frequency-selective processing for sounds or vibrations and 2) determining the spatial overlap in frequency-selective sound and vibration processing. By leveraging across-modality repetition effects, these univariate analyses were aimed at establishing evidence for interactive auditory and tactile frequency processing which would support the notion of shared frequency representations. We supplemented this node-based analysis with region-of-interest (ROI) analyses that specifically focused on putative auditory and somatosensory regions in the temporal and parietal lobes defined using a combination of functional masks (i.e., the localizer results) and anatomical masks (see below). We conducted separate ROI-based analyses to evaluate adaptation effects, to compare representational geometries based on multivariate response patterns, and to relate brain regions based on their patterns of spontaneous BOLD signal fluctuations. Each analysis is described in detail in the following sections.

fMRI data pre-processing
We used AFNI software (Cox, 1996) to perform data pre-processing and univariate analyses. Three-dimensional surface models were created with FreeSurfer (Fischl et al., 1999) (freesurfer-Darwin-lion-stable-pub-v5.3.0) and visualized in SUMA (Saad and Reynolds, 2012). Functional datasets were corrected for slice timing, motion corrected, and despiked. Volumes with head motion exceeding 0.3 mm/TR or a fraction of 0.05 outlier voxels (calculated using the 3dToutcount function) were excluded from analysis.
Volume data were projected into surface space using the 3dVol2Surf function in AFNI with the node-wise signal computed as the mean signal averaged over 10 evenly-spaced points between the smooth white matter and pial matter. Data from the localizer scan and the frequency adaptation scans that were included in whole-brain univariate analyses were normalized to the N27 atlas brain (Mazziotta et al., 2001) and spatially smoothed (4-mm FWHM 2D Gaussian kernel) in standard surface space. All data were expressed in percent signal change with respect to the mean signal in each scan.

fMRI data analyses: Localizer scans
Localizer scans were modeled using GLMs that included 3 regressors of interest corresponding to tactile stimulation, auditory pure-tone stimulation, and auditory band-pass noise stimulation convolved with gamma-variate functions (Boynton et al., 1996). Head motion parameters and drift parameters (linear, quadratic, and cubic) were included as nuisance regressors. In group-level analyses, we identified the surface nodes that exhibited positive BOLD signal changes in each block type. Because we only used these patterns to define an analysis mask, we used a liberal uncorrected threshold (p < 0.05). Activation maps for auditory and tactile stimulation blocks were combined to generate a single analysis mask.

fMRI data analyses: Adaptation scansgeneral
We fitted node-wise multiple linear regression models in first-level analyses. The simplest GLM included 8 stimulus regressors (4 adaptorprobe modality conditions: AA, TT, AT, TA; 2 adaptation conditions: CHANGE, REPEAT) and catch events (CE) as separate regressors. This model (9 conditions of interest) assumed that there were no response differences related to adapting frequency. We also fitted a more complex GLM which allowed for effects of adapting frequency; this model included 16 adaptation conditions (4 modality conditions; 2 adaptation conditions; 2 adapting frequencies: 100 Hz and 300 Hz) along with CE for a total of 17 regressors of interest. Each regression model also included motion and drift parameters as nuisance regressors. Gamma-variate functions were used to model stimulus responses in the main analyses. Beta coefficients from the first-level analyses were submitted to grouplevel univariate analyses designed to identify nodes that exhibited significant response modulation related to the adaptation conditions.
All group-level univariate test results were constrained by the mask defined by the localizer results and statistically thresholded at a false discovery rate (FDR) corrected q < 0.05. We used the 'SurfClust' AFNI function to summarize the clustered activation patterns identified as significant after FDR correction (comprising 2 or more nodes). Data from all clusters comprising 2 or more surface nodes were included in subsequent analyses. In a separate validation procedure, we additionally analyzed the data using finite impulse response deconvolution in order to visualize response time courses for each stimulus condition (8 time points; 0-19.25s post-event onset).

Group-level univariate analysis: Within-modality adaptation events
To identify brain regions whose responses were consistent with frequency-selective processing for sounds or vibrations, we performed a single analysis on the within-modality adaptation events (AA and TT) by conducting a two-way repeated measures ANOVA (rmANOVA) with modality (auditory or tactile) and adaptation condition (CHANGE or REPEAT) as the factors. This analysis ignores potential modulating effects of adapting frequency. By evaluating the main effects of modality and adaptation condition as well as the modality Â adaptation interaction, this analysis could reveal regions that responded similarly or differentially to the auditory and tactile adaptation events. For each node showing a significant main effect of adaptation condition, we defined an adaptation index (AI) as the response on CHANGE events minus the response on REPEAT events. Large AI absolute magnitudes indicate greater frequency-sensitivity and positive AI values reflect repetition suppression. For each participant, we sorted the data into 11 anatomically-defined ROIs (Destrieux et al., 2010) spanning the parietal and temporal lobes in standard surface space: central sulcus (cS), postcentral gyrus (pcG), postcentral sulcus (pcS), supramarginal gyrus (SMG), subcentral gyrus and sulcus (subcG/S), posterior Sylvian fissure (pSF), planum temporale (PT), transverse temporal sulcus (tTS), superior temporal gyrus (STG), lateral superior temporal gyrus (latSTG), and the inferior Insular cortex (infIns). To determine if AI values differed according to modality and ROI, we conducted a two-way rmANOVA with modality and ROI as factors.

Group-level univariate analysis: Across-modality adaptation events
To identify brain regions whose responses were consistent with shared and interactive processing of sound and vibration frequency, we evaluated repetition effects in the AT and TA data. We first tested for adaptation effects that were invariant to adapting frequency using the two-way rmANOVA with modality and adaptation condition as factors. Because no regions showed significant main or interaction effects of adaptation condition, we then tested for adaptation effects that depended on the adapting frequency. For the AT and TA conditions separately, we conducted two-way rmANOVA with adapting frequency (100 Hz or 300 Hz) and adaptation condition (REPEAT or CHANGE) as factors.

Modeling across-modality adaptation responses
For activation clusters identified as significant for the adapting frequency Â adaptation condition interaction, we evaluated a set of competing models that related β Obsthe observed BOLD response profile over the frequency conditions (adaptor-probe: 100-100, 100-300, 300-100, 300-300)to linear combinations of 4 potential predictor variables: adaptor frequency (f a ), probe frequency (f p ), response pattern to tactile-only events (β TT ), and response pattern to auditory-only events (β AA ). These terms enabled us to quantify how the adapting and probe frequencies contributed to the across-modality response pattern. We also reasoned that a component of a cluster's across-modality responses may reflect its response to the within-modality events. Thus, we tested 4 hypotheses by fitting 4 alternative models: where, w a , w p , w TT , and w AA are the weights to the model terms based on adapting frequency, probe frequency, the TT response pattern, and the AA response pattern. For each participant, the models were trained and tested on data over all of the significant clusters using a 10-fold crossvalidation procedure: On each fold, the model parameters were estimated using 90% of the data and model performance was determined using the held-out data. The cross-validation procedure was repeated 1000 times to generate distributions of parameter estimates and goodness-of-fit values. We performed model comparison using Akaike information criteria (AIC) and Bayesian information criteria (BIC) to identify the preferred model.

Representational similarity analysis (RSA)
We compared the spatial activation patterns associated with each event to define the representational space occupied by the adaptation stimulus sequences for each ROI. These cross-validated analyses included any node identified as significant in the univariate analyses. Each participant's scans were divided into 2 datasets on which separate GLMs were fitted to estimate 2 sets of multivariate patterns for each of the 17 event types (2 modalities Â 2 adapting frequencies Â 2 adaptation conditions þ 1 catch event). We then computed correlations between the activation patterns associated with all possible pairs of events. We calculated a distance metric for each pair of events (1 minus the correlation) where distances ranged between 0 (identical pattern) and 2 (perfectly anti-correlated pattern). The full distance matrix (DM) thus defined the representational space for a given ROI. Note that separate analyses in which response pattern similarity was defined by Mahalanobis distance instead of correlation yielded similar overall results (Supplemental Material).
To quantify how the similarity of activation patterns in each ROI related to the modality of the stimuli comprising the events, we defined a tactile modality similarity index (MSI T ) for each ROI using the 1st order DMs: where D across is the mean distance between the tactile-only events and the across-modality (AT and TA) events and D within is the mean distance between just the tactile-only events. Larger MSI T values indicate greater relative response pattern similarity among tactile-only events, which we expected for conventionally-defined somatosensory areas. We computed an analogous auditory modality similarity index (MSI A ) for each ROI as well.
To establish the similarity of representational spaces across ROIs, we compared the 1st order DMs established for the different ROIs and computed the correlations between the DMs of ROI pairs. We calculated a distance metric for each pair of ROIs (1 minus the correlation) to construct a 2nd order DM where distances ranged from 0 (identical 1st order DMs) to 2 (perfectly anti-correlated 1st order DMs). The full 2nd order DM thus defined the ROI landscape for each participant.
In group-level analyses, we averaged 2nd order DMs over participants to establish a mean 2nd order DM. We applied single-linkage hierarchical clustering on the group DM to quantitatively characterize the group-level ROI landscape. To better visualize the ROI landscape, we used classical multidimensional scaling (MDS) and projected the ROI distance matrix into a 2-dimensional subspace that maximally preserved the distances between ROI pairs. MDS was applied to each participant's 2nd order DM and the subspace representations were averaged after Procrustes alignment (Ejaz et al., 2015).

ROI-based resting-state analysis
This analysis was aimed at characterizing the similarity of spontaneous BOLD signal fluctuations across the 11 parietal and temporal ROIs. After standard preprocessing, BOLD signal time series were estimated for each node by taking the residuals after regressing out head motion and physiological noise signals (Barnes et al., 2014;Jo et al., 2010). Time series were extracted for all of the nodes included in the RSA. For each participant, a mean time series was calculated for each ROI separately. The correlation between the mean time series of pairs of ROIs was calculated and the correlations were converted into distances (1 minus correlation) ranging from 0 (perfectly correlated time series) to 2 (perfectly anti-correlated time series). This distance matrix thus provided an ROI landscape defined according to temporal correlations in participants' resting state data. Distance matrices were averaged across participants for a group DM which was directly compared to the group-averaged 2nd order DM calculated from the RSA. We computed the correlation between the group DMs to determine the similarity of the ROI landscapes determined according to the spatial activation profiles in the task-based data and the temporal correlations present in the resting state data. We tested the significance of this correlation using a randomization test which compared the observed correlation against a null distribution of correlations expected by chance. The null distribution was computed by repeating the correlation calculation 1000 times with the node-labels randomly assigned to the data before computing ROI averages on each iteration.

Data and code availability
All data and analysis code are available upon direct request.

Results
3.1. Behavioral results: performance levels are consistent with participants attending to stimulus frequency Participants performed a frequency monitoring task as they underwent fMRI. Each event comprised a series of 4 brief stimuli delivered to the participant's right hand, and participants covertly judged whether the frequency of the probe stimulus matched the frequency of the preceding three adaptor stimuli. To confirm that participants performed this covert frequency monitoring task, participants were required to report their frequency judgment by button press when visually cued randomly throughout each scan. Participants generally maintained high performance on this task (Fig. 1C). Although performance varied significantly according to the modality conditions (F ¼ 12.4, P ¼ 2.3e-6), performance levels on average exceeded chance (0.5) on both the within-modality trials (AA: 0.93 AE 0.03, t ¼ 15.3, P ¼ 3.8e-12; TT: 0.90 AE 0.03, t ¼ 13.8, P ¼ 2.3e-11) and across-modality trials with auditory adaptation (AT: 0.75 AE 0.06, t¼ 4.5, P ¼ 2.3e-4). Performance on across-modality trials with tactile adaptation was more variable across subjects and the group average was marginally greater than chance (TA: 0.66 AE 0.06, t ¼ 2.7, P ¼ 0.016). These results indicate that most participants were compliant in attending to stimulus frequency during the unimodal and crossmodal events.
3.2. Localizer results: auditory and tactile stimulation are associated with partially-overlapping activity patterns in parietal and temporal cortex Fig. 2A shows the bilateral regions in which BOLD signal changes were associated with blocks of auditory and tactile stimulation. Although sound-related activity tended to be confined to temporal regions and touch-related activity tended to be confined to parietal regions, BOLD signal variations in a number of regions were associated with stimulation in either sensory modality. The analysis mask used to evaluate the adaptation scans included any region that exhibited auditory or tactile responses during the localizer scan.
3.3. Within-modality repetition suppression reveals frequency-selective auditory and tactile processing in overlapping regions We conducted a two-way rmANOVA to identify regions whose response variations were related to stimulus modality or adaptation condition. This analysis revealed a number of regions exhibiting a main effect of stimulus modality (Table 1)those that responded differentially to AA and TT eventsand the modality-dependent responses spanned much of the parietal and temporal regions included in the localizer mask. To identify brain regions whose response patterns were consistent with population-level (voxel-wise) neural selectivity for auditory or tactile frequency, we contrasted responses measured when the probe frequency differed from the adapting frequency (CHANGE events) against responses measured when the probe frequency matched the adapting frequency (REPEAT events). This contrast (main effect of adaptation condition) revealed significant activation clusters (Fig. 3A) that were largely confined to the parietal and temporal lobes in the left hemisphere (Table 1). Notably, these clusters all exhibited greater responses to CHANGE events compared to REPEAT events, revealing a pattern consistent with repetition suppression. This pattern is evident in the GLM coefficients (Fig. 3B) and estimated temporal response profiles (Fig. 3C) for auditory-only and tactile-only events in an example ROI, the posterior Sylvian fissure. In fact, greater responses to CHANGE events compared to REPEAT events were observed in nearly all of the defined temporal and parietal ROIs with the auditory-only events (Fig. 3D) and tactile-only events (Fig. 3E). Although auditory-only and tactile-only events were both associated with distributed responses over putative auditory and somatosensory regions, auditory responses were relatively greater in temporal regions while tactile responses were relatively greater in parietal regions. Importantly, no regions showed significant modality Â adaptation interaction effects, indicating that adaptation response patterns tended to be comparable for the AA and TT events.
We compared the magnitude of response differences between CHANGE and REPEAT unisensory events by calculating an adaptation index for each sensory modality over the 9 ROIs comprising significant clusters (Fig. 3F). AI values differed significantly according to brain region (2-way rmANOVA; ROI main effect: F ¼ 2.12, P ¼ 0.04; ROI Â modality interaction effect: F ¼ 3.63, P ¼ 0.0007), but the main effect of modality failed to achieve statistical significance (F ¼ 0, P ¼ 0.94). Because of the significant interaction effect, we compared adaptation indices for the auditory and tactile events in each ROI and found significant differences only in the supramarginal gyrus (t¼ 4.09, P ¼ 2.2e-4) and the subcentral gyrus/sulcus (t ¼ 3.14, P ¼ 0.003) with larger adaptation effects observed for the tactile events in both ROIs. We additionally determined if there was a general relationship between auditory and tactile adaptation strength across temporal and parietal regions and found that adaptation indices for the two sensory modalities, averaged over participants, tended to be negatively correlated (mean ρ ¼ À0.084 AE 0.026) (one-sample t-test versus 0; t ¼ 3.212, P ¼ 0.0048). When performed on each ROI separately, this analysis revealed only significant correlations in the supramarginal gyrus (mean ρ ¼ À0.122 AE 0.038; t ¼ 3.208, P ¼ 0.0046). These results indicate that the neural populations most sensitive to auditory and tactile temporal modulation effects do not necessarily reside in the same voxels.

Across-modality adaptation events reveal frequency-selective temporal modulation effects in parietal cortex
To identify brain regions whose response patterns were consistent with population-level crossmodal frequency adaptation, we contrasted CHANGE and REPEAT responses in analyses that combined the AT and TA events (i.e., analysis independent of adapting modality) and analyses performed on the AT and TA events separately. These contrasts failed to reveal any significant clusters. We then tested the possibility that acrossmodality interaction effects varied according to the adaptor and probe frequencies. Accordingly, for the AT and TA events separately, we performed contrasts that involved interactions between the adaptation condition (CHANGE or REPEAT) and adapting frequency (100 Hz or 300 Hz) factors. Although no clusters displayed significant effects for this interaction in analysis of AT events, we identified 15 significant clusters in anterior parietal cortex and lateral parietal cortex with TA events (Fig. 4A). Similar spatial patterns were revealed in analyses that excluded participants who achieved <60% performance on the frequency monitoring task, which indicates that the activations shown in Fig. 4A were not dominated by the responses of poorly-performing subjects. The response profile in these clusters revealed a clear interaction pattern of distinct repetition effects which depended on adapting frequency (Fig. 4B). We observed repetition facilitation with 100-Hz adaptors as REPEAT events were associated with greater responses compared to CHANGE events. In contrast, we observed repetition suppression with 300-Hz adaptors as responses to REPEAT events were reduced relative to CHANGE events.
To better understand the TA interaction pattern, we fit a set of simple linear models to the TA responses. Model comparison using both AIC and BIC identified the same preferred model (Fig. 4C), which comprised terms related to frequency-specific effects of the tactile adaptors and auditory probes as well as the TT response patterns. Fig. 4D shows that the preferred model, fitted to responses over all 15 clusters exhibiting significant TA interactions in a cross-validated manner (Materials and Methods), explained greater than 80% of the response variance across clusters. Inspection of the model weights, capturing adaptor (0.000619 AE 0.000048) and probe (À0.00076 AE 0.0000382) frequency effects as well as the impact of the TT response profiles (0.76 AE 0.04), reveals that 1) the tactile adaptors contributed positively to cluster responses with Significant clusters identified in group-level analyses (q < 0.05 FDR corrected). Data are shown for clusters containing 10 or more nodes. A total of 82 clusters (2 node minimum) exhibited a significant main effect of modality (AUDITORY vs TACTILE). The 32 clusters comprising 10 or more nodes are sorted according to their relative responses to AA and TT conditions. A total of 53 clusters (2 node minimum) exhibited a significant main effect of adaptation condition (CHANGE vs REPEAT). The 23 clusters comprising 10 or more nodes all exhibited larger response to CHANGE compared to REPEAT events. larger modulations associated with the 300-Hz adaptors, 2) the auditory probes contributed negatively to cluster responses with larger modulations associated with the 300-Hz probes, and 3) the scaled TT response profile comprised a portion of the TA response profile. These results, along with the fact that no constraints were placed on which clusters contributed model training data and testing data, reveal that TA response patterns were highly consistent across the distributed clusters spanning anterior and lateral parietal cortex. In fact, the same results pattern was observed when data from the large parietal cluster were excluded from the model fitting and testing procedure which indicates that the TA response profile was not driven by the largest activity cluster. In sum, the model results reveal that TA response patterns in parietal regions can be understood as a combination of repetition suppression (represented by the TT pattern), frequency-specific response enhancement from the tactile stimulus component, and frequency-specific auditory response inhibition from the auditory stimulus component.

Spatial activation patterns are consistent with modality-based cortical organization
We performed multivariate analyses to further explore how unimodal and crossmodal adaptation sequences were represented across temporal and parietal regions (Material and Methods). In each ROI separately, we determined the pairwise similarity of spatial activation patterns associated with each event type in a split-half cross-validated manner (Fig. 5A).
Importantly, by performing these analyses using cross-validation, we minimized the possibility that representational structures in the data were trivially attributable to correlated noise. Fig. 5B shows a 1st order distance matrix calculated for the superior temporal gyrus in a representative participant. The elements of this DM describe the response pattern similarity between event pairs and the collective DM can be considered as reflecting the representational space occupied by the unimodal, crossmodal, and catch-trial events in STG neural space. The STG DM clearly shows that the spatial activation patterns associated with events comprising auditory stimuli (i.e., the AA, AT, TA events) tended to be highly similar to each other and dissimilar from TT events. Patterns associated with TT events also exhibited a degree of similarity with each other. While the representational space of STG is unsurprising given the known role of STG in auditory processing, the modest similarity between TT events reinforces the idea that temporal regions also respond systematically to tactile inputs alone.
To establish the relationship between representational spaces in different temporal and parietal regions, we generated 1st order DMs for each ROI and subsequently quantified the similarity of these DMs between ROI pairs. A 2nd order DM for the representative participant is shown in Fig. 5C, with each element describing the similarity of the representational spaces of ROI pairs. The example DM reveals systematic relationships between representational spaces over the ROIs. First, there were greater similarities among the regions contained within the parietal lobe and the regions within the temporal lobe. Second, the  ¼ 20). A, Red nodes indicate regions that exhibited significant response differences between CHANGE and REPEAT events regardless of modality (FDR corrected, q < 0.05; constrained by functional localizer mask; minimum cluster size ¼ 2 nodes). The yellow outline indicates the boundaries of the posterior Sylvian fissure (pSF) region-of-interest whose responses are summarized in B and C. B, Group-averaged response magnitudes estimated for auditory and tactile within-modality adaptation events for the left pSF region-of-interest. This ROI exhibited responses to both auditory-only and tactile-only stimulus sequences. CHANGE events (blue bars) were associated with greater responses compared to REPEAT events (red bars) for auditory and tactile comparisons. Error bars indicate S.E.M. C, Group-averaged temporal response profiles estimated for auditory and tactile within-modality adaptation events for left pSF. Responses to CHANGE events (blue traces) tended to larger than responses to REPEAT events (red traces). Thin traces indicate individual subject data. D, Average responses to auditory-only REPEAT and CHANGE events. Each marker indicates mean response over nodes in different parietal and temporal ROIs (inset) exhibiting a significant main effect of adaptation condition (A). Error bars indicate S.E.M. E, Average responses to tactile-only REPEAT and CHANGE events. Conventions as in D. F, Box and whisker plot depicts adaptation index (Materials and Methods) for auditory-only (white columns) and tactile-only events (gray columns) for different ROIs. Horizontal black lines indicate medians. Box edges indicate interquartile range (IQR). Whiskers extend to most extreme data points within 1.5x IQR and individual points indicate data outside of 1.5x IQR. pcG, postcentral gyrus; pcS, postcentral sulcus; SMG, supramarginal gyrus; subcG/S, subcentral gyrus and sulcus; pSF, posterior silvian fissure; PT, planum temporale; tTS, tranverse temporal sulcus; latSTG, lateral superior temporal gyrus; infIns, inferior insula.
representational spaces of parietal regions tended to be dissimilar from those of temporal regions. This systematic pattern was also obvious in the average 2nd order DM over participants (Fig. 6A) with hierarchical clustering revealing 2 ROI groupings that were clearly organized according to regional affiliation to parietal or temporal cortex. To better visualize the relationship between ROIs, we projected and aligned each participant's 2nd order DM in a two-dimensional "brain region" subspace (Material and Methods). The configuration of the brain regions in this subspace (Fig. 6B) clearly reveals structured relationships over the ROIs with a parietal "stream"spanning anterior parietal cortex (central sulcus, postcentral gyrus, postcentral sulcus), lateral parietal cortex (subcentral gyrus and sulcus), and posterior parietal cortex (supramarginal gyrus)and a separate temporal "stream"spanning STG, tranverse temporal sulcus, lateral STG, planum temporale, and posterior Silvian fissure. Notably, the MDS pattern also indicates that the inferior insular cortex contains a representational space that is distinct from the other ROIs. Analogous cross-validated analysis performed using Mahalanobis distance rather than Pearson correlation produced a highly similar cortical landscape in the 2D subspace (Fig. S1). To establish a better intuition for the range of representational spaces observed over the ROIs, we highlighted the average 1st order DMs for the postcentral gyrus and the superior temporal gyrus (Fig. 6C). We also highlighted the average 1st order DM for the planum temporale ROI, which occupied a location between the parietal and temporal groups in the brain region subspace. As expected, the representational spaces in the parietal and temporal regions reflected a clear modality-dependence: the postcentral gyrus DM showed greater similarities among events comprising tactile cues (TT, TA, AT) while the STG DM showed greater similarities among events comprising auditory cues, as seen in the example participant ( Fig. 5B). In contrast, because sensory events were generally associated with more similar activity patterns in planum temporale, the DM computed for PT was less characterized by modality-dependence. This pattern implies that regions like the PT and pSF, which are traditionally considered to be higher-order sensory association areas, differentiate between the auditory and tactile stimulus components less in their responses to the unimodal and crossmodal events than regions like the pcG and STG. To quantify these ROI differences, for each ROI we defined separate indices that expressed the modality-dependence of response pattern similarity and how response magnitudes varied according to modality (Materials and Methods). Both metrics (Fig. 7) varied significantly over ROIs (MSI: main effect of ROI: F ¼ 6.5, P ¼ 1.3e-8, main effect of modality: F ¼ 21.2, P ¼ 0.0002, ROI Â modality interaction: F ¼ 28.2, P ¼ 2.2e-16; Response magnitude: main effect of ROI: F ¼ 25.8, P ¼ 2.2e-16, main effect of modality: F ¼ 134.5, P ¼ 4.6e-10, ROI Â modality interaction: F ¼ 85.1, P ¼ 2.2e-16). These indices reveal cortical landscapes that are consistent with the conventional view of sensory systems: early sensory areas are more dedicated to a single modality and higherorder association areas are responsive to multiple modalities. Thus, the multivariate results are consistent with the traditional notion of modality-based sensory systems.
3.6. Resting state functional connectivity patterns are also consistent with modality-based cortical organization We next determined whether the structured relationships between the parietal and temporal regions, established from the similarity of representational spaces in each ROI, were intrinsic characteristics of our participants' brains or merely a reflection of task-evoked responses given Fig. 4. Results with tactile adaptation and auditory probes. (n ¼ 20) A, Crossmodal adaptation sequences comprising tactile adaptors and auditory probes (TA events) were associated with response modulation patterns in anterior and lateral parietal cortex (inset) that were characterized by significant interactions between adapting frequency (100 Hz vs 300 Hz) and adaptation condition (CHANGE vs REPEAT). Red nodes indicate clusters comprising a minimum of 2 surface nodes. B, Average activation profile for example cluster. With 100-Hz adaptors, larger responses are observed on REPEAT events (i.e., with the 100-Hz probes compared to the 300-Hz probes). In contrast, with 300-Hz adaptors, larger responses are observed on CHANGE events (i.e., with the 100-Hz probes compared to the 300-Hz probes). Error bars indicate S.E.M. C, TA response profiles in each cluster can be explained as a linear combination of the cluster's TT response profile (green bars), frequencyspecific response modulation attributed to the tactile adaptor stimuli (orange bars), and frequency-specific response modulation attributed to the auditory probe stimuli (purple bars) (Materials and Methods). The weighted sum of the depicted response components accounts for the activity pattern shown in B. D, Observed and model predicted response patterns over 15 clusters (minimum cluster size ¼ 2 nodes) exhibiting significant TA interaction patterns. Each cluster's profile includes 4 responses (gray bar indicates profile for Cluster 1). Cross-validated models fit to response patterns in 90% of the data accounted for an average of 81% of the response variance in the held-out data (Materials and Methods).
our particular paradigm and stimulus conditions. Accordingly, we performed an analysis on each participant's resting state fMRI data that was analogous to the RSA performed on their task-based fMRI data. For each participant, we computed the average BOLD signal time series in each of the 11 temporal and parietal ROIs (Fig. 8A). We then calculated the correlation between the time series in each ROI pair to generate a "connectivity" matrix expressing the degree of similarity in the temporal variations of spontaneous signal fluctuations across the ROIs. The average connectivity matrix over participants was significantly correlated with the group-averaged 2nd order DM defined in the RSA (Fig. 8B; ρ ¼ 0.28, P < 0.001; mean null ρ ¼ 0.16 AE 0.03). In fact, direct comparisons between the connectivity matrix and 2nd order DMs within participants yielded correlations ranging from 0.84 to 0.95 with a significant mean correlation of 0.92 AE 0.03. Indeed, projecting and aligning the brain regions across participants in a two-dimensional subspace (Materials and Methods) resulted in a ROI landscape (Fig. 8C) that closely resembled the landscape generated in the multivariate pattern analysis. The close correspondence between the RSA results (based on spatial variations) and the resting state connectivity analysis results (based on temporal variations) implies that the modality-based organization of regions spanning parietal and temporal cortex is an intrinsic characteristic rather than a structured pattern that emerged as a consequence of our experimental manipulations.

Discussion
We used a fMRI adaptation paradigm to identify brain regions exhibiting selectivity for auditory and tactile temporal frequency. With both unimodal tactile and auditory stimulus sequences, we observed repetition suppression in temporal regions, which we assumed would exhibit frequency selectivity a priori, as well as a number of parietal regions traditionally associated with somatosensory functions. To establish evidence for shared neural representations and multimodal frequency interactions, we tested crossmodal adaptation sequences and only observed significant repetition effectswhich differed for the 100-Hz compared to 300-Hz adaptor eventsin parietal cortex for sequences comprising tactile adaptors followed by an auditory probe. Collectively, our adaptation results imply that both tactile and auditory frequency information modulate parietal and temporal cortex activity. We additionally analyzed the spatial activation patterns associated with the adaptation sequences in the parietal and temporal regions and found relationships between the ROI representational spaces that clearly reflected modality-based organization. The RSA results revealed a cortical landscape characterized by a set of modality-preferring regions flanking regions that were more modality-invariant. This landscape was also observed in the modality-sorted response amplitude profiles across the ROIs as well as the spatiotemporal patterns in spontaneous signal fluctuations measured in these areas at rest. Thus, over parietal and temporal cortex, we found evidence for a traditional modality-based sensory cortex model alongside evidence for frequency-specific auditory and tactile processing.
The results of our block design localizer scans revealed somatosensory and auditory activations in parietal and temporal regions, respectively, along with overlapping activations consistent with earlier reports (Foxe et al., 2002;Nordmark et al., 2012;Schurmann et al., 2006). To test whether response modulation in these regions could be associated with frequency-specific processing, we quantified BOLD signal changes in a fMRI frequency adaptation paradigm using an event-related design. Within-modality adaptation events were associated with BOLD signal changes that were consistent with repetition suppression effects (Barron et al., 2016;Grill-Spector et al., 2006;Krekelberg et al., 2006): Responses to stimulus sequences comprising repeats of a single frequency were significantly reduced relative to responses to sequences in which the probe frequency changed from the adaptor frequency. Many of the regions exhibiting frequency-based repetition suppression have also been shown to exhibit BOLD adaptation effects in other somatosensory Fig. 5. Representational structure of unimodal and multimodal stimulus sequences. A, Activity matrix depicts the spatial activation patterns associated with each event (rows) over the nodes comprising a ROI (columns). Rows are additionally ordered according to frequency conditions (adaptor-probe frequencies) in each modality set (AA, AT, TA, and TT modality conditions). B, First-order distance matrix (DM) indicating the pair-wise similarity between activation patterns associated with different event types in the superior temporal gyrus ROI. Lower distance values indicate more similar activation patterns. Similarity was quantified in a cross-validated manner with a single participant's data acquired in different runs (Materials and Methods). Due to run-specific noise, identical events over different runs were associated with slightly different activation patterns which resulted in non-zero distances along the main diagonal. C, Second-order DM in example participant indicating the similarity of first-order DMs between different ROI pairs. Lower distance values indicate that greater similarities in the representational spaces of two ROIs. (Hegner et al., 2007;Tame et al., 2012) and auditory (Millin et al., 2018) contexts. Our BOLD adaptation effects could relate to perceptual adaptation effects reported for frequency perception by somatosensation Hollins, 1994, 1993;Hollins et al., 1990;Tommerdahl et al., 2005) and audition (Alais et al., 2015;Parra and Pearlmutter, 2007;Zwicker, 1964). Although the magnitude of repetition suppression BOLD effects can be correlated to behavior, particularly in priming or memory paradigms (Dobbins et al., 2004;Voss et al., 2009;Wig et al., 2005), we did not observe any meaningful or significant brain-behavior relationships. The challenges associated with establishing brain-behavior correlations with fMRI data are well-documented (Rousselet and Pernet, 2012) and the absence of such correlations in our study may be due to a number of factors. First, our behavioral paradigm was primarily designed merely to verify participants' compliance on the frequency monitoring task rather than to provide robust performance estimates, measure changes in behavior (as in learning or priming), or even to quantify perceptual sensitivities. The poor performance of some participants, particularly in the TA condition, may have indicated a lack of task compliance and failure to attend to the stimuli, but it is important to note that BOLD repetition effects can occur in the absence of attention (Bentley et al., 2003;Hsu et al., 2014;Naccache and Dehaene, 2001;Vuilleumier et al., 2005). Second, with our coarse performance measures, there was little between-subject variance in performance, particularly with the within-modality events, and this likely influenced the sensitivity of the linear regression models to establish relationships between behavior and brain responses. Alternative behavioral paradigms and fMRI designs are likely required to establish the relationship between BOLD signal changes and perceptual adaptation to temporal frequency information.
A limitation of fMRI adaptation paradigms and the interpretation of repetition suppression effects is the potential confounding effect of attention. Specifically, because stimulus sequences consisting of changes can be more perceptually salient than sequences consisting of repeating stimuli, the larger BOLD signal changes associated with the CHANGE events may be attributable to general attentional modulation effects rather than feature-specific processing. We addressed this critical consideration with two additional analyses. First, to address the concern that the repetition suppression effects we observed only reflected  attentional modulation effects based on change detection, we characterized brain responses in areas previously implicated in automatic change detection and novel event processing (Ranganath and Rainer, 2003). We reasoned that, if our stimuli effectively caused participants to engage in change detection, activity levels in these regions should also differentiate between CHANGE and REPEAT events. We compared CHANGE and REPEAT activations (Supplemental Material) in ROIs centered in the right temporoparietal juncture, the right anterior insula, and the left cingulate motor area because these ROIs have been shown to respond more to changing versus repeating sensory conditions (Downar et al., 2000;Huettel et al., 2002) and these regions did not overlap with the analysis masks defined by our block design localizer scans. Unlike the parietal and temporal regions in the left hemisphere, which we identified as showing within-modality repetition suppression effects (Fig. 3, Table 1), no significant response differences between CHANGE and REPEAT events were observed in these control ROIs (Fig. S2). Although these supplemental results do not allow us to definitively infer participants' cognitive states during our experiment, they at least show that the sensory events in our experiments did not significantly engage brain regions previously implicated in automatic change detection, though they were associated with robust repetition suppression effects in sensory regions spanning parietal and temporal cortex.
In a separate effort to address the possibility that our results reflected attentional modulation rather than feature-specific processing, we performed cross-validated multivariate pattern analysis (MVPA; Supplemental Material) and trained classifiers to decode auditory or tactile frequency information from distributed signals measured in parietal and temporal cortex, specifically the regions that exhibited significant withinmodality adaptation effects (Table 1). We trained and tested the classifiers on only the response patterns associated with the AA and TT REPEAT events. We reasoned that the activation patterns in the parietal and temporal regions should represent sound and vibration frequency information if these regions support frequency-specific processing. For auditory decoding (Fig. S3), significant classifier performance was observed when using only temporal cortex data, but not when using only parietal cortex data. However, when a classifier was trained and tested on data combined over parietal and temporal regions, decoding accuracy was significantly better than that achieved with temporal cortex data alone, implying that some meaningful auditory frequency information resides in parietal cortex, consistent with earlier reports . For tactile decoding (Fig. S3), classifier performance failed to significantly exceed chance performance with only parietal cortex or temporal cortex data. Critically, significant decoding accuracy was achieved when the classifier was trained and tested on the combined data over parietal and temporal regions, implying that tactile frequency information may be represented in a distributed manner over parietal and temporal regions. Moreover, significant decoding accuracy was not achieved using data combined over parietal cortex and a control ROI in visual cortex (Supplemental Material), indicating that merely enlarging the multivariate pattern alone was insufficient to yield improved classifier performance. Thus, decodable tactile frequency information appeared to be specifically available in the activation patterns distributed over parietal and temporal cortex. Admittedly, our experimental paradigm was not designed to produce data well-suited for decoding analyses, but these supplemental MVPA results provide some corroborative evidence for feature-specific processing in the parietal and temporal regions identified by repetition suppression.
Our within-modality adaptation results suggest population-level frequency tuning for auditory and tactile stimulation in both parietal and temporal cortex. For auditory regions, the observation of auditory and tactile frequency adaptation effects is perhaps unsurprising given the established existence of neural populations tuned for frequency (Bendor and Wang, 2005;Wang and Walker, 2012) and the common observation that auditory cortical regions can respond to tactile stimulation alone (Foxe et al., 2002;Kayser et al., 2005;Nordmark et al., 2012;Perez-Bellido et al., 2017;Schurmann et al., 2006), even at the unit level (Fu et al., 2003;Lemus et al., 2010). Our data also revealed repetition suppression effects for both modalities in somatosensory areas. The general observation of auditory responses in parietal regions is consistent with other fMRI studies (Beauchamp and Ro, 2008;Liang et al., 2013;Perez-Bellido et al., 2017). Furthermore, causal manipulation of parietal cortex activity can selectively disrupt auditory frequency perception (Convento et al., 2018). However, the inference that neurons in somatosensory regions are tuned for frequency is more tenuous because there has been scant neurophysiological evidence for explicit rate-based frequency tuning in somatosensory neurons unlike their counterparts in the auditory system. Indeed, despite the suggestion that explicit coding for vibration frequencies exists based on modeling studies (Bensmaia et al., 2005;Crommett et al., 2017;Rahman and Yau, 2019), only frequencies in the flutter range (<65 Hz) appear to be encoded in spike rates (Saal et al., 2016;Salinas et al., 2000); higher frequencies (>100 Hz) like those tested here only appear to be represented in a spike timing code in primate somatosensory cortex (Ferrington and Rowe, 1980;Harvey et al., 2013;Mountcastle et al., 1969). Note that a recent study reported explicit  ). B, Correlations between the time series of ROI pairs defined a second-order "connectivity" matrix that summarized the network architecture in the resting state data of each participant. This architecture was significantly correlated with the second-order distance matrix computed from representational similarity analyses (Fig. 6A): The actual correlation (red) between ROI landscapes generated from task-based data and resting state data far exceeded the correlations in the null distribution (blue). C, Multidimensional scaling of the ROI distances (based on resting state correlations) in two-dimensional space. Ellipses show S.E.M. after Procrustes alignment across participants. Although the resting state analyses were performed independently of the task-based multivariate analyses, these yielded highly similar ROI landscapes.
coding for high frequency vibrations in somatosensory neurons for the first time in the mouse (Prsa et al., 2019). Thus, while our results imply that frequency-tuned neural populations may reside in parietal cortex of primates, this conjecture remains to be tested in neurophysiology experiments.
The recruitment of auditory cortex by touch has been interpreted as evidence for a function-based organization scheme in which putative auditory neural circuits that are specialized to represent temporal frequency information also support tactile frequency processing (Yau et al., 2009b). Supramodal processing systems have been invoked for the processing of shape (Lacey et al., 2009;Pascual-Leone and Hamilton, 2001) and motion (Konkle et al., 2009;Wacker et al., 2011). Despite our finding that auditory and tactile frequency-selective responses overlap at a voxel level in many regions, our collective results would appear to argue against a strong supramodal hypothesis for a number of reasons. First, if the same neural populations represented auditory and somatosensory frequency, we may have expected that the magnitude of adaptation effects would be positively correlated between the two modalities. We instead only observed negative or weak correlations over all ROIs or in individual ROIs. Second, if somatosensory and auditory frequency information were carried in overlapping neural populations that were distributed over frequency-processing cortex, a classifier trained on data associated with stimulation in one modality would be predicted to decode frequency information from data associated with stimulation in the other modality. In supplemental analyses (Supplemental Material), we failed to observe significant crossmodal decoding even when using the combined parietal and temporal cortex data which yielded significant cross-validated decoding for auditory and tactile frequency separately (Fig. S3). Finally, if neurons carried supramodal frequency representations, we may have expected to find repetition suppression effects with the crossmodal events that were comparable to those seen with unisensory events. Despite psychophysical evidence for crossmodal adaptation effects (Badde et al., 2016;Crommett et al., 2017;Levitan et al., 2015) and crossmodal adaptation of fMRI responses with other stimulus features (Doehrmann et al., 2010;Tal and Amedi, 2009), we did not observe significant simple across-modality adaptation effects. With the TA events, significant frequency-dependent patterns were found in anterior parietal cortex and lateral parietal cortex. These patterns were largely explained as a linear combination of 1) frequency-specific excitatory response modulation attributed to the tactile adaptors, 2) frequency-specific suppressive response modulation attributed to the auditory probes, and 3) a response component consistent with the activations associated with the TT events. Importantly, each of these response components reflected frequency-dependent processing, so the TA results provide some evidence for the existence of neural circuitry in parietal cortex that responds selectively to auditory and tactile frequency. In sum, while our collective univariate analysis results imply that a number of parietal and temporal regions exhibit frequency-specific responses to auditory and tactile stimulation, these representations may be carried in separate or minimally shared neural populations at the voxel level (Driver and Noesselt, 2008;Klemen and Chambers, 2012).
Although exposure to auditory tones has been shown to modulate the subsequent perception of tactile frequency in a frequency-specific manner , our fMRI data did not reveal significant response modulation with the AT adaptation events. One possibility is that the response components associated with the tactile probes were too weak relative to the responses associated with the auditory adaptors. Indeed, responses to tactile stimulation tended to be weaker than responses to auditory stimulation (Fig. 3). Likewise, spatial activity patterns were more consistent and robust for auditory stimuli compared to tactile stimuli in the higher-order regions that exhibited prominent multimodal responsiveness (Fig. 7). A strong activity imbalance favoring the auditory modality may have obscured potential frequency-dependent tactile responses or interaction effects. This imbalance may reflect differences in the perceived intensities of the auditory and tactile events, despite our efforts to perform across-modality equating. Alternatively, fMRI responses may simply be more robust for sounds experienced through ear-buds compared to vibrations experienced on the fingertips of one hand. Another possibility is that longer adaptation periods are required to reveal the fMRI correlates of AT adaptation effects. The event-related design used here was only suited to reveal adaptation effects that could emerge within a period of seconds, but perceptual crossmodal adaptation effects have been seen following prolonged adaptation on the order of minutes. Future studies should manipulate the relative strength of the auditory and tactile responses, perhaps by testing auditory band-passed noise stimuli which also interact with touch in a frequency-dependent manner Yau et al., 2009b), and adaptation durations to establish a more comprehensive understanding of audiotactile crossmodal adaptation effects.
Beyond the potential confounding effects of attention on our experimental paradigm addressed above, there are additional caveats to consider. First, there may be concerns that our results reflect in part the response demands of the behavioral task used in the adaptation scans. Specifically, because we required participants to perform covert frequency judgements on each event while overtly reporting these judgements only when cued explicitly, activity associated with motor preparation may have differentiated between the REPEAT and CHANGE events even when no responses were required. However, this may be an unlikely explanation of our repetition suppression results, which were mainly confined to the left hemisphere, because 1) motor preparation likely would have preferentially involved the (right hemisphere) sensorimotor regions contralateral to the (left) hand used for the button box responses and 2) the response demands were identical on within-and across-modality events while we only observed repetition suppression effects with the former. Another critical consideration relates to the interpretation of localized repetition suppression effects: Given the relatively coarse temporal resolution of fMRI and the sluggishness of the BOLD signal, there is ambiguity in relating repetition suppression effects to neural processing in any particular region. As with all fMRI, activation patterns likely reflect the processing of feedback or recurrent signals as much as feedforward responses. Furthermore, adaptation patterns can be inherited from upstream processes (Mur et al., 2010;Tolias et al., 2005). Accordingly, adaptation effects observed in higher-order brain regions may simply reflect the repetition suppression effects in earlier cortical regions or even adaptation in the periphery or subcortical systems. Analogous experiments using imaging and recording methods with higher temporal resolutions will be needed to address these concerns. Similarly, additional experiments that more finely sample the stimulus parameter space (adaptation durations and frequencies) may better enable efforts to link the BOLD adaptation results to neural computations using quantitative encoding models as has been done in studies of compressive temporal summation in the human visual system (Zhou et al., 2017). A more complete sampling of stimulus space may also support pattern component modeling efforts (Diedrichsen et al., 2018) that could link the RSA results to neural encoding models. Such experiments are required to move beyond the more descriptive adaptation results presented here that may only be suited to inform questions regarding the sensory modality preferences of parietal and temporal cortex.
Our multivariate analysis of the task-based data revealed an ROI landscape that was consistent with the traditional view of parietal and temporal cortical systems. Even as we observed multimodal responses in the parietal and temporal regions, there was a general organization seemingly defined according to sensory modality preferences. Crucially, our resting state analyses revealed a highly similar organization pattern between the parietal and temporal regions (Fig. 8). This correspondence was critical because it implies that the ROI landscape cannot be attributed to stimulus-related activity only. In control experiments, we took particular care to confirm that the tactile stimulation used in the scanner were inaudible (Materials and Methods). Moreover, in an earlier experiment , we confirmed that parietal responses to auditory stimulation could not be trivially attributed to mechanical vibrations generated via the in-ear bud headphones by demonstrating that the similar activations can result from sounds heard through circumaural air-conduction headphones. These control experiments alleviate concerns that the overlap between the auditory and tactile responses could be explained by some physical explanation related to poor stimulus control; the resting state results, free of any sensory stimulation effects, directly address this issue. The relationship between intrinsic fluctuations and representational similarity analyses has received recent attention (Henriksson et al., 2015) and our data are consistent with the notion that spatially-patterned intrinsic cortical dynamics could underlie or enhance apparent relationships in the representational geometries over different brain regions. Regardless of how the RSA and functional connectivity results relate, it is indisputable that both reveal a ROI landscape resembling the traditional model of sensory cortex organization. This landscape is also evident in anatomical connectivity profiles determined from invasive tracer studies (Cappe et al., 2009a;Cappe and Barone, 2005;Hackett et al., 2007) and diffusion-tensor imaging (Ro et al., 2013). How information is dynamically routed through these connected systems remains an open question. Attention may play a critical role in gating signal transmission between the parietal and temporal networks (Convento et al., 2018). Future studies are required to investigate the state-dependence of frequency information processing in parietal and temporal cortex.
In summary, our results provide evidence for frequency-selective processing of tactile and auditory stimulation over brain regions spanning parietal and temporal cortex. The frequency-selectivity, which is observable at the scale of individual voxels, is revealed in unimodal repetition suppression effects as well as the frequency-dependent suppression seen with the auditory probes in the TA events. These results contribute to the growing literature showing that sensory regions traditionally thought to be dedicated to a single modality can respond to inputs from multiple sensory modalities. Our data highlight the featurespecificity of the multimodal responses, consistent with a functionbased scheme for brain organization; however, it is critical to note that the sparse examples of across-modality adaptation suggest only weak evidence for supramodal frequency representations. Importantly, we also see evidence of a modality-based scheme in which brain regions are organized according to their preferences for audition or touch. Thus, our results reveal that multimodal frequency responses can be distributed over traditionally-defined sensory cortical systems. Because we only characterized responses to sinusoidal stimulation, our analyses could only address processing differences associated with frequency or modality manipulations. Future studies will need to test more complex and naturalistic stimuli, like frequency sweeps (Crommett et al., 2019) or textures (Lederman, 1979;Manfredi et al., 2014;Yau et al., 2009a), to more thoroughly investigate how auditory and tactile representations are elaborated and maintained over hierarchical processing streams in parietal and temporal cortex.

Declaration of competing interest
MT discloses that he is a cofounder of Cortical Metrics, which built and provided the tactile stimulator used in these experiments, and he receives royalties through this. The remaining authors declare no competing financial interests.