def create_rois_from_clusters(contrast_tmap, mask, threshold=3.09,
                              height_control='brute', cluster_threshold=10,
                              save_path=None):
    if save_path is not None:
        if not os.path.exists(save_path):
            os.makedirs(save_path)

    thresholded = map_threshold(contrast_tmap, mask, threshold,
                                height_control, cluster_threshold)
    cluster_map, n_cluster = label(thresholded.get_data() > 0)

    clusters = []
    masker = NiftiMasker(mask_img=mask)
    masker.fit()
    mask_affine = nib.load(mask).get_affine()
    for label_ in range(1, n_cluster + 1):
        cluster = cluster_map.copy()
        cluster[cluster_map != label_] = 0
        cluster[cluster_map == label_] = 1
        cluster = nib.Nifti1Image(cluster, mask_affine)
        clusters.append(cluster)
        if save_path is not None:
            nib.save(cluster, os.path.join(save_path,
                     'cluster_{0}.nii'.format(label_)))

    return clusters

def preprocess(num, subj, subj_dir, subj_warp_dir, force_warp=False):
    bold_path = 'BOLD/task001_run00%i/bold_dico_bold7Tp1_to_subjbold7Tp1.nii.gz' % (num+1)
    bold_path = os.path.join(DATA_DIR, subj, bold_path)
    template_path = os.path.join(DATA_DIR, 'templates', 'grpbold7Tp1', 'brain.nii.gz')
    warp_path = os.path.join(DATA_DIR, subj, 'templates', 'bold7Tp1', 'in_grpbold7Tp1', 'subj2tmpl_warp.nii.gz')

    output_path = os.path.join(subj_warp_dir, 'run00%i.nii.gz' % num)

    if force_warp or not os.path.exists(output_path):
        print 'Warping image #%i...' % num
        subprocess.call(['fsl5.0-applywarp', '-i', bold_path, '-o', output_path, '-r', template_path, '-w', warp_path, '-d', 'float'])
    else:
        print 'Reusing cached warp image #%i' % num

    print 'Loading image #%i...' % num
    bold = load(output_path)

    masker = NiftiMasker(load(MASK_FILE))
    # masker = niftimasker(load(MASK_FILE), detrend=true, smoothing_fwhm=4.0,
    #                     high_pass=0.01, t_r=2.0, standardize=true)
    masker.fit()
    print 'Removing confounds from image #%i...' % num
    data = masker.transform(bold, confounds(num, subj))
    print 'Detrending image #%i...' % num
    filtered = np.float32(savgol_filter(data, 61, 5, axis=0))
    img = masker.inverse_transform(data-filtered)
    print 'Smoothing image #%i...' % num
    img = image.smooth_img(img, 4.0)
    print 'Saving image #%i...' % num
    save(img, os.path.join(subj_dir, 'run00%i.nii.gz' % num))
    print 'Finished with image #%i' % num

def MaskFlatten(concat_dict, mask, iter_n):
    '''Mask image data, convert to 2D feature matrix'''
    nifti_masker = NiftiMasker(mask_img=mask)
    masked_dict = {}
    for i in range(iter_n):
        masked_dict[i] = nifti_masker.fit_transform(concat_dict[i])
    return masked_dict

def ts(img_path,
	mask=False,
	substitution={},
	):
	"""
	Return the mean and median of a Region of Interest (ROI) time course.

	Parameters
	----------

	img_path : str
		Path to NIfTI file from which the ROI is to be extracted.
	maks : nilearn.NiftiMasker or str, optional
		Nilearn `nifti1.Nifti1Image` object to use for masking the desired ROI, or a string specifying the path of a maskfile.
	substitution : dict, optional
		A dictionary with keys which include 'subject' and 'session'.
	"""
	if substitution:
		img_path = img_path.format(**substitution)
	img_path = path.abspath(path.expanduser(img_path))
	img = nib.load(img_path)
	try:
		masked_data = mask.fit_transform(img)
	except:
		mask = path.abspath(path.expanduser(mask))
		mask = NiftiMasker(mask_img=mask)
		masked_data = mask.fit_transform(img).T
	ts_means = np.mean(masked_data, axis=0)
	ts_medians = np.mean(masked_data, axis=0)
	return ts_means, ts_medians

def nilearn_denoise(in_file, brain_mask, wm_mask, csf_mask,
                      motreg_file, outlier_file,
                      bandpass, tr ):
    """Clean time series using Nilearn high_variance_confounds to extract 
    CompCor regressors and NiftiMasker for regression of all nuissance regressors,
    detrending, normalziation and bandpass filtering.
    """
    import numpy as np
    import nibabel as nb
    import os
    from nilearn.image import high_variance_confounds
    from nilearn.input_data import NiftiMasker
    from nipype.utils.filemanip import split_filename

    # reload niftis to round affines so that nilearn doesn't complain
    wm_nii=nb.Nifti1Image(nb.load(wm_mask).get_data(), np.around(nb.load(wm_mask).get_affine(), 2), nb.load(wm_mask).get_header())
    csf_nii=nb.Nifti1Image(nb.load(csf_mask).get_data(), np.around(nb.load(csf_mask).get_affine(), 2), nb.load(csf_mask).get_header())
    time_nii=nb.Nifti1Image(nb.load(in_file).get_data(),np.around(nb.load(in_file).get_affine(), 2), nb.load(in_file).get_header())
        
    # infer shape of confound array
    # not ideal
    confound_len = nb.load(in_file).get_data().shape[3]
    
    # create outlier regressors
    outlier_regressor = np.empty((confound_len,1))
    try:
        outlier_val = np.genfromtxt(outlier_file)
    except IOError:
        outlier_val = np.empty((0))
    for index in np.atleast_1d(outlier_val):
        outlier_vector = np.zeros((confound_len, 1))
        outlier_vector[index] = 1
        outlier_regressor = np.hstack((outlier_regressor, outlier_vector))
    
    outlier_regressor = outlier_regressor[:,1::]
        
    # load motion regressors
    motion_regressor=np.genfromtxt(motreg_file)
    
    # extract high variance confounds in wm/csf masks from motion corrected data
    wm_regressor=high_variance_confounds(time_nii, mask_img=wm_nii, detrend=True)
    csf_regressor=high_variance_confounds(time_nii, mask_img=csf_nii, detrend=True)
    
    # create Nifti Masker for denoising
    denoiser=NiftiMasker(mask_img=brain_mask, standardize=True, detrend=True, high_pass=bandpass[1], low_pass=bandpass[0], t_r=tr)
    
    # denoise and return denoise data to img
    confounds=np.hstack((outlier_regressor,wm_regressor, csf_regressor, motion_regressor))
    denoised_data=denoiser.fit_transform(in_file, confounds=confounds)
    denoised_img=denoiser.inverse_transform(denoised_data)
        
    # save  
    _, base, _ = split_filename(in_file)
    img_fname = base + '_denoised.nii.gz'
    nb.save(denoised_img, img_fname)
    
    confound_fname = os.path.join(os.getcwd(), "all_confounds.txt")
    np.savetxt(confound_fname, confounds, fmt="%.10f")
    
    return os.path.abspath(img_fname), confound_fname

    def transform(self, imgs, confounds=None):
        """

        Parameters
        ----------
        imgs: list of Niimg-like objects
        """
        self._check_fitted()

        if self.smoothing_fwhm:
            imgs = smooth_img(imgs, self.smoothing_fwhm)

        imgs = [_utils.check_niimg_3d(img) for img in imgs]

        for i, roi in enumerate(self.mask_img_):
            masker = NiftiMasker(mask_img=roi)
            x = masker.fit_transform(imgs)
            if self.extract_funcs is not None:
                x = np.array([FDICT[f][0](x, **FDICT[f][1]) for f in self.extract_funcs])
            if i == 0:
                X = x
            else:
                X = np.concatenate((X, x), axis=0)

        return X.swapaxes(0, 1)

    def apply_mask(self, mask):
        """ Mask Brain_Data instance

        Args:
            mask: mask (Brain_Data or nifti object)
            
        """

        if isinstance(mask,Brain_Data):
            mask = mask.to_nifti() # convert to nibabel
        if not isinstance(mask, nib.Nifti1Image):
            if type(mask) is str:
                if os.path.isfile(mask):
                    mask = nib.load(mask)
               # Check if mask need to be resampled into Brain_Data mask space
                if not ((self.mask.get_affine()==mask.get_affine()).all()) & (self.mask.shape[0:3]==mask.shape[0:3]):
                    mask = resample_img(mask,target_affine=self.mask.get_affine(),target_shape=self.mask.shape)
            else:
                raise ValueError("Mask is not a nibabel instance, Brain_Data instance, or a valid file name.")

        masked = deepcopy(self)
        nifti_masker = NiftiMasker(mask_img=mask)
        masked.data = nifti_masker.fit_transform(self.to_nifti())
        if len(self.data.shape) > 2:
            masked.data = masked.data.squeeze()
        masked.nifti_masker = nifti_masker
        return masked

def extract_brain_rad(db, rad_column, rad_dir, stat, include_chim=False):
    """Replaces radiation presence by stat on whole brain ROI.

    Assumes brain mask and radiation nifti file is in rad_dir."""
    brain_mask_file = 'BrainMask_to_rd.nii.gz'
    extracted_rad_stat = {}  # Memoization of radiation statistic
    for idx, row in db.iterrows():
        if row[rad_column] == 1:
            sub_id = row['patient']
            if sub_id in extracted_rad_stat:
                db.loc[idx, rad_column] = extracted_rad_stat[sub_id]
            else:
                mask_path = os.path.join(rad_dir, sub_id, brain_mask_file)
                mask_check = os.path.isfile(mask_path)
                rad_path = os.path.join(rad_dir, sub_id, sub_id + '.nii')
                rad_check = os.path.isfile(rad_path)
                if mask_check and rad_check:
                    masker = NiftiMasker(mask_path)
                    rad_stat = stat(masker.fit_transform(rad_path))
                    extracted_rad_stat[sub_id] = rad_stat
                    db.loc[idx, rad_column] = rad_stat
                else:
                    db.loc[idx, rad_column] = None
        elif not include_chim:
            db.loc[idx, rad_column] = None

    db = db[db[rad_column].notnull()]
    return db

	def _run_interface(self, runtime):
		from nilearn.input_data import NiftiMasker, NiftiLabelsMasker
		from nipype.utils.filemanip import split_filename
		import nibabel as nib
		import os

		functional_filename = self.inputs.in_file
		atlas_filename = self.inputs.atlas_filename
		mask_filename = self.inputs.mask_filename

		# Extracting the ROI signals
		masker = NiftiLabelsMasker(labels_img=atlas_filename,
                           background_label = 0,
                           standardize=True,
                           detrend = True,
                           verbose = 1
                           )
		time_series = masker.fit_transform(functional_filename)

		# Removing the ROI signal from the time series
		nifti_masker = NiftiMasker(mask_img=mask_filename)
		masked_data = nifti_masker.fit_transform(functional_filename, confounds=time_series[...,0])
		masked_img = nifti_masker.inverse_transform(masked_data)

		# Saving the result to disk
		outputs = self._outputs().get()
		fname = self.inputs.in_file
		_, base, _ = split_filename(fname)
		nib.save(masked_img, os.path.abspath(base + '_regressed.nii.gz'))
		return runtime

def significant_signal(data_path,
	substitution={},
	mask_path='',
	exclude_ones=False,
	):
	"""Return the mean and median inverse logarithm of a p-value map.

	Parameters
	----------

	data_path : str
		Path to a p-value map in NIfTI format.
	mask_path : str
		Path to a region of interest map in NIfTI format.
		THIS IS ALMOST ALWAYS REQUIRED, as NIfTI statistic images populate the whole 3D circumscribed space around your structure of interest,
		and commonly assign null values to the background.
		In an inverse logarithm computation, null corresponds to infinity, which can considerably bias the evaluation.
	substitution : dict
		Dictionary whose keys are format identifiers present in `data_path` and whose values are strings.

	Returns
	-------

	mean : float
	median : float
	"""

	if substitution:
		data_path = data_path.format(**substitution)
	data_path = path.abspath(path.expanduser(data_path))
	try:
		img = nib.load(data_path)
	except FileNotFoundError:
		return float('NaN'), float('NaN')
	if mask_path:
		mask_path = path.abspath(path.expanduser(mask_path))
		masker = NiftiMasker(mask_img=mask_path)
		masked_data = masker.fit_transform(img).T
		data = masked_data[~np.isnan(masked_data)]
	else:
		data = img.get_data()
		data = data[~np.isnan(data)]
	# We interpret zero as the lowest p-value, and conservatively estimate it to be equal to just under half of the smallest value in the defined range
	nonzero = data[np.nonzero(data)]
	data_min = np.min(nonzero)
	data_min = data_min*0.49
	data[data == 0] = data_min
	if exclude_ones:
		data = data[data!=1]
	data = -np.log10(data)
	# We use np.ma.median() because life is complicated:
	# https://github.com/numpy/numpy/issues/7330
	median = np.ma.median(data, axis=None)
	mean = np.mean(data)

	return mean, median

    def similarity(self, image, method='correlation'):
        """ Calculate similarity of Brain_Data() instance with single Brain_Data or Nibabel image

            Args:
                self: Brain_Data instance of data to be applied
                image: Brain_Data or Nibabel instance of weight map

            Returns:
                pexp: Outputs a vector of pattern expression values

        """

        if not isinstance(image, Brain_Data):
            if isinstance(image, nib.Nifti1Image):
                image = Brain_Data(image)
            else:
                raise ValueError("Image is not a Brain_Data or nibabel instance")
        dim = image.shape()

        # Check to make sure masks are the same for each dataset and if not create a union mask
        # This might be handy code for a new Brain_Data method
        if np.sum(self.nifti_masker.mask_img.get_data()==1)!=np.sum(image.nifti_masker.mask_img.get_data()==1):
            new_mask = intersect_masks([self.nifti_masker.mask_img, image.nifti_masker.mask_img], threshold=1, connected=False)
            new_nifti_masker = NiftiMasker(mask_img=new_mask)
            data2 = new_nifti_masker.fit_transform(self.to_nifti())
            image2 = new_nifti_masker.fit_transform(image.to_nifti())
        else:
            data2 = self.data
            image2 = image.data


        # Calculate pattern expression
        if method is 'dot_product':
            if len(image2.shape) > 1:
                if image2.shape[0]>1:
                    pexp = []
                    for i in range(image2.shape[0]):
                        pexp.append(np.dot(data2, image2[i,:]))
                    pexp = np.array(pexp)
                else:
                    pexp = np.dot(data2, image2)
            else:
                pexp = np.dot(data2, image2)
        elif method is 'correlation':
            if len(image2.shape) > 1:
                if image2.shape[0]>1:
                    pexp = []
                    for i in range(image2.shape[0]):
                        pexp.append(pearson(image2[i,:], data2))
                    pexp = np.array(pexp)
                else:
                    pexp = pearson(image2, data2)
            else:
                pexp = pearson(image2, data2)
        return pexp

def _vectorize_nii(in_data_file, mask_file, parcellation_path, fwhm):
    from nilearn.input_data import NiftiMasker, NiftiLabelsMasker
    import nibabel as nib

    if parcellation_path is None:
        masker = NiftiMasker(mask_img=mask_file, smoothing_fwhm=fwhm)
    else:
        masker = NiftiLabelsMasker(labels_img=parcellation_path, smoothing_fwhm=fwhm)

    vectorized_data = masker.fit_transform(in_data_file)
    return vectorized_data, masker

def map_threshold(stat_img, mask_img, threshold, height_control='fpr',
                  cluster_threshold=0):
    """ Threshold the provvided map

    Parameters
    ----------
    stat_img : Niimg-like object,
       statistical image (presumably in z scale)

    mask_img : Niimg-like object,
        mask image

    threshold: float,
        cluster forming threshold (either a p-value or z-scale value)

    height_control: string
        false positive control meaning of cluster forming
        threshold: 'fpr'|'fdr'|'bonferroni'|'none'

    cluster_threshold : float, optional
        cluster size threshold

    Returns
    -------
    thresholded_map : Nifti1Image,
        the stat_map theresholded at the prescribed voxel- and cluster-level
    """
    # Masking
    masker = NiftiMasker(mask_img=mask_img)
    stats = np.ravel(masker.fit_transform(stat_img))
    n_voxels = np.size(stats)

    # Thresholding
    if height_control == 'fpr':
        z_th = norm.isf(threshold)
    elif height_control == 'fdr':
        z_th = fdr_threshold(stats, threshold)
    elif height_control == 'bonferroni':
        z_th = norm.isf(threshold / n_voxels)
    else:  # Brute-force thresholding
        z_th = threshold
    stats *= (stats > z_th)

    stat_map = masker.inverse_transform(stats).get_data()

    # Extract connected components above threshold
    label_map, n_labels = label(stat_map > z_th)
    labels = label_map[(masker.mask_img_.get_data() > 0)]
    for label_ in range(1, n_labels + 1):
        if np.sum(labels == label_) < cluster_threshold:
            stats[labels == label_] = 0

    return masker.inverse_transform(stats)

def make_ttest(reg1, reg2):
    masker = NiftiMasker(nib.load(MASK_FILE), standardize=False)
    masker.fit()

    subjects = [1, 2, 3, 5, 6, 7, 8, 9, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20]

    a = np.arctanh(join_all_subjects(reg1, subjects, masker))
    b = np.arctanh(join_all_subjects(reg2, subjects, masker))
    t, prob = ttest_rel(a, b)

    tt = masker.inverse_transform(t)
    pp = masker.inverse_transform(prob)
    return tt, pp

def load_data():
    with open(expanduser('~/data/HCP_unmasked/data.json'), 'r') as f:
        data = json.load(f)
        for this_data in data:
            this_data['array'] += '.npy'
        mask_img = expanduser('~/data/HCP_mask/mask_img.nii.gz')
    masker = NiftiMasker(mask_img=mask_img, smoothing_fwhm=4,
                         standardize=True)
    masker.fit()
    smith2009 = fetch_atlas_smith_2009()
    init = smith2009.rsn70
    dict_init = masker.transform(init)
    return masker, dict_init, sorted(data, key=lambda t: t['filename'])

def preprocess_varpar(num, subj, subj_dir, **kwargs):
    from nistats.design_matrix import make_design_matrix
    from nistats.first_level_model import run_glm
    bold_path = 'BOLD/task001_run00%i/bold_dico_bold7Tp1_to_subjbold7Tp1.nii.gz' % (num+1)
    bold_path = os.path.join(DATA_DIR, subj, bold_path)
    mask = os.path.join(DATA_DIR, subj, 'templates', 'bold7Tp1', 'brain_mask.nii.gz')
    bold = load(bold_path)
    masker = NiftiMasker(mask)
    data = masker.fit_transform(bold)
    dmat = make_design_matrix(np.arange(data.shape[0])*TR, hrf_model='fir', drift_order=5,
                              **kwargs)
    labels, results = run_glm(data, dmat, noise_model='ols', verbose=1)
    img = masker.inverse_transform(StandardScaler().fit_transform(results[0.0].resid))
#    return StandardScaler().fit_transform(results[0.0].resid)
    save(img, os.path.join(subj_dir, 'run00%i.nii.gz' % num))

class SmoothResampleMasker(BaseMasker):

    def __init__(self, mask_img=None, smoothing_fwhm=None, resampling=None, searchlight=False):

        self.mask_img = mask_img
        self.smoothing_fwhm = smoothing_fwhm
        self.resampling = resampling
        self.searchlight = searchlight

        self.masker = None

    def fit(self):

        if self.resampling is not None:
            self.mask_img = resample_img(self.mask_img, target_affine=np.diag(self.resampling * np.ones(3)))
        self.masker = NiftiMasker(mask_img=self.mask_img)
        self.masker.fit()

        return self

    def transform(self, imgs, confounds=None):

        smooth_prefix = '' if self.smoothing_fwhm is None else 's%g' % self.smoothing_fwhm
        resample_prefix = '' if self.smoothing_fwhm is None else 'r%g' % self.smoothing_fwhm

        if not isinstance(imgs, list):
            imgs = [imgs]

        path_first = imgs[0] if isinstance(imgs[0], str) else imgs[0].get_filename()

        path_first_resampled = os.path.join(os.path.dirname(path_first), resample_prefix + os.path.basename(path_first))
        path_first_smoothed = os.path.join(os.path.dirname(path_first), smooth_prefix + resample_prefix + os.path.basename(path_first))

        if self.resampling is not None and self.smoothing_fwhm is not None:
            if self.resampling is not None:
                if not os.path.exists(path_first_resampled) and not os.path.exists(path_first_smoothed):
                    imgs = resample_img(imgs, target_affine=np.diag(self.resampling * np.ones(3)))
                else:
                    imgs = []
            if self.smoothing_fwhm is not None:
                if not os.path.exists(path_first_smoothed):
                    imgs = smooth_img(imgs, self.smoothing_fwhm)
                else:
                    imgs = []
        else:
            imgs = [check_niimg_3d(img) for img in imgs]

        return self.masker.transform(imgs)

class signal_extractor():

    def __init__(self, dataset = None):
        self.dataset = dataset
        if dataset.has_key('mask'):
            self.masker = NiftiMasker(mask_img = self.dataset.mask,
                                low_pass = .1,
                                high_pass = .01,
                                smoothing_fwhm =6.,
                                t_r = 1.05,
                                detrend = True,
                                standardize = False,
                                memory_level = 0,
                                verbose=5)
        else:
            self.masker = NiftiMasker(
                                low_pass = .1,
                                high_pass = .01,
                                smoothing_fwhm =6.,
                                t_r = 1.05,
                                detrend = True,
                                standardize = False,
                                memory_level = 0,
                                verbose=5)
    def extract(self):
        for idx, func in enumerate([self.dataset.func1]):
            #add mask, smoothing, filter and detrending


            for i in range(len(self.dataset.subjects)):
                tic = time.clock()
                #extract signal to x
                x = self.masker.fit_transform(func[i])
                print "loading time : "+ str(time.clock() - tic)
                yield x, self.masker

def loader(anat, downsample, target_affine, dataroot, subject, maskpath, nrun,
           niifilename, labels, **kwargs):
    ''' 
    All parameters are submitted as cfg dictionary.
    Given parameters in cfg, return masked and concatenated over runs data 
    
    Input
    anat: MNI template
    downsample: 1 or 0
    target_affine: downsampling matrix
    dataroot: element of path to data
    subject: folder in dataroot with subject data
    maskpath: path to mask
    nrun: number of runs
    niifilename: how is the data file called
    labels: labels from load_labels function
    
    Output
    dict(nii_func=nii_func,nii_mean=nii_mean,masker=masker,nii_mask=nii_mask)
    nii_func: 4D data
    nii_mean: mean over 4th dimension
    masker: masker object from nibabel
    nii_mask: 3D mask
    '''
    nii_func = list()
    for r in range(nrun):
        fname = '{0}/{1}/run{2}/{3}'.format(dataroot, subject, r+1, niifilename) # Assumption about file location
        nii_img = load(fname, mmap=False)
        nii_img.set_sform(anat.get_sform())
        # Get mean over 4D
        nii_mean = mean_img(nii_img)
        # Masking
        nii_mask = load(maskpath)
        nii_mask.set_sform(anat.get_sform())
        # Binarize the mask
        nii_mask = check_binary(nii_mask)
        if downsample:
            nii_img = resample_img(nii_img, target_affine=target_affine)
            nii_mask = resample_img(nii_mask, target_affine=target_affine, interpolation='nearest')
        masker = NiftiMasker(nii_mask, standardize=True)
        nii_img = masker.fit_transform(nii_img)
        # Drop zero timepoints, zscore
        nii_img = drop_labels(nii_img, labels.get('to_drop_zeros')[r])
        nii_func.append(stats.zscore(nii_img, axis=0)) # zscore over time
    # throw data together
    nii_func = np.concatenate(nii_func)
    return dict(nii_func=nii_func, nii_mean=nii_mean, masker=masker, nii_mask=nii_mask)

    def __init__(self, data=None, Y=None, X=None, mask=None, output_file=None, **kwargs):
        if mask is not None:
            if not isinstance(mask, nib.Nifti1Image):
                if type(mask) is str:
                    if os.path.isfile(mask):
                        mask = nib.load(mask)
            else:
                raise ValueError("mask is not a nibabel instance")
            self.mask = mask
        else:
            self.mask = nib.load(os.path.join(get_resource_path(),'MNI152_T1_2mm_brain_mask.nii.gz'))
        self.nifti_masker = NiftiMasker(mask_img=self.mask)

        if data is not None:
            if type(data) is str:
                data=nib.load(data)
                self.data = self.nifti_masker.fit_transform(data)
            elif type(data) is list:
                # Load and transform each image in list separately (nib.concat_images(data) can't handle images of different sizes)
                self.data = []
                for i in data:
                    if isinstance(i,six.string_types):
                        self.data.append(self.nifti_masker.fit_transform(nib.load(i)))
                    elif isinstance(i,nib.Nifti1Image):
                        self.data.append(self.nifti_masker.fit_transform(i))
                self.data = np.array(self.data)
            elif not isinstance(data, nib.Nifti1Image):
                raise ValueError("data is not a nibabel instance")

            # Collapse any extra dimension
            if any([x==1 for x in self.data.shape]):
                self.data=self.data.squeeze()
        else:
            self.data = np.array([])

        if Y is not None:
            if type(Y) is str:
                if os.path.isfile(Y):
                    Y=pd.read_csv(Y,header=None,index_col=None)
            if isinstance(Y, pd.DataFrame):
                if self.data.shape[0]!= len(Y):
                    raise ValueError("Y does not match the correct size of data")
                self.Y = Y
            else:
                raise ValueError("Make sure Y is a pandas data frame.")
        else:
            self.Y = pd.DataFrame()

        if X is not None:
            if self.data.shape[0]!= X.shape[0]:
                raise ValueError("X does not match the correct size of data")
            self.X = X
        else:
            self.X = pd.DataFrame()

        if output_file is not None:
            self.file_name = output_file
        else:
            self.file_name = []