#!/usr/bin/env fslpython
# test fslstats with a bunch of different options
#
# Usage: fslstats [preoptions] <input> [options]
#
# preoption -t will give a separate output line for each 3D volume of a 4D timeseries
# preoption -K < indexMask > will generate seperate n submasks from indexMask, for indexvalues 1..n where n is the maximum index value in indexMask, and generate statistics for each submask
# Note - options are applied in order, e.g. -M -l 10 -M will report the non-zero mean, apply a threshold and then report the new nonzero mean
# -l <lthresh> : set lower threshold
# -u <uthresh> : set upper threshold
# -r : output <robust min intensity> <robust max intensity>
# -R : output <min intensity> <max intensity>
# -e : output mean entropy ; mean(-i*ln(i))
# -E : output mean entropy (of nonzero voxels)
# -v : output <voxels> <volume>
# -V : output <voxels> <volume> (for nonzero voxels)
# -m : output mean
# -M : output mean (for nonzero voxels)
# -s : output standard deviation
# -S : output standard deviation (for nonzero voxels)
# -w : output smallest ROI <xmin> <xsize> <ymin> <ysize> <zmin> <zsize> <tmin> <tsize> containing nonzero voxels
# -x : output co-ordinates of maximum voxel
# -X : output co-ordinates of minimum voxel
# -c : output centre-of-gravity (cog) in mm coordinates
# -C : output centre-of-gravity (cog) in voxel coordinates
# -p <n> : output nth percentile (n between 0 and 100)
# -P <n> : output nth percentile (for nonzero voxels)
# -a : use absolute values of all image intensities
# -n : treat NaN or Inf as zero for subsequent stats
# -k <mask> : use the specified image (filename) for masking - overrides lower and upper thresholds
# -d <image> : take the difference between the base image and the image specified here
# -h <nbins> : output a histogram (for the thresholded/masked voxels only) with nbins
# -H <nbins> <min> <max> : output a histogram (for the thresholded/masked voxels only) with nbins and histogram limits of min and max
# Note - thresholds are not inclusive ie lthresh<allowed<uthresh
import sys
import os
import fsl.utils.run as run
import numpy as np
import nibabel as nib
# set the random seed so that the tests are reproducible
np.random.seed(0)
options = [
{"option": "-r", "expected": "0.018533 0.978824"},
{"option": "-R", "expected": "0.000546 0.999809"},
{"option": "-e", "expected": "0.906103"},
{"option": "-E", "expected": "0.906103"},
{"option": "-v", "expected": "1000 1000.000000"},
{"option": "-V", "expected": "1000 1000.000000"},
{"option": "-m", "expected": "0.495922"},
{"option": "-M", "expected": "0.495922"},
{"option": "-s", "expected": "0.290744"},
{"option": "-S", "expected": "0.290744"},
{"option": "-w", "expected": "0 10 0 10 0 10 0 1"},
{"option": "-x", "expected": "9 7 4"},
{"option": "-X", "expected": "8 2 1"},
{"option": "-c", "expected": "4.498706 4.545873 4.613188"},
{"option": "-C", "expected": "4.498706 4.545873 4.613188"},
{"option": "-p 50", "expected": "0.482584"},
{"option": "-P 50", "expected": "0.482584"},
]
# create a list of tests and expected results
tests = []
for option in options:
tests.append(
{
"preoptions": "",
"mask": "",
"options": option["option"],
"expected": option["expected"],
}
)
tests_with_preoptions = [
{
"preoptions": "-K",
"mask": "mask.nii.gz",
"options": "-m",
"expected": "0.948930 \n0.947671 \n1.003258 \n1.010696",
},
{
"preoptions": "-t -K",
"mask": "mask.nii.gz",
"options": "-m",
"expected": "0.459736 \n0.476035 \n0.504080 \n0.549485 \n0.489194 \n0.471636 \n0.499178 \n0.461211",
},
{
"preoptions": "-t",
"mask": "",
"options": "-m",
"expected": "0.496675 \n0.487950",
},
]
# taken from fslchpixdim test
def create_image(shape, pixdim):
pixdim = list(pixdim)
data = np.random.random(shape).astype(np.float32)
hdr = nib.Nifti1Header()
hdr.set_data_dtype(np.float32)
hdr.set_data_shape(shape)
hdr.set_zooms(pixdim[:len(shape)])
return nib.Nifti1Image(data, np.eye(4), hdr)
def create_mask(input_img, n_rois=4):
'''
Create a mask image from the input image.
The mask image will have n_rois different values with random sizes randomly placed in the image.
'''
data = input_img.get_fdata()
mask = np.zeros_like(data)
for i in range(n_rois + 1):
roi_size = np.random.randint(1, data.size)
roi_value = i
roi = np.random.choice(data.size, roi_size, replace=False)
mask.flat[roi] = roi_value
return nib.Nifti1Image(mask, input_img.affine, input_img.header)
def test_fslstats():
imgfile = 'test.nii.gz'
shape = (10,10,10)
img = create_image(shape, [1] * len(shape))
img.to_filename(imgfile)
mask = create_mask(img)
mask.to_filename('mask.nii.gz')
# run the tests witoout preoptions
for test in tests:
cmd = f"fslstats {test['preoptions']} {test['mask']} {imgfile} {test['options']}"
# remove any double spaces from empty test options
cmd = " ".join(cmd.split())
print(cmd)
output = run.runfsl(cmd)
# trim any trailing whitespace from the output (fslstats adds a newline at the end)
output = output.rstrip()
try:
assert output == test["expected"]
except AssertionError:
print(cmd)
print(output)
print(test["expected"])
raise
# run the tests with preoptions
imgfile = 'test_4d.nii.gz'
shape = (10,10,10, 2)
img = create_image(shape, [1] * len(shape))
img.to_filename(imgfile)
for test in tests_with_preoptions:
cmd = f"fslstats {test['preoptions']} {test['mask']} {imgfile} {test['options']}"
# remove any double spaces from empty test options
cmd = " ".join(cmd.split())
print(cmd, flush=True)
output = run.runfsl(cmd)
# trim any trailing whitespace from the output (fslstats adds a newline at the end)
output = output.rstrip()
try:
assert output == test["expected"]
except AssertionError:
print('Failed test')
print(cmd)
print('Output')
print(output)
print('Expected')
print(test["expected"])
raise
if __name__ == "__main__":
if len(sys.argv) > 1:
os.chdir(sys.argv[1])
test_fslstats()