Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

DWI preprocessing #13

Closed
wants to merge 13 commits into from
68 changes: 55 additions & 13 deletions Dockerfile
Original file line number Diff line number Diff line change
@@ -1,12 +1,61 @@
FROM bids/base_freesurfer
# Use Ubuntu 16.04 LTS
FROM ubuntu:xenial-20161213

## Install the validator
RUN apt-get update && \
apt-get install -y curl && \
curl -sL https://deb.nodesource.com/setup_6.x | bash - && \
apt-get remove -y curl && \
apt-get install -y nodejs && \
rm -rf /var/lib/apt/lists/* /tmp/* /var/tmp/*

RUN npm install -g [email protected]

# Download FreeSurfer
RUN apt-get -y update \
&& apt-get install -y wget && \
wget -qO- ftp://surfer.nmr.mgh.harvard.edu/pub/dist/freesurfer/5.3.0-HCP/freesurfer-Linux-centos4_x86_64-stable-pub-v5.3.0-HCP.tar.gz | tar zxv -C /opt \
--exclude='freesurfer/trctrain' \
--exclude='freesurfer/subjects/fsaverage_sym' \
--exclude='freesurfer/subjects/fsaverage3' \
--exclude='freesurfer/subjects/fsaverage4' \
--exclude='freesurfer/subjects/fsaverage5' \
--exclude='freesurfer/subjects/fsaverage6' \
--exclude='freesurfer/subjects/cvs_avg35' \
--exclude='freesurfer/subjects/cvs_avg35_inMNI152' \
--exclude='freesurfer/subjects/bert' \
--exclude='freesurfer/subjects/V1_average' \
--exclude='freesurfer/average/mult-comp-cor' \
--exclude='freesurfer/lib/cuda' \
--exclude='freesurfer/lib/qt' && \
apt-get install -y tcsh bc tar libgomp1 perl-modules curl && \
rm -rf /var/lib/apt/lists/* /tmp/* /var/tmp/*

# Set up the environment
ENV OS Linux
ENV FS_OVERRIDE 0
ENV FIX_VERTEX_AREA=
ENV SUBJECTS_DIR /opt/freesurfer/subjects
ENV FSF_OUTPUT_FORMAT nii.gz
ENV MNI_DIR /opt/freesurfer/mni
ENV LOCAL_DIR /opt/freesurfer/local
ENV FREESURFER_HOME /opt/freesurfer
ENV FSFAST_HOME /opt/freesurfer/fsfast
ENV MINC_BIN_DIR /opt/freesurfer/mni/bin
ENV MINC_LIB_DIR /opt/freesurfer/mni/lib
ENV MNI_DATAPATH /opt/freesurfer/mni/data
ENV FMRI_ANALYSIS_DIR /opt/freesurfer/fsfast
ENV PERL5LIB /opt/freesurfer/mni/lib/perl5/5.8.5
ENV MNI_PERL5LIB /opt/freesurfer/mni/lib/perl5/5.8.5
ENV PATH /opt/freesurfer/bin:/opt/freesurfer/fsfast/bin:/opt/freesurfer/tktools:/opt/freesurfer/mni/bin:/usr/local/sbin:/usr/local/bin:/usr/sbin:/usr/bin:/sbin:/bin:$PATH

# Install FSL 5.0.9
RUN apt-get update && \
apt-get install -y --no-install-recommends curl && \
curl -sSL http://neuro.debian.net/lists/trusty.us-ca.full >> /etc/apt/sources.list.d/neurodebian.sources.list && \
curl -sSL http://neuro.debian.net/lists/xenial.us-ca.full >> /etc/apt/sources.list.d/neurodebian.sources.list && \
apt-key adv --recv-keys --keyserver hkp://pgp.mit.edu:80 0xA5D32F012649A5A9 && \
apt-get update && \
apt-get install -y fsl-core=5.0.9-3~nd14.04+1 && \
apt-get install -y fsl-core=5.0.9-1~nd+1+nd16.04+1 && \
rm -rf /var/lib/apt/lists/* /tmp/* /var/tmp/*

# Configure environment
Expand All @@ -23,7 +72,7 @@ ENV FSLOUTPUTTYPE=NIFTI_GZ
RUN echo "cHJpbnRmICJrcnp5c3p0b2YuZ29yZ29sZXdza2lAZ21haWwuY29tXG41MTcyXG4gKkN2dW12RVYzelRmZ1xuRlM1Si8yYzFhZ2c0RVxuIiA+IC9vcHQvZnJlZXN1cmZlci9saWNlbnNlLnR4dAo=" | base64 -d | sh

# Install Connectome Workbench
RUN apt-get update && apt-get -y install connectome-workbench=1.2.3-1~nd14.04+1
RUN apt-get update && apt-get -y install connectome-workbench=1.2.3-1~nd16.04+1

ENV CARET7DIR=/usr/bin

Expand Down Expand Up @@ -54,18 +103,11 @@ ENV HCPPIPEDIR_Global=${HCPPIPEDIR}/global/scripts
ENV HCPPIPEDIR_tfMRIAnalysis=${HCPPIPEDIR}/TaskfMRIAnalysis/scripts
ENV MSMBin=${HCPPIPEDIR}/MSMBinaries

RUN apt-get update && apt-get install -y --no-install-recommends python-pip python-six python-nibabel && \
RUN apt-get update && apt-get install -y --no-install-recommends python-pip python-six python-nibabel python-setuptools && \
rm -rf /var/lib/apt/lists/* /tmp/* /var/tmp/*
RUN pip install https://github.com/INCF/pybids/archive/158dac2062dc6b5a4ab2f92090108eedc3387575.zip
RUN pip install pybids==0.0.1
ENV PYTHONPATH=""

RUN apt-get update && \
curl -sL https://deb.nodesource.com/setup_4.x | bash - && \
apt-get install -y nodejs && \
rm -rf /var/lib/apt/lists/* /tmp/* /var/tmp/*

RUN npm install -g [email protected]

COPY run.py /run.py
RUN chmod +x /run.py

Expand Down
119 changes: 86 additions & 33 deletions run.py
Original file line number Diff line number Diff line change
Expand Up @@ -138,7 +138,7 @@ def run_generic_fMRI_volume_processsing(**args):
run(cmd, cwd=args["path"], env={"OMP_NUM_THREADS": str(args["n_cpus"])})

def run_generic_fMRI_surface_processsing(**args):
print(args)
# print(args)
args.update(os.environ)
cmd = '{HCPPIPEDIR}/fMRISurface/GenericfMRISurfaceProcessingPipeline.sh ' + \
'--path={path} ' + \
Expand All @@ -153,6 +153,7 @@ def run_generic_fMRI_surface_processsing(**args):
run(cmd, cwd=args["path"], env={"OMP_NUM_THREADS": str(args["n_cpus"])})

def run_diffusion_processsing(**args):
# print(args)
args.update(os.environ)
cmd = '{HCPPIPEDIR}/DiffusionPreprocessing/DiffPreprocPipeline.sh ' + \
'--posData="{posData}" ' +\
Expand All @@ -162,13 +163,16 @@ def run_diffusion_processsing(**args):
'--echospacing="{echospacing}" '+ \
'--PEdir={PEdir} ' + \
'--gdcoeffs="NONE" ' + \
'--dwiname="{dwiname}" ' + \
'--printcom=""'
cmd = cmd.format(**args)
print('\n',cmd, '\n')
run(cmd, cwd=args["path"], env={"OMP_NUM_THREADS": str(args["n_cpus"])})


__version__ = open('/version').read()

parser = argparse.ArgumentParser(description='HCP Pipeliens BIDS App (T1w, T2w, fMRI)')
parser = argparse.ArgumentParser(description='HCP Pipelines BIDS App (T1w, T2w, fMRI)')
parser.add_argument('bids_dir', help='The directory with the input dataset '
'formatted according to the BIDS standard.')
parser.add_argument('output_dir', help='The directory where the output files '
Expand All @@ -193,7 +197,7 @@ def run_diffusion_processsing(**args):
'fMRISurface', 'DiffusionPreprocessing'],
default=['PreFreeSurfer', 'FreeSurfer', 'PostFreeSurfer',
'fMRIVolume', 'fMRISurface',
'DiffusionPreprocessing'])
'DiffusionPreprocessing', 'TaskfMRIAnalysis'])
parser.add_argument('--license_key', help='FreeSurfer license key - letters and numbers after "*" in the email you received after registration. To register (for free) visit https://surfer.nmr.mgh.harvard.edu/registration.html',
required=True)
parser.add_argument('-v', '--version', action='version',
Expand Down Expand Up @@ -343,15 +347,15 @@ def run_diffusion_processsing(**args):
type='bold',
extensions=["nii.gz", "nii"])]
for fmritcs in bolds:
fmriname = "_".join(fmritcs.split("sub-")[-1].split("_")[1:]).split(".")[0]
fmriname = "_".join(fmritcs.split("sub-")[-1].split("_")[1:-1]).split(".")[0]
assert fmriname

fmriscout = fmritcs.replace("_bold", "_sbref")
if not os.path.exists(fmriscout):
fmriscout = "NONE"

fieldmap_set = layout.get_fieldmap(fmritcs)
print(fieldmap_set)
# print(fieldmap_set)
if fieldmap_set and fieldmap_set["type"] == "epi":
SEPhaseNeg = None
SEPhasePos = None
Expand Down Expand Up @@ -410,31 +414,80 @@ def run_diffusion_processsing(**args):
if stage in args.stages:
stage_func()

dwis = layout.get(subject=subject_label, type='dwi',
extensions=["nii.gz", "nii"])

# print(dwis)
# acqs = set(layout.get(target='acquisition', return_type='id',
# subject=subject_label, type='dwi',
# extensions=["nii.gz", "nii"]))
# print(acqs)
# posData = []
# negData = []
# for acq in acqs:
# pos = "EMPTY"
# neg = "EMPTY"
# dwis = layout.get(subject=subject_label,
# type='dwi', acquisition=acq,
# extensions=["nii.gz", "nii"])
# assert len(dwis) <= 2
# for dwi in dwis:
# dwi = dwi.filename
# if "-" in layout.get_metadata(dwi)["PhaseEncodingDirection"]:
# neg = dwi
# else:
# pos = dwi
# posData.append(pos)
# negData.append(neg)
#
# print(negData)
# print(posData)
posData = []
negData = []
PEdir = "None"
dwiname = "Diffusion"
dirnums = []
onerun = False

numruns = set(layout.get(target='run', return_type='id',
subject=subject_label, type='dwi',
extensions=["nii.gz", "nii"]))
# accounts for multiple runs, number of directions, and phase encoding directions
if not numruns:
onerun= True
numruns = {'run-01'}
if numruns:
for session in numruns:
if not onerun:
bvals = [f.filename for f in layout.get(subject=subject_label,
type='dwi', run=session,
extensions=["bval"])]
else:
bvals = [f.filename for f in layout.get(subject=subject_label,
type='dwi', extensions=["bval"])]
dwi_dict = {'bvalFile':[], 'bval':[], 'dwiFile':[], 'direction':[]}
for bvalfile in bvals:
with open(bvalfile) as f: # get number of directions
bvalues = [bvalue for line in f for bvalue in line.split()]
dwi_dict['bvalFile'].append(bvalfile)
dwi_dict['bval'].append(len(bvalues) - 1)
dwiFile = glob(os.path.join(os.path.dirname(bvalfile),'{0}.nii*'.format(os.path.basename(bvalfile).split('.')[0]))) # ensures bval file has same name as dwi file
assert len(dwiFile) == 1
dwi_dict['dwiFile'].append(dwiFile[0])
dwi_dict['direction'].append(layout.get_metadata(dwiFile[0])["PhaseEncodingDirection"][0])

# check if length of lists in dictionary are the same
n = len(dwi_dict['bvalFile'])
assert all(len(dwi_dict[k]) for k,v in dwi_dict.items())

for dirnum in set(dwi_dict['bval']):
idxs = { i for k,v in dwi_dict.iteritems() for i in range(0,len(dwi_dict['bval'])) if v[i] == dirnum }
PEdirNums = set([dwi_dict['direction'][i] for i in idxs])
for PEdirNum in PEdirNums:
dwis = [ dwi_dict['dwiFile'][i] for i in idxs if dwi_dict['direction'][i] == PEdirNum ]
assert len(dwis) <= 2
dwiname = "Diffusion" + "_dir-" + str(dirnum) + "_" + session + "_corr_" + str(PEdirNum)
if "j" in PEdirNum:
PEdir = 2
elif "i" in PEdirNum:
PEdir = 1
else:
RuntimeError("Phase encoding direction not specified for diffusion data.")
pos = "EMPTY"
neg = "EMPTY"
gdcoeffs = "None"
for dwi in dwis:
if "-" in layout.get_metadata(dwi)["PhaseEncodingDirection"]:
neg = dwi
else:
pos = dwi
try:
echospacing = layout.get_metadata(pos)["EffectiveEchoSpacing"] * 1000
dwi_stage_dict = OrderedDict([("DiffusionPreprocessing", partial(run_diffusion_processsing,
posData=pos,
negData=neg,
path=args.output_dir,
subject="sub-%s" % subject_label,
echospacing=echospacing,
PEdir=PEdir,
gdcoeffs="NONE",
dwiname=dwiname,
n_cpus=args.n_cpus))])
for stage, stage_func in dwi_stage_dict.iteritems():
if stage in args.stages:
stage_func()
except:
print("You may have missing diffusion data in the positive phase encoding direction.")