Spatio-temporal profile of brain activity during gentle touch investigated with magnetoencephalography

Positive affective touch plays a central role in social and inter-personal interactions. Low-threshold mechanoreceptive afferents, including slowly-conducting C-tactile (CT) afferents found in hairy skin, transmit such signals from gentle touch to the brain. Tactile signals are processed, in part, by the posterior insula, where it is the thought to be the primary target for CTs. We used magnetoencephalography (MEG) to assess brain activity evoked by gentle, naturalistic stroking touch on the arm delivered by a new MEG-compatible brush robot. We aimed to use high temporal resolution MEG to allow us to distinguish between brain responses from fast-conducting Aβ and slowly-conducting CT afferents. Brush strokes were delivered to the left upper arm and left forearm of 15 healthy participants. We hypothesized that late brain responses, due to slow CT afference, would appear with a time shift between the two different locations on the arm. Our results show that gentle touch rapidly activated somatosensory, motor, and cingulate regions within the first 100 ms of skin contact, which was driven by fast-conducting mechanoreceptive afference, and that these responses were sustained during touch. Peak latencies in the posterior insula were shifted as a function of stimulus location and temporally-separate posterior insula activations were induced by Aβ and CT afference that may modulate the emotional processing of gentle touch on hairy skin. We conclude that the detailed information regarding temporal and spatial brain activity from MEG provides new insights into the central processing of gentle, naturalistic touch, which is thought to underpin affective tactile interactions.


Introduction
Many studies have investigated brain activity to stimulation of the glabrous (non-hairy) skin of the hand, typically used in discriminative touch; however, less is known about the processing of biologicallymeaningful touch to the hairy skin, which is pertinent in the reception of touch and in affective interactions. Although it is a simple stimulus, a gentle stroke on the arm gives rise to a complex pattern of peripheral neural signals, where both discriminative and affective components are important in providing a complete sensory experience. Gentle touch on the hairy (non-glabrous) skin activates both fast-conducting Aβ mechanoreceptive afferents, with conduction velocities in the 35-80 m/s (Macefield et al., 1989;Kakuda, 1992) and slowly-conducting C-tactile (CT) afferents that conduct at~1 m/s (Vallbo et al., 1999;Watkins et al., 2017), thus although the timing of afference is different, it is pertinent to investigate these together.
Aβ afferents signal discriminative aspects of touch, i.e. exact and rapid changes in skin deformation (Vallbo and Johansson, 1984;Vallbo et al., 1995;Johnson, 2001;Abraira and Ginty, 2013), which is required for fast sensorimotor integration, whereas CTs are unlikely to play a key role in such processes due to their slow conduction velocity. CTs have never been found in human glabrous skin and the temporal delay before CT information reaches the central nervous system has implications for its processing. The preferential stimulation for evoking maximal frequency CT firing is a skin-temperature, caress-like stroke moving at a velocity between 1 and 10 cm/s over the skin. CT mean firing frequency correlates with perceived pleasantness of caress-like stimuli Ackerley et al., 2014b), suggesting that CTs are important in mediating positive affective touch (McGlone et al., 2014), e.g., the soothing touch between a parent and child, or in strengthening conspecific bonding.
Centrally, it is known that the posterior insula receives vast input from small diameter primary afferents (Aδ-and C-fibers) and processes affective information (Craig, 2002), whereas the majority of Aβ input appears to project to the S1 and S2 that process more discriminative aspects of touch (Ackerley et al., 2012;Sanchez-Panchuelo et al., 2012Case et al., 2016), which has led to a dichotomy in the processing of this information. The separation of these tactile pathways is useful when interpreting specific aspects of touch, yet it is clear that gentle touch on hairy skin activates both Aβ-and CT afferents and produces both affective and discriminative components of touch, thus tactile perception appears to involve the integration of both sensory pathways. Previous studies with functional magnetic resonance imaging (fMRI) have indicated that CTs project to the contralateral posterior insula, but not to primary (S1) or secondary (S2) somatosensory cortices (Olausson et al., 2002(Olausson et al., , 2008, but also that higher order somatosensory association areas (for review see: McGlone et al., 2014) are activated during gentle touch. However, the role of mechanoreceptive afferents in gentle touch and their subsequent activity in these brain regions is still unclear. The dichotomy between Aβ/somatosensory cortex activation and CT/posterior insula activation has been emphasized in the literature, but it is likely that the co-processing and integration of tactile information is more complex than this, especially as the posterior insula receives many different types of bodily information (Davis et al., 1998;Bj€ ornsdotter et al., 2009;Mazzola et al., 2012).
The majority of studies on the brain's response to gentle touch have used fMRI, which has excellent spatial resolution, but the blood-oxygenlevel-dependent (BOLD) response is too slow to detect rapid changes in neural activity (Hall et al., 2014). Magnetoencephalography (MEG), on the other hand, can track brain dynamics on a millisecond timescale and spatial resolution is generally 1 cm or better (Baillet, 2017). Hence, the present study aimed to define the location and progression of brain activity related to gentle touch on the hairy skin in humans, using MEG in combination with a new MEG-compatible tactile stimulator able to deliver naturalistic, caress-like stimuli. We hypothesized that S1 and S2 should be activated by fast-conducting Aβ afference first (all within 150 ms of skin stimulation (Wegner et al., 2000)), whereas the slowly-conducting CT-driven activity should occur later. Based on previous electroencephalography (EEG) work (Ackerley et al., 2013), we hypothesized that, due to the slow conduction velocity of CT-afferents, the subsequent brain activity would arise, in the posterior insula, with a latency of 600-700 ms when stimulating the hairy skin of the forearm, thus, we specifically aimed to clarify CT-driven activation of the insula. Furthermore, if the stimulation site were shifted more proximally to the upper arm, then the insular activity should show a corresponding reduction in latency.

Participants
A total of 21 healthy volunteers participated in the study. Six participants were excluded from final data analysis due to discontinuing the study (n ¼ 1), lack of an anatomical MRI-scan (n ¼ 2), or noisy MEGrecordings (n ¼ 3), hence the data presented in this study comes from 15 participants (6 males; 2 left handed; age 29.8 AE 5.7 (mean AE SD) years; range 22-45 years). A priori exclusion criteria were any previous or current psychiatric or neurological disorder and use of any psychoactive medication. The local ethics committee of Gothenburg, Sweden, approved the study. Experiments were carried out according to the declaration of Helsinki. All participants received written and oral information about the study before giving written informed consent to participate. Participants were reimbursed at 200 SEK per hour.

Stimuli
Naturalistic caress-like brush stroking with an approximate velocity of 3 cm/s was delivered to the left arm. This was conducted using a novel, custom-made, non-magnetic stimulator that stroked a soft brush over the arm, driven by pneumatic artificial muscles (PAMs; fluidic muscles by Festo AG & Co, Esslingen, Germany) mounted on a sturdy support frame, using glass-fiber reinforced composite tubes (Exelens S30; Exel Composites Oyj, Vantaa, Finland) and injection-molded plastic polypropylene connectors (S30 connectors; Exel Composites Oyj, Vantaa, Finland) (Fig. 1). The fluidic muscle contracted by 20% when compressed air at 6 bar was applied and relaxed when the compressed air was released. Two fluidic muscles were used to move the brush horizontally (muscle specification: DMSP-10-1000 AM-CM, diameter 10 mm, length of the contracting part 1000 mm, max travel 200 mm), and vertically (muscle specification: DMSP-10-140 AM-CM, diameter 10 mm, length of the contracting part 140 mm, max travel 28 mm). The compressed air was controlled by a pneumatic valve (SY5220-6LOU-01F-Q, SMC Corporation, Tokyo, Japan) and controlled by computer-generated TTL trigger pulses to lower the brush, stroke the skin, lift the brush, and pull the brush back to the initial position. The electronics and solenoid valves were placed outside the magnetically shielded room (MSR), and a 3.5-m pneumatic tube (internal diameter 2.5 mm) connected the pneumatic valve to the PAM.
The brush stimulator was tested for magnetic artefacts and the potential for it to produce non-tactile sensory responses (e.g. visual or auditory cues about the movement). First, we ran a stroking protocol similar to that in the main experiment, where the stimulator was located in position (i.e. approximately where it would be if a participant were present) and moved in the air. No MEG artefacts were found from this. Secondly, we tested the stimulator to see whether during the movement, non-tactile MEG responses were found by having the stimulator brushing on a phantom arm, while the participant had their left arm in their lap. The participant's view of the stimulator was obscured by a curtain and they wore earplugs. Again, we saw no correlated activity in the MEG that corresponded to the movement of the stimulator.
The brush stimulator was equipped with an accelerometer (ADXL335, Analog Devices Inc., Norwood, MA), two optoswitches (Omron E3X-DA41-S: Omron Corporation, Tokyo, Japan), and a load cell (Tedea-Huntleigh model 1004; VPG Corporate, Malvern, PA) to measure the acceleration, skin contact, velocity of the brush, and load on the skin, respectively. Similar approaches in MEG have been utilized previously (Jousm€ aki et al., 2007;Bourguignon et al., 2011;Piitulainen et al., 2013). The unique brush stimulator used in the present study was built at Aalto NeuroImaging, Aalto University, Finland.
The brush was a soft 5 cm wide paint brush, made out of goat's hair. A multifilament fiber-optic sensor was attached alongside the bristles of the brush, marking the timing of brush contact with the skin (t ¼ 0), and a load cell was used to measure the pressure applied on the skin. It is important to note that the trigger signal (i.e., t ¼ 0) was defined as the onset of brush contact with the skin. As such, contact between the brush and hair on the skin, and the resulting activation of hair mechanoreceptive afferents, occurred previous to this timepoint (i.e., for t < 0). Presentation ® (Neurobehavioral Systems, Berkley, CA) software was used to operate the brush robot with TTL-pulses.

Experimental design
Two conditions were performed. In the first condition (upper arm), brush strokes were delivered to the left lateral aspect of the upper arm proximal to the elbow. In the second condition (lower arm), brush strokes were delivered to the left dorsolateral aspect of the forearm proximal to the wrist. The upper arm condition was performed prior to the lower arm condition for all participants.
The robot was programmed to brush the skin for 1500 ms. However, the true brush contact duration showed slight variations between each participant due to differences in how the brush travelled over the skin (differences in anatomy, i.e. long vs. short arms and muscle mass affected how the brush was placed). The contact duration for all participants was 1500 AE 123 ms (mean AE SD), with no significant difference between upper and lower arm durations.
A block design containing random oddball trials was used in order to control participants' attention throughout the entire experiment. A total of 10 blocks per condition was performed. Each block consisted of 22 stimuli, with 20 strokes of 1500 ms duration, and two strokes of 600 ms duration, all at~3 cm/s. The 600 ms duration stimulus was the "oddball" trial and participants were instructed to press a button (which was registered as an event in the continuous MEG recording) with their right index finger whenever they felt a brush stroke that deviated from the other strokes. The oddball trials occurred randomly within each block. A total of 200 "long strokes" and 20 oddball strokes were delivered in each condition. The inter-stimulus interval (ISI) was set to 5060 ms þ a random duration of 1-1000 ms in order to avoid anticipatory effects.
Participants were dressed in non-magnetic clothes and seated comfortably in an upright position in the MEG, with their left arm supported by a table and soft cushioning. A curtain blocked the participants from seeing the movements of the brush robot and earplugs were used to suppress possible auditory contamination from the robot's pneumatics. Participants were instructed to rest their gaze on a point in front of them and keep their facial muscles as relaxed as possible throughout the experiment. They were asked to keep their head, shoulders, and body still, eyes open, and to blink as little as possible. Participants were asked to pay attention to the sensation of the brush strokes.

Measurements
Experiments were carried out in the NatMEG laboratory, Karolinska Institutet, Sweden (www.natmeg.se). The 306-channel MEG was recorded with an Elekta Neuromag TRIUX™ system (Elekta Oy, Helsinki, Finland) inside a magnetically shielded room (Model Ak3B, Vacuumschmelze GmbH & Co, Hanau, Germany). The sensor array of the Elekta MEG instrument encompasses 102 sensor elements each comprising one magnetometer and two planar gradiometers. The signals were sampled at 1.0 kHz with a band-pass filter set at 0.1-330 Hz.
A two-lead electrocardiogram (ECG) as well as vertical and horizontal electrooculogram (EOG) with bipolar surface electrodes were sampled for offline artifact rejection. The impedances for the EOG and ECG electrodes were checked to be below 10 kΩ. Continuous head position tracking was carried out using 4 head-position-indicator (HPI) coils. Individual anatomical landmarks (nasion, left and right preauricular points), a minimum of 200 head shape points, and HPI coil positions were digitized with a Polhemus Fastrak (Polhemus, Colchester, Vermont).
The distance (mean AE SD) from the two brush sites (the initial contact point for the brush on the skin on the upper and lower arm respectively) to the 7th cervical spinal process (C7) was 35 AE 3 cm (upper arm) and 61 AE 5 cm (lower arm), averaged over all participants. These measures were obtained in order to get an estimation of the distance that the peripheral nerves, innervating the stimulated skin areas, cover until they reach the spinal cord.
Individual T1-weighted anatomical MRI scans of the head were obtained for all participants included in the analysis.

Offline data processing
The MEG data and code are available on request from the authors. MaxFilter 2.2.10 (Elekta Oy, Helsinki, Finland) was used for spatiotemporal signal-space separation (tSSS: Taulu and Simola, 2006;Taulu and Hari, 2009) to reduce noise, and to allow head movement compensation in the raw continuous MEG-data. The MaxFilter parameters were set to 20 s buffer length and a correlation limit of 0.98. MEG-channels exhibiting noise during recordings were manually marked as bad before MaxFilter was applied.
Further data pre-processing and analysis of the evoked response were carried out in MNE Python (Gramfort et al., 2013(Gramfort et al., , 2014. Continuous data were low-pass filtered with a cut off frequency of 100 Hz, and notch filtered at 50 Hz to remove power-line noise. Subsequently, independent component analysis (ICA) was used on the continuous data to identify and remove eye movements and blinks. The data were then segmented into epochs of À1 s-4 s around the onset of the brush contact with the skin, i.e. the onset of the trigger signal from the fiber-optic sensor. Oddball trials were excluded from further analysis. Epochs with signals exceeding peak-to-peak amplitude of 4 pT/m for gradiometers and 4 pT for magnetometers were automatically rejected. The remaining epochs were visually inspected and those containing movement, muscle artifacts, or MEG sensor jumps were removed. The average number of epochs over all participants that were left after cleaning the data was 177 (range 151-192) for the upper arm condition and 168 (range 128-197) for the lower arm condition. The cleaned epochs were subsequently down-sampled to 250 Hz.
Event-related fields (ERF) and global field power (GFP) were calculated separately for magnetometers and gradiometers. GFP is defined as the absolute (magnetic field in the current study) magnitude of the standard deviation of the sensor signals over all sensors; in this case, the standard deviation is estimated across all magnetometer channels (Lehmann and Skrandies, 1980). The pre-stimulus interval between À1 and À0.2 s was used for baseline correction.

Source analysis
Segmentation of individual T1-weighted MRIs, followed by surface reconstruction, were carried out using the Freesurfer software package Fischl et al., 1999). For each participant, digitized head shape and anatomical landmarks were aligned to the individual MRI in order to position the head and MEG-sensors in a common coordinate system. A single layer boundary element method (BEM) model and a source space with 4098 grid points distributed on the cortical surface per hemisphere were used for computing a forward solution with free dipole orientations.
The depth-weighted minimum norm estimate (MNE; H€ am€ al€ ainen and Ilmoniemi, 1994; Lin et al., 2006) method was used for distributed source modeling of the evoked response. The depth-weighting coefficient was set to 0.8 and the parameter for loose dipole orientation was set to 0.2. Subsequently, dynamic statistical parametric maps (dSPM; Dale et al., 2000) were obtained by normalizing the estimated current amplitudes at each location to their respective standard error of the estimate. The noise covariance matrix used for calculating the inverse operator and noise-normalization was estimated from the non-averaged single trials using the pre-stimulus baseline from À1 to -0.2 s. The obtained dSPM values are F-distributed (3 and 200 degrees of freedom; Dale et al., 2000). The dSPM maps indicate the signal-to-noise ratio of the current estimate at each cortical location as a function of time. The obtained dSPM "movies" of cortical activity over time are useful for visualization of the data since they identify locations where MNE amplitudes are statistically above the noise level (Meeren et al., 2013).
For the grand average calculation of the evoked response, the source activity from each individual was morphed into a common source space, i.e. the average brain (fsaverage) provided by Freesurfer. The dSPM values were averaged across all participants and the grand average maps were used to identify spatio-temporal cortical patterns that were consistent across participants. For visualization, the grand averaged dSPM values were thresholded at the 96th percentile. The fsaverage brain surface was anatomically labeled using the Destrieux atlas (Destrieux et al., 2010). A list of labels that were active at the same time during both conditions was assembled (Table 1). A label was considered active if it contained dSPM values belonging to the 96th percentile of the data. The rationale for creating lists with overlapping activity for the two conditions was that it provides a list of brain areas that were repeatedly activated by the same stimulus type, i.e. gentle brushing on hairy skin, and not dependent on location (upper or lower arm).

Statistical analysis
To correct for multiple comparisons across sources and time points, we calculated thresholds for the dSPM maps using False Discovery Rate according to Benjamini/Yekutieli (Genovese et al., 2002), using the MNE-Python software.
Since MEG data are subject to considerable correlation both in the spatial and temporal domains, the multiple comparison problem is usually overestimated, thus we complemented parametric statistics with a non-parametric method based on clusterwise statistics (Maris and Oostenveld, 2007;Maris, 2012). Here, we used the test for paired contrasts with a spatial prior (i.e. the location of the sources in cortical space) implemented in the MNE-Python software package. Briefly, this test selects samples from the dSPM maps where the values exceed an a priori threshold of p < 0.01. Clusters are formed by significant samples that are adjacent in space and time, and the sum of the values within each cluster is used as a cluster-level statistic. The largest-valued cluster is selected as the test statistic and a reference distribution of the test statistic is obtained by permutation of the participants (i.e., a random permutation between baseline and stimulus), and re-calculation of the test statistic 1024 times. Finally, the null hypothesis (i.e., no difference between stimulus and baseline) is tested by comparing the test statistic from the observed data against the reference distribution at the p < 0.05 level. This test will correct for the multiple-comparison problem across time points, participants, and cortical sources (Maris and Oostenveld, 2007;Meeren et al., 2013). However, the test cannot be used for thresholding of the parametric dSPM maps. Instead, analysis of significant clusters provides a way to underpin the dSPM maps by confirming that observed activations correspond to significant clusters. This was done for the activations listed in Table 1.

Region of interest (ROI) analysis
Based on the hypothesis that CT afferents project to the contralateral posterior insular cortex (Olausson et al., 2002(Olausson et al., , 2008Bj€ ornsdotter et al., 2009), source time series from the contralateral (right) posterior insula were extracted. The anatomical label in the Destrieux atlas, which is most equivalent to the posterior insula, is denoted the "long insular gyrus and central sulcus of the insula", hence this label was used to obtain the activity time series for the posterior insula. The data within this ROI were collapsed into a single activity time course by calculating the average of all vertex (dipole) time series.
For visual comparison, the source time series from a post hoc selection of a region neighboring the posterior insula were also obtained. The selected region is denoted the "posterior ramus of the lateral sulcus" in the Destrieux atlas, which is an area that separates the posterior insula from the inferior parietal lobule, and merges superiorly with the parietal operculum. In the following text and figures, this region will be referred to as the S2 cortex.

Evoked responses at the sensor level
Evoked responses to the~3 cm/s gentle brush strokes to the left forearm and upper arm were calculated with reference to the time when the brush first contacted the skin. Fig. 2 contains the temporal evolution of the grand averages of global field power (GFP; Fig. 2A and B) for the upper arm and lower arm conditions respectively, and GFP plots for individual participants (Fig. 2C and D). The grand averages were very similar for the two conditions with four main peaks: two at short latency after the onset of the stimuli and two at short latency after the offset. The sustained field during the ongoing stimulus contained some between-participant variability and the response to the offset of the stimulus was sharper in some participants. The reason for the offset response being variable was likely due to the epochs being centered around the onset of the stimuli and not the offset, making the offset response non-time-locked over epochs.  Table 1 show early activation of the S1 and S2 areas with widespread activity within the first 100 ms. This early activity included  D) indicate where there is a significant deviation from baseline (p < 0.05, FDR-corrected). Note that the brush makes contact with arm hair before the skin, resulting in activations before the trigger onset (i.e., 0 ms) due to this input from fast-conducting hair afferents. responses in the postcentral sulcus, bilateral insula, and midline cortical areas (the superior frontal gyrus, and the midcingulate cortex (MCC)). Activity was sustained throughout the brush stroke and persisted after the end of the stroke. Fig. 3 shows whole head images of the dSPM values from the source analysis for four time instances separately and for the two stimulus conditions: upper arm (Fig. 3A) and lower arm (Fig. 3B) stroking. False discovery rate correction for multiple comparisons (using p < 0.05), yielded thresholds for these maps at dSPM-values of 3.02 for the lower arm condition, and 2.52 for upper arm. In addition, brain activity was analyzed using a non-parametric clusterwise statistical test, and it confirmed that activations from the dSPM maps in Fig. 3 belonged to significant clusters. Table 1 summarizes the overlapping activity that was present in the grand average at the same time instances for both conditions. We next report the activations in more detail for the lateral and medial aspects of the cortical surface.

Evoked responses in source space
Starting with the lateral aspects of the brain, Fig. 3 and Table 1 show contralateral activity in the arm region of the central sulcus, and the superior part of the precentral sulcus, in addition to activity in the operculo-insular region at 50 ms. At 100 ms, the operculo-insular activity was enhanced and both conditions exhibited activity in the contralateral postcentral sulcus. Ipsilateral activity in the central sulcus, operculoinsular region, and the superior part of the precentral sulcus was also present at 100 ms. At the end of the brush stimulation, around 1500 ms, there was some enhancement of the contralateral activity, mainly in frontal motor areas, that was still present at 2200 ms post stimulus onset, i.e. 700 ms after the end of brush contact with the skin.
For cortical areas in the medial aspect of the brain, Fig. 3 and Table 1 show that there was bilateral medial activity in the superior frontal gyrus (SFG) and the anterior part of the midcingulate cortex (aMCC) at 50 ms. In the contralateral hemisphere, there was also activity in the posterior part of the midcingulate cortex (pMCC) and the marginal branch of the cingulate sulcus, whereas, on the ipsilateral side, there was additional activity in the paracentral lobule. At 100 ms, the bilateral medial hemispheres exhibited a widespread activity pattern, involving the SFG, aMCC, pMCC, posterior-dorsal part of the cingulate gyrus (dPCC), the paracentral lobule, the precuneus, and the subparietal sulcus. The activity in SFG and MCC was sustained throughout the stimulation and after the offset of brush contact with the skin.

Activity in posterior insula and S2
To test the hypothesis that the activity in the posterior insula should occur at a time latency corresponding to the slow conduction velocity of CT afferents, we estimated the time it took for peripheral CT afferent activity to reach the central nervous system. To this end, measurements of the distance from the stimulation sites on the upper and lower arm to the 7th cervical spinal level were carried out. The peripheral conduction velocity of CTs is~1 m/s (Vallbo et al., 1999;Watkins et al., 2017). This means that in this sample, based on the average distance between the brushing sites on the arm and the spinal cord, any brain activity related to CT signaling should occur after at least 350 ms for upper arm condition and 610 ms for the lower arm condition. Fig. 4A shows the right hemisphere of the inflated fsaverage brain, with the Destrieux labels of the S2 (blue) and posterior insula (green) highlighted. Fig. 4B and C shows the extracted time series from these labels in Fig. 4A.
The activity in contralateral (right) S2 (Fig. 4B) showed sharp peaks at short latency after the onset and the offset of brush skin contact. The high amplitude sharp waveforms indicated that the S2 response was time-and phase-locked to the fast-conducting Aβ afferent input. The peaks in the contralateral (right) posterior insula (Fig. 4C) did not show such sharp and high amplitude responses to the rapidly changing events as in S2, although it followed a similar pattern of activity. For both conditions, it is clear that the onset of activity in the posterior insula was driven by the onset of the contact between the brush and the skin; however, the largest amplitude peaks did not occur with short latency after the onset or end of stimulation. For the upper arm condition, there were two peaks at long latency, at 420 ms and 780 ms. For the lower arm condition, there were two long-latency peaks at 630 ms and 740 ms after the onset of stimulation. These peaks in the insula did not coincide with the strong medial activations depicted in Fig. 3A and B. The activation of the aMCC peaked later, at 1340 ms for the upper, and 1520 ms for the lower arm condition, and there were no peaks in the 400-800 ms range.

Discussion
We used MEG and a custom-built brush robot to elucidate the brain responses to naturalistic, caress-like strokes on the lower and upper arms of healthy humans. We show that this gentle touch rapidly evoked activity (within~100 ms) over a network of regions, including somatosensory, motor, and cingulate regions. Our findings emphasize the importance of considering the Aβ afference when interpreting the implications of gentle affective touch, especially in neuroimaging studies that have typically focused on the activity from 'CT optimal' touch. Although Aβ input may define the more discriminative aspects of such touch (e.g., when, where), we show that the initial activation of posterior insula and cingulate regions is driven by Aβ afference. Thus, it is reasonable to conclude that Aβ input also underpins affective touch in correspondence with the CT input. While the role of Aβ afference in positive affective touch is somewhat recognized (L€ oken et al., 2011;Ackerley et al., 2014a), we elucidate the potential areas activated via such signals. Moreover, we indicate that separate temporal posterior insula activations are induced by both Aβ and CT afference, which in turn may have a modulating effect on the emotional processing of gentle touch on the hairy skin.
We identified a number of important spatio-temporally distinct responses during gentle touch on the hairy skin. Firstly,~50 ms after the brush made initial contact with the skin, there was an activation mainly of the contralateral central sulcus, operculum and the superior part of the precentral sulcus, as well as bilateral medial areas. The activity was widespread at 100 ms such that it also involved the ipsilateral central sulcus, bilateral postcentral sulci, frontal motor areas, operculo-insular areas, and MCC. Due to its short latency, this activity was driven solely by Aβ afferents, which evidently activate the posterior insula. The continued activity during the brush stroke until the end of the stimulation involved the same areas that were already active at 100 ms; however, the signal strength fluctuated over time. Activity in frontal motor and cingulate cortices persisted beyond the end of the brush stroke.
The activity fluctuations over time were particularly evident in contralateral S2 and posterior insula. There were clear peaks in response to the onset and offset of the stimulation, likely due to Aβ input. Furthermore, there was also a sustained response between these peaks, i.e. during the ongoing stimulation and after the brush stroke had ended, especially in the insula (cf. Ackerley et al., 2013). The sustained insular response between the onset/offset peaks was likely due to CT afference, as there was a shift in the evoked waves between the upper and lower arm conditions, which corresponded with our predicted timecourses. These mid-brushing peaks were not as well defined as those relating to the Aβ input, which is likely due to the much larger range of conduction velocities in the population of C compared to A afferents (Erlanger and Gasser, 1930), which results in a considerably longer duration of the barrage of CT afference reaching the brain. The response was furthermore weaker in the second, i.e. lower arm, condition: this may be due to participant fatigue. A previous EEG study on gentle brush strokes found a frontal midline ultra-late event-related potential that was attributed to CT input (Ackerley et al., 2013). We did not find an exact corresponding ultra-late event-related field in the current MEG study, potentially due to different spatial sensitivities of EEG and MEG (e.g. superficial vs. deep sources, response orientation; Ahlfors et al., 2010); however, our MEG responses indicate a temporal correlation between CT input and posterior insula activity, which match the latency of response that would be predicted based on the slower CT conduction time.
The posterior insula has been implicated in the processing of CT afference (Olausson et al., 2002(Olausson et al., , 2008Bj€ ornsdotter et al., 2009;Morrison et al., 2011aMorrison et al., , 2011bMorrison, 2016). However, only indirect evidence relates its activity with CT input and studies have shown that early insula responses from touch on hairy skin correspond to Aβ afference, akin to the early-onset responses we present herein. Grandi and Gerbella (2016) showed that the monkey posterior insula contains neurons that are speed-selective to a caress delivered on hairy skin, but their response was near-instantaneous. A MEG study investigating the insular-opercular response to innocuous transient mechanical stimulation of the forearm in humans found a sharp deflection of activity in the posterior insulã 100 ms post-stimulation (Hayamizu et al., 2016). Others have found contralateral posterior insula activity to saltatory moving tactile stimulation of the human glabrous hand (Oh et al., 2017) and that the insula was also activated during neutral, unpleasant, and pleasant, moving tactile stimuli of the lower leg (Hua et al., 2008).
Taking the results of our present study together with previous findings on tactile responses in the posterior insula, it is likely that this region plays a general role in coding moving tactile stimuli (cf. Morrison et al., 2011b). We show the importance of having precise temporal information when understanding such responses, as the posterior insula was clearly activated by both Aβ and CT afference, although at different times, and this has implications for the interpretation of many neuroimaging studies of gentle touch where Aβ and CT afference may not be distinguishable and have typically been considered separately.
Previous attempts to make a distinction between Aβ and CT afference have been based on the their relative activations to slow and fast stroking (Morrison et al., 2011a(Morrison et al., , 2011bMcGlone et al., 2012;Perini et al., 2015). The onset and offset responses in the insula driven by Aβs presented herein indicate that further precise investigations in the temporal domain are needed to clarify this. For example, while presently not achievable, a selective and complete block of Aβ or CT-afferent signals would make it possible to further distinguish between their contributions to brain responses.
Concerning other areas, we report activations in the superior precentral sulcus and the superior frontal sulcus and gyrus, both parts of Brodmann area (BA) 6. The lateral portion of BA6 is the premotor cortex (PMC), whereas the medial side constitutes the supplementary motor area (SMA). Anterior to BA6, the superior frontal sulcus and gyrus extend Fig. 4. Region of interest (ROI) analysis of the activity in the contralateral (right) S2 and posterior insula. A, labels used for the ROI analysis: Blue, Destrieux atlas label "posterior ramus of the lateral sulcus" (here called S2); green, Destrieux "long insular gyrus and central sulcus of the insula" (here called posterior insula). B, C, time series of the activity in S2 and posterior insula; peak latencies are marked with vertical dashed lines. Graphs on the left hand side belong to the upper arm condition and graphs on the right hand side belong to the lower arm condition. Note that the amplitude scales differ for the two conditions. As in Fig. 2, the brush makes contact with arm hair before the skin, resulting in activations before the trigger onset (i.e., 0 ms) due to this input from fast-conducting hair afferents. Gray horizontal bars indicate the duration of the brush movement over the skin. Dashed horizontal lines indicate significant deviation from baseline (p < 0.05, FDR-corrected). D, right (contralateral) hemisphere images at the latencies for the peaks in posterior insular activity. The colorbar is the same as in Fig. 3. into BA8. These areas were activated bilaterally, with an early onset of activity that was sustained during the stimulation and beyond the end of the stimulation. Activity in BA6 and 8 has been reported in fMRI studies investigating gentle touch vs. rest (Olausson et al., 2002;Bj€ ornsdotter et al., 2009;Gordon et al., 2013;Davidovic et al., 2016;Sailer et al., 2016). Furthermore, Olausson et al. (2002) report activity in BA6 in response to gentle brush strokes in patient G.L., who lacks Aβ afferents but has an intact CT system. This suggests that Aβ afferents are not necessary for recruiting activity in motor-related brain areas in a passive task.
The SMA is activated in the initiation and inhibition of motor output (Nachev et al., 2008). Monkey studies indicated the SMA is involved in touch processing when the tactile stimulation dictates a motor response (Romo et al., 1993;Hern andez et al., 2002). Our participants were instructed to attend to the stimuli and press a button when there was a deviant (oddball) in order to maintain their attention throughout the experiment. The use of the oddball task could explain the activity in BA6 and 8 in our study. However, these areas are activated even when participants do not have a task involving a motor response (Olausson et al., 2002;Bj€ ornsdotter et al., 2009;Gordon et al., 2013;Davidovic et al., 2016;Sailer et al., 2016). Thus, these more classical motor areas seem to encode parts of moving touch, implying that they may be activated during general cutaneous movement, and not just from muscle activity (Aimonetti et al., 2007). The co-processing of these signals may then aid in the accurate interpretation of afferent feedback in kinesthesia and bodily awareness.
We also report activity in the anterior and posterior mid-cingulate cortex (aMCC, pMCC). In fMRI, gentle touch to the arm activated the dorsal anterior cingulate cortex (dACC) (Gordon et al., 2013) and activity in the pregenual anterior cingulate cortex (pgACC) was seen during forearm massage (Lindgren et al., 2012). However, it should be noted that it has been recently suggested that aMCC is a better label for the brain region previously designated as dACC (Cieslik et al., 2015;Vogt, 2016) and we have chosen to use this nomenclature. The aMCC is activated by painful stimuli and itch, but is also engaged in reward processes (Hadland et al., 2003;Vogt, 2016). It is active, along with the anterior insula, in a wide range of empathic responses during observation of pleasant or unpleasant scenes (Bernhardt and Singer, 2012). Thus, it has been suggested that the anterior insula and aMCC contribute both to subjective/emotional experiences and to adaptive responses in actual and predicted body states. The pMCC is, on the other hand, involved in the processing and regulation of pain (Vogt, 2016) and bodily awareness (Vogt and Laureys, 2005). Further studies will be necessary to define which of these possibilities provide the better explanation for the activation of the MCC during gentle touch to the arm.

Conclusions
Our MEG study demonstrates that gentle touch on the hairy skin of the arm evokes activity in a set of well-defined but widespread regions, including S1 and S2, bilateral insula, motor, premotor, and cingulate areas. The posterior insula activity was driven by gentle touch that activated Aβ and CT afferents together, thus further work should be conducted into the parallel processing of both types of input, and hence, we argue against a clear dichotomy between these afferent systems in positive affective touch.

Declaration of interest
None.

Funding sources
This work was supported by grants to J.W. from the Swedish Research Council (grant 2017-01717), the Knut and Alice Wallenberg Foundation (project NeuroSQUID, 2014.0102), and the Sahlgrenska University Hospital (ALFGBG grant 3161). The NatMEG facility at Karolinska Institutet, Sweden, is supported by Knut & Alice Wallenberg (grant 2011.0207) and the Swedish Research Council (grant D0589201). This study was also supported by the Karolinska Institutet Strategic Neuroscience Program, StratNeuro.