Best way to handle signal-dependent RDM covariance matrix?

[johnmarktaylor91]

Hi all--so, some RDM comparators (whitened cosine, whitened correlation) involve whitening the RDMs using the RDM covariance matrix. This matrix is the sum of two terms: a signal-independent term, and a signal-dependent term (see attached PDF: Heiko has a writeup from awhile ago deriving the formulas for these). Currently, rsatoolbox only uses the signal-independent term.

It seems like we might at least want to have the option of letting the user include the signal-dependent term as well, no? What I'm wondering is, what would be the best way of handling this in terms of the user interface? The signal-dependent term requires quantities such as 1) the estimated true (whitened) patterns U, and 2) the post-whitening voxel covariance matrix sigma_r. It seems possibly confusing to add these as extra arguments to compare_rdm (since they are only relevant in this specific case), but I'm not sure what the best method might be. No rush for this, it's just on my mind as I'm working on a methods paper that uses the full formula for the covariance matrix.

Cov_of_Crossnobis_RDMs-3 (2).pdf

[JasperVanDenBosch]

Could we make it descriptor? And the relevant method then takes it from there? It could have some naming convention.

[johnmarktaylor91]

That seems very sensible to me.

[HeikoSchuett]

Hi @johnmarktaylor91
I think, the main task for implementation is a function that actually returns the full signal dependent covariance estimate.

In the current implementation sigma_k is passed through already as another argument to the compare function. I think consistency would be great here, i.e. we should either make both descriptors or pass the new one through additionally as well.

I think the really clean version needs to modify the RDMs object such that the matrices are correctly handled for subset and subsample functions. This would be a bit of work, but would make all interactions with bootstrapping and crossvaligdation work correctly.

best
Heiko

[johnmarktaylor91]

I think, the main task for implementation is a function that actually returns the full signal dependent covariance estimate.

I agree—I have a sloppy version for my own tinkering but haven't done it cleanly yet.

In the current implementation sigma_k is passed through already as another argument to the compare function. I think consistency would be great here, i.e. we should either make both descriptors or pass the new one through additionally as well.

Making them descriptors might be a bit easier from the user's standpoint, right? The pre-whitened patterns (which are needed for the distance covariance estimate) could be computed as part of calc_rdm and tacked on as a descriptor, rather than the user computing them in a separate function call.

One quick q as long as we're on this topic—in practice, the sigma_r in your formula for the distance covariance (i.e., the leftover voxelwise covariance that remains after whitening) will just be the identity matrix, right? Since in real data we don't have access to the ground truth, and we've already made our best estimate of the voxelwise covariance when we do the whitening.

I think the really clean version needs to modify the RDMs object such that the matrices are correctly handled for subset and subsample functions. This would be a bit of work, but would make all interactions with bootstrapping and crossvaligdation work correctly.

Good point and I will keep this in mind if I end up doing the code changes for this.

One final q: Niko and I are working on a new method, "framed RSA" that involves inserting the all-zeros and all-c patterns (where c is some constant). These require special handling, since they are noiseless. The method is looking quite promising (improved ROI and DNN layer identification performance) and I suspect we will want to add it to rsatoolbox. Any initial thoughts on what might be the best way of handling these patterns in terms of the plumbing? One thought is to add an argument to calc_rdm to include these patterns, and give them their own dedicated label in each of the pattern_descriptors, so that downstream functions know what they are. No rush on this if it requires more discussion, just curious if there's an obvious answer.

[JasperVanDenBosch]

Adding framed RSA to the toolbox would be great of course. Wrt to the inputs, I'd have a slight preference to keep it more object-oriented to prevent the list of arguments becoming too long, maybe these could be properties of a FramedRSADataset class or optional properties of the base Dataset? The calc methods, if doing framed RSA, would then know where to find these on the input objects.

[jdiedrichsen]

Hi @johnmarktaylor91 :

Signal-depend part of the covariance: Thanks for you note - Before starting, I am a bit confused though. Is your formula 40 not identical to Eq 27 in the Diedrichsen et al. 2022? OR am I missing something? If yes, then you need to tell me what the difference in the simulation is...
For using the signal-dependent part of the covariance matrix for model comparison, I think it is not consistent to use the estimated U's in this expression, as their variance will substantially overestimate the true signal (i.e. it will contain noise). My guess in you simulations you used the true U's? Instead, in our previous paper (Diedrichsen & Kriegeskorte, 2017), we used the signal predicted under the current model. This allows you to estimate a likelihood of the distances under the current model p(D|M), which is useful for subsequent Bayesian model comparison, but it requires you to iteratively update the covariance matrix with each iteration (IRLS-style). Ultimately, for fMRI data, the signal-dependent part tends to be very small, so that we did not see much practical utility in taking the signal dependence into account. Do you have a compelling example where it would make much of a difference?

Frame RSA: This is an interesting idea - we have been using something quite similar for some of our papers. Actually, I don't think that there is anything special that needs to be done to include rest. Resting baseline is of course not noiseless, but it is also measured with noise. If all measures are expressed relative to rest, they carry this variance in them. You can undo this, by adding a row of zeros and then re-centering the activity estimates of each voxel / channel. Sigma_k is then rank-deficient, but that is also the case when you just estimate the condition-wise distances. Adding an explicit rest simply means that you retain the information that is normally lost going from second-moment to distances.
Dealing with the all-one direction is of course more tricky. as distances (and subsequent MDS's) depend on the constant that you use for scaling this vector. In the past, I have not included the all-one into the model fit, but simply projected an all-one unit vector into the MDS and then shown it for visualization. Do you have an application where the all-one would be used in model fit or alignment between RDMs?

Happy to discuss over the break if you are interested.

Joern