Comparing representational geometries using whitened unbiased-distance-matrix similarity β part 4/6
Cleaned-text code points [46745, 62698). Bibliography omitted here; source.md includes it.
The optimal method of course always depends on the data and models that need to be compared (Figure 8). The first decision is whether the models are meant to predict the dissimilarities quantitatively (ratio scale) or only their ranks (ordinal scale). Quantitative predictions can often be derived if we have an explicit model of the shape of the underlying activity profiles. The distribution of activity profiles may also be predicted from activities in an artificial neural network model [ 20 ] , directly from perceptual judgements [ 25 ] , or the statistics of external training data [ 9 ] . In other cases, the model may only predict the rank ordering of the dissimilarities, but not by how much one dissimilarity is larger than another. In such cases, rank-correlations are most appropriate [ 8 ] . While this approach can be statistically less powerful [ 19 ] , it is robust against any possible monotonic transformation of the dissimilarities.
The next decision is whether the activity patterns can be estimated independently and with approximately the same variance across all partitions. If we cannot be sure that this is the case, then crossvalidated, unbiased distance estimates should be used. This is important, because the bias on the standard Euclidean or Mahalanobis distances will be structured if the noise is not i.i.d., such that the model comparison will be biased. Even in situations in which the measured activity-pattern estimates can be assumed to be i.i.d., the unbiased estimation approach can be more powerful than using the biased estimates and Pearsonβs correlation. This is because the meaningful zero point (which indicates that there is no pattern difference) can help distinguish models. Which approach is better depends on the number of partitions, the signal-to-noise ratio, the experimental conditions, and the structure of the models (Fig. 4). Overall, however, the increased robustness of violations of noise assumptions will generally outweigh the cost of increased variance, especially if the number of partitions is large.
Whether biased or unbiased distance estimates are used, RDMs should always be compared using whitened RDM correlations or cosine similarities. These measures perform often substantially better, but never worse than standard approaches. Because we can approximate the true covariance structure well using the covariance structure under the assumption that there is no true signal, the approach can be implemented in a computationally efficient manner. Thus, we do not see any reason to use standard Pearson correlation or cosine similarity for RDM comparison.
### 3.3 Probabilistic versions of RSA
An alternative to using unbiased distance estimates is to rely on an explicit model of the measurement noise. This approach is used in PCM [ 15 , 14 ] , and two other related Bayesian formulations of RSA [ 26 , 10 ] . All three approaches rely on estimating the likelihood of the data under the RDM model, using the assumption of Gaussian signal and noise. While this method can be used to derive maximum-likelihood estimates of the distances [ 10 ] , these are by construction not unbiased. Rather, all three approaches rely on a valid estimate of the trial-by-trial ( πΊ K \bm{\Sigma}_{K} ) and spatial ( πΊ P \bm{\Sigma}_{P} ) noise covariance for valid model inference. Here we show that, if these estimates are valid, PCM provides a theoretical upper bound for the performance of any RSA method. However, if the trial-by-trial noise covariance is not specified correctly (as in the example of Fig. 4a), model selection bias would arise, similar to the one shown for biased RDM correlations (see Fig 4a). In cases in which there is a high uncertainty about ( πΊ K \bm{\Sigma}_{K} ), the use of WUC may therefore be preferable to PCM.
### 3.4 Influence of spatial noise correlations
While we have extensively discussed the influence of the trial-by-trial covariance ( πΊ K \bm{\Sigma}_{K} ) on the distance estimates, many neuroscientific datasets also have substantial spatial noise correlations between channels ( πΊ P \bm{\Sigma}_{P} ). This is especially true for fMRI data, where strong correlations between voxels exist. For this reason, it is common practice to attenuate these correlations by spatial pre-whitening activity estimates before multivariate analysis [ 6 ] . However, given that our estimate of the spatial noise structure is not perfect and often needs to be regularized, we will be left with some residual noise correlations between voxels.
Eq. 29 shows that such correlations simply inflate the (co-)variance of the distance estimates under the null-hypothesis by a constant factor, which is related to the "effective" number of voxels, tr β ( πΊ P ) 2 / tr β‘ ( πΊ P β πΊ P ) \mathrm{tr}(\bm{\Sigma}_{P})^{2}/\mathrm{tr}(\bm{\Sigma}_{P}\bm{\Sigma}_{P}) . Because subsequent RSA inference, for example a t-test on distance estimates across subjects, takes the variability of distances into account directly, the inflation does not affect the validity of such inferences. In Appendix 6.3, we provide a simulation with spatially correlated noise. The results confirm the advantage of using whitened RDM correlations or RDM cosine similarities generalizes to situations with substantial noise correlations between channels.
In contrast, the likelihood estimates from PCM directly rely on valid spatial covariance estimates. This is not a problem if we simply use the these estimates for model comparison using frequentist tests [ 15 ] , as again, the increased variability of the likelihood estimates will be directly taken into account. However, a direct interpretation of the likelihoods, for example for a subsequent Bayesian group analysis [ 27 ] , needs to be done with some caution.
Interestingly, the simulation (see 6.3) also shows that for very high spatial noise correlations, model comparison using WUC is slightly superior to model comparison using PCM. Thus, the likelihood-based test of PCM outperforms RSA only if the spatial noise structure is relatively close to the assumption of independence.
### 3.5 Electrophysiological and other non-Gaussian data
While the derivation of the bias (Eq. 4 ) and the unbiasedness of the crossvalidated distance estimator (Eq. 24 ) holds for any data distribution, the derivation of the covariance matrix of distances (Eq. 29 ) depends on the assumption of a Gaussian distribution of the noise.
One important application, for which this assumption is clearly violated, are neural spiking data obtained from single- or multi-unit electrophysiological recordings. Estimates of instantaneous firing rates show a Poisson-like distribution across repeated measurements, with increasing firing variability for higher rates and a rightward skew. For such data, we recommend to take the square-root of the instantaneous firing rates before calculating distances [ 28 ] . The square-root transform makes the distribution more symmetric the variance approximately constant across firing rates.
To test the approach suggested here against modest violations of the Gaussian noise assumption, we repeated the simulations shown in Fig. 7 using noise distributions with either a modest skew ( Ο 2 \chi^{2} -distribution with df=6) or a modest kurtosis (t-distribution with df=6). Overall, we were able to replicate the increasing advantage of the advantage of the WUC with increasing number of conditions. Thus, model comparison using WUC appears to be robust against modest violations of the Gaussian noise assumption, as we would find in real fMRI or electrophysiological data.
### 3.6 Summary
Taken together, we believe that the WUC provides an important new measure for RDM comparison that should become standard for many applications. The new measure extends the linear centered kernel alignment (CKA) [ 23 , 12 ] and RV coefficient [ 11 ] by removing the biasing influence of measurement noise.
The WUC and whitened RDM Pearson correlation have been implemented in a new Python-based RSA toolbox, released by the team of authors ( https://github.com/rsagroup/rsatoolbox , [ 16 ] ). We hope that the results presented in this paper, together with the accessible implementation, will accelerate the adoption of what we consider to be current best practice in RSA.
## 4 Methods
### 4.1 Extended definitions
To derive the mean and full variance-covariance matrix of the distance estimates, it is useful to make some more general definitions and assumptions. We assume that each measured activity profile (column of π ^ m \hat{\mathbf{B}}_{m} ) has covariance πΊ K \bm{\Sigma}_{K} between conditions. For fMRI, this correlation structure is caused by the sluggish nature of the hemodynamic response, as well from the low-frequency noise inherent to the measurements. A reasonable estimate of πΊ K \bm{\Sigma}_{K} can be derived directly from the first level linear model [ 29 ] . For other modalities, it may be more reasonable to assume independence of measurements.
We also assume that each measured activity pattern (row of π ^ m \hat{\mathbf{B}}_{m} ) has a covariance of πΊ P \bm{\Sigma}_{P} between channels. Variability of fMRI, EEG, MEG measurements clearly shows substantial spatial structure. Again, this noise structure can usually be estimated from the residuals of the first-level linear model (see Mahalanobis distance). To remove the redundancy of πΊ K \bm{\Sigma}_{K} and πΊ P \bm{\Sigma}_{P} in terms of the overall scaling of the noise, we restrict trace β‘ ( πΊ P ) = P \mathrm{trace}(\bm{\Sigma}_{P})=P .
For the derivation of the full covariance matrix of distance estimates, we need to make the slightly more restrictive assumption that π ^ m \hat{\mathbf{B}}_{m} has a matrix-normal distribution across partitions m m . While this assumption is reasonable for fMRI data, it is recommended to apply the square-root transform to neuronal spiking data to make it conform to the normal assumption [ 28 ] .
To derive distances between conditions in a matrix notation, we define a D Γ K D\times K contrast matrix π \mathbf{C} . The k t β h k^{th} row of this matrix contains a 1 1 and a β 1 -1 for the two conditions that are contrasted in the k t β h k^{th} distance, all other entries are 0 0 . The product ππ \mathbf{CB} then results in a D Γ P D\times P matrix that contains the pattern differences πΉ k = π i β π j \bm{\delta}_{k}=\mathbf{b}_{i}-\mathbf{b}_{j} in its rows. We define:
π« β‘ ππ β π T β π T / P \bm{\Delta}\equiv\mathbf{C}\mathbf{B}\,\mathbf{B}^{T}\mathbf{C}^{T}/P (22)
The diagonal of π« \bm{\Delta} contains the squared distances d k d_{k} (divisively normalized by the number of channels). On the basis of πΊ K \bm{\Sigma}_{K} , we can also define the D Γ D D\times D variance-covariance matrix of the pattern-difference estimates ( π β π ^ m \mathbf{C}\hat{\mathbf{B}}_{m} ):
π΅ β‘ Var β‘ ( π β π ^ m ) = π β πΊ K β π T . \bm{\Xi}\equiv\mathrm{Var}(\mathbf{C}\,\hat{\mathbf{B}}_{m})=\mathbf{C}\bm{\Sigma}_{K}\mathbf{C}^{T}. (23)
β’ Size is given in number of rows Γ \times number of columns. For consistency of notation, vectors are defined to be in either row or column orientation.
### 4.2 Bias of distance estimates
Eq. 4 can be derived by expressing the estimated pattern difference ( πΉ ^ k , m \hat{\bm{\delta}}_{k,m} ) as the sum of the true pattern-difference vectors πΉ k \bm{\delta}_{k} and the measurement noise ( π k , m \bm{\xi}_{k,m} ). By substituting this into Eq. 3 and taking the expected value ( E E ), it is straightforward to show that the distance estimator is positively biased:
E β‘ ( d k ^ ) = E β‘ ( 1 M 2 β β m β n ( πΉ k + π k , m ) β ( πΉ k + π k , n ) T / P ) = ( πΉ k β πΉ k T + E β‘ ( β m π k , m β π k , m T / M ) ) / P = d k + π΅ k β k / M . \begin{split}\mathrm{E}(\hat{d_{k}})&=\mathrm{E}\left(\frac{1}{M^{2}}\sum_{m}\sum_{n}(\bm{\delta}_{k}+\bm{\xi}_{k,m})(\bm{\delta}_{k}+\bm{\xi}_{k,n})^{T}/P\right)\\ &=\left(\bm{\delta}_{k}\bm{\delta}_{k}^{T}+\mathrm{E}\left(\sum_{m}\bm{\xi}_{k,m}\bm{\xi}_{k,m}^{T}/M\right)\right)/P\\ &=d_{k}+\bm{\Xi}_{kk}/M.\end{split} (24)
The bias arises by multiplying the noise with itself (i.e. for the cases of m = n m=n ). For all other cases ( m β n m\neq n ), the noise terms are independent, and the expected value of their product is zero.
### 4.3 Variance of distance estimates
An analytical expression for the variance-covariance matrix of the vector of distance estimates can be derived using the following general result (see Appendix B1, B2 for details). If the matrix π \mathbf{A} has matrix normal distribution β³ β π© β ( π , π΅ , πΊ ) \mathcal{MN}(\mathbf{B},\bm{\Xi},\bm{\Sigma}) , then the diagonal of ππ T \mathbf{A}\mathbf{A}^{T} has the expected value and variance:
E β‘ ( diag β‘ ( ππ T ) ) \displaystyle\rm E(\rm diag(\mathbf{A}\mathbf{A}^{T})) = \displaystyle= diag β‘ ( ππ T + tr β‘ ( πΊ ) β π΅ ) \displaystyle\rm diag(\mathbf{B}\mathbf{B}^{T}+\rm tr(\bm{\Sigma})\bm{\Xi}) (25) Var β‘ ( diag β‘ ( ππ T ) ) \displaystyle\rm Var(\rm diag(\mathbf{A}\mathbf{A}^{T})) = \displaystyle= 4 β π β πΊ β π T β π΅ + 2 β t β r β ( πΊ β πΊ ) β ( π΅ β π΅ ) , \displaystyle 4\mathbf{B}\bm{\Sigma}\mathbf{B}^{T}\circ\bm{\Xi}+2\rm tr(\bm{\Sigma\Sigma})(\bm{\Xi}\circ\bm{\Xi}), (26) where β \circ is the element-by-element multiplication of two matrices, and tr \mathrm{tr} () is the trace of a matrix. When setting π \mathbf{A} to the mean of the pattern differences across partitions, we can easily derive the variance of the biased distance estimate (Eq. 3 ).
Var β‘ ( π ^ ) \displaystyle\mathrm{Var}(\hat{\mathbf{d}}) = \displaystyle= 1 P 2 β ( 2 β tr β ( πΊ P β πΊ P ) M 2 β π΅ β π΅ + 4 β P M β π« β β π΅ ) \displaystyle\frac{1}{P^{2}}\left(\frac{2\mathrm{tr}(\bm{\Sigma}_{P}\bm{\Sigma}_{P})}{M^{2}}\bm{\Xi}\circ\bm{\Xi}+\frac{4P}{M}\bm{\Delta}^{*}\circ\bm{\Xi}\right) (27) π« β \displaystyle\bm{\Delta}^{*} = \displaystyle= ππ β πΊ P β π T β π T \displaystyle\mathbf{CB}\bm{\Sigma}_{P}\mathbf{B}^{T}\mathbf{C}^{T} (28)
The expression for the variance of the unbiased estimate of the distance (Eq. 5 ) can be derived by taking the covariance of all pairs of partitions into account (see Appendix B3).
Var β‘ ( π ~ ) = 1 P 2 β ( 2 β t β r β ( πΊ P β πΊ P ) M β‘ ( M β 1 ) β π΅ β π΅ + 4 β P M β π« β β π΅ ) \mathrm{Var}(\tilde{\mathbf{d}})=\frac{1}{P^{2}}\left(\frac{2\mathrm{tr}(\bm{\Sigma}_{P}\bm{\Sigma}_{P})}{M(M-1)}\bm{\Xi}\circ\bm{\Xi}+\frac{4P}{M}\bm{\Delta}^{*}\circ\bm{\Xi}\right) (29)
Intuitively the variances of distance estimates come from the product of signal and noise (averaged over M M partitions) and product of noise with noise (averaged over M 2 M^{2} pairs of partitions for the biased distance estimate and over M β‘ ( M β 1 ) M(M-1) pairs of partitions for the unbiased estimate).
### 4.4 Spatial pre-whitening and Mahalanobis distances
In the result section, we focus on biased and unbiased estimates of the Euclidean distance. Previous work [ 6 ] , however, demonstrates clearly that taking into account the spatial covariance structure of fMRI noise ( πΊ P \bm{\Sigma}_{P} ) can increase the reliability of distance estimates.
In the simplest case, we ignore the correlation between voxels and simply divide the activity estimates for each voxel by the square root of the diagonal elements of πΊ P \bm{\Sigma}_{P} . This step already prevents noisy voxels to influence the distance estimate overly much. Additionally, we can use multivariate pre-whitening, i.e. post-multiplication of πΊ P β 1 / 2 \bm{\Sigma}_{P}^{-1/2} . This step gives less weight to the information contained in two voxels that are highly correlated in their random variability than to information contained in two uncorrelated voxels. Calculating Euclidean distances on multivariate pre-whitened data is equivalent to calculating a Mahalanobis distance.
In practice, we do not have access to the voxel-by-voxel covariance matrix. However, we can use the residuals π m \mathbf{R}_{m} from the first-level general linear model to derive an estimate,