ICA FIX implemented

Dockerfile

@@ -6,6 +6,19 @@ ARG DEBIAN_FRONTEND=noninteractive
 ENV LANG="C.UTF-8" \
     LC_ALL="C.UTF-8"
 
+#adding jessie backport
+#Installing R prerequisites for fix
+
+RUN sh -c 'echo "deb http://ftp.debian.org/debian jessie-backports main" >> /etc/apt/sources.list' && \
+    apt-get -qq update && \
+    apt-get -t jessie-backports -y install r-base r-base-dev mc nano wget
+
+RUN Rscript -e 'install.packages(c("kernlab","ROCR", "class", "party", "e1071", "randomForest"), dependencies=TRUE,repos="http://cran.cnr.berkeley.edu/")'
+#older version of party is required (I do it this way, becasue otherwise I would need all the dependencies of party which was already installed above)
+RUN wget https://cran.r-project.org/src/contrib/Archive/party/party_1.0-25.tar.gz && \
+	R CMD INSTALL party_1.0-25.tar.gz && \
+	rm party_1.0-25.tar.gz
+
 # Download FreeSurfer
 RUN apt-get -qq update && \
     apt-get install -yq --no-install-recommends \
@@ -55,10 +68,6 @@ ENV OS=Linux \
     MNI_PERL5LIB=/opt/freesurfer/mni/lib/perl5/5.8.5 \
     PATH=/opt/freesurfer/bin:/opt/freesurfer/fsfast/bin:/opt/freesurfer/tktools:/opt/freesurfer/mni/bin:$PATH
 
-# Install MCR 2016b
-ENV MATLABCMD="/opt/matlabmcr-2016b/v91/toolbox/matlab" \
-    MATLAB_COMPILER_RUNTIME="/opt/matlabmcr-2016b/v91" \
-    LD_LIBRARY_PATH="/opt/matlabmcr-2016b/v91/runtime/glnxa64:/opt/matlabmcr-2016b/v91/bin/glnxa64:/opt/matlabmcr-2016b/v91/sys/os/glnxa64:$LD_LIBRARY_PATH"
 
 RUN apt-get update -qq \
     && apt-get install -y -q --no-install-recommends \
@@ -67,12 +76,36 @@ RUN apt-get update -qq \
        libxt6 \
        unzip \
     && apt-get clean \
-    && rm -rf /var/lib/apt/lists/* /tmp/* /var/tmp/* \
-    && echo "Downloading MATLAB Compiler Runtime ..." \
-    && curl -fsSL --retry 5 -o /tmp/mcr.zip https://ssd.mathworks.com/supportfiles/downloads/R2016b/deployment_files/R2016b/installers/glnxa64/MCR_R2016b_glnxa64_installer.zip \
-    && unzip -q /tmp/mcr.zip -d /tmp/mcrtmp \
-    && /tmp/mcrtmp/install -destinationFolder /opt/matlabmcr-2016b -mode silent -agreeToLicense yes \
-    && rm -rf /tmp/*
+    && rm -rf /var/lib/apt/lists/* /tmp/* /var/tmp/* 
+
+
+#    \
+#    && echo "Downloading MATLAB Compiler Runtime ..." \
+#    && curl -fsSL --retry 5 -o /tmp/mcr.zip https://ssd.mathworks.com/supportfiles/downloads/R2016b/deployment_files/R2016b/installers/glnxa64/MCR_R2016b_glnxa64_installer.zip \
+#    && unzip -q /tmp/mcr.zip -d /tmp/mcrtmp \
+#    && /tmp/mcrtmp/install -destinationFolder /opt/matlabmcr-2016b -mode silent -agreeToLicense yes \
+#    && rm -rf /tmp/*
+
+# Install MCR 2016b
+#ENV MATLABCMD="/opt/matlabmcr-2016b/v91/toolbox/matlab" 
+#ENV MATLAB_COMPILER_RUNTIME="/opt/matlabmcr-2016b/v91"
+#ENV LD_LIBRARY_PATH="/opt/matlabmcr-2016b/v91/runtime/glnxa64:/opt/matlabmcr-2016b/v91/bin/glnxa64:/opt/matlabmcr-2016b/v91/sys/os/glnxa64:$LD_LIBRARY_PATH"
+
+
+
+  ##download and install fix
+RUN wget http://www.fmrib.ox.ac.uk/~steve/ftp/fix.tar.gz -O fix.tar.gz && \
+	cd /opt/ && \
+	tar zxvf /fix.tar.gz && \
+    mv /opt/fix* /opt/fix && \
+    rm /fix.tar.gz && \
+    cd /opt/fix && \
+    cd /
+
+
+##environmnetal variables mode=0 means MATLAB compiled version , mode=2 would mean Octave   
+ENV FSL_FIX_MATLAB_MODE=0
+ENV FSL_FIXDIR=/opt/fix		
 
 # Install miniconda2
 ENV PATH="/usr/local/miniconda/bin:$PATH"
@@ -154,12 +187,38 @@ RUN echo "deb http://ftp.de.debian.org/debian stretch main" >> /etc/apt/sources.
     apt-get install -y libstdc++6
 
 # overwrite matlab mcr shared object
-RUN rm /opt/matlabmcr-2016b/v91/sys/os/glnxa64/libstdc++.so.6 && \
-    ln -s /usr/lib/x86_64-linux-gnu/libstdc++.so.6 /opt/matlabmcr-2016b/v91/sys/os/glnxa64/libstdc++.so.6
+#RUN rm /opt/matlabmcr-2016b/v91/sys/os/glnxa64/libstdc++.so.6 && \
+#    ln -s /usr/lib/x86_64-linux-gnu/libstdc++.so.6 /opt/matlabmcr-2016b/v91/sys/os/glnxa64/libstdc++.so.6
+
+RUN ln -s $FSLDIR /usr/local/fsl
+RUN mkdir -p /vols/Data/HCP && \
+	cd /vols/Data/HCP && \
+	wget https://ftp.humanconnectome.org/workbench/workbench-linux64-v1.2.3.zip && \
+	unzip workbench-linux64-v1.2.3.zip && \
+	cd /	    
 
 # install gradient_unwarp.py (v1.0.3)
 RUN pip install https://github.com/Washington-University/gradunwarp/archive/v1.0.3.zip
 
+#try MALTAB already compiled
+RUN cd $FSL_FIXDIR/compiled/Linux/x86_64 && \
+	cp MCRInstaller.zip /tmp && \
+	 cd /tmp && \
+	 unzip MCRInstaller.zip && \ 
+	 ./install -mode silent -agreeToLicense yes
+
+COPY settings.sh /opt/fix/settings.sh
+RUN chmod +x /opt/fix/settings.sh
+RUN mv /opt/fix/hcp_fix /opt/fix/hcp_fix.old
+COPY hcp_fix /opt/fix/hcp_fix
+RUN chmod +x /opt/fix/hcp_fix
+
+RUN wget http://argyelan.com/ZHH/ZHH30_hp2000.RData && \
+    mv ZHH30_hp2000.RData /opt/fix/training_files/ZHH30_hp2000.RData
+COPY zhh_fix /opt/fix/zhh_fix
+RUN chmod +x /opt/fix/zhh_fix
+
+#COPY MCR.version /opt/fix/MCR.version
 COPY run.py version /
 RUN chmod +x /run.py
 

hcp_fix

@@ -0,0 +1,112 @@
+#!/bin/sh
+
+#   hcp_fix - wrapper script for FIX, for HCP pipelines
+#
+#   Stephen Smith, FMRIB Analysis Group
+#
+#   Copyright (C) 2012-2013 University of Oxford
+#
+#   SHCOPYRIGHT
+
+FSL_FIXDIR=$( cd $(dirname $0) ; pwd)
+export FSL_FIXDIR
+# All fix settings are held in the settings.sh file - edit this file to suit your setup
+. ${FSL_FIXDIR}/settings.sh
+
+#############################################################
+
+Usage() {
+    cat <<EOF
+
+hcp_fix <4D_FMRI_data> <highpass>
+  with <highpass> being the temporal highpass full-width (2*sigma) in seconds
+
+e.g.   hcp_fix BOLD_REST1_RL/BOLD_REST1_RL.nii.gz 200
+
+for detrending-like behaviour, set <highpass> to 2000
+
+EOF
+    exit 1
+}
+
+[ "$2" = "" ] && Usage
+
+fmri=$1
+cd `dirname $fmri`
+fmri=`basename $fmri`
+fmri=`$FSLDIR/bin/imglob $fmri`
+[ `imtest $fmri` != 1 ] && echo No valid 4D_FMRI input file specified && exit 1
+fmri_orig=$fmri
+
+hp=$2
+
+tr=`$FSLDIR/bin/fslval $fmri pixdim4`
+
+echo "processing FMRI file $fmri with highpass $hp"
+
+if [ $hp -gt 0 ] ; then
+  echo "running highpass"
+  hptr=`echo "10 k $hp 2 / $tr / p" | dc -`
+  ${FSLDIR}/bin/fslmaths $fmri -Tmean ${fmri}_Tmean
+  ${FSLDIR}/bin/fslmaths $fmri -bptf $hptr -1 -add ${fmri}_Tmean ${fmri}_hp$hp
+  fmri=${fmri}_hp$hp
+fi
+
+echo "running MELODIC"
+mkdir -p ${fmri}.ica
+$FSLDIR/bin/melodic -i $fmri -o ${fmri}.ica/filtered_func_data.ica -d -250 --nobet --report --Oall --tr=$tr
+
+cd ${fmri}.ica
+
+$FSLDIR/bin/imln ../$fmri filtered_func_data
+$FSLDIR/bin/imln filtered_func_data.ica/mask mask
+
+if [ `$FSLDIR/bin/imtest ../${fmri_orig}_SBRef` = 1 ] ; then
+  $FSLDIR/bin/imln ../${fmri_orig}_SBRef mean_func
+else
+  $FSLDIR/bin/imln filtered_func_data.ica/mean mean_func
+fi
+
+if [ -f ../${fmri_orig}_Atlas.dtseries.nii ] ; then
+  $FSLDIR/bin/imln ../${fmri_orig}_Atlas.dtseries.nii Atlas.dtseries.nii
+fi
+
+mkdir -p mc
+if [ -f ../Movement_Regressors.txt ] ; then
+  cat ../Movement_Regressors.txt | awk '{ print $4 " " $5 " " $6 " " $1 " " $2 " " $3}' > mc/prefiltered_func_data_mcf.par
+else
+  echo "ERROR: Movement_Regressors.txt not retrieved properly." 
+  exit -1
+fi 
+
+${FSL_FIXDIR}/call_matlab.sh -l .fix.log -f functionmotionconfounds $tr $hp 
+
+mkdir -p reg
+cd reg
+$FSLDIR/bin/imln ../../../../T1w_restore_brain highres
+$FSLDIR/bin/imln ../../../../wmparc wmparc
+$FSLDIR/bin/imln ../mean_func example_func
+$FSLDIR/bin/makerot --theta=0 > highres2example_func.mat
+if [ `$FSLDIR/bin/imtest ../../../../T2w` = 1 ] ; then
+  $FSLDIR/bin/fslmaths ../../../../T1w -div ../../../../T2w veins -odt float
+  $FSLDIR/bin/flirt -in ${FSL_FIXDIR}/mask_files/hcp_0.7mm_brain_mask -ref veins -out veinbrainmask -applyxfm
+  $FSLDIR/bin/fslmaths veinbrainmask -bin veinbrainmask
+  $FSLDIR/bin/fslmaths veins -div `$FSLDIR/bin/fslstats veins -k veinbrainmask -P 50` -mul 2.18 -thr 10 -min 50 -div 50 veins
+  $FSLDIR/bin/flirt -in veins -ref example_func -applyxfm -init highres2example_func.mat -out veins_exf
+  $FSLDIR/bin/fslmaths veins_exf -mas example_func veins_exf
+fi
+cd ../..
+
+echo "running FIX"
+if [ $hp != 2000 ] ; then
+  ${FSL_FIXDIR}/fix ${fmri}.ica ${FSL_FIXDIR}/training_files/HCP_hp200.RData 10 -m -h 200
+else
+  ${FSL_FIXDIR}/fix ${fmri}.ica ${FSL_FIXDIR}/training_files/HCP_hp2000.RData 10 -m -h 2000
+fi
+
+$FSLDIR/bin/immv ${fmri}.ica/filtered_func_data_clean ${fmri}_clean
+
+if [ -f ${fmri}.ica/Atlas_clean.dtseries.nii ] ; then
+  /bin/mv ${fmri}.ica/Atlas_clean.dtseries.nii ${fmri_orig}_Atlas_hp${hp}_clean.dtseries.nii
+fi
+

run.py

@@ -152,7 +152,16 @@ def run_generic_fMRI_surface_processsing(**args):
       '--regname="{regname}"'
     cmd = cmd.format(**args)
     run(cmd, cwd=args["path"], env={"OMP_NUM_THREADS": str(args["n_cpus"])})
-
+
+def run_generic_fMRI_ICAFIX_processsing(**args):
+    args.update(os.environ)
+    cmd = '/opt/fix/hcp_fix ' +  '{path}' + '/' + '{subject}' + '/MNINonLinear/Results/' + '{fmriname}'  + '/' + '{fmriname}' + '.nii.gz 2000'
+    cmd = cmd.format(**args)
+    run(cmd, cwd=args["path"], env={"OMP_NUM_THREADS": str(args["n_cpus"])})
+    cmd = '/opt/fix/zhh_fix ' +  '{path}' + '/' + '{subject}' + '/MNINonLinear/Results/' + '{fmriname}'  + '/' + '{fmriname}' + '.nii.gz 2000'
+    cmd = cmd.format(**args)
+    run(cmd, cwd=args["path"], env={"OMP_NUM_THREADS": str(args["n_cpus"])})
+
 def run_diffusion_processsing(**args):
     args.update(os.environ)
     cmd = '{HCPPIPEDIR}/DiffusionPreprocessing/DiffPreprocPipeline.sh ' + \
@@ -191,9 +200,9 @@ def run_diffusion_processsing(**args):
 parser.add_argument('--stages', help='Which stages to run. Space separated list.',
                    nargs="+", choices=['PreFreeSurfer', 'FreeSurfer',
                                        'PostFreeSurfer', 'fMRIVolume',
-                                       'fMRISurface', 'DiffusionPreprocessing'],
+                                       'fMRISurface', 'ICAFIX', 'DiffusionPreprocessing'],
                    default=['PreFreeSurfer', 'FreeSurfer', 'PostFreeSurfer',
-                            'fMRIVolume', 'fMRISurface',
+                            'fMRIVolume', 'fMRISurface', 'ICAFIX',
                             'DiffusionPreprocessing'])
 parser.add_argument('--coreg', help='Coregistration method to use',
                     choices=['MSMSulc', 'FS'], default='MSMSulc')
@@ -362,11 +371,17 @@ def run_diffusion_processsing(**args):
                 SEPhaseNeg = None
                 SEPhasePos = None
                 for fieldmap in fieldmap_set:
-                    enc_dir = layout.get_metadata(fieldmap["epi"])["PhaseEncodingDirection"]
-                    if "-" in enc_dir:
+                    enc_dir = fieldmap['epi'].split('_dir-')[1].split('_')[0]
+                    if "AP" in enc_dir:
                         SEPhaseNeg = fieldmap['epi']
                     else:
                         SEPhasePos = fieldmap['epi']
+                #for fieldmap in fieldmap_set:
+                #    enc_dir = layout.get_metadata(fieldmap["epi"])["PhaseEncodingDirection"]
+                #    if "-" in enc_dir:
+                #        SEPhaseNeg = fieldmap['epi']
+                #    else:
+                #        SEPhasePos = fieldmap['epi']
                 echospacing = layout.get_metadata(fmritcs)["EffectiveEchoSpacing"]
                 unwarpdir = layout.get_metadata(fmritcs)["PhaseEncodingDirection"]
                 unwarpdir = unwarpdir.replace("i","x").replace("j", "y").replace("k", "z")
@@ -412,7 +427,15 @@ def run_diffusion_processsing(**args):
                                                        n_cpus=args.n_cpus,
                                                        grayordinatesres=grayordinatesres,
                                                        lowresmesh=lowresmesh,
-                                                       regname=args.coreg))
+                                                       regname=args.coreg)),
+                                ("ICAFIX", partial(run_generic_fMRI_ICAFIX_processsing,
+                                                       path=args.output_dir,
+                                                       subject="sub-%s"%subject_label,
+                                                       fmriname=fmriname,
+                                                       fmrires=fmrires,
+                                                       n_cpus=args.n_cpus,
+                                                       grayordinatesres=grayordinatesres,
+                                                       lowresmesh=lowresmesh))                                                                    
                                 ])
             for stage, stage_func in func_stages_dict.iteritems():
                 if stage in args.stages: