Python nibabel.Nifti2Image() Examples

def test_check_bs_len(fnames, lengths, expected_out, tmp_path):
    affine = np.eye(4)
    fpaths = []
    for fname, length in zip(fnames, lengths):
        fpath = tmp_path / fname
        nib.Nifti2Image(np.zeros((1, 1, 1, length)), affine=affine).to_filename(str(fpath))
        fpaths.append(str(fpath))

    if all(i < 3 for i in lengths):
        with pytest.raises(RuntimeError) as rterr:
            _check_bs_len(fpaths)
        assert "None of the beta series" in str(rterr.value)
    else:
        assert [op.basename(f) for f in _check_bs_len(fpaths)] == expected_out 

def _calc_beta_map(model, trial_type, hrf_model, tstat):
    """
    Calculates the beta estimates for every voxel from
    a nistats model

    Parameters
    ----------
    model : nistats.first_level_model.FirstLevelModel
        a fit model of the first level results
    trial_type : str
        the trial to create the beta estimate
    hrf_model : str
        the hemondynamic response function used to fit the model
    tstat : bool
        return the t-statistic for the betas instead of the raw estimates

    Returns
    -------
    beta_map : nibabel.nifti2.Nifti2Image
        nifti image containing voxelwise beta estimates
    """
    import numpy as np

    # make it so we do not divide by zero
    TINY = 1e-50
    raw_beta_map = _estimate_map(model, trial_type, hrf_model, 'effect_size')
    if tstat:
        var_map = _estimate_map(model, trial_type, hrf_model, 'effect_variance')
        tstat_array = raw_beta_map.get_fdata() / np.sqrt(np.maximum(var_map.get_fdata(), TINY))
        return nib.Nifti2Image(tstat_array, raw_beta_map.affine, raw_beta_map.header)
    else:
        return raw_beta_map 

def decouple_volumes(v1, v2, mode, se=None, iterations=1):
    """
    
    mode : {inner-from-outer, outer-from-inner, neighbors}
        inner-from-outer: this changes v1 by removing voxels
        outer-from-inner: this changes v2 by adding voxels
        neighbors: this changes v2 by removing voxels
    
    """
    assert mode in ["inner-from-outer","outer-from-inner","neighbors"]
    
    if isinstance(v1, str) and os.path.isfile(v1):
        v1 = nib.load(v1)
    assert isinstance(v1, nib.Nifti1Image) or isinstance(v1, nib.Nifti2Image)
    d1 = v1.get_data()
    if isinstance(v2, str) and os.path.isfile(v2):
        v2 = nib.load(v2)
    assert isinstance(v2, nib.Nifti1Image) or isinstance(v2, nib.Nifti2Image)
    d2 = v2.get_data()
    
    assert d1.ndim is d2.ndim
    
    
    if se is None:
        se = mrph.generate_binary_structure(d1.ndim,1)
    
    if mode == "inner-from-outer":
        # make v2/d2 the inner volume
        d1, d2 = d2, d1
        v1, v2 = v2, v1        
        d2 = d2 & mrph.binary_erosion(d1, se, iterations)
        
    if mode == "outer-from-inner":
        d2 = d2 | mrph.binary_dilation(d1, se, iterations)
        
    if mode == "neighbors":
        d2 = d2 & ~mrph.binary_dilation(d1, se, iterations)
    
    d2 = nib.Nifti1Image(d2, v2.affine, header=v2.header)
    d2.set_filename(v2.get_filename())
    return d2 

def nifti2_from_filename(self):
        """
        Use nibabel.Nifti2Image.from_filename to load nifti2_image.nii.
        """
        file_name = os.path.join(self.data_dir, "nifti2_image")
        nib.Nifti2Image.from_filename(file_name) 

def nifti2_to_filename(self):
        """
        Use nibabel.Nifti2Image.to_filename to save out_nifti2_image.nii.
        """
        file_name = os.path.join(self.data_dir, "out_nifti2_image")
        img = nib.Nifti2Image(self.get_data(), self.get_affine())
        img.to_filename(file_name)
        os.remove(file_name + ".nii") 

def create_sample_data(self):
        """
        Create sample data files. The files created are:

        * analyze_image.hdr + .img: plain ANALYZE image
        * mgh_image.mgh: MGH image
        * nifti1_image.nii: NIfTI1 image
        * nifti2_image.nii:  NIfTI2 image
        * spm2_image.hdr + .img + .mat: SPM2 ANALYZE image
        * spm99_image.hdr + .img + .mat: SPM99 ANALYZE image
        """
        file_name = os.path.join(self.data_dir, "analyze_image")
        analyze_img = nib.AnalyzeImage(self.get_data(), np.eye(4))
        analyze_img.to_filename(file_name)
        file_name = os.path.join(self.data_dir, "mgh_image")
        mgh_img = nib.freesurfer.mghformat.MGHImage(self.get_data(),
                                                    np.eye(4))
        mgh_img.to_filename(file_name)
        file_name = os.path.join(self.data_dir, "nifti1_image")
        nifti1_img = nib.Nifti1Image(self.get_data(), self.get_affine())
        nifti1_img.to_filename(file_name)
        file_name = os.path.join(self.data_dir, "nifti2_image")
        nifti2_img = nib.Nifti2Image(self.get_data(), self.get_affine())
        nifti2_img.to_filename(file_name)
        file_name = os.path.join(self.data_dir, "spm2_image")
        spm2_img = nib.spm2analyze.Spm2AnalyzeImage(self.get_data(),
                                                    np.eye(4))
        spm2_img.to_filename(file_name)
        file_name = os.path.join(self.data_dir, "spm99_image")
        spm99_img = nib.spm99analyze.Spm99AnalyzeImage(self.get_data(),
                                                       np.eye(4))
        spm99_img.to_filename(file_name) 

def load_nifti(nifti_volume):
    """Load and return a nifti file.

    This will apply path resolution if a filename without extension is given. See the function
    :func:`nifti_filepath_resolution` for details.

    Args:
        nifti_volume (string): The filename of the volume to use.

    Returns:
        :class:`nibabel.nifti2.Nifti2Image`
    """
    path = nifti_filepath_resolution(nifti_volume)
    return nifti_info_decorate_nibabel_image(nib.load(path)) 

def write_nifti(data, output_fname, header=None, affine=None, use_data_dtype=True, **kwargs):
    """Write data to a nifti file.

    This will write the output directory if it does not exist yet.

    Args:
        data (ndarray): the data to write to that nifti file
        output_fname (str): the name of the resulting nifti file, this function will append .nii.gz if no
            suitable extension is given.
        header (nibabel header): the nibabel header to use as header for the nifti file. If None we will use
            a default header.
        affine (ndarray): the affine transformation matrix
        use_data_dtype (boolean): if we want to use the dtype from the data instead of that from the header
            when saving the nifti.
        **kwargs: other arguments to Nifti2Image from NiBabel
    """
    if header is None:
        header = nib.nifti2.Nifti2Header()

    if use_data_dtype:
        header = copy.deepcopy(header)
        dtype = data.dtype
        if data.dtype == np.bool:
            dtype = np.char
        try:
            header.set_data_dtype(dtype)
        except nib.spatialimages.HeaderDataError:
            pass

    if not (output_fname.endswith('.nii.gz') or output_fname.endswith('.nii')):
        output_fname += '.nii.gz'

    if not os.path.exists(os.path.dirname(output_fname)):
        os.makedirs(os.path.dirname(output_fname))

    if isinstance(header, nib.nifti2.Nifti2Header):
        format = nib.Nifti2Image
    else:
        format = nib.Nifti1Image

    format(data, affine, header=header, **kwargs).to_filename(output_fname) 

def _estimate_map(model, trial_type, hrf_model, output_type):
    """
    Calculates model output for every voxel from
    a nistats model

    Parameters
    ----------
    model : nistats.first_level_model.FirstLevelModel
        a fit model of the first level results
    trial_type : str
        the trial to create the beta estimate
    hrf_model : str
        the hemondynamic response function used to fit the model
    output_type : str
        Type of the output map.
        Can be ‘z_score’, ‘stat’, ‘p_value’, ‘effect_size’, or ‘effect_variance’

    Returns
    -------
    map_img : nibabel.nifti2.Nifti2Image
        nifti image containing voxelwise output_type estimates
    """
    import numpy as np

    # calculate the beta map
    map_list = []
    map_base = model.compute_contrast(trial_type, output_type=output_type)
    map_list.append(map_base.get_fdata())
    sign = np.where(map_list[0] < 0, -1, 1)
    if 'derivative' in hrf_model:
        td_contrast = '_'.join([trial_type, 'derivative'])
        map_list.append(
            model.compute_contrast(
                td_contrast, output_type=output_type).get_fdata())
    if 'dispersion' in hrf_model:
        dd_contrast = '_'.join([trial_type, 'dispersion'])
        map_list.append(
            model.compute_contrast(
                dd_contrast, output_type=output_type).get_fdata())

    if len(map_list) == 1:
        map_img = map_base
    else:
        map_array = sign * \
            np.sqrt(
                np.sum(
                    np.array([np.power(c, 2) for c in map_list]), axis=0))

        map_img = nib.Nifti2Image(
            map_array,
            map_base.affine,
            map_base.header)

    return map_img 

def binarize(vols, return_empty=False):
    """Binarize a list of input volumes by finding the maximum posterior 
    probability of each and assigning the voxel to this volume.
    
    PARAMETERS
    ----------
    vols : list
        List of filenames or nibabel image objects (describing probabilities of
        different tissue types).
    return_empty : bool
        If true, return an array containing all voxels which are not assigned
        to either of the other volumes (default: False)

    RETURNS
    ----------
    bin_vols : list
        List of ndarrays describing binarized versions of the input volumes.
    unassign : ndarray
        Array containing any unassigned voxels.
    """    
    # if filenames are provided, load data
    volsi = [None]*len(vols)
    for i in range(len(vols)):
        if isinstance(vols[i], str) and os.path.isfile(vols[i]):
            volsi[i] = nib.load(vols[i]).get_data()
        elif type(vols[i]) in [nib.Nifti1Image, nib.Nifti2Image]:
            volsi[i] = vols[i].get_data()
        else:
            # assume numpy array
            volsi[i] = vols[i]

    # Concatenate arrays/images
    imgs = np.concatenate(tuple([v[...,np.newaxis] for v in volsi]), axis=3)
    imgs = np.concatenate((np.zeros_like(volsi[0])[...,np.newaxis], imgs),
                          axis=3)
    
    # Find max indices
    max_idx = np.argmax(imgs, axis=3)
    
    # Binarize. Here vols_bin[0] contain voxels not assigned to any other
    # volume
    vols_bin=[]
    for i in range(imgs.shape[-1]):
        vols_bin.append(max_idx == i)
    #vols_bin = [vb.astype(np.uint8) for vb in vols_bin]

    if return_empty:
        return vols_bin[1:]+[vols_bin[0]]
    else:
        return vols_bin[1:]