ofgulban / bvbabel Goto Github PK

Hi Faruk:

In the read_fmr function, I have notice your note:

# NOTE(Faruk): I dont like this matrix reader but I don't see a
# more elegant way for now.

From what I have researched and tested (only slice timings' code is shown here):

slice_timings = np.ndarray.tolist(np.fromstring(''.join(lines[j + 1 : j + 1 + int(content[1])]), sep='\n'))
info_multiband["Slice timings"] = slice_timings

This could potentially replace the original method. The logic is, since numpy is a dependency of bvbabel, we can use their fromstring function to convert text data into 1D array. This way we could avoid using two loops for reading affine and slice timings arrays. I have tested on the data I have, and compared to the original method, I printed out their slice timing array, and they seems identical. I think this is an indication that this method is working.

If this catches your eyes, I can make a pull request.

Thank you very much!

Hi Faruk,

I would like to read both lag and correlation value from a vmp file from linear correlation analysis.
Currently, when I read the vmp file, I get only a matrix of lag value.
However, I would like to get correlation map as well to mask my ROI.

Thank you so much in advance.
Kenshu

As a feature request, it would be great if this software could support MTC format and include an example MTC file.

@ju-ec , thanks for your contribution #28 . Would it be possible for you to send me some test data for the ROI format (you can also email if you wish)? So that write, create and test functions can be developed in the future.

I have added an initial roundtrip test. @nausikaa8 , would you have time to have a look and see if sth is missing (see here)? If it looks good to you, I am going to proceed adding other file formats to the roundtrip tests.

Also consider adding an example script for converting TRF information to BIDS task events format: /home/faruk/Documents/test_bvbabel/vtc2/anatomy_tmean.vmr

@nausikaa8 , maybe you would be interested in contributing this. Let me know.

Hi Faruk,

I am trying to convert Brainvoyager's .VTC and .VMP outputs in NIFTI using your read_vtc_export_nifti.py and read_vtc_export_nifti.py scripts from the example folder. Unfortunately, the exported NIFTIs don't overlap when opened in an external software.

My current solution for the VTC:

header, data = bvbabel.vtc.read_vtc(FILE, rearrange_data_axes=False)
new_data = np.transpose(data, (0, 2, 1, 3))
new_data = new_data[::-1, ::-1, ::-1, :]

Alessandra

@ofgulban: would it be useful to have the BrainVoyager DWI formats implemented in bvbabel? Most prominent are:

• *.dmr/*.dwi (slice-based data)
• *.vdw (volume-based data)

Some others are:

• *.ddt (diffusion tensors)
• *.grb (bvec/bval)
• *.fbr (reconstructed fiber tracts).

bvbabel cannot read any PRTs that save onset and offset event times separated with tabs.
To fix this problem I would suggest to replace the following lines in prt.py:

values = lines[i+2+j].split(" ")

with

values = lines[i+2+j].split()

This will make sure that any whitespace can be used as separator between the timings.Or use values = lines[i+2+j].split(maxsplit=1) to make sure that each line is splitted only once.

Probably less common, but also possible would be the TAB separator for the color definition. Therefore, I would also suggest to replace the following line:

values = lines[i+2+n].split(" ")

with

values = lines[i+2+n].split() 

Based on the sequential reading order, I believe there is a mismatch in the header information of XStart, XEnd, YStart, YEnd, ZStart, ZEnd. Since I can only rely on sanity checks (couldn't find any documentation) this would be good to double check (this is also the reason I did not submit a pull request).

What it looks like is that the header data is loaded z, x, y (as is the case in many brainvoyager files) instead of x, y, z what is currently the case. if this is true, coding wise, header['ZStart'] and header['ZEnd'] have to be moved above/before header['XStart'].

Create SDM files function is needed

Hi @ofgulban ,
while trying to save a .vmp file (vmp.write_wmp) with a single map I got the following error:

Traceback (most recent call last):

  File "<ipython-input-11-c03dcbcc2ac7>", line 3, in <module>
    bvbabel.vmp.write_vmp(OUTNAME, new_vmp_header, new_vmp_voi)

  File "C:\Users\apizz\miniconda3\lib\site-packages\bvbabel-0.0.2-py3.8.egg\bvbabel\vmp.py", line 413, in write_vmp
    data_img = data_img[::-1, ::-1, ::-1, :]  # Flip BV axes

IndexError: too many indices for array

By default it expects multiple maps. My work around solution was to create a 4D array and add 'None' to the last dimension.

The NIfTI specification defines cal_min and cal_max as follows:

float cal_max;       /*!< Max display intensity */  /* float cal_max;       */
float cal_min;       /*!< Min display intensity */  /* float cal_min;       */

It is standard convention to set both to zero if they are to be ignored. Unfortunately, BrainVoyager seems to set these to 0 and 255 respectively when creating NIfTI data.

The BIDS format BrainVoyager Getting Started Guide (GSG) files exhibit this, e.g. GSGData.zip. Consider the fMRI image that has a mean intensity of 4022 (with a range of 0..49183). Here the cal_min and cal_max of 0..255 is inappropriate, with even dark air voxels having values ~700.

$ fslinfo sub-01_ses-04_task-blocked_run-1_bold.nii.gz   
data_type	FLOAT32
dim1		100
dim2		100
dim3		64
dim4		291
datatype	16
pixdim1		2.000000
pixdim2		2.000000
pixdim3		2.000000
pixdim4		2.000000
cal_max		255.000000
cal_min		0.000000
file_type	NIFTI-1+
$ fslstats sub-01_ses-04_task-blocked_run-1_bold.nii.gz -m -R
4021.941837 0.000000 49183.000000 

It would be great if BrainVoyager could fix this, but bvlabel could also provide a mechanism for detecting and fixing existing bogus values.

Submitted by Judith:

Since I am working on a script for a customer, I tested your prt implementation and it works well – thanks for that!
There is one minor thing that I noticed:
The timings are imported as float arrays, this will result in problems when saving and re-importing any changed PRTs in BV.
Therefore, I would suggest to use int arrays.

error for me in line 296, in write_fmr:

    data = header["ProtocolFile"]
    KeyError: 'ProtocolFile'

This error came up after using the same (unaltered) header file obtained from fmr.read_fmr for saving.
What seems to be the problem is that the read_fmr - because of its approach - is flexible in writing header info
(e.g: 'elif content[0] == "ProtocolFile" for line parsing). However in the write_fmr function is not, leading to KeyErrors for non existing items.

The timings are imported as float arrays, this will result in problems when saving and re-importing any changed PRTs in BV.
Using int arrays should solve the problem.

when I use the command to read a VTC file that was created by brain voyager it returns 0 for all voxel values. When I look at the VTC data from the same file inside brain voyager the values are not zero. The header has the correct information though.
VTC_header, VTC_data = bvbabel.vtc.read_vtc(VTC_path, rearrange_data_axes=False)

Submitted by Judith:

I don’t know whether you have implemented the PRT -> TSV function already.
Just in case you haven’t, I wanted to send you this script using bvbabel for writing TSV files.

# -*- coding: utf-8 -*-
"""
Save PRT as TSV file
Assuming the data has been saved according to the BIDS definition and that 
the PRTs are saved in the rawdata/func folder of each participant and session
and are named sub-ID_ses-ID_task-name_run-ID.prt
"""


__author__ = "Judith Eck"
__version__ = "0.1.0"
__date__ = "06-12-2022"
__name__ = "ConvertPRTtoTSV.py"


# =============================================================================
# Import required packages
# =============================================================================

import numpy as np
import bvbabel


# Define folder structure and file identifiers
project_path = 'C:/Users/JEck/Documents/BrainVoyager/Projects/GSG/rawdata/'
subjects = [1]
runs = [1]
sessions = [4]
task_name = 'task-blocked'

# specify the repetion time
TR = 2 # specified in seconds

# loop over all subjects, sessions and runs
for sub in range(len(subjects)):
    sub_id = 'sub-' + f'{subjects[sub]:02d}'
    for ses in range(len(sessions)):
        ses_id = 'ses-' + f'{sessions[ses]:02d}'
        
        for run in range(len(runs)):
            run_id = 'run-' + f'{runs[run]:02d}'
            
# Read the PRT file, convert it to _events.tsv and save in the same folder
            prtfile = project_path + sub_id + '/' + ses_id + '/func/' + sub_id + '_' + ses_id +  '_' + task_name + '_' + run_id + '.prt'
            prt_header, prt_data = bvbabel.prt.read_prt(prtfile)
            temp = np.empty(shape=(0,3), dtype = str)
            
            for cond in range(len(prt_data)):
                if prt_header['ResolutionOfTime'] == 'msec':
                    prt_data[cond]['Time start'] = np.int_(prt_data[cond]['Time start'])/1000
                    prt_data[cond]['Time stop'] = np.int_(prt_data[cond]['Time stop'])/1000
                else:
                    prt_data[cond]['Time start'] = ((prt_data[cond]['Time start'])-1)*TR
                    prt_data[cond]['Time stop'] = prt_data[cond]['Time stop']*TR
                duration = np.round(np.subtract(prt_data[cond]['Time stop'],prt_data[cond]['Time start']), decimals = 4)
                onset = prt_data[cond]['Time start']
                event_name = np.repeat(prt_data[cond]['NameOfCondition'],len(prt_data[cond]['Time start']))
                temp  = np.append(temp, np.array([onset.astype(str),duration.astype(str),event_name.astype(str)]).T, axis=0) 
                del (duration, onset, event_name)
            np.savetxt((prtfile.split('.')[0] + '_events.tsv'), temp, fmt='%s', delimiter='\t', header='onset\tduration\ttrial_type',comments='')
            del(prt_header, prt_data, prtfile, temp)
            

Hi Faruk, it would be great to have a read/write .glm function.
Thank you in advance,
Alessandra

Thanks for this terrific work. The permissive license will provides a nice mechanism for developers of open source tools to support data from the proprietary BrainVoyager products.

It would be great to include test images, both to validate and regression test the present code and to allow developers without BrainVoyager licenses to port your work to other languages.

There are two approaches:

1. nibabel includes testing scripts and validation images in their tests folder.

2. dcm2niix includes validation images as submodules (dcm_qa, dcm_qa_uih, dcm_qa_nih). The advantage of this approach is faster/smaller installation for end users, as they get the functionality without the bulky validation datasets.

This an example of a TBV ROI:
TEST_TBV.txt

Here is the same data saved in a BV ROI:
TEST_BV.txt

I am wondering if @nausikaa8 has additional test data for PRT? I have finished implementing read write functions. It would be helpful if these functions are tested with other PRT variants.

@nausikaa8 would you have time to help with this? Nothing urgent but your help is highly appreciated.

For many instances it could be preferred to read only the header information instead of the head + img information.
For example, one might want to know the start / end positions of a vtc file, with this information being in the header there is no need to load gigabytes of data to obtain this.
While for many people there is a easy fix (i.e. copy the read_xxx function and remove the data read portion), in the long run it might be good to have a solution to only read header data.

It would be useful to have an showing how to work with SRF (triangular mesh surface) and POI (patches of interest) files.

Example can be based on this dataset: https://search.kg.ebrains.eu/instances/Dataset/ff71a4d1-ea14-4ed6-898e-b92d95b3c446

*image from https://github.com/ofgulban/cortical-auditory-atlas .

This might require an additional flag to rearrange axes in vmr read write functions

Requested by @ju-ec , initial script provided by @ju-ec . Items to consider:

• Plots the motion estimates (as saved by BV)
• Plot framewise displacement with respect to the PRT files (optional).
• Plot a carpet plot to identify problematic volumes in the data.
• Save an SDM file containing a predictor of spike volumes (*_FDspikes.sdm) as identified based on a FD threshold set by the user in the script.
• Censor FMR volumes based on *_FDspikes.sdm.
• Interpolate in place of the censored volumes using linear op spline interpolation.