Connectivity estimation

This example pipeline, ConnectivityEstimation.pyp, shows how NeuroPype can be used to extract the activity of multiple brain regions from an EEG stream, and how these activity estimates can be used to estimate the directed flow of information between the regions, in real time, using Granger-causal approaches. The result is shown on a connectivity plot and on a heatmap, and sent out over LSL.

The pipeline can be divided into six parts, which we discuss in the following: data acquisition and preprocessing, source estimation, connectivity modeling, postprocessing, plots, and output over LSL.

This pipeline is designed for a live data stream, for which it estimates the connectivity in real time. The pipeline plays back a recording over LSL, so that it can be run without a device (see below).

Data Acquisition and Preprocessing

Includes: LSL Input, Assign Channel Locations, Remove Unlocalized Channels, FIR Filter, Decimate and Re-referencing

The pipeline receives its EEG data from the lab streaming layer through the LSL Input node, whose query property is set to name='MyOutputStream'. The three nodes at the bottom of the pipeline (Import File, Stream Data and LSL Output, which are not shown in the screenshots) play back a recording, VisualArrest.xdf, over LSL under that name, so that the example can be run without a live data source. For use with live data, remove them and change the query property of LSL Input so that it finds the stream of your device. See Data input in the example pipelines for more on this setup.

The preprocessing nodes are:

• Assign Channel Locations and Remove Unlocalized Channels: devices often send out non-EEG channels, such as triggers or accelerometers. These two nodes assign each channel a 3D location on the scalp, based on standard montages (channels in the international 10-20 system, or its finer-grained 10-5 variant, work out of the box), and then remove the channels that are not found in a montage (i.e., the non-EEG channels). This is important so that they do not throw off the source estimation. If you have other kinds of montages, you may need a custom montage file (see the documentation of the node).
• The FIR Filter (titled FIR Bandpass) is a band-pass filter (frequencies of 0.5, 1, 45 and 50 Hz), which removes slow drifts as well as high-frequency activity such as power line noise.

• Decimate is a very fast way of downsampling the data (here by a factor of 2), to increase the processing speed of the pipeline. If you have plenty of CPU power, you can drop it (but if your EEG system has a very high sampling rate, such as 2 kHz or more, you should keep it). It does not use an anti-aliasing filter, which is why it follows the low-pass of the FIR Filter.
• Re-referencing is set up to do common average referencing (its reference_range is :, all channels). This is crucial for source estimation, since the head models of NeuroPype assume that the data is common-average referenced.

Source Estimation

Includes: sLORETA, ROI Activations and Merge ROIs

• The sLORETA node calculates the activity for each source location of the head model (its headmodel property, here the standard model based on the Colin27 brain, Colin-339ch-4495v-scalar-4shell-1.1.hm), which are the over 4000 vertices of a mesh that follows the brain surface. The output is a time series with one channel for each vertex. Its regu property adjusts the smoothness of the reconstruction. You could use the LCMV Beamforming node instead of sLORETA, as described in the region of interest spectral power example, which also describes the properties of sLORETA in more detail.
• The ROI Activations node computes the average activity in the regions of interest (ROIs) defined by a brain atlas, which are a parcellation of the brain surface into named anatomical regions. The head model contains one atlas, the Desikan-Killiany atlas, with 467 regions. With the default settings, the node outputs a time series with one channel for each region.
• We want the activity of only a few, larger regions of the brain, to reduce the dimensionality and make it tractable to compute the information flow between them. So the Merge ROIs node, which is a Combine Channels node (its operation is average), merges the ROIs of the atlas into nine new regions, named AG_R, OC_L, OC_R, PoMC, MPFC, IPCS, AG_L, DLPFC and SMA. Its mapping property lists the ROIs that are averaged for each new region, for instance AG_R, which is the average of inferiorparietal R.5, inferiorparietal R.7, inferiorparietal R.8 and inferiorparietal R.11. If you had only a few target regions, you could use Select Regions of Interest and Mean nodes instead (see the region of interest spectral power example), but with this many regions, the wiring would be cumbersome. When editing the mapping, it can help to copy it to a text editor and break the lines to see what is merged into what. The names of the available ROIs are listed in a file with the same name as the head model, ending in -atlas-labels.txt, in the resources/headmodels folder of NeuroPype. To learn where each region is located, see the companion article of the Desikan-Killiany atlas.

Connectivity Modeling

Includes: Group Lasso MVAR and Direct Directed Transfer Function (dDTF)

Network connectivity is estimated in two stages. First, a (time-varying) multivariate autoregressive (MVAR) model is estimated, which describes the dynamics of the system of the nine regions, and then the directed information flow between the regions is calculated from this model.

• The Group Lasso MVAR node estimates the MVAR model over a sliding window of the data (window_length of 0.5 seconds and a window_step of 0.1 seconds). It assumes that at any time, not every region strongly drives every other one, but that the connectivity is "sparse". It is fast enough for real-time use as long as there are not too many regions, or too high a sampling rate or model order (model_order of 10). The main property to adjust is the degree of sparsity (the regularization strength, lambda_reg, 0.2 here). When applying the pipeline to a different set of regions or another EEG system, check that the parameters are well-chosen by reviewing the output: the estimated connectivity should be neither non-existent (too sparse) nor have everything connected to everything else (not sparse enough). Other properties govern the complexity of the model (model order), the amount of data it is estimated from (window length), the fitting procedure (rho and alpha) and the stopping criteria (max_iter and the tolerances), and trade off running time against accuracy. See the documentation of the node.
• The Direct Directed Transfer Function (dDTF) node computes a measure of information flow from the MVAR model. Other measures are found in the Connectivity category of the nodes (for instance partial directed coherence). Its frequencies property selects the frequencies for which it is computed (here 76), and with temporal_normalization enabled, the short-time dDTF (sdDTF) is computed, which is one of the highest-quality measures. The output, for every time step and every frequency, is a square matrix that describes the degree of information flow from each region to every other one: a four-way tensor indexed by time, frequency, region (source) and region (target).

Postprocessing

Includes: Moving Average, Select Range and Mean

The connectivity solution is a rich representation, which can be used in many ways. Here, we prepare it for visualization. A Moving Average node smooths it over time (30 samples, order and order_unit), which gives a more slowly changing estimate. The Select Range node selects the frequencies (on the frequency axis) that are averaged, and the Mean node (titled Avg Frequency) averages over the frequency axis. Since the selection of the Select Range node is :, which selects everything, the result is the average connectivity across all frequencies. If you set it to a range, for instance 8:12, you get the average connectivity within that frequency band only.

Plots

Includes: Connectivity Plot, Operate Constant and Heatmap Plot

The Connectivity Plot node shows the connectivity between the regions, as links drawn between them on a 3D plot of the brain (title of Network Connectivity, view of left), with the 50 strongest links (max_links), and the size of each region's node showing its outflow (node_measure).

A Heatmap Plot shows the same information as an image (a matrix plot), with the source regions on one axis (the x_axis, space.source), and the target regions on the other (the y_axis, space.target), which shows the flow of information between the regions. The Operate Constant node before it multiplies the values by a constant (100000, the constant property) to get them into the color range of the plot (the color_limits, 0 to 100).

Output over LSL

Includes: Fold Into Axis and LSL Output

The connectivity can be sent out over LSL for another application. Since LSL transmits a two-dimensional matrix only, the Vectorize Channels node (a Fold Into Axis node) merges the two space axes of the connectivity (the source and the target regions) into a single space axis, and the result is sent by an LSL Output node as a stream called GB ROI Connectivity (with a stream_type of CONN).