[ENH] Add conn.plot_circle to Dynamic Connectivity Example (and Estimation of Significance?)

[alexrockhill]

I think it would be nice to actually see the connectivity structure in the dynamic connectivity example. I'm happy to make a quick PR to do plot_circle and make it look nice.

I think a deeper question is whether the connections are strong enough not to be likely to be observed by chance, which I think is a question people doing this analysis might be interested as the main result from this analysis. Perhaps, since there is a decent number of epochs, a cross-validation approach would be reasonable. Perhaps with a bootstrap approach of shuffling the data in time to create a null distribution (maybe also preserving autocorrelations by just sliding the channels around relative to each other as in this paper https://pubmed.ncbi.nlm.nih.gov/34717794/). Curious to know what others think, especially @adam2392. Maybe this is too much computation for an example also...

[adam2392]

I think the plot_circle addition would be nice if you can add it. However, I would also add a stipulation that this really doesn't imply "connectivity" but is simply fitting a linear system model to the data. It "prolly can" reproduce dynamics, etc., which then make analysis of the system properties interesting.

This brings up your second point. Note that our entire discussion assumes that the connections/dynamics are time-invariant, so we could potentially use the multiple epochs to answer the question. I think in general interpreting the connections as significant is up to the user's stated assumptions. E.g. in a fully Granger-Causal model, one might set a threshold of "0.05" and fit the model, and then anything above that implies there is a connection and anything below is not; you just assume that the model fit is all you needed to do to test for connections. On the other hand, maybe they employ a bootstrap approach as you mentioned, which is also testing the same thing, but here you're trying to estimate a null distribution of your data to then test for connections. I guess end-goal is the same, but means are different. Just like someone doing t-test vs doing Wilcoxon rank-sum vs bootstrapping.

The challenge with this second point is unfortunately, it is rather difficult to convey to users, hence the "Note" on the GH README heh.

So tldr: In favor of:

• adding the plot_circle to show it off and add some text interpreting it.
• possibly even adding a function for time-series window shuffling to allow bootstrapping with explicit warnings and notes in the docstring mentioning what this is "intended to be used for"

[alexrockhill]

That all sounds super reasonable but the one thing that wasn't clear was how you would approach estimating significant connections in this case where the question is about resting-state connectivity for the given ECoG grid. It sounds like maybe the bootstrapping approach?

[adam2392]

Yep I would say doing the bootstrap is fine, but needs to come with a warning/caveat because it's a common mistake that ppl don't fully specify their assumptions when doing this type of statistical approach :)

Are you attempting to add a PR?

[alexrockhill]

I haven't yet, but I could give it a shot. I wanted to synchronize on this issue so that we're in agreement of the recommendation for best practice. The goal here (and my goal for other projects hopefully if this works 😄) is to take in epochs as in the dynamic connectivity example and assign significance values to a model of connectivity across those epochs. I'd prefer not to assume the connectivity is quasi-static/stationary but I think a typical use-case assumption is that the connectivity should be the same across epochs. As I'm writing this, I'm thinking maybe this should be a separate example with task-based data because that better satisfies the second assumption. Dynamic connectivity for resting state seems a bit odd, am I missing something as to why other methods of connectivity don't suffice in this case?

What do you think? Is that a reasonable example use-case? And, do you think this is the best connectivity framework for this type of example? Looking at the other examples, something that seems to be lacking is statistical quantification of the connectivity results (i.e. the output is all very nice plots but essentially qualitative) so maybe the real need is for an example that uses any one of these or maybe multiple and quantifies the probability that they would be observed by chance?

[adam2392]

Re resting state: Depends I guess :p. Dynamics might change on a slow time-scale (i.e. completely time-varying), or you assume that dynamics are changing but connectivity is the same (okay so the values of the connectivity matrix change but the "true non-zero" edges are the same throughout), or you assume that the dynamics are not changing and the connectivity is not changing.

Tbh it kind of is subjective at this point, hence it's not as straightforward :(

Re your second paragraph: statistical quantification also for "graphs" is much harder compared to iid unfortunately too. Hence, I don't think we can solve those issues on our end alone, but some ppl I know work on it at https://github.com/microsoft/graspologic, so I think we can make API wraps to them for doing statistical testing.

[alexrockhill]

It makes the most sense to me to assume that there are regions that are and are not anatomically connected, which would eventually be verified with DWI, and that the dynamics are changing. Maybe as a first pass, this could be approximated by ignoring the lack of anatomical connections and assuming that they will be non-significant based on any statistics such that you can assume all-to-all connectivity and just model the dynamics.

Can we have some somewhat reasonable bootstrap approach that avoids graph quantification and just evaluates each connection on its own? I think clearly this would be suboptimal as the graph information incorporates reasonable assumptions that would make the statistical quantification better (like how assuming spatiotemporal clusters makes time and time-frequency stats tests much better) but a first pass could use a simpler FDR-type correction, and possibly then only these significant connections could be made into a graph. Clearly the graph assumption, like the clusters, is something that would be really nice to have though... Maybe it wouldn't be too hard to add grasplogic as a dependency (it looks well-maintained) and just go for it right away. What do you think?

[adam2392]

Sure, I think the time-series shuffling approach is reasonable and pretty common enough that it could be implemented and used.

Down to add graspologic as soft-dependency as well for wraps for certain statistical tests.

[alexrockhill]

Okay, I have a few other things to wrap up, then I'll give a new example a shot.

[alexrockhill]

I started working on it but one main thing is that conn = vector_auto_regression(...) returns an object that when you go to call conn.plot_circle you get

<ipython-input-61-812c3494b704> in <module>
----> 1 conn.plot_circle()

~/projects/mne-connectivity/mne_connectivity/base.py in plot_circle(self, **kwargs)
    775     @copy_function_doc_to_method_doc(plot_connectivity_circle)
    776     def plot_circle(self, **kwargs):
--> 777         plot_connectivity_circle(
    778             self.get_data(),
    779             node_names=self.names,

~/projects/mne-connectivity/mne_connectivity/viz/circle.py in plot_connectivity_circle(con, node_names, indices, n_lines, node_angles, node_width, node_height, node_colors, facecolor, textcolor, node_edgecolor, linewidth, colormap, vmin, vmax, colorbar, title, colorbar_size, colorbar_pos, fontsize_title, fontsize_names, fontsize_colorbar, padding, ax, fig, subplot, interactive, node_linewidth, show)
    142             ax = plt.subplot(*subplot, polar=True)
    143 
--> 144     return _plot_connectivity_circle(
    145         con=con, node_names=node_names, indices=indices, n_lines=n_lines,
    146         node_angles=node_angles, node_width=node_width,

~/projects/mne-python/mne/viz/circle.py in _plot_connectivity_circle(con, node_names, indices, n_lines, node_angles, node_width, node_height, node_colors, facecolor, textcolor, node_edgecolor, linewidth, colormap, vmin, vmax, colorbar, title, colorbar_size, colorbar_pos, fontsize_title, fontsize_names, fontsize_colorbar, padding, ax, interactive, node_linewidth, show)
    185         con = con[indices]
    186     else:
--> 187         raise ValueError('con has to be 1D or a square matrix')
    188 
    189     # get the colormap

ValueError: con has to be 1D or a square matrix

which should at least be more informative, but ideally would animate some kind of connectivity plot over time. I think it would be reasonable to start with the animation in an example and then integrate it into the source code if that's okay with you @adam2392.

[adam2392]

Ah yeah so plot_circle will only work rn for a 2D connectivity array (i.e. static). However, yeah I agree that either i) the error can be more informative, or ii) we could just "average" over time(?), and iii) yeah possibly produce some animation over time instead.

If you're up for PRing, I'm up for reviewing :)

[alexrockhill]

Ok, I think I got it working besides some colorbar proliferation bugs but it's complete chaos...

@adam2392, any thoughts on how to make sense of the data?

EDIT: This only went up to 84 epochs, my bad, should go all the way to the 302 in the example.

Here's the full one: