PCC random noise + measurement reduction

breaker.py

import nibabel as nb
import numpy.ma as ma
import numpy as np
from nilearn.image import resample_img
import pp
from glob import glob
import os
import re

def main():
    workshopDataLoc = '/media/workshop/openData_for_workshop/NITRC/workshop'
    subjects = [x for x in os.listdir(workshopDataLoc) if x.startswith('subj')]

    group = 0
    for subject in subjects:
        subj_REST_in = glob(os.path.join(workshopDataLoc, subject, 'REST/*.nii'))[0]
        if group == 0:
            print subject, 'is group 0'
            subj_REST_out = re.sub('fMRI.nii', 'DMN-broken-fMRI.nii.gz', subj_REST_in)
            group = 1
        elif group == 1:
            print subject, 'is group 1'
            subj_REST_out = re.sub('fMRI.nii', 'DMN-broken-fMRI.nii.gz', subj_REST_in, False)
            group = 0


        DMN_break(subj_REST_in, subj_REST_out)

def DMN_break(inputImg, outputImg, pccCenter=[38, 49, 23], breakpcc=True):
    size = 3

    x_start = pccCenter[0] - size
    x_end = pccCenter[0] + size
    y_start = pccCenter[1] - size
    y_end = pccCenter[1] + size
    z_start = pccCenter[2] - size
    z_end = pccCenter[2] + size

    data = nb.load(inputImg)
    fMRI_data = data.get_data()

    shape = data.get_data().shape

    mask = np.zeros(data.get_data()[:,:,:,:].shape)
    mask[x_start:x_end,
         y_start:y_end,
         z_start:z_end,
         :] = 1

    percentile = 90
    mu = fMRI_data[fMRI_data > np.percentile(fMRI_data, percentile)].mean()
    sigma = fMRI_data[fMRI_data > np.percentile(fMRI_data, percentile)].std()

    s = np.random.normal(mu, sigma, fMRI_data.shape).reshape(fMRI_data.shape)

    if breakpcc:
        fMRI_data[mask==1] = np.absolute(s[mask==1])
    else:
        pass

    modified_fMRI = nb.Nifti1Image(fMRI_data[:,:,:,:150], data.affine, header=data.header)
    modified_fMRI.header.set_zooms(data.header.get_zooms()[:3] + (2,))


    modified_fMRI.to_filename(outputImg)


if __name__=='__main__':
    main()