A Fast, Validated Open Source Toolkit for Dynamic Contrast Enhanced MRI Analysis
This is the repository for the paper published in PeerJ: DCEMRI.jl: a fast, validated, open source toolkit for dynamic contrast enhanced MRI analysis.
From the Julia website,
Julia is a high-level, high-performance dynamic programming language for technical computing, with syntax that is familiar to users of other technical computing environments. It provides a sophisticated compiler, distributed parallel execution, numerical accuracy, and an extensive mathematical function library. The library, largely written in Julia itself, also integrates mature, best-of-breed C and Fortran libraries for linear algebra, random number generation, signal processing, and string processing.
Put simply, it looks like Matlab, which is simple to learn and familiar to most MRI researchers, but it works better and faster and is completely free. In particular, for the problem of DCE MRI, Julia's simple and flexible parallel computing model allows almost perfect parallelization of the nonlinear least squares fitting problem. In my informal testing, the intrinsic speed of Julia coupled to my parallel implementation produced a factor of 20-40 speedup over comparable Matlab and Python.
Installation is simple. First you need Julia. The simplest way to get Julia is to grab the current release version from julialang.org.
Next you need DCEMRI.jl. Open Julia. You should see a terminal window with the julia> prompt. This is analogous to the command line in Matlab.
To install DCEMRI.jl, first type ] at the prompt to get into the package shell then type:
(@v1.5) pkg> add DCEMRI
This might take a minute, because the validation data must be downloaded, as well as a few supporting packages.
If you want to DCEMRI.jl to create plots for you, you also need Python with Matplotlib installed. Most OS X and Linux machines probably already have a version of Python with Matplotlib installed. If you don't have Python with Matplotlib, such as if you are on Windows, you can grab the excellent, free Anaconda python distribution, which comes with Matplotlib pre-installed.
If you have installation problems, check the next section on platform specific issues. If none of that helps, please file an issue at the github repository, and we will try to help you.
If you get an error that xcrun is missing on Mac OS X, follow the instructions to install the Developer Tools Command Line Tools, and then switch back to your Julia window and run Pkg.build("HDF5"). This should fix everything. If you don't get this error, but installation fails, make sure you have Xcode installed, along with the optional Command Line Tools.
There have also been reports of Julia popping up behind other windows on OS X. Look behind your windows if you think Julia isn't starting.
On Windows, you will almost certainly need to install Anaconda to have a Python distribution. Some Windows versions claim to be 64-bit, but aren't quite. If you are having trouble on Windows, and you have installed a 64-bit Anaconda and Julia, you can try uninstalling them and installing 32-bit versions of both.
After installing the Julia and the DCEMRI module, you should run the validations, to make sure the calculations work correctly on your machine. The easiest way to do this is to start Julia and then run
julia> using DCEMRI
julia> validate()
This will run both validations (4 and 6), which could take up to an hour, depending on the number of cores you started Julia with and how fast your hardware is. Examine the results to make sure that the parameters have been recovered accurately. The text output of the scripts will also print quantitative measures of parameter accuracy. By default validation results will be written to subdirectories of the DCEMRI.jl module directory, which is queried with Pkg.dir("DCEMRI"). You can customize this directory for the validation by passing the output path as the second argument. This only works the individual validation functions, e.g. validate(6, "/my/path") or validate(4, "/my/path"). Since both validations write figures with identical names, they can't go to the same output directory.
An example output is shown here:
julia> validate(4)
Running analysis of noise-free QIBA v4 data ...
running models
found multi-flip data
fitting R1 relaxation rate to multi-flip data
fitting 6 x 23 points on each of 4 workers
processed 90 voxels in 2.2 s (41.5 vox/s)
computing signal enhancement ratios
converting DCE signal to effective R1
converting effective R1 to tracer tissue concentration Ct
fitting DCE data
attempting Extended Tofts-Kety model
fitting 661 x 23 points on each of 4 workers
processed 90 voxels in 3.8 s (23.5 vox/s)
saving results to /Users/dss/.julia/v0.3/DCEMRI/test/q4/results/results.mat
Plotting results ...
Kt RMSE (%): 6.97465437361441
Kt max error (%): 23.493640353851994
Kt CCC: 0.9998009845162595
ve RMSE (%): 18.02170557638968
ve max error (%): 99.99999999999996
ve CCC: 0.8904290685710147
vp RMSE (%): 23.770196145538407
vp max error (%): 92.10583127104924
vp CCC: 0.9999200988268792
Running analysis of noisy QIBA v4 data ...
running models
found multi-flip data
fitting R1 relaxation rate to multi-flip data
fitting 6 x 2250 points on each of 4 workers
processed 9000 voxels in 0.5 s (19436.3 vox/s)
computing signal enhancement ratios
converting DCE signal to effective R1
converting effective R1 to tracer tissue concentration Ct
fitting DCE data
attempting Extended Tofts-Kety model
fitting 661 x 2250 points on each of 4 workers
processed 9000 voxels in 341.7 s (26.3 vox/s)
saving results to /Users/dss/.julia/v0.3/DCEMRI/test/q4/results_noisy/results.mat
Plotting results ...
Kt RMSE (%): 11.311615941962662
Kt max error (%): 100.0
Kt CCC: 0.9742179876687028
ve RMSE (%): 18.238054961776477
ve max error (%): 100.0
ve CCC: 0.7026132423939505
vp RMSE (%): 12.654024477709797
vp max error (%): 100.0
vp CCC: 0.9717255972607232
Validation complete. Results can be found in /Users/dss/.julia/v0.3/DCEMRI/test/q4.
To perform the validation on the Quantitative Imaging Biomarkers Alliance phantoms for yourself from the original DICOMS (not necessary, see below), you will need to download the DICOMS from Daniel Barboriak's Lab. Then the scripts in the q4 and q6 folders will help you translate the DICOM data to MAT files suitable for input into the Julia code.
BUT I have already done this step for you and included the MAT files. This also avoids you needing to install Python if you don't have it already. If you want to install Python and run the scripts to convert the DICOM data to MAT files, then I recommend the Anaconda Python distribution. It has everything you need for scientific programming with Python.
You can run the in vivo data demo with the command
demo(). It will save the output by default to Pkg.dir("DCEMRI","demo/results") by default. You can change this location by passing a path string to demo(). After a few seconds to a few minutes, depending on the speed of your machine, you will see the following output text:
julia> demo()
Processing in vivo data ...
running models
found multi-flip data
fitting R1 relaxation rate to multi-flip data
fitting 10 x 4582 points on each of 4 workers
processed 18327 voxels in 1.0 s (19055.6 vox/s)
computing signal enhancement ratios
converting DCE signal to effective R1
converting effective R1 to tracer tissue concentration Ct
fitting DCE data
attempting Standard Tofts-Kety model
fitting 25 x 1694 points on each of 4 workers
processed 6774 voxels in 1.1 s (5928.8 vox/s)
saving results to results/results.mat
Plotting results ...
Demo run complete.
Results can be found in /Users/dss/.julia/v0.3/DCEMRI/demo/results
All units in the code are SI where possible. Sometimes, due to numerical accuracy issues, they have been converted internally. But all data should be supplied to the code in SI units. In particular, time should be in seconds, and relaxation rates in inverse seconds. The two exceptions to this rule are that flip angles should be in degrees and Ktrans is output in min^-1.
In the simplest incarnation, if you already have a MAT file containing your data, you can run the analysis from within Julia using
julia> using DCEMRI
julia> results = fitdata(datafile="/path/to/your/datafile.mat")
DCEMRI.jl will look for parameters in the input MAT file, and if they are found will use them. Anything not found in the MAT file will be initialized from the defaults. These defaults can be viewed with the defaults() command. You may also override both the MAT file and the defaults by passing keyword arguments to fitdata.
DCEMRI.jl has two basic modes of operation. The first is command-line invocation, like an operating system command. To run it as a shell command, first edit the first line of bin/dcefit to point to where you installed your Julia binary, as in
#!/path/to/julia/binary
Next, make sure bin/dcefit is executable. It should already be, but it doesn't hurt to check. Next copy it to a directory that is in your shell's search path. A good place on UNIX systems, such as Mac OS X, is /usr/local/bin.
dcefit can parse arguments passed on the command line to configure the model and point to the input data and output file. To see the available options, run dcefit -h at the terminal prompt, you will get
usage: dcefit [-O OUTFILE] [-R RELAXIVITY] [-r TR] [-d DCEFLIP]
[-c SERCUTOFF] [-t T1FLIP [T1FLIP...]]
[-m MODELS [MODELS...]] [-p] [-w WORKERS] [-v] [-h]
[datafile]
Process DCE-MRI data. Optional arguments can be used to override any
values found in input files. For questions, contact David Smith
<[email protected]>. For bug reports and feature requests, file an
issue at http://github.com/davidssmith/DCEMRI.jl
positional arguments:
datafile path to MAT file containing DCE and T1 data
(default: "input.mat")
optional arguments:
-O, --outfile OUTFILE
path to MAT file to contain the ouput
(default: "results.mat")
-R, --relaxivity RELAXIVITY
contrast agent relaxivity (1/s) (type:
Float64)
-r, --TR TR repetition time (ms) (type: Float64)
-d, --DCEflip DCEFLIP
flip angle of DCE data (type: Float64)
-c, --SERcutoff SERCUTOFF
minimum SER to include in processing mask
(type: Float64)
-t, --T1flip T1FLIP [T1FLIP...]
list of flip angle(s) of T1 data (type:
Float64)
-m, --models MODELS [MODELS...]
list of models: 1=plasma only, 2=Standard,
3=Extended (type: Int64)
-p, --plotting plot intermediate results
-w, --workers WORKERS
number of parallel workers to use (one per CPU
core is good) (type: Int64, default: 4)
-v, --verbose show verbose output
-h, --help show this help message and exit
To process a DCEMRI data set from the command line, the minimum invocation is
dcefit /path/to/my/dce/data.mat.
The input data MAT file must contain the following:
Cp: Arterial input function as a vector, resampled to the DCE time points.DCEdata: DCE data as a 3-D array (1 time by 2 space dimensions).DCEflip: flip angle in deg of DCE datat: time vector representing the dcedata samples.TR: repetition time of DCE scan- R1 information as either
R10andS0, representing pre-calculated R1 relaxation maps, or asT1data, indicating that a multi-flip scan was performed and must be analyzed. IfT1datais supplied, the code also needsT1flip, a vector of flip angles (in degrees) for the multi-flip data.
The results will be saved in the current directory as results.mat. You can override the output file name and location with the --outfile flag.
If you've made it this far, you are ready to run the DCE analysis on your own data. Congratulations! If you have problems or find bugs, please file an issue on the github repository or email us. If you find ways to make it better, please submit your improvements as well. We hope that this can become a community effort that leads to an outstanding, rock solid, trustworthy tool.
To keep your installation of DCEMRI.jl up to date, periodically run Pkg.update() at the julia> prompt.
This project was funded by the National Cancer Institute of the National Institutes of Health, under Award Numbers K25CA176219, U01CA142565, R01CA129961, R25CA092043. The content is solely the responsibility of the authors and does not necessarily represent the official views of the National Institutes of Health.