Back to the main page.
Bug 2753 - Getting started with fNIRS (and fNIRS/EEG)
| Status | CLOSED FIXED |
| Reported | 2014-11-04 11:57:00 +0100 |
| Modified | 2019-08-10 12:33:05 +0200 |
| Product: | FieldTrip |
| Component: | core |
| Version: | unspecified |
| Hardware: | PC |
| Operating System: | Windows |
| Importance: | P5 minor |
| Assigned to: | Jörn M. Horschig |
| URL: | |
| Tags: | |
| Depends on: | |
| Blocks: | |
| See also: | http://bugzilla.fieldtriptoolbox.org/show_bug.cgi?id=2882http://bugzilla.fieldtriptoolbox.org/show_bug.cgi?id=2881 |
Jörn M. Horschig - 2014-11-04 11:57:07 +0100
As discussed in the FT meeting on Wed 29/10/2014: The fNIRS community is getting more and more interested in combining EEG/fNIRS analysis. In a meeting of the NIRS Network Nijmegen (NNN) interest was raised to use FieldTrip to combine fNIRS and EEG analysis in one toolbox. To support this endeavor, additional documentation and code is required. We came up with the following, general ideas: - we need to extend the 'getting started' guide with 'getting started on fNIRS data', explaining how to get data into FieldTrip. This requires settling on how to represent fNIRS data best in FieldTrip and how to incorporate wavelengths, optical density transformations, transmitter/receiver distinction, etc. This needs to be documented properly and already existing functions need to be extended (e.g. ft_datatype should correctly detect fNIRS data). - getting started with specific fNIRS systems requires input from those responsible for or experienced with the specific system. Individual companies (or researchers?) are responsible for writing and maintaining XXX2fieldtrip code, where XXX is the recording system (this is the same for EEG and MEG systems). All XXX2fieldtrip functions should result in the fNIRS FieldTrip data representation, which should be as close as possible to the EEG/MEG representation. - code for fNIRS specific analysis have to be written by those who are experienced in fNIRS analysis. All functions should adhere to FT guidelines, and be maintained by some responsible person (ideally the expert who initially commits the function). The functions should be located in fieldtrip/contrib/fnirs. This requires input from the fNIRS community, as the FT team does lack experience in fNIRS analysis. All functions should be generically written such that they are compatible with any fNIRS system (if possible). - once enough functions are available to conduct a full analysis, a tutorial should be written that illustrates how to conduct fNIRS analysis using FieldTrip. We should not strive to explain all pro's and con's of different (preprocessing) methods but stick to the basics. The tutorial should also explain how to combine fNIRS and EEG data using FieldTrip. - additional tutorial can be rewritten on demand and availability by those who are interested to contribute. Specific action points: - I will take over the initial start by creating the getting started fNIRS guide for the Artinis Oxymon system and think of a good way to represent fNIRS specific properties. I'll create a development/project page on the wiki to track current progress and thoughts. - Mark van Wanrooij (DCN/biophysics department) indicated that he and his group would gladly help in contributing code and/or documentation as they will start EEG/fNIRS analysis from next year onwards and already have experience with fNIRS analysis. Once I finished the getting started guide, I will contact him again.
Robert Oostenveld - 2014-11-04 12:31:48 +0100
could you share a NIRS dataset that was recorded in a task-related experiment in its original (Artinis) data format?
Jörn M. Horschig - 2014-12-18 16:43:46 +0100
I created a project page, which includes a todo list: http://fieldtrip.fcdonders.nl/development/nirs The next days/weeks, I will update the project page with suggestions on how FieldTrip could present NIRS data and peculiarities. These are not too different from how FieldTrip is presenting MEG data (e.g. a NIRS channel is defined by the transmitter-receiver relation, similarly as a channel in MEG can be defined by the two coils of a gradiometer). Therefore, I will try to stick closely to what is defined in MEG. (In reply to Robert Oostenveld from comment #1) I am not sure whether sharing a dataset be of any help. We maintain our own data format, and seem to be reluctant in sharing details about this (though, no magic happening there). Our Matlab import function is obfuscated for that reason. Is that a problem for FieldTrip to not having the binary format of the exported data explained somewhere?
Robert Oostenveld - 2014-12-19 09:44:46 +0100
(In reply to Jörn M. Horschig from comment #2) > Our Matlab import function is obfuscated for that reason. > Is that a problem for FieldTrip to not having the binary > format of the exported data explained somewhere? It is not ideal from a usability perspective, but legally or technically not a problem. I suggest we create a fieldtrip/external/artinis folder. It can (besides the p-code) contain its own README and COPYRIGHTS. You can even add a splash license or disclaimer that is shown once upon the first call of the reader code. I suggest you look into external/openmeeg, which I think is the best example of a non-GPL extension of FT. See openmeeg_license.m for the splash license (or external/mne/mne_license.m). The external/ctf code also has a disclaimer (although there not in a separate function). You could add an artinis_license.m and have that one called from the p-code. Alternative is that we call artinis_license from ft_read_header and ft_read_data just prior to calling the reading function.
Jörn M. Horschig - 2014-12-19 09:51:27 +0100
@Robert: could you check the two ideas below? Short update (and summary of the project page): upon thinking further about how to incorporate, I came to the conclusion that it could be sufficient to just introduce a new type of sensor description, called 'opto' for optodes: opto - contains information about the optodes. opto.tra - contains information about how receiver and transmitter form channels. opto.balance - information about balancing and transformations from ODs to concentrations opto.optopos - contains information about the position of the optodes. opto.optotype - contains information about the type of optode ('receiver' or 'transmitter'). opto.wavelength - NxM matrix, where N is the number of transmitters and M the number of wavelengths per transmitter (note that in theory, some transmitters can have only 2 wavelengths and others more - mostly this is constant in a given system though) opto.chanpos - contains information about the position of the channels (i.e. averages of optopos) opto.chantype - contains information about the channel type ('nirs') As channels measure (at least) two different things simultaneously (oxygenated hemoglobin and deoxygenated hemoglobin), I propose to have this reflected in the channel labels, e.g. data.label = {'Rx1 - Tx1 [HbO]', 'Rx1 - Tx1 [HbR]'}; selection of e.g. only oxygenated data can then be achieved using ft_channelselection.
Jörn M. Horschig - 2014-12-19 09:53:49 +0100
(In reply to Robert Oostenveld from comment #3) I agree that it is not the most user-friendly solution. I can also think that the obfuscation is partly due to make the Matlab import 'foolproof', which the average Matlab user might require but the average FieldTrip user not. I'll think about it and discuss here at Artinis.
Robert Oostenveld - 2014-12-19 16:27:07 +0100
following the phone call, I just noticed this docu % FT_DATATYPE_SENS describes the FieldTrip structure that represents % an EEG, ECoG, or MEG sensor array. This structure is commonly called % "elec" for EEG and "grad" for MEG, or more general "sens" for either % one. which already mentions ECoG. It would have to be changed to % FT_DATATYPE_SENS describes the FieldTrip structure that represents % an MEG, EEG, ECoG or NIRS sensor array. This structure is commonly called % "grad" for MEG, "elec" for EEG/ECoG and "opto" for NIRS, or more % general "sens" for either one of these. It would also have to be mentioned on http://fieldtrip.fcdonders.nl/faq/how_are_electrodes_magnetometers_or_gradiometers_described Furthermore, anywhere in the code where is now something like isfield(data, 'elec') the 'opto' case would have to be dealt with as well. Looking at mac011> grep -l isfield.*elec *.m | wc -l 21 and the functions it returns, opto would not have to be built into 21, but still in a considerable number of functions. Unrelated, but imagine this: you could make a neighbour definition of different chomophores (instead of Rx-Tx pairs), and cluster in the statistics over multiple concentrations being affected at the same time by the experimental manipulation, thereby compile more evidence for a significant effect.
Jörn M. Horschig - 2014-12-19 16:57:01 +0100
It just struck my mind, when having a leadfield-like solution for transforming ODs to concentrations, we are still having a problem for transforming the data. We need to have some function that takes data as input and returns data so that data labels are updated accordingly. (In reply to Robert Oostenveld from comment #6) Okay, I will take care of that as well. About stats: hmm, nice thought. Usually, however, we are looking for inverse effects in oxy- vs. deoxy, i.e. an increase in HbO and a decrease in HbR (otherwise the signal is probably not physiological). Would that work without dirty tricks (such as inversing one of the two concentrations)?
Jörn M. Horschig - 2015-01-07 10:08:08 +0100
Created attachment 689 NIRS header v1 Happy new year ;) I worked on implementing the header reading function into FT by adding relevant parts to ft_read_header and creating an external\artinis folder containing our code. I attached the resulting header of a 24-channel NIRS recording as a .mat-file. The header looks as follows: hdr = Fs: 10 nChans: 48 label: {48x1 cell} nSamples: 71994 nSamplesPre: 0 nTrials: 1 chantype: {48x1 cell} chanunit: {48x1 cell} opto: [1x1 struct] hdr.opto = tra: [24x48 double] optopos: [24x2 double] optotype: {24x1 cell} chanpos: [48x2 double] chantype: {48x1 cell} wavelength: [856 764 856 764 856 764 856 765 858 764 858 764 858 764 858 764] transmits: [24x16 logical] Although in common terms, we speak of a transmitter as the combination of multiple laser, I decided to treat each laser as one individual transmitter in FT. I did this because different channels of the same receiver-transmitter combination are defined by different wavelengths (i.e. laser). An example: Say we have two channels, channel1 is named 'Rx1-Tx1 [850nm]' (Receiver 1 receives light from Transmitter 1 emitting light at 850nm) and channel2 is named 'Rx1-Tx1 [760nm]'. The recorded data is physically based on the same optodes on the scalp, as Transmitter 1 emits both 850nm and 760nm (time multiplexed). However, the channel data is defined by absorption of light of different wavelengths (i.e. laser). I see the now implemented way as the only way to be explicit about this. Defining a transmitter as a combination of wavelengths would render it impossible to have an unambiguous mapping from transmitter-wavelengths to channel (Rx-Tx combination). Next stop: ft_read_data and related functions. Btw, I plan to commit my changes just when one can actually use FT to analyse NIRS data (i.e. when ft_read_data is extended and other functions are compatible with NIRS data).
Robert Oostenveld - 2015-01-08 11:03:05 +0100
(In reply to Jörn M. Horschig from comment #8) Overall it look fine to me, but I have some suggestions. Perhaps it would be more explicit what we are talking about if we would use hdr.opto.fiberpos hdr.opto.fibertype This is more like "coil" in the MEG case. Missing is hdr.opto.label, which should have as many elements as chanpos and chantype. I am a bit surprised by opto.chanpos and opto.optopos being Nx2 rather than Nx3. Why not in 3 dimensions? hdr.opto.tra should be transposed, i.e. the number of rows should be equal to the number of channels. The way that you have coded the wavelengths does not match my initial expectation, i.e. I would have done it differently (more compact), like this >> hdr.opto.wavelength ans = 764 765 856 858 >> hdr.opto.transmits ans = 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 % first transmit fiber, at 856nm 1 0 0 0 0 0 1 0 1 0 0 0 0 0 1 0 1 0 0 0 0 0 1 0 0 1 0 0 0 0 0 1 1 0 0 0 0 0 0 1 1 0 0 0 0 0 0 1 1 0 0 0 0 0 0 1 1 0 0 0 In this example there are no transmit fibers that cary multiple wavelengths, but I guess that would happen, right? Hmm, I now see that hdr.opto.optopos row 9 and 10 are the same. That strongly suggests that it is the same fiber! In fact, looking at octopus/fiberpos, I would now conclude that there are 8 transmit fibers, each with two wavelengths. Would this not better be represented as >> hdr.opto.transmits ans = 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0 % first transmit fiber, at 764 and 856 1 0 1 0 1 0 1 0 0 1 1 0 1 0 0 1 1 0 0 1 1 0 0 1 1 0 0 1 Consequently, optopos/fiberpos and optotype/filbertype would also have length 8+8 rather than 8+2*8. --- I typed the following section prior to identifying the 24/16 confusion, so I still assumed that there were 24 fibers in total, each with only a single wavelength. ---- Looking at imagesc(hdr.opto.tra') I can clearly see the 8 receivers (coded positively) and 16 transmitters. But I cannot see how a channel is constructed, since wavelength info is missing. I would think that for channel 1, which is Rx1a-Tx @ 856nm, it could have hdr.opto.tra(1,:) = zeros(1,24) hdr.opto.tra(1,1) = 3; % receive from fiber 1 at wavelength 3 hdr.opto.tra(1,9) = -3; % transmit from fiber 9 at wavelength 3 I am not sure yet whether tra should contain negative numbers. Whether a column (fiber) is a receiver to transmitter is also present in optotype, and hence does not have to be coded again. How would the matrix be used elsewhere? Would it be used as a lookup table, or for computations? Perhaps an interesting idea to further flesh out how tra should be organised is to consider some unusual applications: say that we average two channels that have the same transmit fiber (also same wavelength) but different receive fibers. How would that be represented in a row of tra? Or if we compute a difference (ratio?) between one transmit-receive pair at wavelength 1 versus wavelength 2? ---- Now considering that opto.transmits is 16x4 rather than than 24x16, I would say that hdr.opto.tra = zeros(48,16) % 48 channels constructed from 16 fibers >> hdr.label{1} 'Rx1a-Tx1 [856nm]' hdr.opto.tra(1,1) = 3; % receive from fiber 1 at wavelength 3 hdr.opto.tra(1,9) = -3; % transmit from fiber 9 at wavelength 3 >> hdr.label{2} 'Rx1a-Tx1 [764nm]' hdr.opto.tra(1,1) = 1; % receive from fiber 1 at wavelength 1 hdr.opto.tra(1,9) = -1; % transmit from fiber 9 at wavelength 1 ----- cheers Robert
Robert Oostenveld - 2015-01-08 11:16:59 +0100
(In reply to Robert Oostenveld from comment #9) It now occurs to me that compared to yours I only made the description more compact. Your version with the split-wavelength optodes did allow to lookup the wavelength and fibers for a specific channel from opto.tra. Please think for some time about my more compact representation, and check whether it is really compatible/convertible to yours (and vice versa). I think I prefer mine because the explicit link to fibers as physical element in the cap (which people can see sticking out) and the more compact representation of wavelengths (there are 4, hence a vector of length 4). So my feeling is that my description can be more easily explained to a non-tech person.
Jörn M. Horschig - 2015-01-12 15:25:44 +0100
Hi Robert, I agree that calling it fibers makes it more intuitive and overall better structured. I went forth and implemented most of your suggestions. However, now I start doubting the compactness of the opto.wavelength field. In the previous version, opto.wavelength was [1xL], where L is the number of physical lasers. You suggested to represent it as [1xW], where W is the number of unique wavelengths. This makes the opto-structure cleaner, and less redundant, at first sight. I am currently weighing off this cleanliness against correctness and being intuitive. Cleanliness: - By representing unique wavelengths, it is much easier to see which fibers and channels used what wavelengths. - It also looks aesthetically more appealing. Well, as appealing as an integer-vector can look like ;) - From a programmer's point of view, this solution is also much more elegant. Correctness: - I am not sure if our lasers are really emitting at 764nm or e.g. at 764.2nm. We are rounding here to integers. We are one of the few companies who are honest by admitting that not all lasers are emitting at the exact same wavelengths (e.g. in a 24ch device, there can be lasers with e.g. 762nm, 763nm and 764nm - we aim to have them as close to each other as possible). Other companies claim that there are strictly two wavelengths in their system, thus for most systems opto.wavelength will be 1x2 or [1x3] (which could be in fact incorrect). In my opinion, the compact view supports such incorrectness, and I can imagine many users wondering why there are multiple wavelengths in some systems and only 2 in others. In fact actually, it is quite essential to know the exact wavelength as that will be used to convert from optical densities to changes in chromophore concentration. Intuition: - As said, in contrast to the previous opto-structure, the current design neglects individual lasers. These are physical devices connected to the mainboard. Each fiber is connected to 'a number' of lasers. No laser is, however, connected to multiple fibers! It is thus a surjective but not an injective mapping (I knew my Bachelor pays off some day by being able to throw around such fancy terms!). Naive users could thus get the impression that two transmitters emitting at 764nm were connected to the same physical laser. On the other hand, by simply looking at the device, they might realize that this impression is wrong. Thus, in total I see Cleanliness as a big plus for the new design, correctness as a small plus for the old design (mathematically they are equivalent) and intuition as a small plus for the old design. --------------------------------------------------------- Concerning combing channels, updating the opto-structure: There are several points to consider: 1) I think opto.tra is cannot be dealt with analogously to any other .tra matrix, simply because a transmitter does not receive any data. So, the vector 1 0 0 -1 0 0 representing 'Rx1 - Tx1 [852nm]' is not obtained by subtracting activity in channel Tx1 from activity in channel Rx1. It's data should be interpreted as 'the light going out of Tx1 at the first wavelength that is received by Rx1'. 2) A combination of channels is no simple addition of data points. For example a bipolar combination of two NIRS-channels does not actually yield activity that is spatially 'in between' these two channels (as it is commonly interpreted for EEG data). For NIRS, such a derivation can spatially only be interpreted as the 'difference in location A vs. location B'. Things such as different wavelengths etc. complicate that matter (the difference between two channels that are based on different wavelengths is hardly interpretable other than 'there is difference'). 3) Although we could update the .tra-matrix according to such computations, I have no idea what the content of the matrix would then represent. For example, say that the raw data is represented in the .tra-matrix as such: 1 0 0 -1 0 0 'Rx1 - Tx1 [764nm]' 0 1 0 0 0 -1 'Rx2 - Tx3 [764nm]' After computing a bipolar derivation "chan1-chan2", how should the .tra-matrix look like? 1 -1 0 -1 0 1 would not make much sense wrt my earlier explanation. Rx2 is not transmitting light at any wavelength, Tx3 is not receiving light. 1 1 0 -1 0 -1 would also not make much sense, as there would be no difference between "chan1-chan2" and "chan1+chan2". Most importantly, the receiver-transmitter mapping gets completely lost. Thus, it might be better to re-think the .tra matrix. I can imagine that we might want an additional matrix similar to .transmits (or change .transmits), which reflects what fiber is transmitting which wavelengths. However, we also need some matrix representing what receiver receives light from which transmitter at what wavelength (thus, the current .tra matrix). At the moment, I am considering whether we should swap the data-representation of .tra and .transmit, i.e. .tra becomes a boolean-matrix, where a '1' simply means 'contributing to this channel', and .transmit represents at which wavelength a fiber is transmitting light. Since a fibers are transmitting light of multiple wavelengths, we would however still need a clean mapping from wavelengths (or laser) to channel. I'll think about this a bit and will drop my thoughts here when done ;)
Robert Oostenveld - 2015-01-12 17:58:39 +0100
(In reply to Jörn M. Horschig from comment #11) Hi Jorn, Thanks for sharing the elaborate considerations, and learning me some new words along the way ;-) Regarding Cleanliness/Intuition and the wavelengths: ideally the structure would allow either representation without it being interpreted differently. Similar like the rows of grad.chanori that map onto columns of grad.tra, which both can be multiplied by -1 without the interpretation changing (just flipping orientation twice). The nuance of the (only) slightly different wavelengths was not obvious to me before, but now I get it. I would be also be perfectly fine with wavelength = 1*Nlasers. I also had not thought about the number of lasers at the same wavelength and guess that I had in mind that a single laser would be shared between channels with some optical switching in between. I am still not happy with the tra. As you say, it cannot be used for the same types of computations. So why call it tra? We have lasers, fibers, tissue, fibers again, and detectors (are these diodes?). So why not map it as laser-fiber-photodetector or transmitter-fiber-receiver? Each transmitter is connected to a fiber, each receiver is connected to a fiber, a channel is a combination of transmitter+fiber and receiver+fiber. The number of transmitters is the number of lasers, and each of them can have its own wavelength. This then leads me to think in the direction of opto.chanlabel = nchan * 1 opto.chanpos = nchan * 3 opto.fiberpos = nfiber * 3 opto.wavelength = 1 * nlaser opto.transmitter = nlaser * nfiber opto.receiver = nreceiver * nfiber opto.sensor = ... Hmm, I am still not happy. A channel requires the specification of the wavelength and the fiber, plus the specification of receiving fiber. Not of the receiving photodetector, right? That means that a single channel needs three numbers/indices/boolvecs. Better perhaps... use http://en.wikipedia.org/wiki/Transceiver (as a mix of both) opto.transmitter = nfiber * nwavelength/nlaser opto.receiver = nfiber * nwavelength/nlaser opto.transceiver = nchan * 2, where the transmitter is in the 1st and the receiver in the 2nd column. Or each row of transceiver is a vector of positive and negative numbers, coding transmitter and receiver indices, respectively. The transmitter and receiver code in a boolean matrix what sending/receiving wavelength is connected to each fiber.
Jörn M. Horschig - 2015-01-15 15:51:48 +0100
Interesting thought, however not possible, unfortunately. There are L Lasers (and thus L wavelengths), T transmitters and R receivers. The total number of fibers is N=T+R. The relation between T and L is dependent on the system. Mostly it is T=L/2, i.e. every transmitter transmits light of 2 different wavelengths. Sometimes it is T=L/3, but in theory (at least for our system), one transmitter can transmit light at X wavelengths, another at Y, the next one at Z, etc.. Let's recap on what we settled on, and what needs to be revised: The header seems to be fixed as: hdr = Fs: scalar nChans: scalar C label: cell-array Cx1 nSamples: scalar nSamplesPre: scalar=0 nTrials: scalar=1 chantype: cell-array Cx1 chanunit: cell-array Cx1 opto: struct The opto-structure needs more work. Approved fields in the opto-structure are: opto = fiberpos: matrix Nx2 (we have only layouts atm, can also be Nx3) fibertype: cell Nx1 chanpos: matrix Cx2 or Cx3 chantype: cell Cx1 wavelength: vector 1xL I think we should also add: fiberlabel: cell Nx1 We are missing information on a) what optodes are receiver and transmitter? b) how are optodes or lasers combined to form channels? We thought about a) previously and thought to solve it with a 'transmits'-field: transmits: matrix NxL, boolean For b) we had the .tra field in the past, which took values akin to three-valued logic (-1, 0 or 1). But .tra is used for fields of different purposes in FieldTrip, therefore we would like to refine this field. I see three main solutions: 1) Very elegant would be to solve a) and b) in one field. Let's settle on the proposal to use 'transceiver' as the name (I was always thinking of it as optoCmb, analogously to chanCmb): transceiver: matrix NxC where the matrix can take values from 1 to L, or -1 if it is a receiver. Then, we wouldn't need to have a 'transmits' field. 2) Alternatively, we could also take a step back and go to your suggestions. As we know that a receiver always sees all wavelengths of a certain transmitter, we need 'transceiver' to be a NxC boolean matrix or, as you proposed, a Cx2-matrix in addition to the 'transmits'. The second 'transceiver' dimension determines the receiver and the wavelength (or laser) to create the channel. If we settle on wavelengths as being 1xL, there is a unique transmitter-laser mapping. by combining 'transmits' and 'transceiver' 3) 'transceiver' enters the third dimension as transceiver: matrix Cx3 with receiver, transmitter and wavelengths. This would abandon the .transmits-field. Btw, thinking about it, 'transmits' seems like an important field, as it shows a basic property of the system independent of hows channels are formed.
Robert Oostenveld - 2015-01-15 17:02:25 +0100
> There are L Lasers (and thus L wavelengths), T transmitters > and R receivers. The total number of fibers is N=T+R. How about this with your option 1: opto = fiberpos: matrix Nx2 or Nx3 fibertype: cell Nx1 chanpos: matrix Cx2 or Cx3 chantype: cell Cx1 tranceiver: matrix CxN wavelength: vector 1xL Note that I made transceiver CxN, to be consistent with grad.tra. We could even again call it opto.tra(nsceiver). The ambiguity in "tra" is that it can also be "transmits" and that a negative number would naturally be interpreted as "receives". Let me try this out.... If we assume 2 (approximate) wavelengths tranceiver [ 1 0 0 0 -1 0 0 0 % fiber 1 sends at wavelength 1, fiber 5 receives 2 0 0 0 -1 0 0 0 % fiber 1 sends at wavelength 2, fiber 5 receives 0 1 0 0 -1 0 0 0 % fiber 2 sends at wavelength 1, fiber 5 receives 0 2 0 0 -1 0 0 0 % fiber 2 sends at wavelength 2, fiber 5 receives ... ] or with 2*T wavelengths tranceiver [ 1 0 0 0 -1 0 0 0 % fiber 1 sends at wavelength 1, fiber 5 receives 2 0 0 0 -1 0 0 0 % fiber 1 sends at wavelength 2, fiber 5 receives 0 3 0 0 -1 0 0 0 % fiber 2 sends at wavelength 3, fiber 5 receives 0 4 0 0 -1 0 0 0 % fiber 2 sends at wavelength 4, fiber 5 receives ... ] This allows for fibertype{any(transceiver>0, 1)} = 'transmitter'; fibertype{any(transceiver<0, 1)} = 'receiver'; and colums/fibers that have both positive and negative number are to be flagged as error (until even more fancy hardware comes along that combines the Tx and Rx in a single fiber). I would be in favour of using similarly sized negative numbers for the receiving end, so tranceiver [ 1 0 0 0 -1 0 0 0 % fiber 1 sends at wl 1, fiber 5 receives at wl 1 2 0 0 0 -2 0 0 0 % fiber 1 sends at wl 2, fiber 5 receives at wl 2 0 3 0 0 -3 0 0 0 % fiber 2 sends at wl 3, fiber 5 receives at wl 3 0 4 0 0 -4 0 0 0 % fiber 2 sends at wl 4, fiber 5 receives at wl 4 ... ] as this gives a more useful representation in terms of mathematical operations. I can also imagine that the lasers are somehow pulsed, and that a certain pulse scheme of the laser is being used to time the receiver and to prevent cross-talk. That would also nicely be coded in here, without actually having to specify it in detail. Splitting it over multiple matrices makes it more complicated. 3 dimensions are difficult to explain. Shall we(*) try it like this and see how it works out in practice? *) that means you ;-)
Jörn M. Horschig - 2015-01-16 08:50:32 +0100
okay, that sounds like a nice plan. I especially like the symmetry of transmitters and receivers being coded as '+/- [wavelengths]', respectively. Then, one can also easily find all receiver's that are connected to a specific laser/transmitter. I am still not sure if we should call it .tra because it might be very confusing. Also, if activity across channels is combined (e.g. averaged, or if ICA is perfoemd), I am not sure if and how this should be reflected in the opto-structure. But maybe we(*) can best iterate our way down there (i.e. face the problems when they occur). *) ;)
Jörn M. Horschig - 2015-02-02 16:12:43 +0100
Created attachment 696 updated nirs header I updated the code according to our discussion. The updated nirs header is attached. It'd be nice if you can check whether it looks good (I added fiberlabel and chanlabel to opto). As already proposed in my last comment, I think it's best to keep the structure like this and see if I run into problems further down in the implementation. If the header is alright, I will move on to implementing the data representation and the associated necessary transformations from optical densities (ODs) to (changes in) chromophore concentration.
Jörn M. Horschig - 2015-02-03 15:12:12 +0100
Since I am running into too much trouble already, I'll stop for now and start a new discussion. There are several suboptimalities with our own data representation (imho) with respect to the FiedlTrip way of reading data in, but there is also a more general, high-level issue we need to resolve first. For every measurement, the so-called 'differential path length factor' (DPF) needs to be set. This factor is the estimated increase in length of photon traveling between the transmitter and the receiver taking the underlying tissue and scattering into account. Light is transmitted perpendicular to the skin into the brain, and it travels along an infamous banana-shaped path until it arrives at the receiver. While the receiver and the transmitter are only a few centimeters away, the actual photon traveling along the banana-shape is longer. This is what the DPF represents. This factor is later needed to transform the change in optical density to changes in chromophore concentration (using simple regression). As we already discussed, Robert, you suggested that the DPF is part of the subject-specific data, and thus should be part of the data-structure. It is not part of the optode-setup, as I could in theory put the same sensor arrangement on top of someone else's head, who has a different DPF. I agree with this idea, but I do not see the DPF as more subject-specific than optode-positions.Optodes are be attached to the subject's skull, hence they are in fact also subject-specific. I could try to fit the exact same sensor-arrangement to your head than I did on my, but I bet that this will not yield a satistfactory data quality ;) Similarly, I could use the same DPF for your data as for mine, but it also will make the data estimate less accurate. Thus, I also consider the DPF as 'meta'-data, which is not directly related to the measurement, thus belonging to the header - similar to optode-positions. The main reason why I would like the DPF in the header (or optode-structure), however, is a mere practical one: the current FieldTrip implementation of ft_read_data needs to only return a data-matrix. There is simply no room to return the DPF as well. As already said, the DPF will be needed later on for transformations, thus it needs to be read out of the data-file somewhere. My suggestion is to make the DPF part of the header. I think you are right that it does not belong to the opto-structure, as it has no influence on any properties of the optodes themselves. Unfortunately this makes the whole data processing pipeline rather complicated, because for the transformation we need three ingredients: 1) the data 2) the DPF (i.e. the header-structure) 3) the distance between optodes (i.e. the opto-structure) I would like the header to contain hdr.DPF - scalar, differential path length factor and the data to contain data.hdr - structure, header of the data incl. optode-structure Then, all parameters to work with are returned to user upon a call to ft_preprocessing. Any counter-suggestions? And sorry for the long tale again ;)
Robert Oostenveld - 2015-02-03 15:21:58 +0100
(In reply to Jörn M. Horschig from comment #17) Previously I had related the "opto" to "elec/grad" and the DPF to the volume conduction model. One way of creating the EEG volume conduction model in absence of an MRI is by fitting spheres (or something else) to the measured electrode positions. This sounds what you are suggesting with the DPFs and optode positions. It is fine with me to add the DPF details to the opode structure. That is better than to add it to data.hdr, since that is not kept consistently with the data.
Jörn M. Horschig - 2015-02-04 12:35:14 +0100
is hdr.orig usually used to dump all of system-specific metainformation that might be useful later? That could save me multiple calls to read in the datafile, e.g. when reading in events. But it looks kinda ugly...
Jörn M. Horschig - 2015-02-05 12:22:59 +0100
Created attachment 698 ft_preprocessing on nirs btw, just struck my mind - I should probably also adapt the BUCN_NIRS functions then, shouldn't I?
Jörn M. Horschig - 2015-02-05 16:56:07 +0100
Do we need to specify manufacturer specific nirs sensor labels? I would generally argue for saying that NIRS data in Fieldtrip needs to either end with '[???nm]', or with '[O2Hb]' or '[HHb]', so the senstype 'nirs' should be sufficient. I will go ahead and implement it like this. Feel free to disagree ;)
Robert Oostenveld - 2015-02-05 17:51:17 +0100
(In reply to Jörn M. Horschig from comment #19) yes. But you should only assume it to exist in ft_read_header, ft_read_data and ft_read_event and not in any other (higher level) fieldtrip function.
Jörn M. Horschig - 2015-02-06 09:22:00 +0100
hmm, for ft_senslabel I run into trouble (it's more a style issue). There are near infinitely possible combinations of sensor labels (any receiver id can go with any transmitter id, and transmitters can in theory have any wavelength). So, I could code ft_senslabel to return all these possibilities from say IDs 1-100 and wavelengths anything between 600nm and 900nm in steps of 1nm. But you will probably easily see that this is a) ugly and b) takes up a lot of memory for such a simple task. I am thinking about also doing a regexp matching here, but then there would be no senslabel for NIRS-data. I could put this regexp-check in an external function if desired. Since you usually have smart ideas about these things, is there anything you can advise? Another issue: NIRS-data can have anything from 1 to >100 channels, and mostly auxiliary channels are also in the data. So I cannot just check in ft_senstype if some ratio of channels is NIRS, I would need to check if there is any NIRS channel. This would need to go at the end of the if-statement so that this check is only done if the data is of no other type.
Jörn M. Horschig - 2015-02-23 11:52:43 +0100
any thoughts on the senslabel issue?
Robert Oostenveld - 2015-02-23 14:40:15 +0100
(In reply to Jörn M. Horschig from comment #24) I think that senslabel should not be needed. Its primary reason for existence is to deal with different MEG systems, allowing ft_senstype to figure out what kind of system it is. But the implemented strategy in general is not very robust, because we can have combined MEG+EEG, just as NIRS+EEG. So the assumption of a single senstype is not very robust.
Jörn M. Horschig - 2015-02-23 15:14:32 +0100
(In reply to Robert Oostenveld from comment #25) the reason I'd like it is to be able to use 'NIRS' in ft_channelselection. But I can also go on without this (until people start requesting this).
Robert Oostenveld - 2015-02-23 17:55:56 +0100
(In reply to Jörn M. Horschig from comment #26) You could also just extend ft_channelselection for that. It presently supports groups such as EOG ECG which are also not dealt with using senslabel.
Jörn M. Horschig - 2015-03-05 10:39:57 +0100
good idea, I just implemented this. The tiny neat-nerd in me dislikes that some meta-channels are in lowercase like 'mua' and 'lfp' and others in uppercase like 'MEG' and 'EEG'. Any thoughts on why that is and if a case insensitive matching might be appropriate? Just thinking about whether to 'NIRS' or 'nirs'... minor issue though ;)
Robert Oostenveld - 2015-03-05 14:40:14 +0100
(In reply to Jörn M. Horschig from comment #28) I think that in channelselection there is a tendency to use upper case, and in ft_chantype and other low-level code a tendency for using lower case. Quite often it does not matter, but for channels themselves the p in Fp1 (Frontal pole) means something else than the P in P1 (Parietal). It would be good if it were made more consistent, and where possible the code would be made case insensitive. In the code itself I am in favour of lower case, hence I would prefer that as default.
Jörn M. Horschig - 2015-03-05 14:50:31 +0100
while I'm on it, should I then also change things like findmeg = find(strcmp(channel, 'MEG')); to findmeg = find(strcmpi(channel, 'MEG')); ?
Jörn M. Horschig - 2015-03-05 16:42:23 +0100
one more (again with more text to express and explain my thoughts): we talked about how transformation from "changes in optical density" to "concentration changes" such as "changes in (de-)oxygenated blood" can be handled. This transformation is simply a remixing of the channels, i.e. a low and a high wavelength measurement get weighted and summed to achieve an estimate of the chromophore concentration. You suggested to implement this in FT similar to how leadfields are handled, i.e. have a function computing a transform-matrix that can then be applied to the data. While this works in theory, this does not resolve the re-labelling of the channels. Or in other words, we need to know to what chromophore concentration estimate the transformation is. While thinking about a solution for this, I came across montages. They would contain a 'tra'nsformation matrix, and an array containing future channel labels. a) Would it be alright if I utilize the montage structure for these transformations? b) Apart from the fact that a leadfield does not require a labelnew-field, why is a montage not used to represent the leadfield? It would contain a labelorg-field for checking consistency with the data. To me it sounds that it would resolve the problems we had the past couple of years with the leadfield.
Robert Oostenveld - 2015-03-13 18:34:11 +0100
(In reply to Jörn M. Horschig from comment #31) Re a) fine with me. I do suggest to make a high level function to have a place where you can document and explain the functionality, otherwise people will never find it. The high level functio can call apply_montage. Re b) never thought of it. The montage is used at a slightly higher level in the code, where bookkeeping is done. The leadfields are used in sections where labels dont exist most of the time. But there is the code in between. Note that leadfields now should have labels (not that i wrote the code, but so i was told).
Jörn M. Horschig - 2015-04-07 10:32:24 +0200
Hi Robert, a short status update: (1) I finished the work on the conversion functions, so the only thing left to do is the automatic layoutting, for which most functionality is already built-in. (2) I discussed licensing here at Artinis. We would like to have our specific read-in function obfuscated, but everything else open source. The only restriction we would like to have is to restrict commercial usage, but that everything may well be modified and re-distributed. Seems that GPL, but CC looks like a good option. Also, it'd be great if there would be a mention that Artinis created and maintains that fNIRS part of FieldTrip. Should I prepare a patch once I am done with all the ft_ main functions that I modified (e.g. ft_preprocessing and ft_channelselection), or should I immediately upload them myself? Best, Jörn
Jörn M. Horschig - 2015-04-24 13:32:53 +0200
I committed the first version of the plugin, and after a little bugfix it seems to work fine for the Biphysics lab (van Opstal's group). So, I will close this one and create a new bug concerning topographic plotting, which does not work, yet. All commits FieldTrip functions: https://code.google.com/p/fieldtrip/source/detail?r=10334 https://code.google.com/p/fieldtrip/source/detail?r=10335 https://code.google.com/p/fieldtrip/source/detail?r=10340 https://code.google.com/p/fieldtrip/source/detail?r=10341 external/artinis: https://code.google.com/p/fieldtrip/source/detail?r=10337 https://code.google.com/p/fieldtrip/source/detail?r=10338 https://code.google.com/p/fieldtrip/source/detail?r=10339 https://code.google.com/p/fieldtrip/source/detail?r=10342
Jörn M. Horschig - 2015-06-22 10:36:22 +0200
Hi Robert, just to let you know that the first version is out and running stable after some initial bugfixes - someone at HBM told me that they heard from you that this is not released yet, so therefore just a reminder that it is ;)
Robert Oostenveld - 2016-05-13 12:58:25 +0200
I made some updates, see https://github.com/fieldtrip/fieldtrip/pull/159
Jörn M. Horschig - 2016-05-13 14:01:47 +0200
looks all good, apart from the documentation of ft_datatype_sens and subsequent handling of NIRS data. I would propose to leave this for now and let me fix this once we get going with this, alright? Of course, you are free to implement this as well, but I would expect that it's quite a time investment, and you maybe want to spend your time with more imminent problems right now ;)