diff --git a/.gitignore b/.gitignore index 5f05b47..c872d38 100644 --- a/.gitignore +++ b/.gitignore @@ -1,5 +1,8 @@ *to*_Atlas_v2/ +HArtMuT_mix_Colin27_small.mat +HArtMuT_NYhead_small.mat mni2atlas/ sa_nyhead.mat sa_nyhead_sereegareduced.mat -*.asv \ No newline at end of file +*.asv +*.afdesign \ No newline at end of file diff --git a/docs/README.md b/docs/README.md index 9502e3d..57b24b7 100644 --- a/docs/README.md +++ b/docs/README.md @@ -18,24 +18,27 @@ Reference: - [Installation](#installation) - [Configuration](#configuration) - [Obtaining a lead field](#obtaining-a-lead-field) - - [Using a leadfield generated in Brainstorm](#using-a-leadfield-generated-in-brainstorm) + - [Using a leadfield generated in Brainstorm](#using-a-leadfield-generated-in-brainstorm) - [Picking a source location](#picking-a-source-location) + - [Using atlases](#using-atlases) - [Orienting a source dipole](#orienting-a-source-dipole) - [Defining a source activation signal](#defining-a-source-activation-signal) - - [Brain components and scalp data](#brain-components-and-scalp-data) + - [Components and scalp data](#components-and-scalp-data) - [Finalising a dataset](#finalising-a-dataset) - [Variability](#variability) - [Noise and multiple signal classes](#noise-and-multiple-signal-classes) - [A note on components and sources](#a-note-on-components-and-sources) - - [Event-related spectral perturbations](#event-related-spectral-perturbations) - - [Frequency burst for event-related synchronization](#frequency-burst-for-event-related-synchronization) - - [Inverse frequency burst for event-related desynchronization](#inverse-frequency-burst-for-event-related-desynchronization) - - [Amplitude modulation](#amplitude-modulation) - - [Pre-generated data as activation signal](#pre-generated-data-as-activation-signal) - - [Autoregressive models](#autoregressive-models) + - [More signal types](#more-signal-types) + - [Event-related spectral perturbations](#event-related-spectral-perturbations) + - [Frequency burst for event-related synchronization](#frequency-burst-for-event-related-synchronization) + - [Inverse frequency burst for event-related desynchronization](#inverse-frequency-burst-for-event-related-desynchronization) + - [Amplitude modulation](#amplitude-modulation) + - [Pre-generated data as activation signal](#pre-generated-data-as-activation-signal) + - [Autoregressive models](#autoregressive-models) - [Random class generation and multiple components](#random-class-generation-and-multiple-components) - [Source identification](#source-identification) - [Component variability](#component-variability) + - [Relativity of the lead field units](#relativity-of-the-lead-field-units) - [Using the GUI](#using-the-gui) - [Extending SEREEGA](#extending-sereega) - [Lead fields](#lead-fields) @@ -45,15 +48,14 @@ Reference: - [Anonymous class creation function](#anonymous-class-creation-function) - [Connectivity benchmarking framework](#connectivity-benchmarking-framework) - [Simulating data with a specific signal-to-noise ratio](#simulating-data-with-a-specific-signal-to-noise-ratio) -- [Relativity of the Units of Measure](#relativity-of-the-units-of-measure) - [Contact](#contact) - [Special thanks and acknowledgements](#special-thanks-and-acknowledgements) ## Introduction -SEREEGA is an open-source MATLAB-based toolbox to generate simulated, event-related EEG toy data. Starting with a forward model obtained from a head model or pre-generated lead field, dipolar _brain components_ can be defined. Each component has a specified position and orientation in the brain. Different activation signals can be assigned to these components. EEG data is simulated by projecting all activation signals from all components onto the scalp and summing them together. +SEREEGA is an open-source MATLAB-based toolbox to generate simulated, event-related EEG toy data. Starting with a forward model obtained from a head model or pre-generated lead field, dipolar _components_ can be defined. Each component has a specified position and orientation in the head model. Different activation signals can be assigned to these components. EEG data is simulated by projecting all activation signals from all components onto the scalp and summing them together. -SEREEGA is modular in that different head models and lead fields are supported, as well as different activation signals. Currently, SEREEGA supports the [New York Head (ICBM-NY)](https://www.parralab.org/nyhead) pre-generated lead field, the [Pediatric Head Atlases](https://www.pedeheadmod.net) pre-generated lead fields, and [FieldTrip](http://www.fieldtriptoolbox.org) custom lead field generation. Four types of activation signals are provided: _ERP_ (event-related potential, i.e. systematic deflections in the time domain), _ERSP_ (event-related spectral perturbation, i.e. systematic modulations of oscillatory activations), _noise_ (stochastic processes in the time domain), and signals based on an _autoregressive model_. A final _data_ class allows the inclusion of any already existing time series as an activation signal. +SEREEGA is modular in that different head models and lead fields are supported, as well as different activation signals. Currently, SEREEGA supports the [New York Head (ICBM-NY)](https://www.parralab.org/nyhead) pre-generated lead field, the [Pediatric Head Atlases](https://www.pedeheadmod.net) pre-generated lead fields, the [HArtMuT (Head Artefact Model using Tripoles)](https://github.com/harmening/HArtMuT) pre-generated lead fields, and [FieldTrip](http://www.fieldtriptoolbox.org) custom lead field generation. Four types of activation signals are provided: _ERP_ (event-related potential, i.e. systematic deflections in the time domain), _ERSP_ (event-related spectral perturbation, i.e. systematic modulations of oscillatory activations), _noise_ (stochastic processes in the time domain), and signals based on an _autoregressive model_. A final _data_ class allows the inclusion of any already existing time series as an activation signal. SEREEGA is intended as a tool to generate data with a known ground truth in order to evaluate neuroscientific and signal processing methods, e.g. blind source separation, source localisation, connectivity measures, brain-computer interface classifier accuracy, derivative EEG measures, et cetera. @@ -99,9 +101,9 @@ Download [SEREEGA](https://github.com/lrkrol/SEREEGA) to your computer and add a It is recommended to have MATLAB 2014b or higher with the DSP and Parallel Computing toolboxes installed. SEREEGA's core functionality requires R2013b. This is primarily due to the use of `addParameter` rather than `addParamValue`; if your MATLAB does not support `addParameter`, try exchanging all those references with `addParamValue`. That should restore basic functionality. Some plotting functions rely on the `boundary` function introduced in R2014b. Some signal generation functions depend on the DSP toolbox version 8.6 (R2014a) or higher. Data simulation tries to use parallel processing from the Parallel Computing Toolbox where available. -[EEGLAB](https://sccn.ucsd.edu/eeglab) is used for a number of functions, and should be started before generating a lead field as EEGLAB's `readlocs` is used to add channel location information. SEREEGA was tested with EEGLAB 13.6.5b. +[EEGLAB](https://sccn.ucsd.edu/eeglab) is used for a number of functions, and should be started before generating a lead field as EEGLAB's `readlocs` is used to add channel location information. SEREEGA was tested with EEGLAB 13.6.5b and 2021.1. -When using the New York Head, as in this tutorial, make sure you have the [New York Head (ICBM-NY) lead field in MATLAB format](http://www.parralab.org/nyhead) in your path. Similarly, the [Pediatric Head Atlases](https://www.pedeheadmod.net) must also first be obtained separately, if you intend to use those. When using FieldTrip to generate a custom lead field, the file `/fileio/ft_read_sens` from the FieldTrip directory will be necessary. FieldTrip can be installed as an EEGLAB plug-in. +When using the New York Head, as in this tutorial, make sure you have the [New York Head (ICBM-NY) lead field in MATLAB format](http://www.parralab.org/nyhead) in your path. Similarly, the [Pediatric Head Atlases](https://www.pedeheadmod.net) and [HArtMuT head models](https://github.com/harmening/HArtMuT) must also first be obtained separately, if you intend to use those. When using FieldTrip to generate a custom lead field, the file `/fileio/ft_read_sens` from the FieldTrip directory will be necessary. FieldTrip can be installed as an EEGLAB plug-in. To obtain a lead field from [Brainstorm](https://neuroimage.usc.edu/brainstorm/), Brainstorm needs to be in MATLAB's path. When some functionality in SEREEGA relies on external files which you don't have, a warning will appear with instructions on where to find them. Reasons for these files not having been included out of the box include file size and license conflicts. ### Configuration @@ -117,7 +119,7 @@ epochs.srate = 1000; % their sampling rate in Hz epochs.length = 1000; % their length in ms ``` -Additional, general parameters can be added here as well, such as `prestim` to indicate a pre-stimulus time period and shift the x-axis, or `marker` to indicate the event marker accompanying each epoch at t=0. +Additional, general parameters can be added here as well, such as `prestim` to indicate a pre-stimulus time period (in ms) and shift the x-axis accordingly, or `marker` to indicate the event marker accompanying each epoch at t=0. ### Obtaining a lead field @@ -136,6 +138,8 @@ lf = lf_generate_fromnyhead('montage', 'S64'); lf2 = lf_generate_fromnyhead('labels', {'Fz', 'Cz', 'Pz', 'Oz'}) ``` +When not indicating any montage or labels, most functions will default to using all available electrode positions. + When generating the lead field using FieldTrip, we can also indicate the resolution of the source model (i.e. the density of the source grid): ```matlab @@ -148,61 +152,61 @@ We can inspect the channel locations in the obtained lead field and their relati plot_headmodel(lf); plot_headmodel(lf, 'style', 'boundary', 'labels', 0); ``` + + #### Using a leadfield generated in Brainstorm -Brainstorm allows you to create personalised leadfield (aka headmodels in Brainstorm). The `lf_generate_frombrainstorm` enables you to convert a leadfield generated in Brainstorm to a format usable by SEREEGA. +Most lead fields supported by SEREEGA are more or less pre-packaged. A separate software called [Brainstorm](https://neuroimage.usc.edu/brainstorm/) allows you to create personalised lead fields (aka headmodels in Brainstorm) based on MRI scans. The `lf_generate_frombrainstorm` enables you to convert a lead field generated in Brainstorm to a format usable by SEREEGA. + +To use this function, Brainstorm needs to be installed and in MATLAB's path. + +To use your own, personalised lead fields, you need to provide three file paths (plus one optional) to the following files generated by Brainstorm: +* The headmodel *mat* file generated in Brainstorm. +* The channel location *mat* file used for the generation of the headmodel in Brainstorm. +* The structural MRI *mat* file to which the electrodes have been coregistered. +* (Optional) the atlas file if used. + +A description of how to obtain these files is included in the help section of `lf_generate_frombrainstorm`. The command then becomes: ```matlab -% Run with default files +% use personal files +lf = lf_generate_frombrainstorm( ... + 'headmodel', 'PATH_TO\my_headmodel.mat', ... + 'chanloc', 'PATH_TO\my_chanloc.mat', ... + 't1', 'PATH_TO\my_T1_MRI.mat', ... + 'atlas', 'PATH_TO\my_atlas.mat'); +``` + +The function also comes with a precomputed head model and associated files, which it uses by default if no other files are supplied, i.e., when you run: + +```matlab +% use defaults lf = lf_generate_frombrainstorm(); ``` -To use this function you need to provide three file paths (plus one optional) to the following files generated by Brainstorm: -* The headmodel *mat* file generated in Brainstorm -* The channel location *mat* file used for the generation of the headmodel in Brainstorm -* The structural MRI *mat* file to which the electrodes have been coregistered. -* (Optional) the atlas file if used +Note that the atlas remains optional. Therefore, to use this, you would run: -Description of how to obtain these files is included in the help section of `lf_generate_frombrainstorm` +```matlab +lf = lf_generate_frombrainstorm('atlas', 'scout_Mindboggle_62.mat'); +``` -The current function is provided with a precomputed headmodel and associated files. For instance, if you have created a headmodel in Brainstorm from a personal MRI scan. The default model and files were obtained as follows: +These default files were obtained as follows: * *EGI HydroCel 256* channel locations already coregistered by Brainstorm to the default T1 image (next point) * *ICBM152 T1* image used in the default anatomy by Brainstorm * Surface *Headmodel* generated from the above files by applying OpenMEEG to the BEM surfaces and costraining the dipoles to the cortex surface * *Mindboggle 62* atlas -If you have created a personal headmodel, add the required files to the path, then use: - -```matlab -% Use personal files -lf = lf_generate_frombrainstorm('chanloc', 'PATH_TO\my_chanloc.mat', 't1', 'PATH_TO\my_T1_MRI.mat', 'headmodel', 'PATH_TO\my_headmodel.mat', 'atlas', 'PATH_TO\my_atlas.mat'); -``` Brainstorm headmodels are expressed in $\frac{V}{A-m}$. `lf_generate_frombrainstorm` automatically converts these into $\frac{\mu V}{nA-m}$ according to the following formula: $$ \frac{V}{A-m} \cdot \frac{10^6 \mu V}{V} \cdot \frac{A}{10^9 nA} = 10^{-3} \frac{10^{-3} \mu V}{nA-m} $$ -The conversion allows working with realistic units, which would help define and interpret the results, especially for ERPs. However, if you desire to work with the International System, you can set the *scaleUnits* argument to 0 to turn the conversion off. The help file of the function provides more details about this. +The conversion provides more realistic units, which helps when defining and interpreting results, especially for ERPs. However, if you desire to work with the International System, you can set the `scaleUnits` argument to `0` to turn the conversion off. The help file of the function provides more details about this. -```matlab -% Use SI units (V/A-m) -lf = lf_generate_frombrainstorm('chanloc', 'PATH_TO\my_chanloc.mat', 't1', 'PATH_TO\my_T1_MRI.mat', 'headmodel', 'PATH_TO\my_headmodel.mat', 'atlas', 'PATH_TO\my_atlas.mat', 'scaleUnits', 0); -``` +Furthermore, note that this function converts the coordinate system from the *Subject-Coordinate-System* (SCS) used in Brainstorm, to the MNI system used in SEREEGA. You can turn this option off by setting the parameter `useMNI` to `0`. -Furthermore, note that this function converts the coordinate system from the *Subject-Coordinate-System* (SCS) used in Brainstorm, to the MNI system used in SEREEGA. You can tunr this -option off by setting the parameter `useMNI` to 0. +Finally, as Brainstorm uses the SCS, the channel and dipole coordinates are expressed in meters. `lf_generate_frombrainstorm` converts these measures to mm, for consistency with the other leadfields. The argument `useMm` can be used to override this. -```matlab -% Use SI units (V/A-m) -lf = lf_generate_frombrainstorm('chanloc', 'PATH_TO\my_chanloc.mat', 't1', 'PATH_TO\my_T1_MRI.mat', 'headmodel', 'PATH_TO\my_headmodel.mat', 'atlas', 'PATH_TO\my_atlas.mat', 'scaleUnits', 0, 'useMNI', 0); -``` -Finally, as Brainstorm uses the SCS, the channel and dipole coordinates are expressed in meters. `lf_generate_frombrainstorm' converts these measures in mm, for consistency with the other leadfields. -However, if you wish to use meters, then set the parameter `useMm` to 0. Note that the electrode names will not be displayed when calling `plot_headmodel()` with measures expressed in meters. - -```matlab -% Use SI units (V/A-m) -lf = lf_generate_frombrainstorm('chanloc', 'PATH_TO\my_chanloc.mat', 't1', 'PATH_TO\my_T1_MRI.mat', 'headmodel', 'PATH_TO\my_headmodel.mat', 'atlas', 'PATH_TO\my_atlas.mat', 'scaleUnits', 0, 'useMNI', 0, 'useMm', 0); -``` ### Picking a source location @@ -225,6 +229,27 @@ plot_source_location(source, lf, 'mode', '3d'); Other options to obtain source locations are `lf_get_source_inradius` to get all sources in a specific radius from either given location coordinates, or from another source in the lead field. If you want more than one source picked at random but with a specific spacing between them, use `lf_get_source_spaced`. +#### Using atlases + +Starting with SEREEGA v1.4.0, lead fields can have an associated atlas. This means that each individual source in the lead field can be associated with a specific region of the brain. Depending on the lead field and atlas used, some sources may be in the 'Prefrontal cortex', others in the 'Superior temporal gyrus', and others in 'Brodmann area 32', for example. This allows you to pick a source based on the names of these regions, as opposed to using numerical coordinates. + +As of December 2022, the New York Head also comes with such an atlas. (If you have downloaded this head model before that time, you may want to download it again.) The following code uses `plot_headmodel` to visualise the region taken up by one particular region, defined as the right insular cortex. It then uses `lf_get_source_random` we saw before to obtain a random source from within this region. + +```matlab +plot_headmodel(lf, 'region', {'Brain_Right_Insular_Cortex'}, 'labels', 0); +source_insularrnd = lf_get_source_random(lf, 'region', {'Brain_Right_Insular_Cortex'}); +``` + +Other functions allow you to e.g. obtain a source closest to some calculated 'middle' of a region. All `lf_get_source_*` functions can be constrained to specific regions. + +```matlab +source_insularmid = lf_get_source_middle(lf, 'region', {'Brain_Right_Insular_Cortex'}); +plot_source_location(source_insularmid, lf, 'mode', '3d'); +``` + +Note that few lead fields come with an appropriate atlas. The latest version of the New York Head does, and the HArtMuT head model has a detailed atlas for all its non-brain sources, using those obtained from the New York Head for its cortical sources. Currently, no other supported third-party lead field comes with an atlas. A utility called [mni2atlas](https://github.com/dmascali/mni2atlas) has been made compatible with SEREEGA through the `utl_add_atlas_frommni2atlas` function, allowing source coordinates from a lead field to be mapped according to a number of available atlases. However, this may result in a poorly-fitting solution, so use with caution. + + ### Orienting a source dipole The sources in the lead field are represented by dipoles at specific locations. These dipoles can be oriented in space into any direction: they can be pointed towards the nose, or the left ear, or anywhere else. Their orientation is indicated using an [x, y, z] orientation array. Along with its location, a dipole's orientation determines how it projects onto the scalp. @@ -236,7 +261,7 @@ plot_source_projection(source, lf, 'orientation', [0 1 0], 'orientedonly', 1); plot_source_projection(source, lf, 'orientation', [0 0 1], 'orientedonly', 1); ``` -These projections can be linearly combined to get any simulated dipole orientation. In the ICBM-NY lead field, default orientations are included, which orient the dipole perpendicular to the cortical surface. If no orientation is provided, this default orientation is used. FieldTrip-generated lead fields and the Pediatric Head Atlases have no meaningful default orientation and thus revert to all-zero orientations. It is thus necessary to explicitly indicate an orientation for these lead fields to work. For that, you can use `utl_get_orientation_random` to get a random orientation, `utl_get_orientation_pseudoperpendicular` for an orientation pointing outward toward the scalp surface, or `utl_get_orientation_pseudotangential` for an orientation that attempts to approximate a tangent to the cortical surface. +These projections can be linearly combined to get any simulated dipole orientation. In the ICBM-NY lead field, default orientations are included, which orient the dipole perpendicular to the cortical surface. If no orientation is provided, this default orientation is used. FieldTrip-generated lead fields and the Pediatric Head Atlases have no meaningful default orientation and thus revert to all-zero orientations. It is thus necessary to explicitly indicate an orientation for these lead fields to work. For that, you can use `utl_get_orientation_random` to get a random orientation, `utl_get_orientation_pseudoperpendicular` for an orientation pointing outward toward the scalp surface, or `utl_get_orientation_pseudotangential` for an orientation that attempts to approximate a tangent to the surface. ```matlab plot_source_projection(source, lf, 'orientation', [1, 1, 0]); @@ -245,7 +270,7 @@ plot_source_projection(source, lf); When looking at source localisation results in EEG literature, authors usually report the _location_ of a source, but not its _orientation_ per se. What they do report, however, is its projection pattern, usually in the form of topoplots like the ones we plotted in this section. To mimic known effects, you could thus either find a souce from your region of interest whose default orientation results in the projection pattern you are seeking to mimic, or you can pick a source and adjust its orientation to match the desired projection. -Note that the `plot_source_projection` function used here merely plots the projection pattern at the indicated orientation; it does not by itself change the orientation of any source in the simulation. The actual source orientations to be used in the simulation are indicated at the level of [brain components](#brain-components-and-scalp-data), described further below. +Note that the `plot_source_projection` function used here merely plots the projection pattern at the indicated orientation; it does not by itself change the orientation of any source in the simulation. The actual source orientations to be used in the simulation are indicated at the level of [components](#components-and-scalp-data), described further below. ### Defining a source activation signal @@ -274,7 +299,7 @@ erp = utl_check_class(erp, 'type', 'erp') (You may notice one of the parameters of the final class is its `type`, which is set to `erp`. This is how SEREEGA knows how to interpret this class. If we had manually added this field to our class definition, we no longer would have needed to pass this information to `utl_check_class`, and instead could have simple called `erp = utl_check_class(erp)`.) -We can now plot what this ERP would look like. +We can now plot what this ERP would look like. For this, we also need to know the length of the epoch and the sampling rate, which we defined earlier in our `epochs` configuration struct. ```matlab plot_signal_fromclass(erp, epochs); @@ -283,9 +308,9 @@ plot_signal_fromclass(erp, epochs); An ERP activation class can have any number of peaks. For _n_ peaks, you should then also indicate _n_ latencies, _n_ widths, _n_ amplitudes, et cetera--one for each peak. For example, `erp = struct('peakLatency', [450, 500], 'peakWidth', [200, 200], 'peakAmplitude', [-1, 1])` produces an ERP with two peaks. -### Brain components and scalp data +### Components and scalp data -Having defined both a signal (the ERP) and a source location plus projection pattern for this signal, we can now combine these into a single component. _Brain components_ again are represented as structure arrays in this toolbox, with separate fields for the component's source location, orientation, and its activation signal. +Having defined both a signal (the ERP) and a source location plus projection pattern for this signal, we can now combine these into a single component. _Components_ again are represented as structure arrays in this toolbox, with separate fields for the component's source location, orientation, and its activation signal. ```matlab c = struct(); @@ -295,7 +320,7 @@ c.signal = {erp}; % ERP class, defined above c = utl_check_component(c, lf); ``` -Note that, just as for classes, there is a function, `utl_check_component`, to validate the component structure and fill in any missing parameters. For example, if none is indicated, this function reverts the source's orientation to a default value obtained from the lead field. +Note that, just as for classes, there is a function, `utl_check_component`, to validate the component structure and fill in any missing parameters. For example, if no orientation is indicated, this function reverts the source's orientation to a default value obtained from the lead field. (Also have a look at `utl_create_component` in the examples further below for a shorthand function to replace manually assigning structure arrays with `signal` and `source` fields. It also automatically verifies the resulting components.) @@ -317,7 +342,7 @@ scalpdata = generate_scalpdata(c, lf, epochs); After the above step, `scalpdata` contains the channels-by-samples-by-epochs simulated data, but no additional information, such as time stamps, channel names, et cetera. -We can turn this into a dataset according to the EEGLAB format which has a standard structure to save such information. Doing so will also allow us to use EEGLAB's analysis tools, for example, to see the scalp projection time course of the ERP we just simulated. At this point, two optional parameters in the configuration array `epochs` are taken into account. +We can turn this into a dataset according to the EEGLAB format, which has a standard structure to save such information. Doing so will also allow us to use EEGLAB's analysis tools, for example, to see the scalp projection time course of the ERP we just simulated. At this point, two optional parameters in the configuration array `epochs` are taken into account. ```matlab epochs.marker = 'event 1'; % the epochs' time-locking event marker @@ -347,7 +372,7 @@ Keep in mind that EEGLAB uses the variable `EEG` internally to refer to the curr If we scroll through the data, we see that all 100 epochs we generated are exactly the same, having a peak at exactly the indicated latency with a width of exactly 200. ```matlab -pop_eegplot( EEG, 1, 1, 1); +pop_eegplot(EEG, 1, 1, 1); ``` This is of course unrealistic and defeats the entire purpose of even simulating multiple epochs in the first place. @@ -364,7 +389,7 @@ EEG = utl_create_eeglabdataset(generate_scalpdata(c, lf, epochs), ... pop_eegplot(EEG, 1, 1, 1); ``` -Note that `peakLatencyDv` applies to each peak latency value separately. Another parameter, `peakLatencyShift`, works the same way but applies to all values equally. +When multiple peaks are defined, note that `peakLatencyDv` applies to each peak latency value separately. Another parameter, `peakLatencyShift`, works the same way but applies to all values equally. Indicating a slope results in a consistent change over time. An amplitude of 1 and an amplitude slope of -.75, for example, results in the signal having a peak amplitude of 1 in the first epoch, and .25 in the last. @@ -403,15 +428,17 @@ noise = struct( ... noise = utl_check_class(noise); c.signal = {noise}; + EEG = utl_create_eeglabdataset(generate_scalpdata(c, lf, epochs), ... epochs, lf); pop_eegplot(EEG, 1, 1, 1); ``` -It does not have to be either/or. We can also add this noise and the already-existing ERP activation together. We can do this by simply adding both the ERP class variable, and the noise class variable to the component's signal field. +It does not have to be either/or. We can also add this noise and the already-existing ERP activation together. We can do this by simply adding both the ERP class variable, and the noise class variable to the component's signal field. ```matlab c.signal = {erp, noise}; + EEG = utl_create_eeglabdataset(generate_scalpdata(c, lf, epochs), ... epochs, lf); pop_eegplot(EEG, 1, 1, 1); @@ -434,8 +461,21 @@ You can define as many components as you want using any combination of sources a Keep in mind that adding multiple sources to one and the same component is completely different from adding the same signal(s) to multiple different components that each have one source. In the former case, the _exact same_ signal activation is projected simultaneously from all sources. In the latter case, separately simulated instances of the signal activation are projected from each source. +This, then, results in the following simplified workflow using SEREEGA: + +![SEREEGA workflow](/docs/workflow.png) -### Event-related spectral perturbations + +### More signal types + +We have seen the _event-related potential_ (ERP) and _noise_ types of signal activations above. SEREEGA supports additional signal types, and [can easily be extended]((#activation-signals)) to support more. The next few sections will address the other available signals. + +Although not all parameters can easily be visualised, the below figure provides an illustration of three base signal types (ERP, noise, and ERSP, described next) along with its main parameters. + +![Selected SEREEGA signals and parameters](/docs/signals.png) + + +#### Event-related spectral perturbations ERP and noise classes are examples of two types of signal activations. A third type concerns oscillatory activations. In its basic form, it is merely a sine wave of a given frequency, amplitude and, optionally, phase. @@ -467,7 +507,7 @@ plot_signal_fromclass(ersp, epochs); Note that the above two examples contained a `modulation` field that was set to `none`. When this field is set to different values, these signals serve as base oscillatory signals that are then modulated in the indicated way. -#### Frequency burst for event-related synchronization +##### Frequency burst for event-related synchronization The base oscillatory signal can be modulated such that it appears only as a single frequency burst, with a given peak (centre) latency, width, and taper. @@ -483,7 +523,7 @@ plot_signal_fromclass(ersp, epochs); ``` -#### Inverse frequency burst for event-related desynchronization +##### Inverse frequency burst for event-related desynchronization The inverse of the above is also possible. It results in an attenuation in the given window. In both cases, it is also possible to set a minimum amplitude, in order to restrict the attenuation, which is otherwise 100%. @@ -496,7 +536,7 @@ plot_signal_fromclass(ersp, epochs); ``` -#### Amplitude modulation +##### Amplitude modulation This allows the base signal's amplitude to be modulated according to the phase of another. In the example below, a 20 Hz base frequency is modulated using a 2 Hz wave. @@ -533,7 +573,7 @@ pop_eegplot(EEG, 1, 1, 1); Note that variability parameters (deviation, slope) can be added to almost all parameters defining ERSP signals, as above with ERP parameters. -### Pre-generated data as activation signal +#### Pre-generated data as activation signal Most SEREEGA activation classes are intended for the procedural generation of an activation signal based on given parameters. To include time series that were externally obtained, generated, or modulated, a `data` class is provided. This allows a signal activation to be extracted from a matrix containing time series, based on given indices. @@ -558,7 +598,7 @@ plot_signal_fromclass(data, epochs); This class projects, for each simulated epoch `e`, the `e`th row of data from the matrix `randomdata`. `plot_signal` simply plots the first epoch. -### Autoregressive models +#### Autoregressive models When used as a single signal, an autoregressive model (ARM) class generates a time series where each sample depends linearly on its preceding samples plus a stochastic term. The order of the model, i.e. the number preceding samples the current sample depends on, can be configured. The exact model is generated using random parameters and fixed for that class upon verification. @@ -607,7 +647,7 @@ sourcelocs = lf_get_source_spaced(leadfield, 10, 50); erp = erp_get_class_random([1:3], [200:1000], [25:200], ... [-1:.1:-.5, .5:.1:1], 'numClasses', 10); -% combining into brain components +% combining into components c = utl_create_component(sourcelocs, erp, lf); % plotting each component's projection and activation @@ -645,6 +685,13 @@ Just as activation signals can have an epoch-to-epoch variability, so too can th Indicating more than one source for a component, as mentioned [above](#a-note-on-components-and-sources), results in a signal projected from all those sources. However, in versions of SEREEGA older than v1.1.0, this was interpreted as location variability of the source. For each epoch, it randomly selected one of the indicated sources to project the signal activation from. This behaviour can be reinstated using the `legacy_rndmultsrc` argument when calling `generate_scalpdata`. As such, you could use `lf_get_source_inradius` to get all sources around a particular point, and add all of these to a single component to simulate activity coming from variable locations in the same brain region. +### Relativity of the lead field units + +It is relevant to highlight that the various lead fields supported by SEREEGA can have different units of measure. This is the case, for example, for all three of the originally supported lead fields: NYHead, FieldTrip, and Pediatric Head Atlas. Unfortunately, although knowing the units is important to interpret the results correctly, these have not been reported in the corresponding literature. Consequently, the units of the result are relative to the leadfield employed. As a workaround, the `generate_scalpdata` function contains a `normaliseLeadfield` argument, which normalises the leadfield values to at least maintain comparability. + +An exception to this is the lead field converted from Brainstorm. Brainstorm explicitly utilises the International System; thus, every leadfield generated with it is expressed in $\frac{V}{A-m}$ (more information is provided [here](https://neuroimage.usc.edu/forums/t/eeg-units/1499)). The function `lf_generate_frombrainstorm` automatically converts this to $\frac{\mu V}{nA-m}$, unless otherwise requested. + + ### Using the GUI ![SEREEGA](/docs/SEREEGA-GUI.png) @@ -655,13 +702,13 @@ When using SEREEGA as an EEGLAB plug-in, it will appear as a separate sub-menu i The option **Configure epochs** allows you to indicate how much data is to be simulated. The various options in the **Configure lead field** sub-menu allow you to add one of these lead fields to the dataset, thus making it the basis for the simulation. Note that these options can all be changed at any time, but they must be set before the following functions can be used. -The **Configure components** sub-menu contains three options which can be used to define the brain components that will underly the simulated data. +The **Configure components** sub-menu contains three options which can be used to define the components that will underly the simulated data. First, **Select source locations** provides a dialog window that allows you to find random or specific sources in the brain (i.e. in the configured lead field), determine their desired orientation, and add them to the simulation. Two plots can support you in this task: one for the source's location, and one for its projection. If you have found a source you wish to keep, click *Add source(s)* to add it to the list at the left, and *OK* to save the list to the dataset. **Define signal activations**, like the previous window to select source locations, provides a list of signals currently stored in the dataset, and options to add additional signal classes. Each signal type has its own button, which pops up a second window to input the parameters that define each signal. Values that are indicated with an asterisk (*) in both their row and column are required. Click *OK* to save the list. -Now it must be decided which of the defined signals will project from which of the selected source locations. The **Assign signals to sources** dialog shows a list of the sources selected earlier, and allows you to assign the defined signal classes to these, thus completing the definition of brain components. +Now it must be decided which of the defined signals will project from which of the selected source locations. The **Assign signals to sources** dialog shows a list of the sources selected earlier, and allows you to assign the defined signal classes to these, thus completing the definition of components. Finally, **Simulate data** simulates the data and populates the EEGLAB dataset with the corresponding values. @@ -692,9 +739,14 @@ To add a lead field to SEREEGA, create a `lf_generate_.m` script % zeros(nsources, 3) if no meaning default is available % .pos - nsources x 3 xyz MNI coordinates of each source % .chanlocs - channel information in EEGLAB format +% .atlas - nsources x 1 cell of strings containing the atlas, +% i.e. indciating the name of the corresponding region +% for each source ``` -If the lead field does not come with a meaning default orientation, set the defaults to all zeros to force a manual orientation. +If the lead field does not come with a meaningful default orientation, set the defaults to all zeros to force a manual orientation. + +Each entry in the atlas must start with one of the three generic categories that SEREEGA recognises: Brain, Eye, or Muscle (not case sensitive). If the lead field does not come with a meaningful atlas, simply only indicate one of these three. Use `utl_sanitize_atlas` to make sure the final atlas format is compatible with other SEREEGA functions. Note that SEREEGA's standard plotting functions calculate a brain boundary based on dipole locations. The accuracy of these plots thus depend on the resolution of the lead field. Alternatively, there are two `_dipplot` plotting functions functions which use a standard head model as backdrop. @@ -704,14 +756,14 @@ Some head models, such as the Pediatric Head Atlas of ages up to 2 years old, wi The electrode coordinates require special attention. These should be centred around a similar zero-point as the standard MNI head midpoint, in order for EEGLAB's `topoplot` to be able to properly plot the lead field's projection patterns. Note that EEGLAB's `readlocs`, which is used by SEREEGA, reads the electrode positions with X=Y and Y=-X. -Or perhaps you came to this section because you wanted to add a new source to an existing lead field in your work space. This can be done manually using `lf_add_source`. +(If you came to this section of the documentation because you wanted to add a new source to an existing lead field in your work space: This can be done manually using `lf_add_source`.) ### Activation signals Files and scripts related to signal activation classes are in class-specific subdirectories of the `./signal` directory. To add new classes of signal activations to SEREEGA, the following files, containing functions of the same name, must be supplied in a new subdirectory. denotes the name of the new class. -`_check_class` - Takes a class definition structure, verifies/completes it to pass the requirements of the new class, and returns a class variable that is compatible with all other functions. The class must have a 'type' field that is equal to , allowing `utl_check_class` to find this file. This file is also where the class documentation should be provided. If any deviation and slope fields are provided, add 'Dv' or 'Slope' to the base field name. A slope for a 'frequency' field should thus be called 'frequencySlope', not e.g. 'freqSlope'. This allows e.g. `utl_apply_dvslope` to function properly. Similarly, if any latencies are to be indicated, use a field name that ends in `Latency' and describes what latency it is, as e.g. with `peakLatency` for ERPs and `modLatency` for ERSP modulations. This allows e.g. `utl_shift_latency` to function properly. +`_check_class` - Takes a class definition structure, verifies/completes it to pass the requirements of the new class, and returns a class variable that is compatible with all other functions. The class must have a 'type' field that is equal to , allowing `utl_check_class` to find this file. This file is also where the class documentation should be provided. If any deviation and slope fields are provided, add 'Dv' or 'Slope' to the base field name. A slope for a 'frequency' field should thus be called 'frequencySlope', not e.g. 'freqSlope'. This allows e.g. `utl_apply_dvslope` to function properly. Similarly, if any latencies are to be indicated, use a field name that ends in `Latency` and describes what latency it is, as e.g. with `peakLatency` for ERPs and `modLatency` for ERSP modulations. This allows e.g. `utl_shift_latency` to function properly. `_generate_signal_fromclass` - Takes a (verified) class structure, an epochs configuration structure, and (at least accepts) an epochNumber argument (indicating the number of the currently generated epoch) and 'baseonly' (whether or not to ignore the deviations, slopes, et cetera). Returns a 1-by-nsamples signal activation time course. @@ -829,11 +881,6 @@ data = utl_mix_data(sigdata, noisedata, 1/3); Also note that an additional source of noise, sensor noise, can be added to generated scalp data using the `sensorNoise` argument of `generate_scalpdata`. This noise has no dependencies across channels or samples. -## Relativity of the Units of Measure - -It is relevant to highlight that the three leadfield provided by default (NYhead, Fieldtrip, and Pediatric Head Atlas) have different units of measure. Unfortunately, although knowing the units is important to interpret the results correctly, these have not been reported in the literature. Consequently, the units of the result are relative to the leadfield employed. As a workaround, the `generate_scalpdata` function contains a *normaliseLeadfield* argument, which normalises the leadfield values. - -An exception to this is the leadfield converted from Brainstorm. Brainstorm explicitly utlises the International System; thus, every leadfield generated with it is expressed in $\frac{V}{A-m}$ (for more information: https://neuroimage.usc.edu/forums/t/eeg-units/1499). The function `lf_generate_frombrainstorm` automatically converts this in $\frac{\mu V}{nA-m}$ (unless otherwise requested). This should make the results easier to interpret. ## Contact @@ -842,4 +889,4 @@ Feel free to contact me at lrkrol@gmail.com. ## Special thanks and acknowledgements -Fabien Lotte provided two early drafts of what are now `lf_generate_fromfieldtrip` and `erp_generate_signal`. I'd like to thank Mahta Mousavi for the band-pass filter design, Stefan Haufe for the autoregressive signal generation code and accompanying support, and Ramón Martínez-Cancino for his assistance with the GUI development. Part of this work was supported by the Deutsche Forschungsgemeinschaft (grant number ZA 821/3-1), Idex Bordeaux and LabEX CPU/SysNum, and the European Research Council with the BrainConquest project (grant number ERC-2016-STG-714567). +Fabien Lotte provided two early drafts of what are now `lf_generate_fromfieldtrip` and `erp_generate_signal`. I'd like to thank Mahta Mousavi for the band-pass filter design, Stefan Haufe for the autoregressive signal generation code and accompanying support, and Ramón Martínez-Cancino for his assistance with the GUI development. Part of this work was supported by the Deutsche Forschungsgemeinschaft (grant number ZA 821/3-1), Idex Bordeaux and LabEX CPU/SysNum, and the European Research Council with the BrainConquest project (grant number ERC-2016-STG-714567), and Volkswagen Foundation. \ No newline at end of file diff --git a/docs/changelog.txt b/docs/changelog.txt index b6af022..e99ae76 100644 --- a/docs/changelog.txt +++ b/docs/changelog.txt @@ -1,3 +1,13 @@ +v1.5.0 +- Added utl_call_gui, allowing GUI functions to be called without requiring the EEGLAB GUI itself to be used or a SEREEGA dataset to be loaded +- Added source ID to pop_sereega_sources GUI and console output +- Various documentation updates +- Added support for new HArtMuT (Head Artefact Model using Tripoles) head models containing ocular and muscular sources alongside brain sources +- Added atlas support, allowing sources to be selected from specified named regions. SEREEGA leadfields now have an additional atlas field containing, for each source, the name of the region it belongs to. lf_get_source_* functions have been added/adjusted to allow sources to be obtained from specific regions. lf_add_atlas_frommni2atlas can be used to add atlases to existing lead fields using the mni2atlas function from the FMRIB Software Library (not included). +- Added electrodes argument to plot_headmodel +- Added plot_source_moment to plot source orientation as arrow +- Added support for lead fields generated in Brainstorm + v1.2.2 - Fixed forgotten version increment in eegplugin_sereega.m - Fixed forgotten documentation update reflecting 2018-06-04 utl_apply_dvslopeshift change diff --git a/docs/workflow.png b/docs/workflow.png new file mode 100644 index 0000000..8cd84d5 Binary files /dev/null and b/docs/workflow.png differ diff --git a/eegplugin_sereega.m b/eegplugin_sereega.m index 825adf5..98f5f4c 100644 --- a/eegplugin_sereega.m +++ b/eegplugin_sereega.m @@ -16,6 +16,8 @@ % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% 2022-11-17 lrk +% - Added HArtMuT lead field % 2018 through 2021 lrk % - Version number increments % 2018-04-23 First version @@ -37,7 +39,7 @@ function version = eegplugin_sereega(fig, try_strings, catch_strings) - version = '1.2.0'; + version = '1.5.0'; % find tools menu % --------------- @@ -84,9 +86,10 @@ 'eeglab redraw;' ... ]; cb_epochs = 'EEG = pop_sereega_epochs(EEG);'; + cb_lf_fieldtrip = 'EEG = pop_sereega_lf_generate_fromfieldtrip(EEG);'; + cb_lf_hartmut = 'EEG = pop_sereega_lf_generate_fromhartmut(EEG);'; cb_lf_nyhead = 'EEG = pop_sereega_lf_generate_fromnyhead(EEG);'; cb_lf_pha = 'EEG = pop_sereega_lf_generate_frompha(EEG);'; - cb_lf_fieldtrip = 'EEG = pop_sereega_lf_generate_fromfieldtrip(EEG);'; cb_comp_sources = 'EEG = pop_sereega_sources(EEG);'; cb_comp_signals = 'EEG = pop_sereega_signals(EEG);'; cb_comp_components = 'EEG = pop_sereega_components(EEG);'; @@ -138,9 +141,10 @@ menu_new = uimenu(menu_root, 'label', 'New empty dataset', 'userdata', userdata, 'callback', cb_new); menu_epochs = uimenu(menu_root, 'label', 'Configure epochs', 'userdata', userdata, 'callback', cb_epochs); menu_lf = uimenu(menu_root, 'label', 'Configure lead field'); + menu_lf_fieldtrip = uimenu(menu_lf, 'label', 'FieldTrip', 'userdata', userdata, 'callback', cb_lf_fieldtrip); + menu_lf_hartmut = uimenu(menu_lf, 'label', 'HArtMuT', 'userdata', userdata, 'callback', cb_lf_hartmut); menu_lf_nyhead = uimenu(menu_lf, 'label', 'New York Head', 'userdata', userdata, 'callback', cb_lf_nyhead); menu_lf_pha = uimenu(menu_lf, 'label', 'Pediatric Head Atlas', 'userdata', userdata, 'callback', cb_lf_pha); - menu_lf_fieldtrip = uimenu(menu_lf, 'label', 'FieldTrip', 'userdata', userdata, 'callback', cb_lf_fieldtrip); menu_comp = uimenu(menu_root, 'label', 'Configure components', 'userdata', userdata); menu_comp_sources = uimenu(menu_comp, 'label', 'Select source locations', 'userdata', userdata, 'callback', cb_comp_sources); menu_comp_signals = uimenu(menu_comp, 'label', 'Define signal activations', 'userdata', userdata, 'callback', cb_comp_signals); diff --git a/leadfield/fieldtrip/lf_generate_fromfieldtrip.m b/leadfield/fieldtrip/lf_generate_fromfieldtrip.m index 201ae07..4dfe289 100644 --- a/leadfield/fieldtrip/lf_generate_fromfieldtrip.m +++ b/leadfield/fieldtrip/lf_generate_fromfieldtrip.m @@ -26,6 +26,9 @@ % for each source. % .pos - xyz MNI coordinates of each source % .chanlocs - channel information in EEGLAB format +% .atlas - atlas (region) indication for each source; this +% is simply 'Brain' until more specific +% functionality is added % % Usage example: % >> lf = lf_generate_fromfieldtrip('labels', {'Fz', 'Cz', 'Pz'}, ... @@ -40,6 +43,8 @@ % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% 2021-01-06 lrk +% - Added lf.atlas % 2018-11-08 lrk % - Now continues with empty labels if montage not found % 2018-03-26 lrk @@ -163,5 +168,6 @@ lf.pos = ftlf.pos(ftlf.inside,:); lf.chanlocs = readlocs('chanlocs-standard_1005.elp'); lf.chanlocs = lf.chanlocs(chanidx); +lf.atlas = repmat({'Brain'}, size(lf.pos, 1), 1); end \ No newline at end of file diff --git a/leadfield/hartmut/chanlocs-hartmut-NYhead_small.xyz b/leadfield/hartmut/chanlocs-hartmut-NYhead_small.xyz new file mode 100644 index 0000000..c587b8d --- /dev/null +++ b/leadfield/hartmut/chanlocs-hartmut-NYhead_small.xyz @@ -0,0 +1,231 @@ +1 -36 75.5 25.5 AF3 +2 -28.5 77.5 30 AF3h +3 36 75.5 25.5 AF4 +4 28.5 77.5 30 AF4h +5 -43 72.5 19.5 AF5h +6 43 72.5 19.5 AF6h +7 -55.5 64 -3 AF7 +8 55.5 64 -3 AF8 +9 -12.5 70.5 52 AFF1h +10 12.5 70.5 52 AFF2h +11 -52 60.5 28 AFF5h +12 52 60.5 28 AFF6h +13 -62 54 9 AFF7h +14 62 54 9 AFF8h +15 -14 84.5 18.5 AFp1 +16 14 84.5 18.5 AFp2 +17 0 80 36.5 AFz +18 -37.5 -11 89.5 C1 +19 37.5 -11 89.5 C2 +20 -67.5 -12 65.5 C3 +21 67.5 -12 65.5 C4 +22 -82.5 -13 33.5 C5 +23 82.5 -13 33.5 C6 +24 -38 -29.5 92 CCP1 +25 -20 -30 99 CCP1h +26 38 -29.5 92 CCP2 +27 20 -30 99 CCP2h +28 -68.5 -29.5 68.5 CCP3 +29 -54.5 -29.5 81.5 CCP3h +30 68.5 -29.5 68.5 CCP4 +31 54.5 -29.5 81.5 CCP4h +32 -83.5 -28.5 35.5 CCP5 +33 -78.5 -29 52.5 CCP5h +34 83.5 -28.5 35.5 CCP6 +35 78.5 -29 52.5 CCP6h +36 -37.5 -49.5 91 CP1 +37 37.5 -49.5 91 CP2 +38 -67 -47 67 CP3 +39 67 -47 67 CP4 +40 -81.5 -44.5 35 CP5 +41 81.5 -44.5 35 CP6 +42 -18.5 -69.5 90.5 CPP1h +43 18.5 -69.5 90.5 CPP2h +44 -50 -67 76 CPP3h +45 50 -67 76 CPP4h +46 -72 -62.5 48.5 CPP5h +47 72 -62.5 48.5 CPP6h +48 0 -49.5 100 CPz +49 0 -9 99 Cz +50 -14.5 -111 -54.5 Ex1 +51 60.5 -81.5 -59 Ex10 +52 -68.5 -69 -59.5 Ex11 +53 68.5 -69 -59.5 Ex12 +54 -74.5 -56.5 -60.5 Ex13 +55 74.5 -56.5 -60.5 Ex14 +56 -82.5 -18.5 -62 Ex19 +57 14.5 -111 -54.5 Ex2 +58 82.5 -18.5 -62 Ex20 +59 -81.5 -4 -63 Ex21 +60 81.5 -4 -63 Ex22 +61 -79.5 9 -64.5 Ex23 +62 79.5 9 -64.5 Ex24 +63 -76.5 23 -66.5 Ex25 +64 76.5 23 -66.5 Ex26 +65 -72 36.5 -69 Ex27 +66 72 36.5 -69 Ex28 +67 -65.5 49 -71 Ex29 +68 -28 -108 -57 Ex3 +69 65.5 49 -71 Ex30 +70 -55 59.5 -71 Ex31 +71 55 59.5 -71 Ex32 +72 -41.5 66 -68 Ex33 +73 41.5 66 -68 Ex34 +74 28 -108 -57 Ex4 +75 -40.5 -101.5 -58.5 Ex5 +76 40.5 -101.5 -58.5 Ex6 +77 -51.5 -92.5 -58.5 Ex7 +78 51.5 -92.5 -58.5 Ex8 +79 -60.5 -81.5 -59 Ex9 +80 -13 -106 -81.5 Exx1 +81 55.5 -76.5 -85 Exx10 +82 -62 -65.5 -84.5 Exx11 +83 62 -65.5 -84.5 Exx12 +84 -66.5 -54.5 -84.5 Exx13 +85 66.5 -54.5 -84.5 Exx14 +86 -69 -41.5 -82.5 Exx15 +87 69 -41.5 -82.5 Exx16 +88 -76.5 -20 -90 Exx19 +89 13 -106 -81.5 Exx2 +90 76.5 -20 -90 Exx20 +91 -76.5 -6 -89.5 Exx21 +92 76.5 -6 -89.5 Exx22 +93 -75 5.5 -91 Exx23 +94 75 5.5 -91 Exx24 +95 -72 18.5 -93.5 Exx25 +96 72 18.5 -93.5 Exx26 +97 -68 31 -96.5 Exx27 +98 68 31 -96.5 Exx28 +99 -61 44 -99 Exx29 +100 -25.5 -102.5 -83.5 Exx3 +101 61 44 -99 Exx30 +102 -52 55 -99.5 Exx31 +103 52 55 -99.5 Exx32 +104 -41.5 62.5 -98.5 Exx33 +105 41.5 62.5 -98.5 Exx34 +106 25.5 -102.5 -83.5 Exx4 +107 -37 -96 -85 Exx5 +108 37 -96 -85 Exx6 +109 -46.5 -87 -84.5 Exx7 +110 46.5 -87 -84.5 Exx8 +111 -55.5 -76.5 -85 Exx9 +112 0 -107 -78.5 Exxz +113 0 -112 -52.5 Exz +114 -28.5 56 61.5 F1 +115 28.5 56 61.5 F2 +116 -51.5 50.5 45 F3 +117 51.5 50.5 45 F4 +118 -65.5 44 22.5 F5 +119 65.5 44 22.5 F6 +120 -71.5 39 -2.5 F7 +121 71.5 39 -2.5 F8 +122 -34 25 79 FC1 +123 34 25 79 FC2 +124 -60.5 21 58.5 FC3 +125 60.5 21 58.5 FC4 +126 -76.5 16.5 29 FC5 +127 76.5 16.5 29 FC6 +128 -18 8.5 91.5 FCC1h +129 18 8.5 91.5 FCC2h +130 -51.5 5.5 75 FCC3h +131 51.5 5.5 75 FCC4h +132 -74 3 48 FCC5h +133 74 3 48 FCC6h +134 0 26 87.5 FCz +135 -31.5 41 71.5 FFC1 +136 -16 43 77 FFC1h +137 31.5 41 71.5 FFC2 +138 16 43 77 FFC2h +139 -57 36.5 52 FFC3 +140 -45 39 63.5 FFC3h +141 57 36.5 52 FFC4 +142 45 39 63.5 FFC4h +143 -72 30 26.5 FFC5 +144 -65.5 33 39.5 FFC5h +145 72 30 26.5 FFC6 +146 65.5 33 39.5 FFC6h +147 76 24 -20.5 FFT10h +148 -75 28 12.5 FFT7h +149 75 28 12.5 FFT8h +150 -76 24 -20.5 FFT9h +151 79.5 10.5 -40 FT10 +152 -79.5 12 -2 FT7 +153 79.5 12 -2 FT8 +154 -79.5 10.5 -40 FT9 +155 -82.5 -0.5 15 FTT7h +156 82.5 -0.5 15 FTT8h +157 -31.5 81.5 -2.5 Fp1 +158 31.5 81.5 -2.5 Fp2 +159 0 86.5 -0.5 Fpz +160 0 58 67 Fz +161 -29.5 -111.5 -34 I1 +162 29.5 -111.5 -34 I2 +163 -83.5 -17.5 -38.5 LPA +164 0 32 -179 Nk1 +165 0 -113 -154.5 Nk2 +166 62.5 -41 -168 Nk3 +167 -62.5 -41 -168 Nk4 +168 0 83.5 -41 Nz +169 -31 -114.5 7.5 O1 +170 31 -114.5 7.5 O2 +171 -30.5 -114.5 -12.5 OI1 +172 -15.5 -119 -11.5 OI1h +173 30.5 -114.5 -12.5 OI2 +174 15.5 -119 -11.5 OI2h +175 0 -119 12 Oz +176 -31.5 -83.5 75.5 P1 +177 71.5 -71 -37.5 P10 +178 31.5 -83.5 75.5 P2 +179 -56 -80.5 57.5 P3 +180 56 -80.5 57.5 P4 +181 -70.5 -75 31 P5 +182 70.5 -75 31 P6 +183 -74.5 -70.5 2.5 P7 +184 74.5 -70.5 2.5 P8 +185 -71.5 -71 -37.5 P9 +186 -20 -107 44.5 PO1 +187 54 -95.5 -36 PO10 +188 20 -107 44.5 PO2 +189 -38 -105 35.5 PO3 +190 -29 -106 41.5 PO3h +191 38 -105 35.5 PO4 +192 29 -106 41.5 PO4h +193 -44 -103.5 30 PO5h +194 44 -103.5 30 PO6h +195 -57.5 -96.5 5 PO7 +196 57.5 -96.5 5 PO8 +197 -54 -95.5 -36 PO9 +198 -13.5 -114.5 28.5 POO1 +199 44.5 -107 -14 POO10h +200 13.5 -114.5 28.5 POO2 +201 -26.5 -113 24 POO3 +202 26.5 -113 24 POO4 +203 -37.5 -110.5 15.5 POO5 +204 37.5 -110.5 15.5 POO6 +205 -44.5 -107 -14 POO9h +206 0 -108 47.5 POz +207 -24.5 -97.5 61 PPO1 +208 67 -84.5 -16 PPO10h +209 -12 -97.5 65 PPO1h +210 24.5 -97.5 61 PPO2 +211 12 -97.5 65 PPO2h +212 -36.5 -96 56 PPO3h +213 36.5 -96 56 PPO4h +214 -61.5 -89 27.5 PPO5 +215 61.5 -89 27.5 PPO6 +216 -68 -83.5 3 PPO7 +217 -65.5 -86.5 15.5 PPO7h +218 68 -83.5 3 PPO8 +219 65.5 -86.5 15.5 PPO8h +220 -67 -84.5 -16 PPO9h +221 0 -84.5 81 Pz +222 83.5 -17.5 -38.5 RPA +223 -84.5 -14.5 -1.5 T7 +224 84.5 -14.5 -1.5 T8 +225 -84.5 -43 0 TP7 +226 84.5 -43 0 TP8 +227 -79.5 -57.5 18 TPP7h +228 79.5 -57.5 18 TPP8h +229 -85.5 -29 17.5 TTP7h +230 85.5 -29 17.5 TTP8h +231 0 -117 -30.5 Iz diff --git a/leadfield/hartmut/chanlocs-hartmut-mix_Colin27_small.xyz b/leadfield/hartmut/chanlocs-hartmut-mix_Colin27_small.xyz new file mode 100644 index 0000000..ab3b70e --- /dev/null +++ b/leadfield/hartmut/chanlocs-hartmut-mix_Colin27_small.xyz @@ -0,0 +1,346 @@ +1 -86.076 -19.99 -47.986 LPA +2 85.794 -20.009 -48.031 RPA +3 0.0083 86.811 -39.983 Nz +4 -29.437 83.917 -6.99 Fp1 +5 0.1123 88.247 -1.713 Fpz +6 29.872 84.896 -7.08 Fp2 +7 -48.971 64.087 -47.683 AF9 +8 -54.84 68.572 -10.59 AF7 +9 -45.431 72.862 5.978 AF5 +10 -33.701 76.837 21.227 AF3 +11 -18.472 79.904 32.752 AF1 +12 0.2313 80.771 35.417 AFz +13 19.82 80.302 32.764 AF2 +14 35.712 77.726 21.956 AF4 +15 46.584 73.808 6.034 AF6 +16 55.743 69.657 -10.755 AF8 +17 50.435 63.87 -48.005 AF10 +18 -70.102 41.652 -49.952 F9 +19 -70.263 42.474 -11.42 F7 +20 -64.466 48.035 16.921 F5 +21 -50.244 53.111 42.192 F3 +22 -27.496 56.931 60.342 F1 +23 0.3122 58.512 66.462 Fz +24 29.514 57.602 59.54 F2 +25 51.836 54.305 40.814 F4 +26 67.914 49.83 16.367 F6 +27 73.043 44.422 -12 F8 +28 72.114 42.067 -50.452 F10 +29 -84.076 14.567 -50.429 FT9 +30 -80.775 14.12 -11.135 FT7 +31 -77.215 18.643 24.46 FC5 +32 -60.182 22.716 55.544 FC3 +33 -34.062 26.011 79.987 FC1 +34 0.3761 27.39 88.668 FCz +35 34.784 26.438 78.808 FC2 +36 62.293 23.723 55.63 FC4 +37 79.534 19.936 24.438 FC6 +38 81.815 15.417 -11.33 FT8 +39 84.113 14.365 -50.538 FT10 +40 -85.894 -15.829 -48.283 T9 +41 -84.161 -16.019 -9.346 T7 +42 -80.28 -13.76 29.16 C5 +43 -65.358 -11.632 64.358 C3 +44 -36.158 -9.9839 89.752 C1 +45 0.4009 -9.167 100.24 Cz +46 37.672 -9.6241 88.412 C2 +47 67.118 -10.9 63.58 C4 +48 83.456 -12.776 29.208 C6 +49 85.08 -15.02 -9.49 T8 +50 85.56 -16.361 -48.271 T10 +51 -85.619 -46.515 -45.707 TP9 +52 -84.83 -46.022 -7.056 TP7 +53 -79.592 -46.551 30.949 CP5 +54 -63.556 -47.009 65.624 CP3 +55 -35.513 -47.292 91.315 CP1 +56 0.3858 -47.318 99.432 CPz +57 38.384 -47.073 90.695 CP2 +58 66.612 -46.637 65.58 CP4 +59 83.322 -46.101 31.206 CP6 +60 85.549 -45.545 -7.13 TP8 +61 86.162 -47.035 -45.869 TP10 +62 -73.009 -73.766 -40.998 P9 +63 -72.434 -73.453 -2.487 P7 +64 -67.272 -76.291 28.382 P5 +65 -53.007 -78.788 55.94 P3 +66 -28.62 -80.525 75.436 P1 +67 0.3247 -81.115 82.615 Pz +68 31.92 -80.487 76.716 P2 +69 55.667 -78.56 56.561 P4 +70 67.888 -75.904 28.091 P6 +71 73.056 -73.068 -2.54 P8 +72 73.895 -74.39 -41.22 P10 +73 -54.91 -98.045 -35.465 PO9 +74 -54.84 -97.528 2.792 PO7 +75 -48.424 -99.341 21.599 PO5 +76 -36.511 -100.85 37.167 PO3 +77 -18.972 -101.77 46.536 PO1 +78 0.2156 -102.18 50.608 POz +79 19.878 -101.79 46.393 PO2 +80 36.782 -100.85 36.397 PO4 +81 49.82 -99.446 21.727 PO6 +82 55.667 -97.625 2.73 PO8 +83 54.988 -98.091 -35.541 PO10 +84 -29.413 -112.45 8.839 O1 +85 0.1076 -114.89 14.657 Oz +86 29.843 -112.16 8.8 O2 +87 -29.818 -114.57 -29.216 I1 +88 0.0045 -118.56 -23.078 Iz +89 29.742 -114.26 -29.256 I2 +90 -43.29 75.855 -28.244 AFp9h +91 -38.552 79.953 -4.995 AFp7h +92 -27.986 82.459 2.702 AFp5h +93 -17.195 84.849 10.027 AFp3h +94 -5.9317 86.878 16.2 AFp1h +95 7.1053 87.074 16.469 AFp2h +96 18.923 85.597 11.443 AFp4h +97 28.644 82.976 2.828 AFp6h +98 39.32 80.687 -4.725 AFp8h +99 43.822 76.542 -28.307 AFp10h +100 -63.254 53.857 -30.316 AFF9h +101 -61.351 58.799 0.897 AFF7h +102 -50.8 64.041 23.089 AFF5h +103 -34.316 68.393 41.188 AFF3h +104 -11.436 70.756 50.348 AFF1h +105 13.479 71.201 51.175 AFF2h +106 36.183 69.151 41.254 AFF4h +107 52.397 65.071 22.862 AFF6h +108 62.915 60.045 0.63 AFF8h +109 64.334 54.6 -30.444 AFF10h +110 -79.067 28.081 -31.253 FFT9h +111 -74.5 31.3 4.846 FFT7h +112 -65.238 36.428 36.144 FFC5h +113 -44.41 40.762 61.69 FFC3h +114 -15.424 43.66 77.682 FFC1h +115 17.592 44.054 77.788 FFC2h +116 45.853 41.623 60.647 FFC4h +117 67.128 37.8 35.296 FFC6h +118 78.053 32.982 4.483 FFT8h +119 80.097 28.514 -31.338 FFT10h +120 -84.125 -1.8467 -29.794 FTT9h +121 -82.355 0.8263 8.579 FTT7h +122 -74.692 4.3033 45.307 FCC5h +123 -51.051 7.1772 74.377 FCC3h +124 -18.219 9.0941 92.529 FCC1h +125 18.787 9.2479 91.562 FCC2h +126 51.885 7.7978 73.507 FCC4h +127 77.002 5.3357 45.35 FCC6h +128 83.888 1.9457 8.501 FTT8h +129 84.123 -1.8083 -29.638 FTT10h +130 -86.973 -32.216 -27.848 TTP9h +131 -85.565 -30.629 11.153 TTP7h +132 -76.407 -29.731 49.217 CCP5h +133 -52.928 -28.906 80.304 CCP3h +134 -18.354 -28.322 98.22 CCP1h +135 20.22 -28.148 98.172 CCP2h +136 55.114 -28.386 80.474 CCP4h +137 79.006 -28.986 49.628 CCP6h +138 86 -29.82 11.248 TTP8h +139 88.625 -32.272 -28 TTP10h +140 -78.16 -60.757 -23.824 TPP9h +141 -76.68 -60.832 12.88 TPP7h +142 -68.115 -62.975 47.252 CPP5h +143 -46.914 -64.691 75.296 CPP3h +144 -15.82 -65.6 91.164 CPP1h +145 19.42 -65.595 92.405 CPP2h +146 50.674 -64.482 76.13 CPP4h +147 71.096 -62.624 47.328 CPP6h +148 78.52 -60.432 12.902 TPP8h +149 78.903 -60.955 -23.805 TPP10h +150 -64.597 -87.656 -19.014 PPO9h +151 -62.959 -87.503 12.952 PPO7h +152 -54.01 -89.899 37.332 PPO5h +153 -35.887 -91.667 55.504 PPO3h +154 -12.047 -92.607 65.508 PPO1h +155 13.923 -92.694 66.958 PPO2h +156 37.799 -91.629 56.733 PPO4h +157 54.609 -89.64 37.035 PPO6h +158 63.112 -87.228 12.856 PPO8h +159 65.014 -87.806 -18.952 PPO10h +160 -42.862 -108.07 -13.151 POO9h +161 -40.12 -107.13 12.061 POO7h +162 -31.951 -108.25 23.047 POO5h +163 -19.862 -108.94 29.76 POO3h +164 -6.9194 -109.26 32.71 POO1h +165 6.8036 -109.16 31.582 POO2h +166 20.294 -108.91 28.944 POO4h +167 32.176 -108.25 22.255 POO6h +168 41.098 -107.25 12.138 POO8h +169 43.895 -109.13 -13.17 POO10h +170 -14.85 -117.99 -6.92 OI1h +171 15.095 -118.02 -6.933 OI2h +172 -14.811 87.235 -4.477 Fp1h +173 15.162 88.091 -4.551 Fp2h +174 -54.83 66.413 -29.704 AF9h +175 -51.176 70.836 -1.755 AF7h +176 -39.641 74.867 13.678 AF5h +177 -27.219 78.709 28.375 AF3h +178 -9.1977 80.605 35.133 AF1h +179 10.482 80.865 35.359 AF2h +180 28.58 79.303 28.47 AF4h +181 40.94 75.74 13.86 AF6h +182 52.029 71.847 -1.92 AF8h +183 55.754 67.17 -29.824 AF10h +184 -71.508 41.119 -30.854 F9h +185 -68.556 45.284 3.002 F7h +186 -58.488 50.672 30.192 F5h +187 -39.98 55.26 52.6 F3h +188 -13.384 57.902 64.332 F1h +189 15.834 58.456 64.992 F2h +190 41.794 56.226 51.499 F4h +191 60.052 52.086 28.708 F6h +192 71.959 47.192 2.475 F8h +193 72.798 41.822 -31.026 F10h +194 -82.956 13.32 -30.808 FT9h +195 -80.114 16.39 6.85 FT7h +196 -71.21 20.82 41.324 FC5h +197 -48.512 24.529 69.136 FC3h +198 -17.344 27.024 86.923 FC1h +199 18.418 27.271 86.437 FC2h +200 49.548 25.238 68.43 FC4h +201 73.219 22.007 41.297 FC6h +202 81.58 17.684 6.564 FT8h +203 83.371 13.548 -30.749 FT10h +204 -85.132 -17.056 -28.731 T9h +205 -82.946 -14.883 10.009 T7h +206 -75.294 -12.64 47.904 C5h +207 -51.581 -10.755 78.035 C3h +208 -18.279 -9.4319 97.356 C1h +209 19.678 -9.3041 95.706 C2h +210 53.806 -10.144 77.73 C4h +211 78.125 -11.735 47.84 C6h +212 85.137 -13.906 9.89 T8h +213 86.1 -17.088 -28.756 T10h +214 -84.81 -47.246 -26.22 TP9h +215 -82.704 -46.298 11.974 TP7h +216 -73.301 -46.792 49.109 CP5h +217 -51.049 -47.176 80.016 CP3h +218 -17.354 -47.342 97.41 CP1h +219 20.68 -47.232 98.072 CP2h +220 53.997 -46.89 80.077 CP4h +221 76.55 -46.373 49.14 CP6h +222 85.2 -45.807 12.102 TP8h +223 85.443 -47.221 -26.176 TP10h +224 -72.177 -74.628 -21.536 P9h +225 -70.113 -74.868 12.999 P7h +226 -61.728 -77.624 43.028 P5h +227 -41.673 -79.753 66.715 P3h +228 -13.961 -81.003 81.003 P1h +229 17.298 -80.981 81.641 P2h +230 44.748 -79.611 67.655 P4h +231 63.627 -77.302 43.119 P6h +232 72.104 -74.499 13.025 P8h +233 73.282 -75.077 -21.576 P10h +234 -54.775 -98.977 -16.193 PO9h +235 -51.928 -98.444 12.304 PO7h +236 -43.342 -100.16 30.009 PO5h +237 -28.007 -101.36 42.379 PO3h +238 -9.5034 -102.06 49.418 PO1h +239 10.236 -102.03 48.942 PO2h +240 28.648 -101.39 42.138 PO4h +241 44.221 -100.22 29.808 PO6h +242 52.839 -98.536 12.25 PO8h +243 55.86 -99.894 -16.208 PO10h +244 -14.805 -115.1 11.829 O1h +245 15.146 -115.19 11.833 O2h +246 -15.158 -118.24 -26.048 I1h +247 15.129 -118.15 -26.081 I2h +248 -36.125 72.38 -45.852 AFp9 +249 -43.512 78.58 -9.24 AFp7 +250 -33.285 81.207 -1.14 AFp5 +251 -22.352 83.562 6.071 AFp3 +252 -12.242 86.194 14.188 AFp1 +253 0.1703 87.322 17.442 AFpz +254 13.622 86.758 15.302 AFp2 +255 24.101 84.377 7.433 AFp4 +256 33.913 81.812 -1.035 AFp6 +257 43.948 79.296 -9.3 AFp8 +258 37.712 72.168 -46.197 AFp10 +259 -59.34 52.68 -48.77 AFF9 +260 -63.262 55.992 -11.173 AFF7 +261 -55.82 61.396 11.884 AFF5 +262 -43.382 66.367 32.811 AFF3 +263 -23.582 69.917 47.293 AFF1 +264 0.2763 71.28 52.092 AFFz +265 25.558 70.556 47.827 AFF2 +266 45.152 67.275 32.731 AFF4 +267 58 62.6 11.9 AFF6 +268 64.673 57.274 -11.46 AFF8 +269 60.601 52.267 -49.038 AFF10 +270 -78.484 28.77 -50.522 FFT9 +271 -76.615 28.653 -11.508 FFT7 +272 -71.506 33.926 20.993 FFC5 +273 -55.94 38.716 49.788 FFC3 +274 -30.655 42.415 71.04 FFC1 +275 0.3512 44.074 79.141 FFCz +276 32.645 43.101 70.795 FFC2 +277 57.504 39.852 48.811 FFC4 +278 74.25 35.5 20.38 FFC6 +279 79.034 30.344 -11.997 FFT8 +280 79.92 28.942 -50.914 FFT10 +281 -87.362 -0.5147 -49.837 FTT9 +282 -82.668 -0.9417 -10.284 FTT7 +283 -80.133 2.5853 27.312 FCC5 +284 -64.161 5.8313 60.885 FCC3 +285 -35.749 8.3091 85.459 FCC1 +286 0.3911 9.508 95.56 FCCz +287 36.07 8.6519 83.832 FCC2 +288 65.164 6.6198 60.052 FCC4 +289 81.544 3.6637 27.201 FCC6 +290 83.168 0.1817 -10.364 FTT8 +291 85.393 -0.9523 -49.52 FTT10 +292 -86.632 -31.238 -47.178 TTP9 +293 -85.933 -31.093 -8.474 TTP7 +294 -81.543 -30.173 30.273 CCP5 +295 -66.128 -29.296 65.898 CCP3 +296 -36.93 -28.57 91.734 CCP1 +297 0.3959 -28.163 101.27 CCPz +298 38.54 -28.225 90.976 CCP2 +299 68.854 -28.64 66.41 CCP4 +300 84.553 -29.378 30.878 CCP6 +301 86 -30.28 -8.435 TTP8 +302 86.762 -31.731 -47.253 TTP10 +303 -80.715 -60.646 -43.594 TPP9 +304 -78.599 -59.724 -4.758 TPP7 +305 -73.664 -61.923 30.38 CPP5 +306 -59.411 -63.925 62.672 CPP3 +307 -32.728 -65.32 85.944 CPP1 +308 0.3658 -65.75 94.058 CPPz +309 35.892 -65.138 85.98 CPP2 +310 62.256 -63.615 62.719 CPP4 +311 76.671 -61.548 30.543 CPP6 +312 79.319 -59.303 -4.84 TPP8 +313 81.56 -61.215 -43.8 TPP10 +314 -64.57 -86.432 -38.324 PPO9 +315 -64.583 -86.222 0.033 PPO7 +316 -58.712 -88.705 25.193 PPO5 +317 -46.16 -90.888 47.446 PPO3 +318 -24.648 -92.292 62.076 PPO1 +319 0.2727 -92.758 67.342 PPOz +320 26.437 -92.295 63.199 PPO2 +321 47.144 -90.712 47.678 PPO4 +322 60.813 -88.504 25.662 PPO6 +323 65.152 -85.943 -0.009 PPO8 +324 65.038 -86.718 -38.448 PPO10 +325 -43.128 -107.52 -32.387 POO9 +326 -42.976 -106.49 5.773 POO7 +327 -36.234 -107.72 17.75 POO5 +328 -25.984 -108.62 26.544 POO3 +329 -13.664 -109.27 32.856 POO1 +330 0.1676 -109.28 32.79 POOz +331 13.651 -109.11 30.936 POO2 +332 26.664 -108.67 26.415 POO4 +333 37.701 -107.84 18.069 POO6 +334 43.67 -106.6 5.726 POO8 +335 43.177 -107.44 -32.463 POO10 +336 -29.391 -114.51 -10.02 OI1 +337 0.0525 -119.34 -3.936 OIz +338 29.553 -113.64 -10.051 OI2 +339 -84.161 -16.019 -9.346 T3 +340 -72.434 -73.453 -2.487 T5 +341 85.08 -15.02 -9.49 T4 +342 73.056 -73.068 -2.54 T6 +343 -86.076 -44.99 -67.986 M1 +344 85.794 -45.009 -68.031 M2 +345 -86.076 -24.99 -67.986 A1 +346 85.794 -25.009 -68.031 A2 diff --git a/leadfield/hartmut/lf_generate_fromhartmut.m b/leadfield/hartmut/lf_generate_fromhartmut.m new file mode 100644 index 0000000..63a38d7 --- /dev/null +++ b/leadfield/hartmut/lf_generate_fromhartmut.m @@ -0,0 +1,144 @@ +% lf = lf_generate_fromhartmut(model[, varargin]) +% +% Generates a leadfield based on one of the HArtMuT head models, +% which includes muscular and ocular sources as well as cortical +% sources. +% +% Assumes you have the required HArtMuT leadfield files in MATLAB +% format in the path. As of 2022-11-16, these are available at +% https://github.com/harmening/HArtMuT +% +% HArtMuT publication: +% Harmening, N., Klug, M., Gramann, K., & Miklody, D. (2022). +% HArtMuT - Modeling eye and muscle contributors in neuroelectric +% imaging. Journal of Neural Engineering. In press. doi: +% 10.1088/1741-2552/aca8ce +% +% HArtMuT is copyright 2022 Nils Harmening and licensed under GNU GPL +% 3. +% +% In: +% model - string indicating which HArtMuT model to use. currently +% supported are: 'mix_Colin27_small', 'NYhead_small' +% +% Optional inputs (key-value pairs): +% labels - cell of electrode labels to be used. default uses all 231 +% available channels (including fiducials). +% montage - name of predefined channel montage. see utl_get_montage. +% +% Out: +% lf - the leadfield containing the following fields +% .leadfield - the leadfield, containing projections in three +% directions (xyz) for each source, in a +% nchannels x nsources x 3 matrix +% .orientation - a default orientation for each soure. for the +% New York Head, this gives dipole orientations +% perpendicular to the cortical surface +% .pos - xyz MNI coordinates of each source +% .chanlocs - channel information in EEGLAB format +% .atlas - nsources x 1 cell with atlas (region) +% indication for each source +% +% Usage example: +% >> lf = lf_generate_fromhartmut('NYhead_small'); +% >> plot_headmodel(lf); +% +% Copyright 2022 Laurens R. Krol +% Neuroadaptive Human-Computer Interaction +% Brandenburg University of Technology + +% 2022-11-16 First version + +% This file is part of Simulating Event-Related EEG Activity (SEREEGA). + +% SEREEGA is free software: you can redistribute it and/or modify +% it under the terms of the GNU General Public License as published by +% the Free Software Foundation, either version 3 of the License, or +% (at your option) any later version. + +% SEREEGA is distributed in the hope that it will be useful, +% but WITHOUT ANY WARRANTY; without even the implied warranty of +% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +% GNU General Public License for more details. + +% You should have received a copy of the GNU General Public License +% along with SEREEGA. If not, see . + +function lf = lf_generate_fromhartmut(model, varargin) + +% parsing input +p = inputParser; + +addRequired(p, 'model', @ischar); + +addParameter(p, 'labels', {}, @iscell); +addParameter(p, 'montage', '', @ischar); + +parse(p, model, varargin{:}) + +model = p.Results.model; +labels = p.Results.labels; +montage = p.Results.montage; + +% setting files to load +if strcmpi(model, 'NYhead_small') + hartmutfile = 'HArtMuT_NYhead_small.mat'; +elseif strcmpi(model, 'mix_Colin27_small') + hartmutfile = 'HArtMuT_mix_Colin27_small.mat'; +else + error('SEREEGA:lf_generate_fromhartmut:invalidFunctionArguments', 'Unknown HArtMuT model ''%s''', model); +end + +% loading leadfields +if exist(hartmutfile, 'file') ~= 2 + error('SEREEGA:lf_generate_fromhartmut:fileNotFound', ['Could not find HArtMuT leadfield file: ' hartmutfile '\nMake sure you have obtained the file, and that MATLAB can find it.\nIt should be available at https://github.com/harmening/HArtMuT']); +else + load(hartmutfile, 'HArtMuT'); +end + +if ~isempty(montage) + % taking channel labels from indicated montage + labels = utl_get_montage(montage); +end + +if isempty(labels) + % taking all available EEG electrodes + labels = HArtMuT.electrodes.label; + + if strcmp(model, 'mix_Colin27_small') + % removing channels with duplicate coordinates + labels = setdiff(labels, {'T3', 'T4', 'T5', 'T6'}); + end +end + +% obtaining indices of indicated electrodes +[~, chanidx] = ismember(labels, HArtMuT.electrodes.label); +if any(~chanidx) + missingchans = find(~chanidx); + warning('\nElectrode %s not available', labels{missingchans}); + chanidx(missingchans) = []; +end + +% preparing output +lf = struct(); +lf.leadfield = cat(2, HArtMuT.cortexmodel.leadfield(chanidx,:,:), HArtMuT.artefactmodel.leadfield(chanidx,:,:)); +lf.orientation = cat(1, HArtMuT.cortexmodel.orientation, HArtMuT.artefactmodel.orientation); +lf.pos = cat(1, HArtMuT.cortexmodel.pos, HArtMuT.artefactmodel.pos); +lf.chanlocs = struct( ... + 'type', HArtMuT.electrodes.chantype', ... + 'labels', HArtMuT.electrodes.label', ... + 'X', num2cell( HArtMuT.electrodes.chanpos(:,2)'), ... + 'Y', num2cell(-HArtMuT.electrodes.chanpos(:,1)'), ... + 'Z', num2cell( HArtMuT.electrodes.chanpos(:,3)')); +lf.chanlocs = lf.chanlocs(chanidx); +lf.atlas = cat(1, strcat('Brain_', HArtMuT.cortexmodel.labels), HArtMuT.artefactmodel.labels); +lf.atlas = utl_sanitize_atlas(lf.atlas); + +% converting chanlocs to EEGLAB format +lf.nbchan = numel(lf.chanlocs); +lf = pop_chanedit(lf, 'convert', {'cart2all'}); +lf = rmfield(lf, {'nbchan', 'chaninfo'}); + +fprintf('Loaded HArtMuT %s leadfield.\n', model); + +end \ No newline at end of file diff --git a/leadfield/lf_add_atlas_frommni2atlas.m b/leadfield/lf_add_atlas_frommni2atlas.m new file mode 100644 index 0000000..550f25c --- /dev/null +++ b/leadfield/lf_add_atlas_frommni2atlas.m @@ -0,0 +1,152 @@ +% leadfield = lf_add_atlas_frommni2atlas(leadfield, atlasselector, varargin) +% +% Returns the given lead field with regions assigned to each source +% according to a selected atlas, using mni2atlas from FSL. +% +% WARNING: Region labels will be approximate and possibly entirely +% inaccurate depending on the lead field and atlas used. +% Make sure to what extent your lead field matches the +% atlases used. +% +% Requires the mni2atlas function from the FMRIB Software Library. +% As of 2021-01-08, this can be downloaded from +% https://github.com/dmascali/mni2atlas +% +% Note that mni2atlas is covered by a different license which, i.a., +% does not allow commercial use. +% +% mni2atlas in turn requires some functions from NIfTI/ANALYZE, +% which, as of 2021-01-08, can be downloaded from +% https://mathworks.com/matlabcentral/fileexchange/8797 +% In particular, the required files are load_nii.m, load_nii_hdr.m, +% load_nii_img.m, and xform_nii.m. +% +% In: +% leadfield - the leadfield to which region labels should be added +% atlasselector - the atlas to be used by mni2atlas. the options are: +% 1 - Juelich Histological Atlas +% 2 - Harvard-Oxford Cortical Structural Atlas +% 3 - Harvard-Oxford Subcortical Structural Atlas +% 4 - JHU ICBM-DTI-81 White Matter labels +% 5 - JHU White Matter tractography Atlas +% 6 - Oxford Thalamic Connectivity Atlas +% 7 - Cerebellar Atlas in MNI152 after FLIRT +% 8 - Cerebellar Atlas in MNI152 after FNIRT +% 9 - MNI Structural Atlas +% +% Optional (key-value pairs): +% region - cell of strings indicating the current region +% labels to be mapped and overwritten +% unknownlabel - string label for sources that cannot be found in +% the atlas, default 'Brain Unknown'. should start +% with one of the generic region labels indicated in +% lf_get_regions. +% +% Out: +% leadfield - the leadfield with region labels added to its atlas +% +% Usage example: +% >> lf = lf_generate_fromnyhead(); +% >> lf = lf_add_atlas_frommni2atlas(lf, 9); % not fully accurate +% +% Copyright 2021, 2022 Laurens R. Krol +% Neuroadaptive Human-Computer Interaction +% Brandenburg University of Technology + +% 2022-11-25 lrk +% - Preallocated variables for efficiency +% 2022-11-17 lrk +% - Switched to inpurParser to parse arguments +% - Added optional region argument +% 2021-01-08 First version + +% This file is part of Simulating Event-Related EEG Activity (SEREEGA). + +% SEREEGA is free software: you can redistribute it and/or modify +% it under the terms of the GNU General Public License as published by +% the Free Software Foundation, either version 3 of the License, or +% (at your option) any later version. + +% SEREEGA is distributed in the hope that it will be useful, +% but WITHOUT ANY WARRANTY; without even the implied warranty of +% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +% GNU General Public License for more details. + +% You should have received a copy of the GNU General Public License +% along with SEREEGA. If not, see . + +function leadfield = lf_add_atlas_frommni2atlas(leadfield, atlasselector, varargin) + +% parsing input +p = inputParser; + +addRequired(p, 'leadfield', @isstruct); +addRequired(p, 'atlasselector', @isnumeric); + +addParameter(p, 'region', {'.*'}, @iscell); +addParameter(p, 'unknownlabel', 'Brain Unknown', @ischar); + +parse(p, leadfield, atlasselector, varargin{:}); + +leadfield = p.Results.leadfield; +atlasselector = p.Results.atlasselector; +currentregion = p.Results.region; +unknownlabel = p.Results.unknownlabel; + +if ~exist('mni2atlas', 'file') + error('SEREEGA:lf_add_atlas_frommni2atlas:fileNotFound', 'Could not find the mni2atlas function file (mni2atlas.m) in the path.\nMake sure you have obtained the mni2atlas software, and that MATLAB can find it.\nIt should be available at https://github.com/dmascali/mni2atlas') +end + +currentregionidx = lf_get_source_all(leadfield, 'region', currentregion); + +atlasnames = { ... + 'Juelich Histological Atlas' ... + 'Harvard-Oxford Cortical Structural Atlas' ... + 'Harvard-Oxford Subcortical Structural Atlas' ... + 'JHU ICBM-DTI-81 White Matter labels' ... + 'JHU White Matter tractography Atlas' ... + 'Oxford Thalamic Connectivity Atlas' ... + 'Cerebellar Atlas in MNI152 after FLIRT' ... + 'Cerebellar Atlas in MNI152 after FNIRT' ... + 'MNI Structural Atlas' ... + }; + +numsources = numel(currentregionidx); +atlas = cell(numsources, 1); +unknowns = 0; +confidence = nan(numsources, 1); + +w = waitbar(0, {sprintf('Atlas: %s', atlasnames{atlasselector}), sprintf('Source 0 of %d', numsources), 'Average confidence: 0.0%%', 'Unkowns: 0 (0.0%%)'}, 'Name', 'Looking up sources'); + +for s = 1:numsources + + idx = currentregionidx(s); + + if mod(s, 250) == 0 || s == numsources + % updating waitbar + waitbar(s/numsources, w, {sprintf('Atlas: %s', atlasnames{atlasselector}), sprintf('Source %d of %d', s, numsources), sprintf('Average confidence: %.1f%%', nanmean(confidence)), sprintf('%d unkowns (%3.1f%%)', unknowns, unknowns/s*100)}); + end + + % getting source region + sr = mni2atlas(leadfield.pos(idx,:), atlasselector); + + if isempty(sr.label) + % no result; using default + atlas(s) = {unknownlabel}; + unknowns = unknowns + 1; + else + % parsing string; updating confidence and atlas + sr = strsplit(sr.label{1}, '% '); + confidence(s) = str2double(sr{1}); + atlas(s) = {['Brain ' sr{2}]}; + end + +end + +atlas = utl_sanitize_atlas(atlas); +leadfield.atlas(currentregionidx) = atlas; +delete(w); + +fprintf('Done.\nAtlas: %s\nAverage confidence: %.1f%%\n%d unkowns (%3.1f%%)\n', atlasnames{atlasselector}, mean(confidence), unknowns, unknowns/numsources*100); + +end \ No newline at end of file diff --git a/leadfield/lf_add_source.m b/leadfield/lf_add_source.m index 3c1256b..cc9a076 100644 --- a/leadfield/lf_add_source.m +++ b/leadfield/lf_add_source.m @@ -11,6 +11,7 @@ % patterns, where n equals the number of channels, and % the three columns represent the x, y, and z projection % patterns, respectively +% region - string naming the source's corresponding atlas region % % Optional inputs (key-value pairs): % orientation - 1-by-3 matrix containing default xyz source @@ -20,13 +21,17 @@ % lf - the leadfield with the new source added % % Usage example: -% >> lf = lf_generate_fromnyhead('labels', {'Fz', 'Cz', 'Pz'}) -% >> lf = lf_add_source(lf, [0 0 0], rand(3,3)) +% >> lf = lf_generate_fromnyhead('labels', {'Fz', 'Cz', 'Pz', 'Oz'}); +% >> lf = lf_add_source(lf, [0 0 0], rand(4,3), 'Brain_NewRegion') % -% Copyright 2017 Laurens R Krol +% Copyright 2017, 2022 Laurens R. Krol % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% Neuroadaptive Human-Computer Interaction +% Brandenburg University of Technology +% 2022-12-02 lrk +% - Added region argument % 2017-09-29 First version % This file is part of Simulating Event-Related EEG Activity (SEREEGA). @@ -44,43 +49,61 @@ % You should have received a copy of the GNU General Public License % along with SEREEGA. If not, see . -function lf = lf_add_source(leadfield, position, projection, varargin) - -% parsing input -p = inputParser; - -addRequired(p, 'leadfield', @isstruct); -addRequired(p, 'position', @isnumeric); -addRequired(p, 'projection', @isnumeric); - -addOptional(p, 'orientation', [0 0 0], @isnumeric); - -parse(p, leadfield, position, projection, varargin{:}) - -lf = p.Results.leadfield; -position = p.Results.position; -projection = p.Results.projection; -orientation = p.Results.orientation; - -if iscolumn(position), position = position'; end -if ~all(size(position) == [1, 3]) - error('SEREEGA:lf_add_source:error', 'position variable should be 1-by-3 matrix'); -end +function lf = lf_add_source(leadfield, position, projection, region, varargin) + +if strcmp(region, 'orientation') + + % maintaining backwards compatibility with SEREEGA <= v1.2.2, where + % the region argument did not exist and could now be confused with the + % optional orientation argument + + warning('No atlas region indicated; defaulting to ''Unknown'''); + + lf = lf_add_source(leadfield, position, projection, 'Unknown', 'orientation', varargin{1}); + +else + + % parsing input + p = inputParser; + + addRequired(p, 'leadfield', @isstruct); + addRequired(p, 'position', @isnumeric); + addRequired(p, 'projection', @isnumeric); + addRequired(p, 'region', @ischar); + + addOptional(p, 'orientation', [0 0 0], @isnumeric); + + parse(p, leadfield, position, projection, region, varargin{:}) + + lf = p.Results.leadfield; + position = p.Results.position; + projection = p.Results.projection; + region = p.Results.region; + orientation = p.Results.orientation; + + if iscolumn(position), position = position'; end + if ~all(size(position) == [1, 3]) + error('SEREEGA:lf_add_source:error', 'position variable should be 1-by-3 matrix'); + end -if ~all(size(projection) == [size(lf.leadfield, 1), 3]) - error('SEREEGA:lf_add_source:error', 'projection variable should be %d-by-3 matrix', size(lf.leadfield, 1)); -end + if ~all(size(projection) == [size(lf.leadfield, 1), 3]) + error('SEREEGA:lf_add_source:error', 'projection variable should be %d-by-3 matrix', size(lf.leadfield, 1)); + end -if ~isempty(orientation) - if iscolumn(orientation), orientation = orientation'; end - if ~all(size(orientation) == [1, 3]) - error('SEREEGA:lf_add_source:error', 'orientation variable should be 1-by-3 matrix'); + if ~isempty(orientation) + if iscolumn(orientation), orientation = orientation'; end + if ~all(size(orientation) == [1, 3]) + error('SEREEGA:lf_add_source:error', 'orientation variable should be 1-by-3 matrix'); + end end -end -% adding source to leadfield -lf.pos(end+1,:) = position; -lf.orientation(end+1,:) = orientation; -lf.leadfield(:,end+1,:) = projection; + % adding source to leadfield + lf.pos(end+1,:) = position; + lf.orientation(end+1,:) = orientation; + lf.leadfield(:,end+1,:) = projection; + lf.atlas = [lf.atlas; {region}]; + lf.atlas = utl_sanitize_atlas(lf.atlas); +end + end \ No newline at end of file diff --git a/leadfield/lf_get_regions.m b/leadfield/lf_get_regions.m new file mode 100644 index 0000000..84cd9b5 --- /dev/null +++ b/leadfield/lf_get_regions.m @@ -0,0 +1,75 @@ +% [allregions, numall, generic, numgeneric] = lf_get_regions(leadfield) +% +% Returns a cell of generic region categories present in the given +% leadfield's atlas, as well as a list of all unique regions. +% +% In: +% leadfield - the leadfield from which to get the regions +% +% Out: +% allregions - cell of strings listing all unique regions present in +% the atlas +% numall - numeric list indicating how many sources are present in +% each of the regions in the allregions list +% generic - cell of strings indicating generic region categories that +% are present in the atlas. currently recognised are: +% - brain +% - muscle +% - eye +% numgeneric - numeric list indicating how many sources are present +% in each of the generic region categories +% +% Usage example: +% >> lf = lf_generate_fromnyhead(); +% >> generic = lf_get_regions_fromatlas(lf); +% +% Copyright 2021 Laurens R. Krol +% Neuroadaptive Human-Computer Interaction +% Brandenburg University of Technology + +% 2022-11-25 lrk +% - Renamed from utl_get_regions_fromatlas to lf_get_regions +% 2022-11-18 lrk +% - Changed output order +% 2021-09-07 lrk +% - Improved efficiency for large cell arrays +% 2021-01-07 First version + +% This file is part of Simulating Event-Related EEG Activity (SEREEGA). + +% SEREEGA is free software: you can redistribute it and/or modify +% it under the terms of the GNU General Public License as published by +% the Free Software Foundation, either version 3 of the License, or +% (at your option) any later version. + +% SEREEGA is distributed in the hope that it will be useful, +% but WITHOUT ANY WARRANTY; without even the implied warranty of +% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +% GNU General Public License for more details. + +% You should have received a copy of the GNU General Public License +% along with SEREEGA. If not, see . + +function [allregions, numall, generic, numgeneric] = lf_get_regions(leadfield) + +% getting all unique regions +[allregions, ~, ic] = unique(leadfield.atlas); +allregions = allregions'; +numall = histcounts(ic, numel(allregions)); + +% looking for recognised region categories +cats = {'brain', 'eye', 'muscle'}; + +generic = {}; +numgeneric = []; +for c = cats + % finding any generic regions in list of unique regions + idx = cellfun(@(x) ~isempty(x), regexpi(allregions, ['^' c{1} '.*'], 'once')); + if any(idx) + % counting all sources in generic regions + generic = [generic, c]; + numgeneric = [numgeneric, sum(~cellfun(@isempty, regexpi(leadfield.atlas, ['^' c{1} '.*'], 'once')))]; + end +end + +end \ No newline at end of file diff --git a/leadfield/lf_get_source_all.m b/leadfield/lf_get_source_all.m new file mode 100644 index 0000000..c74b949 --- /dev/null +++ b/leadfield/lf_get_source_all.m @@ -0,0 +1,82 @@ +% sourceIdx = lf_get_source_all(leadfield, varargin) +% +% Returns all source indices in the lead field, optionally +% corresponding to the regions indicated using non-case-sensitive +% regular expressions. +% +% In: +% leadfield - the leadfield from which to get the random source +% +% Optional (key-value pairs): +% region - cell containing strings and/or regex patterns representing +% leadfield.atlas entries. not case sensitive. default: '.*' +% +% Out: +% sourceIdx - the selected source index/indices +% +% Usage example: +% >> lf = lf_generate_fromnyhead(); +% >> sourceIdx = lf_get_source_all(lf, 'region', {'Brain.*'}); +% >> plot_source_location(sourceIdx, lf); +% +% Copyright 2021 Laurens R. Krol +% Neuroadaptive Human-Computer Interaction +% Brandenburg University of Technology + +% 2022-11-30 lrk +% - Escaped special characters when exact match found +% - Improved efficiency for large atlases +% 2022-11-17 lrk +% - Shortened name to lf_get_source_all +% - Switched to inputParser with 'region' as optional parameter +% 2021-01-06 First version + +% This file is part of Simulating Event-Related EEG Activity (SEREEGA). + +% SEREEGA is free software: you can redistribute it and/or modify +% it under the terms of the GNU General Public License as published by +% the Free Software Foundation, either version 3 of the License, or +% (at your option) any later version. + +% SEREEGA is distributed in the hope that it will be useful, +% but WITHOUT ANY WARRANTY; without even the implied warranty of +% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +% GNU General Public License for more details. + +% You should have received a copy of the GNU General Public License +% along with SEREEGA. If not, see . + +function sourceIdx = lf_get_source_all(leadfield, varargin) + +% parsing input +p = inputParser; + +addRequired(p, 'leadfield', @isstruct); + +addParameter(p, 'region', {'.*'}, @iscell); + +parse(p, leadfield,varargin{:}) + +leadfield = p.Results.leadfield; +region = p.Results.region; + +[allregions, ~, ~, ~] = lf_get_regions(leadfield); + +idx = []; +for r = region + if any(strcmp(r{1}, allregions)) + % since regions are assessed using regular expressions, + % 'EyeRetina_Choroid_Sclera_left', for example, also matches + % 'EyeRetina_Choroid_Sclera_leftright'; therefore, when an exact + % match is detected, this is enforced. it can be circumvented by + % using regex wildcards + r{1} = regexprep(r{1}, '([\\\(\)\[\]\{\}\,\.\+\?\^\<\>\!\=])', '\\$1'); + r{1} = ['^' r{1} '$']; + fprintf('Assuming exact match: %s\n', r{1}); + end + + idx = [idx, ~cellfun(@isempty, regexpi(leadfield.atlas, r{1}))]; +end +sourceIdx = find(sum(idx, 2) > 0)'; + +end \ No newline at end of file diff --git a/leadfield/lf_get_source_inradius.m b/leadfield/lf_get_source_inradius.m index 6c6cd63..5f895c9 100644 --- a/leadfield/lf_get_source_inradius.m +++ b/leadfield/lf_get_source_inradius.m @@ -1,4 +1,4 @@ -% sourceIdx = lf_get_source_inradius(leadfield, centre, radius) +% sourceIdx = lf_get_source_inradius(leadfield, centre, radius, varargin) % % Returns the source(s) in the leadfield within a certain radius of % an indicated source or coordinate. @@ -17,10 +17,15 @@ % >> sourceIdx = lf_get_source_inradius(lf, [0 0 0], 10); % >> plot_source_location(sourceIdx, lf); % -% Copyright 2017 Laurens R Krol +% Copyright 2017, 2022 Laurens R. Krol % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% Neuroadaptive Human-Computer Interaction +% Brandenburg University of Technology +% 2022-11-17 lrk +% - Switched to inpurParser to handle arguments +% - Added optional 'region' argument % 2017-04-27 First version % This file is part of Simulating Event-Related EEG Activity (SEREEGA). @@ -38,12 +43,31 @@ % You should have received a copy of the GNU General Public License % along with SEREEGA. If not, see . -function sourceIdx = lf_get_source_inradius(leadfield, centre, radius) +function sourceIdx = lf_get_source_inradius(leadfield, centre, radius, varargin) + +% parsing input +p = inputParser; + +addRequired(p, 'leadfield', @isstruct); +addRequired(p, 'centre', @isnumeric); +addRequired(p, 'radius', @isnumeric); + +addParameter(p, 'region', {'.*'}, @iscell); + +parse(p, leadfield, centre, radius, varargin{:}) + +leadfield = p.Results.leadfield; +centre = p.Results.centre; +radius = p.Results.radius; +region = p.Results.region; + + +regionIdx = lf_get_source_all(leadfield, 'region', region); if numel(centre) == 1, centre = leadfield.pos(centre,:); end sourceIdx = []; -for p = 1:size(leadfield.leadfield, 2) +for p = regionIdx if sqrt( ... (leadfield.pos(p,1) - centre(1))^2 + ... (leadfield.pos(p,2) - centre(2))^2 + ... diff --git a/leadfield/lf_get_source_middle.m b/leadfield/lf_get_source_middle.m new file mode 100644 index 0000000..888bdbd --- /dev/null +++ b/leadfield/lf_get_source_middle.m @@ -0,0 +1,80 @@ +% sourceIdx = lf_get_source_middle(leadfield, region) +% +% Returns the source nearest to either the average coordinates of all +% sources in the the indicated region(s), or to the average of their +% boundaries. +% +% Note that regardless of method, this may or may not be an actual +% 'middle'. +% +% In: +% leadfield - the leadfield from which to get the source +% +% Optional (key-value pairs): +% region - cell containing strings and/or regex patterns representing +% leadfield.atlas entries. not case sensitive. +% method - ('mean'|'minmax') indicating the method how to +% calculate the 'middle'. 'average' takes the mean of all +% coordinate. 'minmax' takes the mean of the minimum and +% maximum values along each axis. default: average +% +% Out: +% sourceIdx - the nearest-to-average source index +% +% Usage example: +% >> lf = lf_generate_fromnyhead(); +% >> sourceIdx = lf_get_source_middle(lf, 'region', {'Brain.*'}); +% >> plot_source_location(sourceIdx, lf); +% +% Copyright 2021 Laurens R. Krol +% Neuroadaptive Human-Computer Interaction +% Brandenburg University of Technology + +% 2022-11-17 First version + +% This file is part of Simulating Event-Related EEG Activity (SEREEGA). + +% SEREEGA is free software: you can redistribute it and/or modify +% it under the terms of the GNU General Public License as published by +% the Free Software Foundation, either version 3 of the License, or +% (at your option) any later version. + +% SEREEGA is distributed in the hope that it will be useful, +% but WITHOUT ANY WARRANTY; without even the implied warranty of +% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +% GNU General Public License for more details. + +% You should have received a copy of the GNU General Public License +% along with SEREEGA. If not, see . + +function sourceIdx = lf_get_source_middle(leadfield, varargin) + +% parsing input +p = inputParser; + +addRequired(p, 'leadfield', @isstruct); + +addParameter(p, 'region', {'.*'}, @iscell); +addParameter(p, 'method', 'average', @ischar); + +parse(p, leadfield, varargin{:}) + +leadfield = p.Results.leadfield; +region = p.Results.region; +method = p.Results.method; + +sourceIdx = lf_get_source_all(leadfield, 'region', region); + +if numel(sourceIdx) < 3 + error('SEREEGA:lf_get_source_middle:error', 'no ''middle'' can be determined from < 3 sources'); +end + +if strcmp(method, 'average') + sourceIdx = lf_get_source_nearest(leadfield, mean(leadfield.pos(sourceIdx, :)), 'region', region); +elseif strcmp(method, 'minmax') + sourceIdx = lf_get_source_nearest(leadfield, mean([max(leadfield.pos(sourceIdx, :)); min(leadfield.pos(sourceIdx, :))]), 'region', region); +else + error('SEREEGA:lf_get_source_middle:invalidFunctionArguments', 'unknown method ''%s''', method); +end + +end diff --git a/leadfield/lf_get_source_nearest.m b/leadfield/lf_get_source_nearest.m index 21a6135..6a2656a 100644 --- a/leadfield/lf_get_source_nearest.m +++ b/leadfield/lf_get_source_nearest.m @@ -1,12 +1,18 @@ -% [sourceIdx, dist] = lf_get_source_nearest(leadfield, pos) +% [sourceIdx, dist] = lf_get_source_nearest(leadfield, pos, varargin) % % Returns the source in the leadfield nearest to the given position, -% and its distance from that position. +% and its distance from that position. The returned source can +% optionally be constrained to indicated region(s), using +% non-case-sensitive regular expressions. % % In: % leadfield - the leadfield from which to get the source % pos - 1-by-3 matrix of x, y, z coordinates % +% Optional (key-value pairs): +% region - cell containing strings and/or regex patterns representing +% leadfield.atlas entries. not case sensitive. default: .* +% % Out: % sourceIdx - the nearest source index % dist - the distance of the found source to the indicated position @@ -16,10 +22,15 @@ % >> sourceIdx = lf_get_source_nearest(lf, [0 0 0]); % >> plot_source_location(sourceIdx, lf); % -% Copyright 2017 Laurens R Krol +% Copyright 2017, 2022 Laurens R. Krol % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% Neuroadaptive Human-Computer Interaction +% Brandenburg University of Technology +% 2022-11-17 lrk +% - Switched to inpurParser to handle arguments +% - Added optional 'region' argument % 2017-04-27 First version % This file is part of Simulating Event-Related EEG Activity (SEREEGA). @@ -37,13 +48,30 @@ % You should have received a copy of the GNU General Public License % along with SEREEGA. If not, see . -function [sourceIdx, dist] = lf_get_source_nearest(leadfield, pos) +function [sourceIdx, dist] = lf_get_source_nearest(leadfield, pos, varargin) + +% parsing input +p = inputParser; + +addRequired(p, 'leadfield', @isstruct); +addRequired(p, 'pos', @isnumeric); + +addParameter(p, 'region', {'.*'}, @iscell); + +parse(p, leadfield, pos, varargin{:}) + +leadfield = p.Results.leadfield; +pos = p.Results.pos; +region = p.Results.region; + +regionIdx = lf_get_source_all(leadfield, 'region', region); distances = sqrt( ... - (leadfield.pos(:,1) - pos(1)).^2 + ... - (leadfield.pos(:,2) - pos(2)).^2 + ... - (leadfield.pos(:,3) - pos(3)).^2); + (leadfield.pos(regionIdx,1) - pos(1)).^2 + ... + (leadfield.pos(regionIdx,2) - pos(2)).^2 + ... + (leadfield.pos(regionIdx,3) - pos(3)).^2); -[dist, sourceIdx] = min(distances); +[dist, i] = min(distances); +sourceIdx = regionIdx(i); end \ No newline at end of file diff --git a/leadfield/lf_get_source_random.m b/leadfield/lf_get_source_random.m index b6d6e6b..6ef6077 100644 --- a/leadfield/lf_get_source_random.m +++ b/leadfield/lf_get_source_random.m @@ -1,25 +1,33 @@ % sourceIdx = lf_get_source_random(leadfield[, number]) % -% Returns (a) random source index/indices. +% Returns (a) random source index/indices, optionally constrained to +% indicated regions. % % In: % leadfield - the leadfield from which to get the random source % -% Optional: +% Optional (key-value pairs): % number - number of sources to return. default: 1 +% region - cell containing strings and/or regex patterns representing +% leadfield.atlas entries. not case sensitive. default: .* % % Out: % sourceIdx - the random source index/indices % % Usage example: % >> lf = lf_generate_fromnyhead(); -% >> sourceIdx = lf_get_source_random(lf, 5); +% >> sourceIdx = lf_get_source_random(lf, 'number', 5); % >> plot_source_location(sourceIdx, lf); % -% Copyright 2017 Laurens R Krol +% Copyright 2017, 2022 Laurens R. Krol % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% Neuroadaptive Human-Computer Interaction +% Brandenburg University of Technology +% 2022-11-17 lrk +% - Switched to inpurParser to handle arguments +% - Added optional 'region' argument % 2017-04-27 First version % This file is part of Simulating Event-Related EEG Activity (SEREEGA). @@ -37,11 +45,38 @@ % You should have received a copy of the GNU General Public License % along with SEREEGA. If not, see . -function sourceIdx = lf_get_source_random(leadfield, num) +function sourceIdx = lf_get_source_random(leadfield, varargin) -if ~exist('num', 'var'), num = 1; end + +if ~isempty(varargin) && isnumeric(varargin{1}) + + % maintaining backwards compatibility with SEREEGA <= v1.2.2, where + % 'number' was a simple optional argument and could be called as e.g. + % lf_get_source_random(leadfield, 5) -sourceIdx = randperm(size(leadfield.leadfield, 2)); -sourceIdx = sourceIdx(1:num); + sourceIdx = lf_get_source_random(leadfield, 'number', varargin{:}); + +else + % parsing input + p = inputParser; + + addRequired(p, 'leadfield', @isstruct); + + addParameter(p, 'number', 1, @isnumeric); + addParameter(p, 'region', {'.*'}, @iscell); + + parse(p, leadfield, varargin{:}) + + leadfield = p.Results.leadfield; + number = p.Results.number; + region = p.Results.region; + + regionIdx = lf_get_source_all(leadfield, 'region', region); + + sourceIdx = regionIdx(randperm(numel(regionIdx))); + sourceIdx = sourceIdx(1:number); + +end + end \ No newline at end of file diff --git a/leadfield/lf_get_source_spaced.m b/leadfield/lf_get_source_spaced.m index 28b3a64..8caf4b5 100644 --- a/leadfield/lf_get_source_spaced.m +++ b/leadfield/lf_get_source_spaced.m @@ -17,6 +17,8 @@ % these must not satisfy the spacing criterion among each % other, but all randomly added sources will be at least % the indicated distance away from these and each other. +% region - cell containing strings and/or regex patterns representing +% leadfield.atlas entries. not case sensitive. default: .* % % Out: % sourceIdx - 1-by-number array containing the source indices of the @@ -27,10 +29,14 @@ % >> sourceIdx = lf_get_source_spaced(lf, 10, 50); % >> plot_source_location(sourceIdx, lf); % -% Copyright 2017 Laurens R Krol +% Copyright 2017, 2022 Laurens R. Krol % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% Neuroadaptive Human-Computer Interaction +% Brandenburg University of Technology +% 2022-11-17 lrk +% - Added optional 'region' argument % 2017-08-01 lrk % - Added option to start with given source array % 2017-06-12 First version @@ -60,35 +66,28 @@ addRequired(p, 'spacing', @(x) isnumeric(x) && ~isempty(x)); addParameter(p, 'sourceIdx', [], @isnumeric); +addParameter(p, 'region', {'.*'}, @iscell); parse(p, leadfield, number, spacing, varargin{:}) leadfield = p.Results.leadfield; number = p.Results.number; spacing = p.Results.spacing; -sourceIdx = p.Results.sourceIdx; +initialSourceIdx = p.Results.sourceIdx; +region = p.Results.region; -if ~all(sourceIdx <= size(leadfield.pos, 1)) +if ~all(initialSourceIdx <= size(leadfield.pos, 1)) warning('not all indicated sources are in the leadfield.'); - sourceIdx = sourceIdx(sourceIdx <= size(leadfield.pos, 1)); + initialSourceIdx = initialSourceIdx(initialSourceIdx <= size(leadfield.pos, 1)); end -% randomising sources -sources = randperm(size(leadfield.pos, 1)); +% getting selected region sources +regionIdx = lf_get_source_all(leadfield, 'region', region); -if isempty(sourceIdx) - % starting with first random source, beginning search at second - sourceIdx = sources(1); - rndsource = 2; -else - % removing given sources from source list, beginning search at first - idx = ismember(sources, sourceIdx); - sources(idx) = []; - rndsource = 1; -end - -% running through remaining sources +% starting search for resulting source indices iteration = 1; +[sources, rndsource] = local_initialise_search(regionIdx, initialSourceIdx); +sourceIdx = initialSourceIdx; while numel(sourceIdx) < number % calculating distance of next source to previously selected sources distances = zeros(numel(sourceIdx), 1); @@ -111,11 +110,30 @@ % the required number of spaced sources if rndsource > numel(sources) warning('could not find %d sources with %d mm spacing; re-randomising (%d)...', number, spacing, iteration); - sources = randperm(size(leadfield.pos, 1)); - sourceIdx = sources(1); - rndsource = 2; + [sources, rndsource] = local_initialise_search(regionIdx, initialSourceIdx); iteration = iteration + 1; end end +end + + +function [sources, rndsource] = local_initialise_search(regionIdx, initialSourceIdx) + +% randomising all eligible sources +sources = regionIdx(randperm(numel(regionIdx))); + +if isempty(initialSourceIdx) + % starting with first random source, beginning search at second + initialSourceIdx = sources(1); + rndsource = 2; +else + % moving given sources from source to the beginning of the list, then + % beginning search at first random source + idx = ismember(sources, initialSourceIdx); + sources(idx) = []; + sources = [initialSourceIdx, sources]; + rndsource = numel(initialSourceIdx) + 1; +end + end \ No newline at end of file diff --git a/leadfield/nyhead/lf_generate_fromnyhead.m b/leadfield/nyhead/lf_generate_fromnyhead.m index dc6696c..76211f7 100644 --- a/leadfield/nyhead/lf_generate_fromnyhead.m +++ b/leadfield/nyhead/lf_generate_fromnyhead.m @@ -30,6 +30,8 @@ % perpendicular to the cortical surface % .pos - xyz MNI coordinates of each source % .chanlocs - channel information in EEGLAB format +% .atlas - nsources x 1 cell with atlas (region) +% indication for each source % % Usage example: % >> lf = lf_generate_fromnyhead('labels', {'Fz', 'Cz', 'Pz'}); @@ -39,6 +41,11 @@ % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% 2022-11-24 lrk +% - Added atlas support from the updated sa_nyhead.mat +% 2021-01-06 lrk +% - Added support for output of lf_prune_nyheadfile +% - Added lf.atlas % 2018-11-08 lrk % - Now continues with empty labels if montage not found % 2017-08-10 lrk @@ -76,7 +83,10 @@ montage = p.Results.montage; % loading the NY Head leadfield -if exist('sa_nyhead.mat') ~= 2 +if exist('sa_nyhead_sereegareduced.mat') + % loading reduced-size version if available, see lf_prune_nyheadfile + load('sa_nyhead_sereegareduced.mat', 'sa'); +elseif exist('sa_nyhead.mat') ~= 2 error('SEREEGA:lf_generate_fromnyhead:fileNotFound', 'Could not find ICBM-NY leadfield file (sa_nyhead.mat) in the path.\nMake sure you have obtained the file, and that MATLAB can find it.\nIt should be available at https://parralab.org/nyhead') else load('sa_nyhead.mat', 'sa'); @@ -107,5 +117,12 @@ lf.pos = sa.cortex75K.vc; lf.chanlocs = readlocs('chanlocs-nyhead231.elp'); lf.chanlocs = lf.chanlocs(chanidx); +if isfield(sa, 'HO_labels') && isfield(sa.cortex75K, 'in_HO') + lf.atlas = strcat('Brain', {' '}, sa.HO_labels(sa.cortex75K.in_HO)); + lf.atlas = utl_sanitize_atlas(lf.atlas); +else + warning('No atlas found; defaulting to Brain_CorticalSurface.\nThere is a newer version of the NY Head that includes an atlas.\nIt should be available at https://parralab.org/nyhead', ''); + lf.atlas = repmat({'Brain_CorticalSurface'}, size(lf.pos, 1), 1); +end end \ No newline at end of file diff --git a/leadfield/nyhead/lf_prune_nyheadfile.m b/leadfield/nyhead/lf_prune_nyheadfile.m index 12819c8..d6f6292 100644 --- a/leadfield/nyhead/lf_prune_nyheadfile.m +++ b/leadfield/nyhead/lf_prune_nyheadfile.m @@ -2,7 +2,7 @@ % % Generates a smaller leadfield file from the original New York Head % (ICBM-NY) file, containing only the information used by SEREEGA. -% This speeds up loading times. The new file will be +% This speeds up loading times. The new file will be called % 'sa_nyhead_sereegareduced.mat', written to the same directory where % 'sa_nyhead.mat' is located. lf_generate_fromnyhead will try % to find the reduced file first. @@ -24,12 +24,14 @@ % takes the first instance it finds. % % Usage example: -% >> lf_nyhead_shrinkfile('C:\Data'); +% >> lf_prune_nyheadfile('C:\SEREEGA\leadfield\nyhead'); % -% Copyright 2021 Laurens R Krol +% Copyright 2021, 2022 Laurens R Krol % Neuroadaptive Human-Computer Interaction % Brandenburg University of Technology +% 2022-11-24 lrk +% - Added atlas compatibility % 2021-01-05 First version % This file is part of Simulating Event-Related EEG Activity (SEREEGA). @@ -65,6 +67,11 @@ function lf_prune_nyheadfile(directory) sa.cortex75K.V_fem = orig.sa.cortex75K.V_fem; sa.cortex75K.normals = orig.sa.cortex75K.normals; sa.cortex75K.vc = orig.sa.cortex75K.vc; + + if isfield(orig.sa, 'HO_labels') && isfield(orig.sa.cortex75K, 'in_HO') + sa.HO_labels = orig.sa.HO_labels; + sa.cortex75K.in_HO = orig.sa.cortex75K.in_HO; + end newfilepath = fullfile(fileparts(filepath), 'sa_nyhead_sereegareduced.mat'); save(newfilepath, 'sa'); diff --git a/leadfield/pha/lf_generate_frompha.m b/leadfield/pha/lf_generate_frompha.m index 2a04ea0..05dfd20 100644 --- a/leadfield/pha/lf_generate_frompha.m +++ b/leadfield/pha/lf_generate_frompha.m @@ -37,7 +37,7 @@ % % In: % atlas - string indicating the atlas to use: '0to2', '4to8', or -% '8to12'. +% '8to18'. % layout - string indicating which layout (i.e. number of channels) % to use: '128', '256', or '2562'. note: '256' is not % available for atlas '0to2'. @@ -57,6 +57,9 @@ % .pos - xyz coordinates of each source (not scaled to % MNI but may be approximate, depending on atlas) % .chanlocs - channel information in EEGLAB format +% .atlas - atlas (region) indication for each source; this +% is simply 'Brain' until more specific +% information is added % % Usage example: % >> lf = lf_generate_frompha('8to18', '2562'); @@ -66,6 +69,9 @@ % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% 2021-01-11 lrk +% - Added lf.atlas +% - Added warning about coordinate system changes % 2018-11-08 lrk % - Now continues with empty labels if montage not found % 2018-03-22 jpaw @@ -183,11 +189,14 @@ pos(i,:) = [-pos(i,2), pos(i,1), pos(i,3)]; end +warning(sprintf('PHA does not use a standardised coordinate system. SEREEGA has\nshifted and rotated all coordinates into a more convenient, but still nonstandard system')); + % preparing output lf = struct(); lf.leadfield = leadfield; lf.orientation = zeros(size(leadfield, 2), 3); lf.pos = pos; lf.chanlocs = chanlocs; +lf.atlas = repmat({'Brain'}, size(lf.pos, 1), 1); end diff --git a/plot/plot_headmodel.m b/plot/plot_headmodel.m index 16f0aaa..f9464bb 100644 --- a/plot/plot_headmodel.m +++ b/plot/plot_headmodel.m @@ -1,14 +1,17 @@ % h = plot_headmodel(leadfield, varargin) % % Plots both the sources and the electrodes from the given lead -% field, illustrating the head model as a whole. The sources can be -% plotted either all individually as points in a scatter plot, or as -% a 3D boundary estimated from these points. +% field, illustrating the model as a whole or selected regions. The +% sources can be plotted either all individually as points in a +% scatter plot, or as a 3D boundary estimated from these points. % % In: % leadfield - the lead field to visualise % % Optional inputs (key-value pairs): +% newfig - (0|1) whether or not to open a new figure window. +% default: 1 +% electrodes - whether or not to plot electrodes (0|1, default: 1) % labels - whether or not to plot electrode labels (0|1, default: 1) % style - 'scatter' plots all brain sources individually; % 'boundary' plots the boundary surface of these sources @@ -20,18 +23,28 @@ % [ 0, 90] = axial % [90, 0] = sagittal % [ 0, 0] = coronal (default: [120, 20] +% region - cell of strings representing leadfield.atlas entries +% indicating which regions should be plotted. default {} +% plots all sources. % % Out: % h - the handle of the generated figure % % Usage example: -% >> lf = lf_generate_frompha('8to18', '2562'); +% >> lf = lf_generate_fromnyhead(); % >> plot_headmodel(lf, 'labels', 0, 'style', 'boundary'); % % Copyright 2018 Laurens R Krol % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% 2022-12-06 lrk +% - Changed length(color) to size(color,1) for when there's only 1 source +% 2022-12-02 lrk +% - Added electrodes argument +% 2021-01-06 lrk +% - Added newfig argument +% - Added region argument % 2018-03-23 First version % This file is part of Simulating Event-Related EEG Activity (SEREEGA). @@ -56,21 +69,31 @@ addRequired(p, 'leadfield', @isstruct); +addParameter(p, 'newfig', 1, @isnumeric); +addParameter(p, 'electrodes', 1, @isnumeric); addParameter(p, 'labels', 1, @isnumeric); addParameter(p, 'style', 'scatter', @ischar); addParameter(p, 'shrink', 1, @isnumeric); addParameter(p, 'view', [120 20], @isnumeric); +addParameter(p, 'region', {}, @iscell); parse(p, leadfield, varargin{:}) leadfield = p.Results.leadfield; +newfig = p.Results.newfig; +electrodes = p.Results.electrodes; labels = p.Results.labels; style = p.Results.style; shrink = p.Results.shrink; viewpoint = p.Results.view; +region = p.Results.region; % drawing figure -h = figure('name', 'Head model', 'NumberTitle', 'off'); +if newfig + h = figure('name', 'Head model', 'NumberTitle', 'off'); +else + h = NaN; +end axis equal; xlabel('X'); ylabel('Y'); @@ -78,21 +101,40 @@ view(viewpoint); hold; -% plotting brain +% selecting sources based on atlas +regionIdx = 1:size(leadfield.pos, 1); +if ~isempty(region) + if isfield(leadfield, 'atlas') && ~isempty(leadfield.atlas) + idx = lf_get_source_all(leadfield, 'region', region); + if any(idx) + regionIdx = idx; + else + warning('indicated region(s) not found in the lead field''s atlas; plotting all sources'); + end + else + warning('no atlas information present in the lead field; plotting all sources'); + end +end + +% plotting regions if strcmp(style, 'scatter') - color = [ones(1, size(leadfield.pos, 1)); linspace(.3, .7, size(leadfield.pos,1)); linspace(.3, .7, size(leadfield.pos,1))]'; - color = color(randperm(length(color)),:); - scatter3(leadfield.pos(:,1), leadfield.pos(:,2), leadfield.pos(:,3), 10, color, 'filled'); + color = [ones(1, numel(regionIdx)); linspace(.3, .7, numel(regionIdx)); linspace(.3, .7, numel(regionIdx))]'; + color = color(randperm(size(color, 1)),:); + scatter3(leadfield.pos(regionIdx,1), leadfield.pos(regionIdx,2), leadfield.pos(regionIdx,3), 10, color, 'filled'); elseif strcmp(style, 'boundary') - k = boundary(leadfield.pos(:,1), leadfield.pos(:,2), leadfield.pos(:,3), shrink); - trisurf(k, leadfield.pos(:,1), leadfield.pos(:,2), leadfield.pos(:,3), 'FaceColor', [1 .75, .75], 'EdgeColor', 'none'); + k = boundary(leadfield.pos(regionIdx,1), leadfield.pos(regionIdx,2), leadfield.pos(regionIdx,3), shrink); + trisurf(k, leadfield.pos(regionIdx,1), leadfield.pos(regionIdx,2), leadfield.pos(regionIdx,3), 'FaceColor', [1 .75, .75], 'EdgeColor', 'none'); light('Position',[1 1 1],'Style','infinite', 'Color', [1 .75 .75]); light('Position',[-1 -1 -1],'Style','infinite', 'Color', [.5 .25 .25]); material dull; end % plotting electrodes -scatter3(-[leadfield.chanlocs.Y]', [leadfield.chanlocs.X]', [leadfield.chanlocs.Z]', 10, [.2 .2 .2], 'filled'); -if labels, text(-[leadfield.chanlocs.Y]'+ 2.5, [leadfield.chanlocs.X]', [leadfield.chanlocs.Z]', {leadfield.chanlocs.labels}, 'Color', [.2 .2 .2], 'FontSize', 8); end +if electrodes + scatter3(-[leadfield.chanlocs.Y]', [leadfield.chanlocs.X]', [leadfield.chanlocs.Z]', 10, [.2 .2 .2], 'filled'); + if labels + text(-[leadfield.chanlocs.Y]'+ 2.5, [leadfield.chanlocs.X]', [leadfield.chanlocs.Z]', {leadfield.chanlocs.labels}, 'Color', [.2 .2 .2], 'FontSize', 8); + end +end end diff --git a/plot/plot_source_location.m b/plot/plot_source_location.m index df68e88..7b7ba2b 100644 --- a/plot/plot_source_location.m +++ b/plot/plot_source_location.m @@ -24,6 +24,9 @@ % [ 0, 90] = axial % [90, 0] = sagittal % [ 0, 0] = coronal (default: [120, 20] +% region - cell of strings representing leadfield.atlas entries +% indicating which regions should be plotted to provide +% context. default {} plots all sources. % % Out: % h - handle of the generated figure @@ -39,6 +42,8 @@ % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% 2021-01-07 lrk +% - Added region argument % 2018-10-04 lrk % - Put the toolbar back in 3d mode % 2018-03-23 First version @@ -70,6 +75,7 @@ addParameter(p, 'shrink', .5, @isnumeric); addParameter(p, 'mode', '2d', @ischar); addParameter(p, 'view', [120, 20], @isnumeric); +addParameter(p, 'region', {}, @iscell); parse(p, sourceIdx, leadfield, varargin{:}) @@ -79,6 +85,22 @@ shrink = p.Results.shrink; mode = p.Results.mode; viewpoint = p.Results.view; +region = p.Results.region; + +% selecting sources based on atlas +regionIdx = 1:size(leadfield.pos, 1); +if ~isempty(region) + if isfield(leadfield, 'atlas') && ~isempty(leadfield.atlas) + idx = lf_get_source_all(leadfield, 'region', region); + if any(idx) + regionIdx = idx; + else + warning('indicated region(s) not found in the lead field''s atlas; plotting all sources'); + end + else + warning('no atlas information present in the lead field; plotting all sources'); + end +end braincolour = [.85 .85 .85]; sourcecolour = [ones(numel(sourceIdx),1), linspace(.6, .3, numel(sourceIdx))', linspace(.6, .3, numel(sourceIdx))']; @@ -88,8 +110,8 @@ if strcmp(mode, '3d') if newfig, h = figure('name', 'Source location', 'NumberTitle', 'off'); end hold on; - k = boundary(leadfield.pos(:,1), leadfield.pos(:,2), leadfield.pos(:,3), shrink); - s = trisurf(k, leadfield.pos(:,1), leadfield.pos(:,2), leadfield.pos(:,3), 'FaceColor', [1 .75, .75], 'EdgeColor', 'none'); + k = boundary(leadfield.pos(regionIdx,1), leadfield.pos(regionIdx,2), leadfield.pos(regionIdx,3), shrink); + s = trisurf(k, leadfield.pos(regionIdx,1), leadfield.pos(regionIdx,2), leadfield.pos(regionIdx,3), 'FaceColor', [1 .75, .75], 'EdgeColor', 'none'); light('Position',[1 1 1],'Style','infinite', 'Color', braincolour); light('Position',[-1 -1 -1],'Style','infinite', 'Color', [.5 .25 .25]); scatter3(leadfield.pos(sourceIdx,1), leadfield.pos(sourceIdx,2), leadfield.pos(sourceIdx,3), markersize, sourcecolour, 'fill'); @@ -101,20 +123,21 @@ elseif strcmp(mode, '2d') if newfig, h = figure('name', 'Source location', 'NumberTitle', 'off', 'ToolBar', 'none'); end - xmin = min(leadfield.pos(:,1)); - xmax = max(leadfield.pos(:,1)); - ymin = min(leadfield.pos(:,2)); - ymax = max(leadfield.pos(:,2)); - zmin = min(leadfield.pos(:,3)); - zmax = max(leadfield.pos(:,3)); + xmin = min(leadfield.pos([regionIdx, sourceIdx],1)); + xmax = max(leadfield.pos([regionIdx, sourceIdx],1)); + ymin = min(leadfield.pos([regionIdx, sourceIdx],2)); + ymax = max(leadfield.pos([regionIdx, sourceIdx],2)); + zmin = min(leadfield.pos([regionIdx, sourceIdx],3)); + zmax = max(leadfield.pos([regionIdx, sourceIdx],3)); % coronal view subplot(1,3,1); pos = get(gca, 'Position'); set(gca, 'Position', [0, pos(2), 1/3, pos(4)]); hold on; - k = boundary(leadfield.pos(:,1), leadfield.pos(:,3), shrink); - fill(leadfield.pos(k,1), leadfield.pos(k,3), braincolour, 'EdgeColor', 'none'); + xz = [leadfield.pos(regionIdx,1), leadfield.pos(regionIdx,3)]; + k = boundary(xz, shrink); + fill(xz(k,1), xz(k,2), braincolour, 'EdgeColor', 'none'); hsxz = scatter(leadfield.pos(sourceIdx, 1), leadfield.pos(sourceIdx, 3), markersize, sourcecolour, 'fill'); xlabel('X'); ylabel('Z'); xlim([xmin xmax]); @@ -126,8 +149,9 @@ pos = get(gca, 'Position'); set(gca, 'Position', [1/2.9, pos(2), 1/3, pos(4)]); hold on; - k = boundary(leadfield.pos(:,2), leadfield.pos(:,3), shrink); - fill(leadfield.pos(k,2), leadfield.pos(k,3), braincolour, 'EdgeColor', 'none'); + yz = [leadfield.pos(regionIdx,2), leadfield.pos(regionIdx,3)]; + k = boundary(yz, shrink); + fill(yz(k,1), yz(k,2), braincolour, 'EdgeColor', 'none'); hsyz = scatter(leadfield.pos(sourceIdx, 2), leadfield.pos(sourceIdx, 3), markersize, sourcecolour, 'fill'); xlabel('Y'); ylabel('Z'); xlim([ymin ymax]); @@ -139,8 +163,9 @@ pos = get(gca, 'Position'); set(gca, 'Position', [2/3, pos(2), 1/3, pos(4)]); hold on; - k = boundary(leadfield.pos(:,1), leadfield.pos(:,2), shrink); - fill(leadfield.pos(k,1), leadfield.pos(k,2), braincolour, 'EdgeColor', 'none'); + xy = [leadfield.pos(regionIdx,1), leadfield.pos(regionIdx,2)]; + k = boundary(xy, shrink); + fill(xy(k,1), xy(k,2), braincolour, 'EdgeColor', 'none'); hsxy = scatter(leadfield.pos(sourceIdx, 1), leadfield.pos(sourceIdx, 2), markersize, sourcecolour, 'fill'); xlabel('X'); ylabel('Y'); xlim([xmin xmax]); diff --git a/plot/plot_source_moment.m b/plot/plot_source_moment.m new file mode 100644 index 0000000..0faa9d1 --- /dev/null +++ b/plot/plot_source_moment.m @@ -0,0 +1,80 @@ +% plot_source_moment(sourceIdx, leadfield, varargin) +% +% Plots source orientations as arrows (dipolar moments), generally in +% an existing plot, e.g. to be used after plot_headmodel. +% +% In: +% sourceIdx - the index of the source in the leadfield to be plotted +% leadfield - the leadfield from which to plot the source +% +% Optional (key-value pairs): +% scale - numeric scaling value to scale the length of the arrow +% orientation - n-by-3 array of xyz source orientation, n being the +% number of sources indicated at sourceIdx. default +% uses the source's default orientation from the lead +% field. +% properties - quiver properties to adjust the drawing style. +% +% Usage example: +% >> lf = lf_generate_fromnyhead; +% >> rgn = {'Brain_Right_Temporal_Pole'}; +% >> plot_headmodel(lf, 'region', rgn, 'electrodes', 0); +% >> plot_source_moment(lf_get_source_middle(lf, 'region', rgn), lf); +% +% Copyright 2022 Laurens R. Krol +% Neuroadaptive Human-Computer Interaction +% Brandenburg University of Technology + +% 2022-12-02 First version + +% This file is part of Simulating Event-Related EEG Activity (SEREEGA). + +% SEREEGA is free software: you can redistribute it and/or modify +% it under the terms of the GNU General Public License as published by +% the Free Software Foundation, either version 3 of the License, or +% (at your option) any later version. + +% SEREEGA is distributed in the hope that it will be useful, +% but WITHOUT ANY WARRANTY; without even the implied warranty of +% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +% GNU General Public License for more details. + +% You should have received a copy of the GNU General Public License +% along with SEREEGA. If not, see . + +function plot_source_moment(sourceIdx, leadfield, varargin) + +% parsing input +p = inputParser; + +addRequired(p, 'sourceIdx', @isnumeric); +addRequired(p, 'leadfield', @isstruct); + +addParameter(p, 'scale', 50, @isnumeric); +addParameter(p, 'orientation', [], @isnumeric); +addParameter(p, 'properties', {'MaxHeadSize', 0.5, 'Color', 'k', 'LineWidth', 3, 'AutoScale', 'off'}, @iscell); + +parse(p, sourceIdx, leadfield, varargin{:}) + +sourceIdx = p.Results.sourceIdx; +lf = p.Results.leadfield; +scale = p.Results.scale; +orientation = p.Results.orientation; +props = p.Results.properties; + +if isempty(orientation) + orientation = lf.orientation(sourceIdx,:); +end + +% drawing quiver arrow(s) +quiver3( ... + lf.pos(sourceIdx,1), ... + lf.pos(sourceIdx,2), ... + lf.pos(sourceIdx,3), ... + orientation(:,1) * scale, ... + orientation(:,2) * scale, ... + orientation(:,3) * scale, ... + props{:} ... + ); + +end \ No newline at end of file diff --git a/pop/pop_sereega_lf_generate_fromhartmut.m b/pop/pop_sereega_lf_generate_fromhartmut.m new file mode 100644 index 0000000..efb2217 --- /dev/null +++ b/pop/pop_sereega_lf_generate_fromhartmut.m @@ -0,0 +1,149 @@ +% EEG = pop_sereega_lf_generate_fromhartmut(EEG) +% +% Pops up a dialog to serve as GUI for lf_generate_fromhartmut, and +% adds a lead field to the given EEGLAB dataset. +% +% The pop_ functions serve only to provide a GUI for some of +% SEREEGA's functions and are not intended to be used in scripts. +% +% In: +% EEG - an EEGLAB dataset +% +% Out: +% EEG - the EEGLAB dataset with added lead field depending on the +% actions taken in the dialog, at EEG.etc.sereega.leadfield +% +% Copyright 2022 Laurens R Krol +% Team PhyPA, Biological Psychology and Neuroergonomics, +% Berlin Institute of Technology + +% 2022-11-17 First version + +% This file is part of Simulating Event-Related EEG Activity (SEREEGA). + +% SEREEGA is free software: you can redistribute it and/or modify +% it under the terms of the GNU General Public License as published by +% the Free Software Foundation, either version 3 of the License, or +% (at your option) any later version. + +% SEREEGA is distributed in the hope that it will be useful, +% but WITHOUT ANY WARRANTY; without even the implied warranty of +% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +% GNU General Public License for more details. + +% You should have received a copy of the GNU General Public License +% along with SEREEGA. If not, see . + +function EEG = pop_sereega_lf_generate_fromhartmut(EEG) + +% strjoin is problematically overloaded by a number of EEGLAB toolboxes; +% making sure that we use the correct one +cdorig = cd(fileparts(which('strfun/strjoin'))); +joinstring = @strjoin; +cd(cdorig); +userdata.joinstring = joinstring; + +models = {'mix_Colin27_small', 'NYhead_small'}; + +% getting channel labels +labels{1} = {}; +fid = fopen('chanlocs-hartmut-mix_Colin27_small.xyz', 'r'); +while ~feof(fid) + line = fgets(fid); + [i, j] = regexp(line, '(?<=\s+)\w+(?=\s*$)'); + if ~isempty(i), labels{1} = [labels{1}, {line(i:j)}]; end +end +fclose(fid); +labels{1} = sort(labels{1}); + +labels{2} = {}; +fid = fopen('chanlocs-hartmut-NYhead_small.xyz', 'r'); +while ~feof(fid) + line = fgets(fid); + [i, j] = regexp(line, '(?<=\s+)\w+(?=\s*$)'); + if ~isempty(i), labels{2} = [labels{2}, {line(i:j)}]; end +end +fclose(fid); +labels{2} = sort(labels{2}); + +% setting defaults +userdata.model = 1; +userdata.labels = labels{userdata.model}; + +% getting montages +montages = {sprintf('All %d channels', length(userdata.labels)), 'Custom selection'}; +montages = [montages, utl_get_montage('?')]; + +if isstr(EEG) + userdata = get(gcbf, 'userdata'); + set(findobj(gcbf, 'tag', 'model'), 'value', userdata.model); + userdata.labels = labels{userdata.model}; + montages{1} = sprintf('All %d channels', length(userdata.labels)); + set(findobj(gcbf, 'tag', 'montage'), 'string', montages); + set(findobj(gcbf, 'tag', 'montage'), 'value', 1); + set(findobj(gcbf, 'tag', 'labeledit'), 'string', joinstring(userdata.labels)); + set(gcf, 'userdata', userdata); + return +end + +% callbacks +cb_model = [ ... + 'userdata = get(gcf, ''userdata'');', ... + 'userdata.model = get(gcbo, ''value'');', ... + 'set(gcf, ''userdata'', userdata);', ... + 'pop_sereega_lf_generate_fromhartmut(''setuserdata'');']; +cb_montage = [ ... + 'userdata = get(gcf, ''userdata'');' ... + 'if get(gcbo, ''value'') == 1', ... + ' set(findobj(''parent'', gcbf, ''tag'', ''labeledit''), ''string'', ''' joinstring(labels{userdata.model}) ''');' ... + 'elseif get(gcbo, ''value'') == 2', ... + ' set(findobj(''parent'', gcbf, ''tag'', ''labeledit''), ''string'', '''');' ... + 'else', ... + ' montage = get(gcbo, ''String'');', ... + ' montage = montage(get(gcbo, ''value''));', ... + ' set(findobj(''parent'', gcbf, ''tag'', ''labeledit''), ''string'', userdata.joinstring(utl_get_montage(montage{1})));' ... + 'end']; +cb_select = [ ... + 'userdata = get(gcf, ''userdata'');', ... + '[~, labelselection] = pop_chansel(userdata.labels);' ... + 'if ~isempty(labelselection)' ... + ' set(findobj(''parent'', gcbf, ''tag'', ''labeledit''), ''string'', labelselection);' ... + 'end']; + +[~, ~, ~, structout] = inputgui( ... + 'geometry', { 1 1 1 1 1 1 1 1 1 [4 1] }, ... + 'geomvert', [1 1 1 3 1 1 1 6 1 1], ... + 'uilist', { ... + { 'style', 'text', 'string', 'Obtain a lead field from HArtMuT', 'fontweight', 'bold' }, ... + { }, ... + { 'style', 'text', 'string', 'Indicate which model to use' }, ... + { 'style', 'listbox', 'string', models, 'tag', 'model', 'callback', cb_model }, ... + { }, ... + { 'style', 'text', 'string', 'Indicate which channels to include in the simulation,' }, ... + { 'style', 'text', 'string', 'using a montage or custom channel selection' }, ... + { 'style', 'listbox', 'string', montages, 'tag', 'montage', 'callback', cb_montage }, ... + { }, ... + { 'style', 'edit', 'string', joinstring(userdata.labels), 'tag', 'labeledit' }, ... + { 'style', 'pushbutton', 'string', '...', 'callback', cb_select } }, ... + 'helpcom', 'pophelp(''lf_generate_fromhartmut'');', ... + 'title', 'Lead field: HArtMuT', ... + 'userdata', userdata); + +if ~isempty(structout) + % user pressed OK, getting label selection in cell format + labelselection = textscan(structout.labeledit, '%s'); + labelselection = labelselection{1}'; + + model = models{structout.model}; + + % adding lead field to EEG structure + disp('Generating lead field...'); + EEG.etc.sereega.leadfield = lf_generate_fromhartmut(model, 'labels', labelselection); + + % also adding code to generate lead field + EEG.etc.sereega.leadfieldfunction = sprintf('lf_generate_fromhartmut(''%s'', ''labels'', %s);', model, ['{''' joinstring(labelselection, ''',''') '''}']); + disp('Done.'); +end + +end + diff --git a/pop/pop_sereega_plot_headmodel.m b/pop/pop_sereega_plot_headmodel.m index 73c6273..391d7b7 100644 --- a/pop/pop_sereega_plot_headmodel.m +++ b/pop/pop_sereega_plot_headmodel.m @@ -17,6 +17,10 @@ % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% 2022-12-02 lrk +% - Added electrodes argument +% 2021-01-07 lrk +% - Added atlas support % 2018-07-16 lrk % - Fixed issue where view was not properly converted % 2018-04-23 First version @@ -38,6 +42,13 @@ function EEG = pop_sereega_plot_headmodel(EEG) +% strjoin is problematically overloaded by a number of EEGLAB toolboxes; +% making sure that we use the correct one +cdorig = cd(fileparts(which('strfun/strjoin'))); +joinstring = @strjoin; +cd(cdorig); +userdata.joinstring = joinstring; + % testing if lead field is present if ~isfield(EEG.etc, 'sereega') || ~isfield(EEG.etc.sereega, 'leadfield') ... || isempty(EEG.etc.sereega.leadfield) @@ -49,14 +60,31 @@ return end +% obtaining list of regions, renaming generic ones +[allregions, ~, generic, ~] = lf_get_regions(EEG.etc.sereega.leadfield); +genreglist = cellfun(@(x) sprintf('All_%s*', x), generic, 'UniformOutput', false); +reglist = [{'All'}, genreglist, allregions]; + +cb_selectregion = [ ... + '[~, regionselection] = pop_chansel({''' joinstring(reglist, ''',''') '''});' ... + 'if ~isempty(regionselection)' ... + ' set(findobj(''parent'', gcbf, ''tag'', ''regionedit''), ''string'', regionselection);' ... + 'end']; + [~, ~, ~, structout] = inputgui( ... - 'geometry', { 1 1 [1 2] 1 [1 2] [1 2] 1 [1 2]}, ... - 'geomvert', [1 1 1 1 1 1 1 1], ... + 'geometry', { 1 1 [2 1] 1 1 [1 2] [1 2] 1 [1 2] [1 2] 1 [1 2]}, ... + 'geomvert', [1 1 1 1 1 1 1 1 1 1 1 1], ... 'uilist', { ... { 'style', 'text', 'string', 'Plot head model', 'fontweight', 'bold' }, ... { }, ... + { 'style', 'text', 'string', 'Select the region(s) to plot' }, ... + { 'style', 'pushbutton', 'string', '...', 'callback', cb_selectregion }, ... + { 'style', 'edit', 'string', joinstring(reglist(1)), 'tag', 'regionedit' }, ... + { }, ... + { 'style', 'text', 'string', 'Electrodes' }, ... + { 'style', 'checkbox', 'string', 'Plot electrodes', 'value', 1, 'tag', 'electrodes' }, ... { 'style', 'text', 'string', 'Labels' }, ... - { 'style', 'checkbox', 'string', 'Plot channel labels', 'value', 1, 'tag', 'labels' }, ... + { 'style', 'checkbox', 'string', 'Plot electrode labels', 'value', 1, 'tag', 'labels' }, ... { }, ... { 'style', 'text', 'string', 'Style' }, ... { 'style', 'popupmenu', 'string', 'scatter|boundary', 'tag', 'style' }, ... @@ -70,8 +98,19 @@ if ~isempty(structout) % user pressed OK + + % getting region selection in cell format + regionselection = textscan(structout.regionedit, '%s'); + regionselection = regionselection{1}'; + + % adding regex to original labels and undoing the renaming + generic = cellfun(@(x) sprintf('^%s.*$', x), generic, 'UniformOutput', false); + allregions = cellfun(@(x) sprintf('^%s$', x), allregions, 'UniformOutput', false); + regexlist = [{'.*'}, generic, allregions]; + region = cellfun(@(x) regexlist{find(strcmp(x, reglist))}, regionselection, 'UniformOutput', false); + styles = {'scatter', 'boundary'}; - plot_headmodel(EEG.etc.sereega.leadfield, 'labels', structout.labels, 'style', styles{structout.style}, 'shrink', str2double(structout.shrink), 'view', str2num(structout.view)); + plot_headmodel(EEG.etc.sereega.leadfield, 'electrodes', structout.electrodes, 'labels', structout.labels, 'style', styles{structout.style}, 'shrink', str2double(structout.shrink), 'view', str2num(structout.view), 'region', region); end end diff --git a/pop/pop_sereega_plot_source_location.m b/pop/pop_sereega_plot_source_location.m index b8cef7b..ab85acb 100644 --- a/pop/pop_sereega_plot_source_location.m +++ b/pop/pop_sereega_plot_source_location.m @@ -1,4 +1,4 @@ -% EEG = pop_sereega_plot_sources(EEG) +% EEG = pop_sereega_plot_source_location(EEG) % % Pops up a dialog to plot the sources associated with the given EEG % dataset. @@ -18,6 +18,8 @@ % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% 2021-01-07 lrk +% - Added atlas support % 2018-04-30 First version % This file is part of Simulating Event-Related EEG Activity (SEREEGA). @@ -37,6 +39,13 @@ function EEG = pop_sereega_plot_source_location(EEG) +% strjoin is problematically overloaded by a number of EEGLAB toolboxes; +% making sure that we use the correct one +cdorig = cd(fileparts(which('strfun/strjoin'))); +joinstring = @strjoin; +cd(cdorig); +userdata.joinstring = joinstring; + % testing if lead field is present if ~isfield(EEG.etc, 'sereega') || ~isfield(EEG.etc.sereega, 'leadfield') ... || isempty(EEG.etc.sereega.leadfield) @@ -55,21 +64,28 @@ return end +% obtaining list of regions, renaming generic ones +[allregions, ~, generic, ~] = lf_get_regions(EEG.etc.sereega.leadfield); +genreglist = cellfun(@(x) sprintf('All_%s*', x), generic, 'UniformOutput', false); +reglist = [{'All'}, genreglist, allregions]; + % general callback functions cbf_get_value = @(tag,property) sprintf('get(findobj(''parent'', gcbf, ''tag'', ''%s''), ''%s'')', tag, property); +cbf_set_value = @(tag,property,value) sprintf('set(findobj(''parent'', gcbf, ''tag'', ''%s''), ''%s'', %s);', tag, property, value); % callbacks -cb_2d = 'plot_source_location([EEG.etc.sereega.components.source], EEG.etc.sereega.leadfield);'; -cb_3d= [ ... - 'shrink = str2num(' cbf_get_value('shrink', 'string') ');' ... - 'view = str2num(' cbf_get_value('view', 'string') ');' ... - 'plot_source_location([EEG.etc.sereega.components.source], EEG.etc.sereega.leadfield, ''mode'', ''3d'', ''shrink'', shrink, ''view'', view);' - ]; +cb_check2d = cbf_set_value('check3d', 'value', ['~' cbf_get_value('check3d', 'value')]); +cb_check3d = cbf_set_value('check2d', 'value', ['~' cbf_get_value('check2d', 'value')]); +cb_selectregion = [ ... + '[~, regionselection] = pop_chansel({''' joinstring(reglist, ''',''') '''});' ... + 'if ~isempty(regionselection)' ... + ' set(findobj(''parent'', gcbf, ''tag'', ''regionedit''), ''string'', regionselection);' ... + 'end']; % building gui -inputgui( ... - 'geometry', { 1 1 1 1 1 1 [1 1] 1 [1 1] [1 1] }, ... - 'geomvert', [1 1 1 1 1 1 1 1 1 1], ... +[~, ~, ~, structout] = inputgui( ... + 'geometry', { 1 1 1 1 1 1 [1 1] 1 [3 1] 1 1 [1 1] [1 1] }, ... + 'geomvert', [1 1 1 1 1 1 1 1 1 1 1 1 1], ... 'uilist', { ... { 'style', 'text', 'string', 'Plot source locations', 'fontweight', 'bold' }, ... { }, ... @@ -77,8 +93,12 @@ { 'style', 'text', 'string', 'To plot individual sources, go to "Select source' }, ... { 'style', 'text', 'string', 'locations" under "Configure components".' }, ... { }, ... - { 'style', 'pushbutton', 'string', 'Plot in 2D', 'callback', cb_2d }, ... - { 'style', 'pushbutton', 'string', 'Plot in 3D', 'callback', cb_3d }, ... + { 'style', 'checkbox', 'string', 'Plot in 2D', 'enable', 'on', 'tag', 'check2d', 'value', 1, 'callback', cb_check2d }, ... + { 'style', 'checkbox', 'string', 'Plot in 3D', 'enable', 'on', 'tag', 'check3d', 'value', 0, 'callback', cb_check3d }, ... + { }, ... + { 'style', 'text', 'string', 'Region(s) to include as context' }, ... + { 'style', 'pushbutton', 'string', '...', 'callback', cb_selectregion }, ... + { 'style', 'edit', 'string', joinstring(reglist(1)), 'tag', 'regionedit' }, ... { }, ... { 'style', 'text', 'string', 'Shrink factor' }, ... { 'style', 'edit', 'string', '0.5', 'tag', 'shrink' }, ... @@ -87,5 +107,25 @@ }, ... 'helpcom', 'pophelp(''plot_source_location'');', 'title', 'Plot source locations'); +if ~isempty(structout) + % user pressed OK + + % getting region selection in cell format + regionselection = textscan(structout.regionedit, '%s'); + regionselection = regionselection{1}'; + + % adding regex to original labels and undoing the renaming + generic = cellfun(@(x) sprintf('^%s.*$', x), generic, 'UniformOutput', false); + allregions = cellfun(@(x) sprintf('^%s$', x), allregions, 'UniformOutput', false); + regexlist = [{'.*'}, generic, allregions]; + region = cellfun(@(x) regexlist{find(strcmp(x, reglist))}, regionselection, 'UniformOutput', false); + + % plotting + if structout.check2d + plot_source_location([EEG.etc.sereega.components.source], EEG.etc.sereega.leadfield, 'region', region); + else + plot_source_location([EEG.etc.sereega.components.source], EEG.etc.sereega.leadfield, 'mode', '3d', 'shrink', str2num(structout.shrink), 'view', str2num(structout.view), 'region', region); + end end +end \ No newline at end of file diff --git a/pop/pop_sereega_sources.m b/pop/pop_sereega_sources.m index 0be4f8c..b9cdafe 100644 --- a/pop/pop_sereega_sources.m +++ b/pop/pop_sereega_sources.m @@ -6,19 +6,28 @@ % The pop_ functions serve only to provide a GUI for some of % SEREEGA's functions and are not intended to be used in scripts. % -% To add sources to the simulation, first use one of the four options -% provided to find sources in the lead field: +% To add sources to the simulation, first select if you want to +% restrict the search to certain regions. By default, all sources +% from the leadfield can be searched for. Click "select" to select +% specific regions from the leadfield's atlas, if available. +% +% Then, use one of the four options provided to find sources in the +% lead field: % - Next to "random", indicate the number of randomly selected % sources you wish to inspect; % - Next to "nearest to", indicate the coordinates in the brain where -% you wish to find a source; +% you wish to find a source. Alternatively, click "middle" to +% obtain the arithmetic mean of all source coordinates in the +% currently selected regions. % - Next to "spaced", indicate the number of sources you wish to % find, and the minimum distance in mm between them; % - Next to "in radius", indicate the location in the brain (or the % ID of the source) and the size of the radius in mm. -% Click "find" to find the source(s). They will show up below under +% Click the "find" button to find the source(s) using the +% corresponding method described above. They will show up below under % "found location", along with their default orientation. Click -% "plot" to inspect these values graphically. +% "plot" to inspect these values graphically. Both plots update +% automatically when new sources are found. % % The default orientation can be changed. Next to "orientation", a % custom orientation can be indicated and applied by clicking @@ -54,10 +63,14 @@ % EEG - the EEGLAB dataset with sources according to the actions % taken in the dialog % -% Copyright 2018 Laurens R Krol +% Copyright 2018 Laurens R. Krol % Team PhyPA, Biological Psychology and Neuroergonomics, % Berlin Institute of Technology +% 2022-11-17 lrk +% - Added source IDs to list of sources, final console output +% 2021-01-08 lrk +% - Added atlas support % 2018-04-25 First version % This file is part of Simulating Event-Related EEG Activity (SEREEGA). @@ -93,7 +106,24 @@ EEG.etc.sereega.components = struct('source', {}, 'signal', {}, 'orientation', {}, 'orientationDv', {}); end +% strjoin is problematically overloaded by a number of EEGLAB toolboxes; +% making sure that we use the correct one +cdorig = cd(fileparts(which('strfun/strjoin'))); +joinstring = @strjoin; +cd(cdorig); +userdata.joinstring = joinstring; + +% obtaining list of regions, renaming generic ones +[allregions, ~, generic, ~] = lf_get_regions(EEG.etc.sereega.leadfield); +genreglist = cellfun(@(x) sprintf('All_%s*', x), generic, 'UniformOutput', false); +reglist = [{'All'}, genreglist, allregions]; + % setting userdata +userdata.region = {'.*'}; +userdata.regiongeneric = generic; +userdata.regionall = allregions; +userdata.regionlist = reglist; +userdata.regioninterpret = @interpret_regionedit; userdata.newsourceidx = []; userdata.newsourceorientation = []; userdata.plotlocationhandle = []; @@ -104,7 +134,7 @@ % generating list of current sources currentsourcelist = {}; ... for i = 1:numel(userdata.currentsourceidx), ... - currentsourcelist = [currentsourcelist, {sprintf('%d at ( %d, %d, %d )', i, round(EEG.etc.sereega.leadfield.pos(EEG.etc.sereega.components(i).source,:)))}]; ... + currentsourcelist = [currentsourcelist, {sprintf('%d at ( %d, %d, %d ), ID %d', i, round(EEG.etc.sereega.leadfield.pos(EEG.etc.sereega.components(i).source,:)), EEG.etc.sereega.components(i).source)}]; ... end % general callback functions @@ -128,7 +158,7 @@ 'end;' ... 'currentsourcelist = {};' ... 'for i = 1:numel(userdata.currentsourceidx),' ... - 'currentsourcelist = [currentsourcelist, {sprintf(''%d at ( %d, %d, %d )'', i, round(EEG.etc.sereega.leadfield.pos(userdata.currentsourceidx(i),:)))}];' ... + 'currentsourcelist = [currentsourcelist, {sprintf(''%d at ( %d, %d, %d ), ID %d'', i, round(EEG.etc.sereega.leadfield.pos(userdata.currentsourceidx(i),:)), userdata.currentsourceidx(i))}];' ... 'end;' ... cbf_set_value('currentsources', 'string', 'currentsourcelist'); ... ]; @@ -150,11 +180,23 @@ ]; % callbacks +cb_selectregion = [ ... + cbf_get_userdata ... + '[~, regionselection] = pop_chansel({''' joinstring(reglist, ''',''') '''});' ... + 'if ~isempty(regionselection)' ... + 'set(findobj(''parent'', gcbf, ''tag'', ''regionedit''), ''string'', regionselection);' ... + 'userdata.region = userdata.regioninterpret(regionselection, userdata.regiongeneric, userdata.regionall, userdata.regionlist);' ... + 'if ~isempty(userdata.plotlocationhandle) && all(ishandle(userdata.plotlocationhandle)),' ... + 'display(''Close and re-plot location to update plotted regions'');' ... + 'end;' ... + 'end;' ... + cbf_set_userdata ... + ]; cb_find_source_random = [ ... cbf_get_userdata ... 'numsources = str2num(' cbf_get_value('fnd_src_random', 'string') '); ' ... 'if ~isempty(numsources),'... - 'userdata.newsourceidx = lf_get_source_random(EEG.etc.sereega.leadfield, numsources);' ... + 'userdata.newsourceidx = lf_get_source_random(EEG.etc.sereega.leadfield, ''number'', numsources, ''region'', userdata.region);' ... 'userdata.newsourceorientation = EEG.etc.sereega.leadfield.orientation(userdata.newsourceidx,:);' ... 'end;' ... cbf_set_userdata ... @@ -162,11 +204,16 @@ cbf_update_plot_location ... cbf_update_plot_projection ... ]; +cb_find_source_middle = [ ... + cbf_get_userdata ... + 'coords = mean(EEG.etc.sereega.leadfield.pos(lf_get_source_all(EEG.etc.sereega.leadfield, ''region'', userdata.region),:));' ... + cbf_set_value('fnd_src_nearest', 'string', 'sprintf(''%.0f '', coords)'); + ]; cb_find_source_nearest = [ ... cbf_get_userdata ... 'xyz = str2num(' cbf_get_value('fnd_src_nearest', 'string') '); ' ... 'if ~isempty(xyz),'... - 'userdata.newsourceidx = lf_get_source_nearest(EEG.etc.sereega.leadfield, xyz);' ... + 'userdata.newsourceidx = lf_get_source_nearest(EEG.etc.sereega.leadfield, xyz, ''region'', userdata.region);' ... 'userdata.newsourceorientation = EEG.etc.sereega.leadfield.orientation(userdata.newsourceidx,:);' ... 'end;'... cbf_set_userdata ... @@ -179,7 +226,7 @@ 'num = str2num(' cbf_get_value('fnd_src_spaced_num', 'string') '); ' ... 'mm = str2num(' cbf_get_value('fnd_src_spaced_mm', 'string') '); ' ... 'if ~isempty(num) && ~isempty(mm),' ... - 'userdata.newsourceidx = lf_get_source_spaced(EEG.etc.sereega.leadfield, num, mm);' ... + 'userdata.newsourceidx = lf_get_source_spaced(EEG.etc.sereega.leadfield, num, mm, ''region'', userdata.region);' ... 'userdata.newsourceorientation = EEG.etc.sereega.leadfield.orientation(userdata.newsourceidx,:);' ... 'end;' ... cbf_set_userdata ... @@ -192,7 +239,7 @@ 'pos = str2num(' cbf_get_value('fnd_src_radius_pos', 'string') '); ' ... 'mm = str2num(' cbf_get_value('fnd_src_radius_mm', 'string') '); ' ... 'if ~isempty(pos) && ~isempty(mm),'... - 'userdata.newsourceidx = lf_get_source_inradius(EEG.etc.sereega.leadfield, pos, mm);' ... + 'userdata.newsourceidx = lf_get_source_inradius(EEG.etc.sereega.leadfield, pos, mm, ''region'', userdata.region);' ... 'userdata.newsourceorientation = EEG.etc.sereega.leadfield.orientation(userdata.newsourceidx,:);' ... 'end;' ... cbf_set_userdata ... @@ -253,7 +300,7 @@ cb_plot_location = [ ... cbf_get_userdata ... 'if isempty(userdata.plotlocationhandle) || ~all(ishandle(userdata.plotlocationhandle)),' ... - '[h, hsxz, hsyz, hsxy] = plot_source_location(userdata.newsourceidx, EEG.etc.sereega.leadfield);' ... + '[h, hsxz, hsyz, hsxy] = plot_source_location(userdata.newsourceidx, EEG.etc.sereega.leadfield, ''region'', userdata.region);' ... 'userdata.plotlocationhandle = [h, hsxz, hsyz, hsxy];' ... cbf_set_userdata ... 'else,' ... @@ -276,6 +323,7 @@ 'userdata.currentsourceidx = [userdata.currentsourceidx, userdata.newsourceidx];' ... 'userdata.currentsourceorientation = [userdata.currentsourceorientation; userdata.newsourceorientation];' ... cbf_set_userdata ... + cbf_set_value('currentsources', 'value', '1') ... cbf_update_fields ... ]; cb_select_source = [ ... @@ -295,78 +343,90 @@ 'userdata.currentsourceorientation(src,:) = [];' ... cbf_set_userdata ... 'if src > numel(userdata.currentsourceidx),' ... - cbf_set_value('currentsources', 'value', '1') ... + 'if numel(userdata.currentsourceidx) == 0,' ... + cbf_set_value('currentsources', 'value', '[]') ... + 'else,' ... + cbf_set_value('currentsources', 'value', '1') ... + 'end;' ... 'end;' ... cbf_update_fields ... cb_select_source ... ]; % geometry -nc = 6; nr = 14; +nc = 6; nr = 16; geom = { ... { nc nr [ 0 0] ... % current sources text [ 2 1] }, ... { nc nr [ 0 1] ... % listbox - [ 2 11] }, ... - { nc nr [ 0 13] ... % remove source + [ 2 13] }, ... + { nc nr [ 0 15] ... % remove source [ 2 1] }, ... { nc nr [ 2 0] ... % add sources text [ 6 1] }, ... - { nc nr [ 2 1] ... % add sources: random text + { nc nr [ 2 1] ... % add sources: from region [ 1 1] }, ... { nc nr [ 3 1] ... % edit [ 2 1] }, ... - { nc nr [ 5 1] ... % find + { nc nr [ 5 1] ... % select [ 1 1] }, ... - { nc nr [ 2 2] ... % add sources: nearest text + { nc nr [ 2 3] ... % add sources: random text [ 1 1] }, ... - { nc nr [ 3 2] ... % edit + { nc nr [ 3 3] ... % edit [ 2 1] }, ... - { nc nr [ 5 2] ... % find + { nc nr [ 5 3] ... % find [ 1 1] }, ... - { nc nr [ 2 3] ... % add sources: spaced text + { nc nr [ 2 4] ... % add sources: nearest text [ 1 1] }, ... - { nc nr [ 3 3] ... % edit # + { nc nr [ 3 4] ... % edit [ 1 1] }, ... - { nc nr [ 4 3] ... % edit mm + { nc nr [ 4 4] ... % middle [ 1 1] }, ... - { nc nr [ 5 3] ... % find + { nc nr [ 5 4] ... % find [ 1 1] }, ... - { nc nr [ 2 4] ... % add sources: radius text + { nc nr [ 2 5] ... % add sources: spaced text [ 1 1] }, ... - { nc nr [ 3 4] ... % edit xyz + { nc nr [ 3 5] ... % edit # [ 1 1] }, ... - { nc nr [ 4 4] ... % edit mm + { nc nr [ 4 5] ... % edit mm [ 1 1] }, ... - { nc nr [ 5 4] ... % find + { nc nr [ 5 5] ... % find + [ 1 1] }, ... + { nc nr [ 2 6] ... % add sources: radius text + [ 1 1] }, ... + { nc nr [ 3 6] ... % edit xyz [ 1 1] }, ... - { nc nr [ 2 6] ... % orientation text + { nc nr [ 4 6] ... % edit mm + [ 1 1] }, ... + { nc nr [ 5 6] ... % find + [ 1 1] }, ... + { nc nr [ 2 8] ... % orientation text [ 1 1] }, ... - { nc nr [ 3 6] ... % orientation edit + { nc nr [ 3 8] ... % orientation edit [ 2 1] }, ... - { nc nr [ 3 7] ... % orientation default + { nc nr [ 3 9] ... % orientation default [ 1 1] }, ... - { nc nr [ 4 7] ... % orientation random + { nc nr [ 4 9] ... % orientation random [ 1 1] }, ... - { nc nr [ 3 8] ... % orientation perpendicular + { nc nr [ 3 10] ... % orientation perpendicular [ 1 1] }, ... - { nc nr [ 4 8] ... % orientation tangential + { nc nr [ 4 10] ... % orientation tangential [ 1 1] }, ... - { nc nr [ 5 6] ... % orientation apply + { nc nr [ 5 8] ... % orientation apply [ 1 1] }, ... - { nc nr [ 2 10] ... % found location text + { nc nr [ 2 12] ... % found location text [ 1 1] }, ... - { nc nr [ 3 10] ... % xyz + { nc nr [ 3 12] ... % xyz [ 2 1] }, ... - { nc nr [ 5 10] ... % plot + { nc nr [ 5 12] ... % plot [ 1 1] }, ... - { nc nr [ 2 11] ... % found orientation text + { nc nr [ 2 13] ... % found orientation text [ 1 1] }, ... - { nc nr [ 3 11] ... % xyz + { nc nr [ 3 13] ... % xyz [ 2 1] }, ... - { nc nr [ 5 11] ... % plot + { nc nr [ 5 13] ... % plot [ 1 1] }, ... - { nc nr [ 3 13] ... % add sources + { nc nr [ 3 15] ... % add sources [ 2 1] }, ... }; @@ -377,11 +437,15 @@ { 'style' 'listbox' 'string' currentsourcelist, 'tag', 'currentsources', 'callback', cb_select_source } ... { 'style' 'pushbutton' 'string' 'Remove source', 'callback', cb_remove_source } ... { 'style' 'text' 'string' 'Add new sources', 'fontweight', 'bold' } ... + { 'style' 'text' 'string' 'From region' } ... + { 'style' 'edit' 'string' reglist{1}, 'tag', 'regionedit' } ... + { 'style' 'pushbutton' 'string' 'Select', 'callback', cb_selectregion } ... { 'style' 'text' 'string' 'Random' } ... { 'style' 'edit' 'string' '#', 'tag', 'fnd_src_random' } ... { 'style' 'pushbutton' 'string' 'Find', 'callback', cb_find_source_random } ... { 'style' 'text' 'string' 'Nearest to' } ... { 'style' 'edit' 'string' 'x y z', 'tag', 'fnd_src_nearest' } ... + { 'style' 'pushbutton' 'string' 'Middle', 'callback', cb_find_source_middle } ... { 'style' 'pushbutton' 'string' 'Find', 'callback', cb_find_source_nearest } ... { 'style' 'text' 'string' 'Spaced' } ... { 'style' 'edit' 'string' '#', 'tag', 'fnd_src_spaced_num' } ... @@ -428,7 +492,23 @@ if ~isempty(userdata.plotprojectionhandle) && ishandle(userdata.plotprojectionhandle) close(userdata.plotprojectionhandle); end + fprintf('Source IDs:%s\n', sprintf(' %d', EEG.etc.sereega.components.source)); fprintf('Number of source locations: %d\n', numel([EEG.etc.sereega.components.source])); end end + + +function region = interpret_regionedit(regioneditvalue, generic, allregions, combinedlist) + + % getting region selection in cell format + regionselection = textscan(regioneditvalue, '%s'); + regionselection = regionselection{1}'; + + % adding regex to original labels and undoing the renaming + generic = cellfun(@(x) sprintf('^%s.*$', x), generic, 'UniformOutput', false); + allregions = cellfun(@(x) sprintf('^%s$', x), allregions, 'UniformOutput', false); + regexlist = [{'.*'}, generic, allregions]; + region = cellfun(@(x) regexlist{find(strcmp(x, combinedlist))}, regionselection, 'UniformOutput', false); + +end diff --git a/utils/utl_call_gui.m b/utils/utl_call_gui.m index 5a0969b..d22d3af 100644 --- a/utils/utl_call_gui.m +++ b/utils/utl_call_gui.m @@ -17,10 +17,13 @@ % >> utl_call_gui('sources', lf); % >> utl_call_gui('plot_headmodel', lf); % -% Copyright 2021 Laurens R Krol +% Copyright 2021, 2022 Laurens R Krol % Neuroadaptive Human-Computer Interaction % Brandenburg University of Technology +% 2022-11-16 lrk +% - Fixed issue where some GUI dialogs relied on variables in the base +% workspace % 2021-10-01 First version % This file is part of Simulating Event-Related EEG Activity (SEREEGA). @@ -40,18 +43,24 @@ function utl_call_gui(popfunction, leadfield) -% creating temporary EEG dataset -EEG.etc.sereega.leadfield = leadfield; %#ok - -% amending popfunction to full pop_sereega_*(EEG) +% amending popfunction to full pop_sereega_* if not already the case if ~startsWith(popfunction, 'pop_sereega_') popfunction = ['pop_sereega_' popfunction]; end -if ~endsWith(popfunction, '(EEG)') - popfunction = [popfunction '(EEG)']; -end + +% some functions need access to the EEG variable in the base workspace: +% copying base EEG to temporary variable and replacing with new local EEG +% that only contains the requested leadfield +localEEG.etc.sereega.leadfield = leadfield; +localEEG.etc.sereega.components = struct('source', {}, 'signal', {}, 'orientation', {}, 'orientationDv', {}); +evalin('base', 'tempEEG_for_utl_call_gui = EEG;'); +assignin('base', 'EEG', localEEG); % calling GUI -eval(popfunction); +evalin('base', [popfunction '(EEG);']); + +% restoring base variables +evalin('base', 'EEG = tempEEG_for_utl_call_gui;'); +evalin('base', 'clear tempEEG_for_utl_call_gui;'); end diff --git a/utils/utl_sanitize_atlas.m b/utils/utl_sanitize_atlas.m new file mode 100644 index 0000000..e41e9f3 --- /dev/null +++ b/utils/utl_sanitize_atlas.m @@ -0,0 +1,67 @@ +% atlas = utl_sanitize_atlas(atlas) +% +% Takes the cell of strings that is a lead field's atlas, and +% sanitizes it for compatibility with other functions. The +% limitations are primarily due to EEGLAB's GUI functions. Currently, +% the following substitutions happen: +% +% - Space (' ') becomes underscore ('_') +% - Apostrophe (') is removed +% +% Also forces the atlas to be a column, and checks whether each atlas +% element belongs to (i.e. each string starts with) one of the known +% categories, i.e. brain, muscle, etc. +% +% In: +% atlas - cell of strings +% +% Out: +% atlas - cell of sanitized strings +% +% Usage example: +% >> atlas = {'Brain Wernicke''s Area', 'Brain Heschl''s Gyrus'}; +% >> atlas = utl_sanitize_atlas(atlas); +% +% Copyright 2022 Laurens R. Krol +% Neuroadaptive Human-Computer Interaction +% Brandenburg University of Technology + +% 2022-11-25 First version + +% This file is part of Simulating Event-Related EEG Activity (SEREEGA). + +% SEREEGA is free software: you can redistribute it and/or modify +% it under the terms of the GNU General Public License as published by +% the Free Software Foundation, either version 3 of the License, or +% (at your option) any later version. + +% SEREEGA is distributed in the hope that it will be useful, +% but WITHOUT ANY WARRANTY; without even the implied warranty of +% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +% GNU General Public License for more details. + +% You should have received a copy of the GNU General Public License +% along with SEREEGA. If not, see . + +function atlas = utl_sanitize_atlas(atlas) + + % forcing column + if isrow(atlas) + atlas = atlas'; + end + + % replacing known problematic characters + atlas = cellfun(@(x) strrep(x, ' ', '_'), atlas, 'UniformOutput', false); + atlas = cellfun(@(x) strrep(x, '''', ''), atlas, 'UniformOutput', false); + + % checking for generic categories + brainidx = strncmpi(atlas, 'brain', 5); + muscleidx = strncmpi(atlas, 'muscle', 6); + eyeidx = strncmpi(atlas, 'eye', 3); + + unknowns = numel(atlas) - sum([brainidx; muscleidx; eyeidx]); + if unknowns > 0 + warning('%d atlas element(s) do not start with generic region categories', unknowns) + end + +end \ No newline at end of file