Copying headers from input images with math_img

This example shows how to copy the header information from one of the input images to the result image when using the function math_img.

The header information contains metadata about the image, such as the dimensions, the voxel sizes, the affine matrix, repetition time (TR), etc. Some of this information might be important for downstream analyses depending on the software used.

Let’s fetch an example fMRI dataset

from nilearn import datasets

dataset = datasets.fetch_adhd(n_subjects=2)

[get_dataset_dir] Dataset found in /home/runner/nilearn_data/adhd
[fetch_single_file] Downloading data from
https://www.nitrc.org/frs/download.php/7783/adhd40_0010064.tgz ...
[_chunk_report_] Downloaded 17940480 of 45583539 bytes (39.4%%,    1.5s
remaining)
[fetch_single_file]  ...done. (3 seconds, 0 min)

[uncompress_file] Extracting data from
/home/runner/nilearn_data/adhd/31769c9cee5cd55f045e62633d651f3d/adhd40_0010064.t
gz...
[uncompress_file] .. done.

Now let’s look at the header of one of these images

from nilearn.image import load_img

subj1_img = load_img(dataset.func[0])
subj2_img = load_img(dataset.func[1])

print(f"Subject 1 header:\n{subj1_img.header}")

Subject 1 header:
<class 'nibabel.nifti1.Nifti1Header'> object, endian='<'
sizeof_hdr      : 348
data_type       : b''
db_name         : b''
extents         : 0
session_error   : 0
regular         : b'r'
dim_info        : 48
dim             : [  4  61  73  61 176   1   1   1]
intent_p1       : 0.0
intent_p2       : 0.0
intent_p3       : 0.0
intent_code     : none
datatype        : float32
bitpix          : 32
slice_start     : 0
pixdim          : [-1.  3.  3.  3.  2.  0.  0.  0.]
vox_offset      : 0.0
scl_slope       : nan
scl_inter       : nan
slice_end       : 60
slice_code      : unknown
xyzt_units      : 10
cal_max         : 0.0
cal_min         : 0.0
slice_duration  : 0.0
toffset         : 0.0
glmax           : 0
glmin           : 0
descrip         : b''
aux_file        : b''
qform_code      : scanner
sform_code      : scanner
quatern_b       : -0.0
quatern_c       : 1.0
quatern_d       : 0.0
qoffset_x       : 90.0
qoffset_y       : -126.0
qoffset_z       : -72.0
srow_x          : [-3. -0. -0. 90.]
srow_y          : [  -0.    3.   -0. -126.]
srow_z          : [  0.   0.   3. -72.]
intent_name     : b''
magic           : b'n+1'

Let’s apply a simple operation using math_img

from nilearn.image import math_img

result_img = math_img("img1 * 1", img1=subj1_img)

By default, math_img simply resets result image’s header to the default Nifti1Header.

This means that it will contain different information as compared to the input image.

We can check that as follows:

print("Following header fields do not match:")
for key in result_img.header:
    if not (subj1_img.header[key] == result_img.header[key]).all():
        print(
            f"For '{key}'\n",
            "\tinput image:",
            subj1_img.header[key],
            "\n\tresult image:",
            result_img.header[key],
        )

Following header fields do not match:
For 'regular'
        input image: b'r'
        result image: b''
For 'dim_info'
        input image: 48
        result image: 0
For 'pixdim'
        input image: [-1.  3.  3.  3.  2.  0.  0.  0.]
        result image: [-1.  3.  3.  3.  1.  1.  1.  1.]
For 'scl_slope'
        input image: nan
        result image: nan
For 'scl_inter'
        input image: nan
        result image: nan
For 'slice_end'
        input image: 60
        result image: 0
For 'xyzt_units'
        input image: 10
        result image: 0
For 'qform_code'
        input image: 1
        result image: 0
For 'sform_code'
        input image: 1
        result image: 2

This could affect some downstream analyses.

For example, here the TR (given as fifth element in pixdim) is changed from 2 in subj1_img to 1 in result_img.

To fix this, we can copy the header of the input images to the result image, like this:

result_img_with_header = math_img(
    "img1 * 1", img1=subj1_img, copy_header_from="img1"
)

Let’s compare the header fields again.

print("Following header fields do not match:")
for key in result_img_with_header.header:
    if not (subj1_img.header[key] == result_img_with_header.header[key]).all():
        print(
            f"For '{key}'\n",
            "\tinput image:",
            subj1_img.header[key],
            "\n\tresult image:",
            result_img_with_header.header[key],
        )

Following header fields do not match:
For 'scl_slope'
        input image: nan
        result image: nan
For 'scl_inter'
        input image: nan
        result image: nan
For 'cal_max'
        input image: 0.0
        result image: 24359.932

We can safely ignore the fields that are still different – scl_scope and scl_inter are just nan and cal_max is supposed to have the maximum data value that is updated automatically by nilearn.

Modifying dimensions in the formula

Now let’s say we have a formula that changes the dimensions of the input images. For example, by taking the mean of the images along the time axis.

Copying the header with the copy_header_from parameter will not work in this case.

So, in such cases we could just use math_img without specifying copy_header_from and then explicitly copy the header from one of the images using new_img_like

result_img = math_img(
    "np.mean(img1, axis=-1) - np.mean(img2, axis=-1)",
    img1=subj1_img,
    img2=subj2_img,
)

Several of the header fields are different:

print("Following header fields do not match:")
for key in result_img.header:
    if not (subj1_img.header[key] == result_img.header[key]).all():
        print(
            f"For '{key}'\n",
            "\tinput image:",
            subj1_img.header[key],
            "\n\tresult image:",
            result_img.header[key],
        )

Following header fields do not match:
For 'regular'
        input image: b'r'
        result image: b''
For 'dim_info'
        input image: 48
        result image: 0
For 'dim'
        input image: [  4  61  73  61 176   1   1   1]
        result image: [ 3 61 73 61  1  1  1  1]
For 'pixdim'
        input image: [-1.  3.  3.  3.  2.  0.  0.  0.]
        result image: [-1.  3.  3.  3.  1.  1.  1.  1.]
For 'scl_slope'
        input image: nan
        result image: nan
For 'scl_inter'
        input image: nan
        result image: nan
For 'slice_end'
        input image: 60
        result image: 0
For 'xyzt_units'
        input image: 10
        result image: 0
For 'qform_code'
        input image: 1
        result image: 0
For 'sform_code'
        input image: 1
        result image: 2

Now we can copy the header explicitly like this:

from nilearn.image import new_img_like

result_img_with_header = new_img_like(
    ref_niimg=subj1_img,
    data=result_img.get_fdata(),
    affine=result_img.affine,
    copy_header=True,
)

Now, only a few not-so-important fields are different.

The modified fields can vary depending upon the formula passed into the function.

In this case, dim and pixdim are different because we took a mean over the time dimension.

And again, cal_min and cal_max are set to minimum and maximum data values respectively, by Nilearn.

print("Following header fields do not match:")
for key in result_img_with_header.header:
    if not (subj1_img.header[key] == result_img_with_header.header[key]).all():
        print(
            f"For '{key}'\n",
            "\tinput image:",
            subj1_img.header[key],
            "\n\tresult image:",
            result_img_with_header.header[key],
        )

Following header fields do not match:
For 'dim'
        input image: [  4  61  73  61 176   1   1   1]
        result image: [ 3 61 73 61  1  1  1  1]
For 'pixdim'
        input image: [-1.  3.  3.  3.  2.  0.  0.  0.]
        result image: [-1.  3.  3.  3.  1.  1.  1.  1.]
For 'scl_slope'
        input image: nan
        result image: nan
For 'scl_inter'
        input image: nan
        result image: nan
For 'cal_max'
        input image: 0.0
        result image: 17736.469
For 'cal_min'
        input image: 0.0
        result image: -15583.777

Estimated memory usage: 837 MB