Calculating Intersubject Correlations

Author

Jason Geller

Intersubject correlation (ISC) is a measure used to assess the degree of synchrony among participants—that is, how similarly their responses change over time. ISC has been widely applied in fMRI and EEG research, particularly in naturalistic contexts, to quantify shared neural activity. More recently, ISC has been extended to webcam-based eye-tracking, allowing researchers to examine shared gaze patterns during dynamic stimuli.

This vignette demonstrates how to calculate ISC using two commonly used approaches:

Pairwise ISC Approach

In this method, each participant’s data is correlated with every other participant’s data, generating a correlation matrix across all pairs.
Correlations are transformed using Fisher’s z to r transformation and averaged and then transformed back into r resulting in one ISC for each participant.

Leave-One-Out ISC Approach - Here, each participant’s data is correlated with the mean time series of all other participants, excluding their own data. - This results in one correlation value for each participant.

To compute ISC, your data should be structured as a p × t matrix, where:

p represents the number of participants
t represents time points

Before highlighting the function calculate_isc() I want to break down each step that goes into calculating the ISC.

Pairwise ISC

Below, we generate a 10 × 1000 matrix

library(tidyverse)
install.packages('webgazeR', repos = c('https://jgeller112.r-universe.dev', 'https://cloud.r-project.org'))
library(webgazeR)

set.seed(123)

data_matrix <- matrix(rnorm(10 * 1000), nrow = 1000, ncol = 10)

data_matrix %>% 
  head()

            [,1]        [,2]       [,3]        [,4]       [,5]       [,6]
[1,] -0.56047565 -0.99579872 -0.5116037 -0.15030748  0.1965498 -0.4941739
[2,] -0.23017749 -1.03995504  0.2369379 -0.32775713  0.6501132  1.1275935
[3,]  1.55870831 -0.01798024 -0.5415892 -1.44816529  0.6710042 -1.1469495
[4,]  0.07050839 -0.13217513  1.2192276 -0.69728458 -1.2841578  1.4810186
[5,]  0.12928774 -2.54934277  0.1741359  2.59849023 -2.0261096  0.9161912
[6,]  1.71506499  1.04057346 -0.6152683 -0.03741501  2.2053261  0.3351310
           [,7]       [,8]       [,9]      [,10]
[1,] -0.6992281 -1.6180367  0.5110004  1.9315759
[2,]  0.9964515  0.3791812  1.8079928 -0.6164747
[3,] -0.6927454  1.9022505 -1.7026150 -0.5625675
[4,] -0.1034830  0.6018743  0.2874488 -0.9899631
[5,]  0.6038661  1.7323497 -0.2691142  2.7312276
[6,] -0.6080450 -0.1468696 -0.3795247 -0.7216662

We use the cor() function to compute pairwise correlations between each participant’s time-binned data and that of every other participant.

We specify use = “pairwise.complete.obs” to ensure that:

If any missing values (NA) exist for a participant, those time points are excluded from the correlation for all comparisons involving that participant.
This allows us to retain all participants in the dataset while ensuring that correlations are computed using only available data.

correlation_matrix <- cor(data_matrix, use = "pairwise.complete.obs")

Since the diagonal of the correlation matrix represents self-correlations (i.e., each participant correlated with themselves), we replace those values with NA, as they are not meaningful for intersubject correlation (ISC):

diag(correlation_matrix) <- NA

correlation_matrix

              [,1]         [,2]         [,3]          [,4]          [,5]
 [1,]           NA  0.086479441 -0.019329537 -0.0029947104  0.0322973377
 [2,]  0.086479441           NA  0.026503334 -0.0070290755  0.0924646624
 [3,] -0.019329537  0.026503334           NA  0.0505608499 -0.0150195347
 [4,] -0.002994710 -0.007029076  0.050560850            NA -0.0004005193
 [5,]  0.032297338  0.092464662 -0.015019535 -0.0004005193            NA
 [6,] -0.018267991 -0.055292983 -0.033003346  0.0438457353  0.0138042312
 [7,] -0.034340849  0.024316458  0.005282974 -0.0200855749 -0.0090534146
 [8,] -0.030156049  0.031077180 -0.025946806 -0.0055278737 -0.0777495779
 [9,] -0.009370688  0.022155213  0.058417535 -0.0064371559  0.0470798891
[10,]  0.022344953  0.022481121  0.029148322  0.0459001272 -0.0047210318
               [,6]         [,7]         [,8]          [,9]        [,10]
 [1,] -0.0182679912 -0.034340849 -0.030156049 -0.0093706879  0.022344953
 [2,] -0.0552929834  0.024316458  0.031077180  0.0221552128  0.022481121
 [3,] -0.0330033458  0.005282974 -0.025946806  0.0584175347  0.029148322
 [4,]  0.0438457353 -0.020085575 -0.005527874 -0.0064371559  0.045900127
 [5,]  0.0138042312 -0.009053415 -0.077749578  0.0470798891 -0.004721032
 [6,]            NA -0.016497452 -0.016411441 -0.0006214018 -0.010177844
 [7,] -0.0164974516           NA  0.037173511  0.0791431861  0.002523005
 [8,] -0.0164114415  0.037173511           NA -0.0422433679  0.010503101
 [9,] -0.0006214018  0.079143186 -0.042243368            NA -0.006174781
[10,] -0.0101778437  0.002523005  0.010503101 -0.0061747810           NA

Next, we apply Fisher’s Z-transformation using atanh() to convert the correlation coefficients into Z-scores before averaging.

Averaging correlation scores directly is bad practice because:

Correlations are bounded between -1 and 1, which skews their distribution.
Fisher’s Z-transformation normalizes correlation values, making them more suitable for averaging.

You can also choose to use the median instead (if this is the case you do not need to apply Fisher transformation)

z_values <- atanh(correlation_matrix)  # Fisher's Z transformation

Once transformed, we compute the average ISC per participant and transform back to r.

isc_values <- tanh(rowMeans(z_values, na.rm = TRUE))  # Convert back to r

the calculate_isc() function with method=“pairwise” does all this for us.

isc_pairwise <- calculate_isc(data_matrix, method = "pairwise", summary_statistic = "mean")

isc_pairwise

 [1]  0.002985107  0.027061018  0.008523871  0.010880947  0.008761466
 [6] -0.010296001  0.007625601 -0.013271741  0.015797990  0.012429909

Leave one out ISC

In the leave-one-out approach, we compute the correlation between each participant’s time series and the average time series of all other participants, excluding their own data.

First, we extract the time series for a single participant:

# Extract the full time series for the current subject
subject_time_series <- data_matrix[, 1]

and calculate the mean of the other subjects data points, excluding the current subject

leave_one_out_mean <- rowMeans(data_matrix[, -1, drop = FALSE], na.rm = TRUE)

We then compute the correlation between the participant’s time series and the leave-one-out mean:

leave_one_out_isc <- cor(subject_time_series, 
                                          leave_one_out_mean, 
                                          use = "pairwise.complete.obs")

Since this needs to be done for every participant, we use the calculate_isc() function to automate the process:

isc_pairwise <- calculate_isc(data_matrix, method = "leave-one-out")

isc_pairwise

 [1]  0.008891827  0.080509746  0.025626567  0.031021518  0.026597776
 [6] -0.029482014  0.023592539 -0.038321505  0.046317774  0.036036027