ntustison/a_background.md

Last active November 17, 2024 20:13

Star () You must be signed in to star a gist
Fork () You must be signed in to fork a gist

Learn more about clone URLs
Clone this repository at <script src="https://gist.github.com/ntustison/12a656a5fc2f6f9c4494c88dc09c5621.js"></script>
Save ntustison/12a656a5fc2f6f9c4494c88dc09c5621 to your computer and use it in GitHub Desktop.

Download ZIP

ANTsX Tutorial

Raw

a_background.md

ANTsX ecosystem

ANTs intro

ANTs --- Set of C++ programs built on top of the Insight Toolkit
ANTsR --- R-based ANTs interface (dependencies include ITKR, ANTsRCore)
ANTsPy --- Python-based ANTs interface
ANTsXNet --- deep learning library built on Keras
- ANTsRNet
- ANTsPyNet

ANTsR / ANTsPy

ANTsXNet

Raw

b_0_ants_interfaces_basic.md

ANTsR / ANTsPy

Installation

Python: git clone https://github.com/ANTsX/ANTsPy; cd ANTsPy; python3 setup.py install
R: devtools::install_github( "ANTsX/ANTsR" )

Basics

Read, write, and visualize images

Python

>>> import ants
>>>
>>> r16_filename = ants.get_ants_data("r16")
>>> r16_image = ants.image_read(r16_filename)
>>>
>>> ants.write_image(r16_image, "r16.nii.gz")
>>>
>>> ants.plot( r16_image )

R

> library( ANTsR )
>
> r16Filename <- getANTsRData( "r16" )
> r16Image <- antsImageRead( r16Filename )
>
> antsImageWrite( r16Image, "r16.nii.gz" )
>
> plot( r16Image )

Image info

Python

>>> import ants
>>>
>>> r16 = ants.image_read(ants.get_ants_data('r16'))
>>> r16
ANTsImage
	 Pixel Type : float (float32)
	 Components : 1
	 Dimensions : (256, 256)
	 Spacing    : (1.0, 1.0)
	 Origin     : (0.0, 0.0)
	 Direction  : [1. 0. 0. 1.]
>>> r16.shape
(256, 256)
>>> r16.origin
(0.0, 0.0)
>>> ants.get_origin(r16)
(0.0, 0.0)
>>> r16.direction
array([[1., 0.],
       [0., 1.]])
>>> ants.get_direction(r16)
array([[1., 0.],
       [0., 1.]])
>>> r16.spacing
(1.0, 1.0)              
>>> ants.get_spacing(r16)
(1.0, 1.0)

R

> library( ANTsR )
>
> r16 <- antsImageRead( getANTsRData( 'r16' ) )
> r16
antsImage
  Pixel Type          : float 
  Components Per Pixel: 1 
  Dimensions          : 256x256 
  Voxel Spacing       : 1x1 
  Origin              : 0 0 
  Direction           : 1 0 0 1 
  Filename           : /Library/Frameworks/R.framework/Versions/4.0/Resources/library/ANTsRCore/extdata/r16slice.jpg 
> dim( r16 )
[1] 256 256
> antsGetOrigin( r16 )
[1] 0 0
> antsGetDirection( r16 )
     [,1] [,2]
[1,]    1    0
[2,]    0    1
> antsGetSpacing( r16 )
[1] 1 1

Conversion to/from native types

Python

>>> import ants
>>>
>>> r16 = ants.image_read(ants.get_ants_data('r16'))
>>> r16_numpy = r16.numpy()
>>> r16_numpy
array([[0., 0., 0., ..., 0., 0., 0.],
       [0., 0., 0., ..., 0., 0., 0.],
       [0., 0., 0., ..., 0., 0., 0.],
       ...,
       [0., 0., 0., ..., 0., 0., 0.],
       [0., 0., 0., ..., 0., 0., 0.],
       [0., 0., 0., ..., 0., 0., 0.]], dtype=float32)
>>> r16_image = ants.from_numpy(r16_numpy, origin=r16.origin, spacing=r16.spacing, direction=r16.direction)

R

> library( ANTsR )
>
> r16 <- antsImageRead( getANTsRData( 'r16' ) )
> r16Array <- as.array( r16 )
> r16Image <- as.antsImage( r16Array, spacing = antsGetSpacing( r16 ), origin = antsGetOrigin( r16 ), direction = antsGetDirection( r16 ) )
> r16Image <- as.antsImage( r16Array, reference = r16 )

Arithmetic operations

Python

>>> import ants
>>>
>>> r16 = ants.image_read(ants.get_ants_data('r16'))
>>> r64 = ants.image_read(ants.get_ants_data('r64'))
>>>
>>> r_x = r16 + r64
>>> r_x = r16 - r64
>>> r_x = r16 * r64
>>> r_x = r16 / r64
>>>
>>> r_x = r16 + 3.14
>>> r_x = r16 - 3.14
>>> r_x = r16 * 3.14
>>> r_x = r16 / 3.14
>>> r_x = r16 ** 3.14

R

> library( ANTsR )
>
> r16 <- antsImageRead( getANTsRData( "r16" ) )
> r64 <- antsImageRead( getANTsRData( "r64" ) )
>
> rx <- r16 + r64
> rx <- r16 - r64
> rx <- r16 * r64
> rx <- r16 / r64
>
> rx <- r16 + 3.14
> rx <- r16 - 3.14
> rx <- r16 * 3.14
> rx <- r16 / 3.14
> rx <- r16 ^ 3.14

Other operations

Python

>>> import ants
>>>
>>> r16 = ants.image_read(ants.get_ants_data('r16'))
>>>
>>> r16.min()
0.0
>>> r16.max()
254.0
>>> r16.std()
84.8173
>>> r16.mean()
50.461365

R

> library( ANTsR )
>
> r16 <- antsImageRead( getANTsRData( "r16" ) )
>
> min( r16 )
[1] 0
> max( r16 )
[1] 254
> sd( r16 )
[1] 84.81795
> mean( r16 )
[1] 50.46136

Raw

b_1_ants_interfaces_image_registration.md

ANTsR / ANTsPy

Image registration

Python

>>> import ants
>>>
>>> r16 = ants.image_read(ants.get_ants_data('r16'))
>>> r16_seg = ants.threshold_image(r16, "Kmeans", 3)
 Kmeans with 3 thresholds
>>> ants.plot(r16, overlay=r16_seg, overlay_alpha=0.5)
>>>
>>> r64 = ants.image_read(ants.get_ants_data('r64'))
>>> r64_seg = ants.threshold_image(r64, "Kmeans", 3)
 Kmeans with 3 thresholds
>>> ants.plot(r64, overlay=r64_seg, overlay_alpha=0.5)
>>>
>>> # Affine only
>>> reg = ants.registration(fixed=r16, moving=r64, type_of_transform="antsRegistrationSyNQuick[a]", verbose=True)
>>> r16xr64_seg = ants.apply_transforms(fixed=r16, moving=r64_seg, transformlist=reg['fwdtransforms'], interpolator='nearestNeighbor')
>>> ants.plot(r16, overlay=r16xr64_seg, overlay_alpha=0.5)
>>> r64xr16_seg = ants.apply_transforms(fixed=r64, moving=r16_seg, transformlist=reg['invtransforms'], whichtoinvert=[True], interpolator='nearestNeighbor')
>>> ants.plot(r64, overlay=r64xr16_seg, overlay_alpha=0.5)
>>>
>>> # Linear and deformable
>>> reg = ants.registration(fixed=r16, moving=r64, type_of_transform="antsRegistrationSyNQuick[b]", verbose=True)
>>> r16xr64_seg = ants.apply_transforms(fixed=r16, moving=r64_seg, transformlist=reg['fwdtransforms'], interpolator='nearestNeighbor')
>>> ants.plot(r16, overlay=r16xr64_seg, overlay_alpha=0.5)
>>> r64xr16_seg = ants.apply_transforms(fixed=r64, moving=r16_seg, transformlist=reg['invtransforms'], whichtoinvert=[True, False], interpolator='nearestNeighbor')
>>> ants.plot(r64, overlay=r64xr16_seg, overlay_alpha=0.5)
>>>
>>> # Multivariate linear and deformable using intensity with wm and gm segmentations
>>> r16_gm = ants.threshold_image(r16_seg, 3, 3, 1, 0)
>>> r64_gm = ants.threshold_image(r64_seg, 3, 3, 1, 0)
>>> r16_wm = ants.threshold_image(r16_seg, 4, 4, 1, 0)
>>> r64_wm = ants.threshold_image(r64_seg, 4, 4, 1, 0)
>>> multivariate_extras = list()
>>> multivariate_extras.append(['MSQ', r16_gm, r64_gm, 2, 1])
>>> multivariate_extras.append(['MSQ', r16_wm, r64_wm, 2, 1])
>>> reg = ants.registration(fixed=r16, moving=r64, type_of_transform="antsRegistrationSyNQuick[b]", multivariate_extras = multivariate_extras, verbose=True)
>>> r16xr64_seg = ants.apply_transforms(fixed=r16, moving=r64_seg, transformlist=reg['fwdtransforms'], interpolator='nearestNeighbor')
>>> ants.plot(r16, overlay=r16xr64_seg, overlay_alpha=0.5)
>>> r64xr16_seg = ants.apply_transforms(fixed=r64, moving=r16_seg, transformlist=reg['invtransforms'], whichtoinvert=[True, False], interpolator='nearestNeighbor')
>>> ants.plot(r64, overlay=r64xr16_seg, overlay_alpha=0.5)

R

> library( ANTsR )
>
> r16 <- antsImageRead( getANTsRData( 'r16' ) )
> r16Seg <- thresholdImage( r16, "Kmeans", 3 )
 Kmeans with 3 thresholds
> plot( r16, r16Seg, alpha = 0.5 )
>
> r64 <- antsImageRead( getANTsRData( 'r64' ) )
> r64Seg <- thresholdImage( r64, "Kmeans", 3 )
 Kmeans with 3 thresholds
> plot( r64, r64Seg, alpha = 0.5 )
>
> # Affine only
> reg <- antsRegistration( fixed = r16, moving = r64, typeofTransform = "antsRegistrationSyNQuick[a]", verbose = TRUE )
> r16xr64Seg <- antsApplyTransforms( fixed = r16, moving = r64Seg, transformlist = reg$fwdtransforms, interpolator = 'nearestNeighbor' )
> plot( r16, r16xr64Seg, alpha = 0.5 ) 
> r64xr16Seg <- antsApplyTransforms( fixed = r64, moving = r16Seg, transformlist = reg$invtransforms, whichtoinvert = c( TRUE ), interpolator = 'nearestNeighbor' )
> plot( r64, r64xr16Seg, alpha = 0.5 ) 
>
> # Linear and deformable
> reg <- antsRegistration( fixed = r16, moving = r64, typeofTransform = "antsRegistrationSyNQuick[b]", verbose = TRUE )
> r16xr64Seg <- antsApplyTransforms( fixed = r16, moving = r64Seg, transformlist = reg$fwdtransforms, interpolator = 'nearestNeighbor' )
> plot( r16, r16xr64Seg, alpha = 0.5 ) 
> r64xr16Seg <- antsApplyTransforms( fixed = r64, moving = r16Seg, transformlist = reg$invtransforms, whichtoinvert = c( TRUE, FALSE ), interpolator = 'nearestNeighbor' )
> plot( r64, r64xr16Seg, alpha = 0.5 ) 
>
> # Linear and deformable using intensity with wm and gm segmentations
> r16Gm <- thresholdImage( r16Seg, 3, 3, 1, 0 )
> r64Gm <- thresholdImage( r64Seg, 3, 3, 1, 0 )
> r16Wm <- thresholdImage( r16Seg, 4, 4, 1, 0 )
> r64Wm <- thresholdImage( r64Seg, 4, 4, 1, 0 ) 
> multivariateExtras <- list()
> multivariateExtras[[1]] <- list( 'MSQ', r16Gm, r64Gm, 2, 1 )
> multivariateExtras[[2]] <- list( 'MSQ', r16Wm, r64Wm, 2, 1 )
> reg <- antsRegistration( fixed = r16, moving = r64, typeofTransform = "antsRegistrationSyNQuick[b]", multivariateExtras = multivariateExtras, verbose = TRUE )
> r16xr64Seg <- antsApplyTransforms( fixed = r16, moving = r64Seg, transformlist = reg$fwdtransforms, interpolator = 'nearestNeighbor' )
> plot( r16, r16xr64Seg, alpha = 0.5 ) 
> r64xr16Seg <- antsApplyTransforms( fixed = r64, moving = r16Seg, transformlist = reg$invtransforms, whichtoinvert = c( TRUE, FALSE ), interpolator = 'nearestNeighbor' )
> plot( r64, r64xr16Seg, alpha = 0.5 )

Template building

Python

>>> import ants
>>>
>>> r16 = ants.image_read(ants.get_ants_data('r16'))
>>> r27 = ants.image_read(ants.get_ants_data('r27'))
>>> r30 = ants.image_read(ants.get_ants_data('r30'))
>>> r62 = ants.image_read(ants.get_ants_data('r62'))
>>> r64 = ants.image_read(ants.get_ants_data('r64'))
>>> r85 = ants.image_read(ants.get_ants_data('r85'))
>>>
>>> rimage_list = [r16, r27, r30, r62, r64, r85]
>>> rtemplate = ants.build_template(initial_template=None, image_list=rimage_list, type_of_transform="antsRegistrationSyNQuick[s]")
>>> ants.plot(rtemplate)

R

> library( ANTsR )
> 
> r16 <- antsImageRead( getANTsRData( 'r16' ) )
> r27 <- antsImageRead( getANTsRData( 'r27' ) )
> r30 <- antsImageRead( getANTsRData( 'r30' ) )
> r62 <- antsImageRead( getANTsRData( 'r62' ) )
> r64 <- antsImageRead( getANTsRData( 'r64' ) )
> r85 <- antsImageRead( getANTsRData( 'r85' ) )
> 
> rimageList <- list( r16, r27, r30, r62, r64, r85 )
> initialAverage <- antsAverageImages( rimageList )
> rtemplate <- buildTemplate( initialTemplate = initialAverage, imgList = rimageList, typeofTransform = "antsRegistrationSyNQuick[s]", verbose = TRUE )
> plot( rtemplate )

Motion correction

Python

>>> import ants
>>> 
>>> # Reduce the number of times points
>>> image = ants.image_read(ants.get_ants_data("pcasl"))
>>> image_array = image.numpy()
>>> image_array = image_array[:,:,:,:10]
>>> image = ants.from_numpy(image_array, spacing=image.spacing, direction=image.direction, origin=image.origin)
>>> 
>>> moco = ants.motion_correction(image, type_of_transform="BOLDRigid", verbose=True)
>>> moco['motion_corrected']
ANTsImage
	 Pixel Type : float (float32)
	 Components : 1
	 Dimensions : (64, 64, 18, 10)
	 Spacing    : (3.4375, 3.4375, 7.2, 4.0)
	 Origin     : (-116.0013, 57.1198, -38.3576, 0.0)
	 Direction  : [ 1.  0.  0.  0.  0. -1.  0.  0.  0.  0.  1.  0.  0.  0.  0.  1.]
>>> moco['motion_parameters']
[['/var/folders/5h/62f7bql14zvgx34nqwd7xvnc0000gn/T/tmpy49zfa890GenericAffine.mat'], ['/var/folders/5h/62f7bql14zvgx34nqwd7xvnc0000gn/T/tmp_qncwnum0GenericAffine.mat'], ['/var/folders/5h/62f7bql14zvgx34nqwd7xvnc0000gn/T/tmpwlnjbm_30GenericAffine.mat'], ['/var/folders/5h/62f7bql14zvgx34nqwd7xvnc0000gn/T/tmppfj9lmyv0GenericAffine.mat'], ['/var/folders/5h/62f7bql14zvgx34nqwd7xvnc0000gn/T/tmpw5y0z0ot0GenericAffine.mat'], ['/var/folders/5h/62f7bql14zvgx34nqwd7xvnc0000gn/T/tmpztsq1ete0GenericAffine.mat'], ['/var/folders/5h/62f7bql14zvgx34nqwd7xvnc0000gn/T/tmptw81m6lx0GenericAffine.mat'], ['/var/folders/5h/62f7bql14zvgx34nqwd7xvnc0000gn/T/tmptxnjmkc_0GenericAffine.mat'], ['/var/folders/5h/62f7bql14zvgx34nqwd7xvnc0000gn/T/tmpwuv83l_a0GenericAffine.mat'], ['/var/folders/5h/62f7bql14zvgx34nqwd7xvnc0000gn/T/tmps8_zjfy60GenericAffine.mat']]
>>> moco['FD']
array([0.        , 0.09816545, 0.18367779, 0.12037325, 0.13567556,
       0.09293036, 0.09830562, 0.18831979, 0.04102131, 0.06663223])

R

> library( ANTsR )
> 
> # Reduce the number of times points
> image <- antsImageRead( getANTsRData( "pcasl" ) )
> imageArray <- as.array( image )
> imageArray <- imageArray[,,,1:10]
> image <- as.antsImage( imageArray, reference = image )
> 
> moco <- antsMotionCalculation( image, txtype = "Rigid", verbose = TRUE )
> moco$moco_img
antsImage
  Pixel Type          : float 
  Components Per Pixel: 1 
  Dimensions          : 64x64x18x10 
  Voxel Spacing       : 3.4375x3.4375x7.19999980926514x4 
  Origin              : -116.0013 57.11981 -38.35765 0 
  Direction           : 1 0 0 0 0 -1 0 0 0 0 1 0 0 0 0 1 
> moco$moco_params
   MetricPre MetricPost     MOCOparam1     MOCOparam2     MOCOparam3
1          0  -1.005069  0.00023146855 -0.00009585432  0.00004607494
2          0  -1.005325 -0.00107804559 -0.00067556275  0.00025181432
3          0  -1.009195 -0.00018273871  0.00001386437 -0.00002405675
4          0  -1.010406 -0.00015799195 -0.00019400520  0.00046687210
5          0  -1.007293  0.00011430427 -0.00004679069 -0.00001049459
6          0  -1.005072  0.00022220720 -0.00007070929 -0.00003112626
7          0  -1.008751  0.00009608059 -0.00001348405  0.00001652053
8          0  -1.003125 -0.00132131224  0.00013327676 -0.00000959006
9          0  -1.005963 -0.00009936649  0.00001233287 -0.00002474547
10         0  -1.007560 -0.00061640351 -0.00064586266 -0.00052936993
     MOCOparam4   MOCOparam5   MOCOparam6
1   0.010914629 -0.008945767 -0.111323305
2   0.001277686 -0.079050128 -0.091623102
3  -0.021036213  0.014788211  0.098308572
4  -0.003656350  0.037733644  0.004530848
5   0.007828627  0.008543149 -0.075201220
6  -0.018174798  0.003812138 -0.146378820
7   0.015950238 -0.030493634 -0.063736015
8  -0.010916564  0.004372850 -0.006914493
9  -0.002273881  0.022307026  0.068820824
10  0.027577054 -0.036157711  0.055110048
> moco$fd
   MeanDisplacement MaxDisplacement
1        0.09661564       0.1930724
2        0.18165174       0.2670051
3        0.10522685       0.1257340
4        0.09715838       0.1322308
5        0.08049656       0.0901812
6        0.10052631       0.1125070
7        0.13801552       0.2513173
8        0.05978751       0.1279247
9        0.06027606       0.1109254
10       0.00000000       0.0000000

Raw

b_2_ants_interfaces_landmark_registration.md

Landmark registration

Python

Two simple examples

>>> import ants
>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> 
>>> examples = ("Big square to little square",
>>>             "Brian's canonical example")
>>> 
>>> moving_points_list = list()
>>> fixed_points_list = list()
>>> 
>>> # Big square to little square
>>> moving_points = np.zeros((4, 2))
>>> moving_points[0, 0] = 64
>>> moving_points[0, 1] = 64
>>> moving_points[1, 0] = 192
>>> moving_points[1, 1] = 64
>>> moving_points[2, 0] = 64
>>> moving_points[2, 1] = 192
>>> moving_points[3, 0] = 192
>>> moving_points[3, 1] = 192
>>> moving_points_list.append(moving_points)
>>> 
>>> fixed_points = np.zeros((4, 2))
>>> fixed_points[0, 0] = 96
>>> fixed_points[0, 1] = 96
>>> fixed_points[1, 0] = 160
>>> fixed_points[1, 1] = 96
>>> fixed_points[2, 0] = 96
>>> fixed_points[2, 1] = 160
>>> fixed_points[3, 0] = 160
>>> fixed_points[3, 1] = 160
>>> fixed_points_list.append(fixed_points)
>>> 
>>> # Brian's example
>>> moving_points = np.zeros((3, 2))
>>> moving_points[0, 0] = 64
>>> moving_points[0, 1] = 64
>>> moving_points[1, 0] = 192
>>> moving_points[1, 1] = 64
>>> moving_points[2, 0] = 64
>>> moving_points[2, 1] = 192
>>> moving_points_list.append(moving_points)
>>> 
>>> fixed_points = np.zeros((3, 2))
>>> fixed_points[0, 0] = 192
>>> fixed_points[0, 1] = 192
>>> fixed_points[1, 0] = 192
>>> fixed_points[1, 1] = 64
>>> fixed_points[2, 0] = 64
>>> fixed_points[2, 1] = 192
>>> fixed_points_list.append(fixed_points)
>>> 
>>> # domain image
>>> r16 = ants.image_read(ants.get_ants_data("r16"))
>>> 
>>> transform_types = ('rigid', 'similarity', 'affine', 'bspline', 'tps', 'diffeo', 'syn', 'time-varying')
>>> 
>>> for e in range(len(examples)):
>>>     warped_grid_arrays = list()
>>>     for i in range(len(transform_types)):
>>>         print("Fitting with transform ", transform_types[i])
>>>         xfrm = ants.fit_transform_to_paired_points(moving_points_list[e], fixed_points_list[e], 
>>>                            transform_type=transform_types[i],
>>>                           domain_image=r16, 
>>>                           number_of_fitting_levels=4, 
>>>                           mesh_size=(4, 4), 
>>>                           number_of_compositions=10000,
>>>                           number_of_integration_steps=10,
>>>                           number_of_time_steps=2,
>>>                           composition_step_size=0.2,
>>>                           convergence_threshold=1e-6,
>>>                           enforce_stationary_boundary=True, 
>>>                           verbose=True)
>>>        grid = ants.create_warped_grid(r16)
>>>        if isinstance(xfrm, dict):
>>>            warped_grid = xfrm['forward_transform'].apply_to_image(grid, r16)
>>>        else:    
>>>            warped_grid = xfrm.apply_to_image(grid, r16)
>>>        warped_grid_arrays.append(warped_grid.numpy()) 
>>>    fig, axes = plt.subplots(2, 4, figsize=(10, 10))
>>>    for i, ax in enumerate(axes.flat):
>>>        ax.imshow(warped_grid_arrays[i], cmap="gray")
>>>        if transform_types[i] == "bspline":
>>>            ax.set_title("Transform: bspline (clamped)")
>>>        else:
>>>            ax.set_title("Transform: " + transform_types[i])
>>>        ax.xaxis.set_visible(False)
>>>        ax.yaxis.set_visible(False)
>>>        for j in range(moving_points_list[e].shape[0]):
>>>            ax.arrow(moving_points_list[e][j,0], moving_points_list[e][j,1],
>>>                     fixed_points_list[e][j,0]-moving_points_list[e][j,0],
>>>                     fixed_points_list[e][j,1]-moving_points_list[e][j,1],
>>>                     head_width=10.0,
>>>                     width=2.0,
>>>                     color='red',
>>>                     ec='red')
>>>
>>>    fig.suptitle(examples[e])
>>>    fig.tight_layout()
>>>    plt.show()
>>>    plt.close()

More complicated example

>>> import ants
>>> import numpy as np
>>> import os
>>> import matplotlib.pyplot as plt
>>> 
>>> os.environ["ITK_GLOBAL_DEFAULT_NUMBER_OF_THREADS"] = "2"
>>> 
>>> x = np.linspace(-1, 1, num=20)
>>> 
>>> moving_x = np.empty(x.shape)
>>> moving_y = np.empty(x.shape)
>>> fixed_x = np.empty(x.shape)
>>> fixed_y = np.empty(x.shape)
>>> 
>>> center_x = 128
>>> center_y = 128
>>> offset_x = 50
>>> offset_y = 50
>>> scale = 64
>>> 
>>> for i in range(len(x)):
>>>     fixed_x[i] = x[i] * scale + center_x + offset_x
>>>     fixed_y[i] = x[i] ** 2 * scale + center_y + offset_y
>>>     moving_x[i] = x[i] * scale + center_x
>>>     moving_y[i] = -x[i] ** 2 * scale + center_y
>>> 
>>> fixed_points = np.vstack((fixed_x, fixed_y)).transpose()
>>> moving_points = np.vstack((moving_x, moving_y)).transpose()
>>> 
>>> #########   rigid transform    ############
>>> 
>>> rigid_xfrm = ants.fit_transform_to_paired_points(moving_points, fixed_points, transform_type="rigid")
>>> 
>>> fixed_points_rigid_warped = np.empty(fixed_points.shape)
>>> for j in range(fixed_points.shape[0]):
>>>     fixed_points_rigid_warped[j,:] = rigid_xfrm.apply_to_point(tuple(fixed_points[j,:]))
>>> 
>>> #########   syn transform    ############
>>> 
>>> r16 = ants.image_read(ants.get_ants_data("r16"))
>>> 
>>> syn_xfrm = ants.fit_transform_to_paired_points(moving_points, fixed_points_rigid_warped, 
>>>                       transform_type="syn",
>>>                       domain_image=r16, 
>>>                       number_of_fitting_levels=4, 
>>>                       mesh_size=(9, 9), 
>>>                       number_of_compositions=100,
>>>                       composition_step_size=0.2,
>>>                       enforce_stationary_boundary=True, 
>>>                       verbose=True)
>>> 
>>> total_forward_xfrm = ants.compose_ants_transforms([rigid_xfrm, syn_xfrm['forward_transform']])
>>> total_fixed_to_middle_xfrm = ants.compose_ants_transforms([rigid_xfrm, syn_xfrm['fixed_to_middle_transform']])
>>> 
>>> fixed_points_syn_warped = np.empty(fixed_points.shape)
>>> fixed_points_syn_mid_warped = np.empty(fixed_points.shape)
>>> moving_points_syn_warped = np.empty(fixed_points.shape)
>>> moving_points_syn_mid_warped = np.empty(fixed_points.shape)
>>> 
>>> for j in range(fixed_points.shape[0]):
>>>     fixed_points_syn_warped[j,:] = total_forward_xfrm.apply_to_point(tuple(fixed_points[j,:]))
>>>     fixed_points_syn_mid_warped[j,:] = total_fixed_to_middle_xfrm.apply_to_point(tuple(fixed_points[j,:]))
>>>     moving_points_syn_warped[j,:] = syn_xfrm['inverse_transform'].apply_to_point(tuple(moving_points[j,:]))
>>>     moving_points_syn_mid_warped[j,:] = syn_xfrm['moving_to_middle_transform'].apply_to_point(tuple(moving_points[j,:]))
>>> 
>>> plt.plot(fixed_points[:,0], fixed_points[:,1], 'g^', label="Fixed points") 
>>> plt.plot(moving_points[:,0], moving_points[:,1], 'go', label="Moving points")
>>> plt.plot(fixed_points_rigid_warped[:,0], fixed_points_rigid_warped[:,1], 'r^', label="Fixed (rigid)")
>>> plt.plot(fixed_points_syn_warped[:,0], fixed_points_syn_warped[:,1], 'r^--', alpha=0.33, label="Fixed (SyN-forward)") 
>>> plt.plot(moving_points_syn_warped[:,0], moving_points_syn_warped[:,1], 'ro--', alpha=0.33, label="Moving (SyN-inverse)") 
>>> plt.plot(fixed_points_syn_mid_warped[:,0], fixed_points_syn_mid_warped[:,1], 'b^', alpha=1.0, label="Fixed (SyN-mid)") 
>>> plt.plot(moving_points_syn_mid_warped[:,0], moving_points_syn_mid_warped[:,1], 'bo--', alpha=0.33, label="Moving (SyN-mid)")
>>> 
>>> plt.legend()
>>> plt.title('Rigid + SyN landmark registration')
>>> plt.show()
>>>
>>> grid = ants.create_warped_grid(r16)
>>> warped_grid = syn_xfrm['forward_transform'].apply_to_image(grid, r16)
>>> ants.plot(warped_grid)

R

Two simple examples

> library( ANTsR )
> library( ggplot2 )
> library( raster )
> 
> examples <- c("Big square to little square",
>               "Brian's canonical example")
> 
> movingPointsList <- list()
> fixedPointsList <- list()
> 
> # Big square to little square
> movingPoints <- array( data = 0, dim = c( 4, 2 ) )
> movingPoints[1, 1] <- 64
> movingPoints[1, 2] <- 64
> movingPoints[2, 1] <- 192
> movingPoints[2, 2] <- 64
> movingPoints[3, 1] <- 64
> movingPoints[3, 2] <- 192
> movingPoints[4, 1] <- 192
> movingPoints[4, 2] <- 192
> movingPointsList[[1]] <- movingPoints
> 
> fixedPoints <- array( data = 0, dim = c( 4, 2 ) )
> fixedPoints[1, 1] <- 96
> fixedPoints[1, 2] <- 96
> fixedPoints[2, 1] <- 160
> fixedPoints[2, 2] <- 96
> fixedPoints[3, 1] <- 96
> fixedPoints[3, 2] <- 160
> fixedPoints[4, 1] <- 160
> fixedPoints[4, 2] <- 160
> fixedPointsList[[1]] <- fixedPoints
> 
> # Brian's example
> movingPoints <- array( data = 0, dim = c( 3, 2 ) )
> movingPoints[1, 1] <- 64
> movingPoints[2, 1] <- 192
> movingPoints[2, 2] <- 64
> movingPoints[3, 1] <- 64
> movingPoints[3, 2] <- 192
> movingPointsList[[2]] <- movingPoints
>
> fixedPoints[1, 1] <- 192
> fixedPoints[1, 2] <- 192
> fixedPoints[2, 1] <- 192
> fixedPoints[2, 2] <- 64
> fixedPoints[3, 1] <- 64
> fixedPoints[3, 2] <- 192
> fixedPointsList[[2]] <- fixedPoints
> 
> # domain image
> r16 <- antsImageRead( getANTsRData( "r16" ) )
> 
> transformTypes <- c( 'rigid', 'similarity', 'affine', 'bspline', 'tps', 'diffeo', 'syn', 'time-varying' )
> 
> for( e in seq.int( length( examples ) ) )
>   {
>   warpedGridArrays <- list()
>   for( i in seq.int( length( transformTypes ) ) )
>     {
>     cat( paste0( "Fitting with transform ", transformTypes[i], "\n" ) )
>     xfrm <- fitTransformToPairedPoints( movingPointsList[[e]], fixedPointsList[[e]], 
>                        transformType = transformTypes[i],
>                        domainImage = r16, 
>                        numberOfFittingLevels = 4, 
>                        meshSize = c( 4, 4 ), 
>                        numberOfCompositions = 10000,
>                        numberOfIntegrationSteps = 10,
>                        numberOfTimeSteps = 2,
>                        compositionStepSize = 0.2,
>                        convergenceThreshold = 1e-6,
>                        enforceStationaryBoundary = TRUE, 
>                        verbose = TRUE )
>     grid <- createWarpedGrid( r16 )
>     if( is.list( xfrm ) )
>       {
>       warpedGrid <- applyAntsrTransformToImage( xfrm$forwardTransform, grid, r16 )
>       } else {
>       warpedGrid <- applyAntsrTransformToImage( xfrm, grid, r16 )
>       }
>     warpedGridArrays[[i]] <- as.array( warpedGrid )
>     }
> 
>   plotDataFrame <- data.frame()
>   for( i in seq.int( length( warpedGridArrays ) ) )  
>     {
>     rasterDataFrame <- as.data.frame( raster( warpedGridArrays[[i]] ), xy = TRUE )
>     rasterDataFrame$TransformType <- rep( transformTypes[i], nrow( rasterDataFrame ) )
>     if( i == 1 )
>       {
>       plotDataFrame <- rasterDataFrame
>       } else {
>       plotDataFrame <- rbind( plotDataFrame, rasterDataFrame )
>       }
>     }
> 
>   plotDataFrame$TransformType[which( plotDataFrame$TransformType == "bspline" )] <- "bspline (clamped)"
>   plotDataFrame$TransformType <- factor( plotDataFrame$TransformType, levels = unique( plotDataFrame$TransformType ) )  
> 
>   p <- ggplot( data = plotDataFrame, aes( x = x, y = y, fill = layer ), colour = "gray" ) + 
>        geom_raster() + 
>        facet_wrap( ~TransformType ) + 
>        scale_fill_gradientn( colours = c( "black", "white" ) ) +
>        theme( legend.position = "none",
>               axis.title.x = element_blank(), 
>               axis.title.y = element_blank(), 
>               axis.text.x = element_blank(), 
>               axis.text.y = element_blank(), 
>               axis.ticks.x = element_blank(), 
>               axis.ticks.y = element_blank()
>               ) +
>        ggtitle( examples[e] )      
> 
>   # ggsave( paste0( "~/Desktop/foo", e, ".pdf" ), plot = p )            
>   }

More complicated example

> library( ANTsR )
> library( ggplot2 )
> 
> Sys.setenv( ITK_GLOBAL_DEFAULT_NUMBER_OF_THREADS = "2" )
> 
> x <- seq( -1, 1, length.out = 20 )
> 
> movingX <- c()
> movingY <- c()
> fixedX <- c()
> fixedY <- c()
> 
> centerX <- 128
> centerY <- 128
> offsetY <- 50
> offsetX <- 50
> offsetY <- 50
> scale <- 64
> 
> for( i in seq.int( length( x ) ) )
>   {
>   fixedX[i] <- x[i] * scale + centerX + offsetX 
>   fixedY[i] <- x[i]^2 * scale + centerY + offsetY
>   movingX[i] <- x[i] * scale + centerX
>   movingY[i] <- -x[i]^2 * scale + centerY
>   }
> 
> fixedPoints <- cbind( fixedX, fixedY )
> movingPoints <- cbind( movingX, movingY )
> 
> #########   rigid transform    ############
> 
> rigidXfrm <- fitTransformToPairedPoints( movingPoints, fixedPoints, transformType = "rigid" )
> 
> fixedPointsRigidWarped <- applyAntsrTransform( rigidXfrm, fixedPoints )
> 
> #########   syn transform    ############
> 
> r16 <- antsImageRead( getANTsRData( "r16" ) )
> synXfrm <- fitTransformToPairedPoints( movingPoints, fixedPointsRigidWarped,
>                 transformType = "syn",
>                 domainImage = r16,
>                 numberOfFittingLevels = 4,
>                 meshSize = c( 9, 9 ),
>                 numberOfCompositions = 100,
>                 compositionStepSize = 0.2,
>                 enforceStationaryBoundary = TRUE,
>                 verbose = TRUE )
> 
> 
> totalForwardXfrm <- composeAntsrTransforms( list( rigidXfrm, synXfrm$forwardTransform ) )
> totalFixedToMiddleXfrm <- composeAntsrTransforms( list( rigidXfrm, synXfrm$fixedToMiddleTransform ) )
> 
> fixedPointsSynWarped <- applyAntsrTransform( totalForwardXfrm, fixedPoints )
> fixedPointsSynMidWarped <- applyAntsrTransform( totalFixedToMiddleXfrm, fixedPoints )
> movingPointsSynWarped <- applyAntsrTransform( synXfrm$inverseTransform, movingPoints )
> movingPointsSynMidWarped <- applyAntsrTransform( synXfrm$movingToMiddleTransform, movingPoints )
> 
> numberOfPoints <- dim( fixedPoints )[1]
> allPointsX <- c( fixedPoints[,1], movingPoints[,1], fixedPointsRigidWarped[,1],
>                  fixedPointsSynWarped[,1], movingPointsSynWarped[,1],
>                  fixedPointsSynMidWarped[,1], movingPointsSynMidWarped[,1] )
> allPointsY <- c( fixedPoints[,2], movingPoints[,2], fixedPointsRigidWarped[,2],
>                  fixedPointsSynWarped[,2], movingPointsSynWarped[,2],
>                  fixedPointsSynMidWarped[,2], movingPointsSynMidWarped[,2] )
> pointsType <- c( rep( "Fixed points", numberOfPoints ), rep( "Moving points", numberOfPoints ), rep( "Fixed (rigid)", numberOfPoints ), 
>                  rep( "Fixed (SyN-forward)", numberOfPoints ), rep( "Moving (SyN-inverse))", numberOfPoints ),
>                  rep( "Fixed (SyN-mid)", numberOfPoints ), rep( "Moving (SyN-mid))", numberOfPoints )
>                )
> 
> pointPlotDf <- data.frame( X = allPointsX, Y = allPointsY, Type = pointsType )
> pointPlot <- ggplot( data = pointPlotDf ) + 
>              geom_point( aes( x = X, y = Y, colour = Type, shape = Type ) ) + 
>              geom_line( aes( x = X, y = Y, colour = Type, group = Type ), alpha = 0.25 )
> pointPlot
> 
> grid <- createWarpedGrid( r16 )
> warpedGrid <- applyAntsrTransform( synXfrm$forwardTransform, grid )
> plot( warpedGrid )

Raw

b_3_ants_velocity_flows.md

Velocity flows across point sets

A more complicated example.

Python

>>> import ants
>>> 
>>> import math
>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> 
>>> domain_array = np.zeros((100, 100))
>>> domain_image = ants.from_numpy(domain_array)
>>> 
>>> number_of_points = 8
>>> dimensionality = 2
>>> 
>>> rectangle = np.zeros((number_of_points, dimensionality))
>>> rectangle[0, :] = [35.0, 15.0]
>>> rectangle[1, :] = [35.0, 50.0]
>>> rectangle[2, :] = [35.0, 85.0]
>>> rectangle[3, :] = [50.0, 85.0]
>>> rectangle[4, :] = [65.0, 85.0]
>>> rectangle[5, :] = [65.0, 50.0]
>>> rectangle[6, :] = [65.0, 15.0]
>>> rectangle[7, :] = [50.0, 15.0]
>>> 
>>> square = np.zeros((number_of_points, dimensionality))
>>> square[0, :] = [25.0, 25.0]
>>> square[1, :] = [25.0, 50.0]
>>> square[2, :] = [25.0, 75.0]
>>> square[3, :] = [50.0, 75.0]
>>> square[4, :] = [75.0, 75.0]
>>> square[5, :] = [75.0, 50.0]
>>> square[6, :] = [75.0, 25.0]
>>> square[7, :] = [50.0, 25.0]
>>> 
>>> circle = np.zeros((number_of_points, dimensionality))
>>> for i in range(number_of_points):
>>>     circle[i, 0] = 50.0 + 25.0 * math.sin(2 * math.pi * i / 8 - 3 * math.pi / 4)
>>>     circle[i, 1] = 50.0 + 25.0 * math.cos(2 * math.pi * i / 8 - 3 * math.pi / 4)
>>> 
>>> point_sets = list()
>>> point_sets.append(rectangle)
>>> point_sets.append(square)
>>> point_sets.append(circle)
>>> 
>>> plt.plot(rectangle[:,0], rectangle[:,1], 'r-', label="Rectangle")
>>> plt.plot(square[:,0], square[:,1], 'g-', label="Square")
>>> plt.plot(circle[:,0], circle[:,1], 'b-', label="Circle")
>>> 
>>> plt.legend()
>>> plt.title('Original data')
>>> plt.show()
>>> 
>>> tv = ants.fit_time_varying_transform_to_point_sets(point_sets, time_points=None,
>>>   number_of_time_steps=None, domain_image=domain_image,
>>>   number_of_fitting_levels=4, mesh_size=2, spline_order=3, displacement_weights=None,
>>>   number_of_compositions=100, composition_step_size=0.5, number_of_integration_steps=10,
>>>   sigma=0.0, convergence_threshold=0.0, rasterize_points=False, verbose=True)
>>> 
>>> # Warp between the endpoints
>>> 
>>> warped_rectangle_to_circle = np.empty_like(rectangle)
>>> warped_rectangle_to_circle[:] = rectangle
>>> for j in range(rectangle.shape[0]):
>>>     warped_rectangle_to_circle[j,:] = tv['forward_transform'].apply_to_point(tuple(rectangle[j,:]))
>>> 
>>> warped_circle_to_rectangle = np.empty_like(circle)
>>> warped_circle_to_rectangle[:] = circle
>>> for j in range(circle.shape[0]):
>>>     warped_circle_to_rectangle[j,:] = tv['inverse_transform'].apply_to_point(tuple(circle[j,:]))
>>> 
>>> plt.plot(warped_rectangle_to_circle[:,0], warped_rectangle_to_circle[:,1], 'r--', label="RectangleToCircle")
>>> plt.plot(circle[:,0], circle[:,1], 'r-', label="Circle")
>>> plt.plot(warped_circle_to_rectangle[:,0], warped_circle_to_rectangle[:,1], 'b--', label="CircleToRectangle")
>>> plt.plot(rectangle[:,0], rectangle[:,1], 'b-', label="Rectangle")
>>> 
>>> plt.legend()
>>> plt.title('Warping between endpoints')
>>> plt.show()
>>> 
>>> # Warp the endpoints (i.e., rectangle and circle) to the middle (i.e., square)
>>> 
>>> forward_xfrm = ants.transform_from_displacement_field(ants.integrate_velocity_field(tv['velocity_field'], 0, 0.5))
>>> warped_rectangle_to_square = np.empty_like(rectangle)
>>> warped_rectangle_to_square[:] = rectangle
>>> for j in range(rectangle.shape[0]):
>>>     warped_rectangle_to_square[j,:] = forward_xfrm.apply_to_point(tuple(rectangle[j,:]))
>>> 
>>> inverse_xfrm = ants.transform_from_displacement_field(ants.integrate_velocity_field(tv['velocity_field'], 1.0, 0.5))
>>> warped_circle_to_square = np.empty_like(circle)
>>> warped_circle_to_square[:] = circle
>>> for j in range(circle.shape[0]):
>>>     warped_circle_to_square[j,:] = inverse_xfrm.apply_to_point(tuple(circle[j,:]))
>>> 
>>> plt.plot(warped_rectangle_to_square[:,0], warped_rectangle_to_square[:,1], 'r--', label="RectangleToSquare")
>>> plt.plot(warped_circle_to_square[:,0], warped_circle_to_square[:,1], 'b--', label="CircleToSquare")
>>> plt.plot(rectangle[:,0], rectangle[:,1], 'r-', label="Rectangle")
>>> plt.plot(square[:,0], square[:,1], 'g-', label="Square")
>>> plt.plot(circle[:,0], circle[:,1], 'b-', label="Circle")
>>> 
>>> plt.legend()
>>> plt.title('Warping to middle')
>>> plt.show()

R

> library( ANTsR )
> library( ggplot2 )
> 
> domainArray <- array( data = 0, dim = c( 100, 100 ) )
> domainImage <- as.antsImage( domainArray )
> 
> numberOfPoints <- 8
> dimensionality <- 2
> 
> rectangle <- array( data = 0, dim = c( numberOfPoints, dimensionality ) )
> rectangle[1,] <- c( 35.0, 15.0 )
> rectangle[2,] <- c( 35.0, 50.0 )
> rectangle[3,] <- c( 35.0, 85.0 )
> rectangle[4,] <- c( 50.0, 85.0 )
> rectangle[5,] <- c( 65.0, 85.0 )
> rectangle[6,] <- c( 65.0, 50.0 )
> rectangle[7,] <- c( 65.0, 15.0 )
> rectangle[8,] <- c( 50.0, 15.0 )
> 
> square <- array( data = 0, dim = c( numberOfPoints, dimensionality ) )
> square[1,] <- c( 25.0, 25.0 )
> square[2,] <- c( 25.0, 50.0 )
> square[3,] <- c( 25.0, 75.0 )
> square[4,] <- c( 50.0, 75.0 )
> square[5,] <- c( 75.0, 75.0 )
> square[6,] <- c( 75.0, 50.0 )
> square[7,] <- c( 75.0, 25.0 )
> square[8,] <- c( 50.0, 25.0 )
> 
> circle <- array( data = 0, dim = c( numberOfPoints, dimensionality ) )
> for( i in seq_len( numberOfPoints ) )
>   {
>   circle[i, 1] <- 50.0 + 25.0 * sin( 2 * pi * ( i - 1 ) / 8 - 3 * pi / 4 )
>   circle[i, 2] <- 50.0 + 25.0 * cos( 2 * pi * ( i - 1 ) / 8 - 3 * pi / 4 )
>   }
> 
> pointSets <- list()
> pointSets[[1]] <- rectangle
> pointSets[[2]] <- square
> pointSets[[3]] <- circle
> 
> allPointsX <- c( rectangle[,1], square[,1], circle[,1] )
> allPointsY <- c( rectangle[,2], square[,2], circle[,2] )
> pointsType <- c( rep( "Rectangle", numberOfPoints ), rep( "Square", numberOfPoints ), rep( "Circle", numberOfPoints ) )
> 
> pointPlotDf <- data.frame( X = allPointsX, Y = allPointsY, Type = pointsType )
> pointPlot <- ggplot( data = pointPlotDf ) +
>              geom_point( aes( x = X, y = Y, colour = Type, shape = Type ) ) +
>              geom_path( aes( x = X, y = Y, colour = Type ), alpha = 0.5 ) +
>              ggtitle( "Original data" )
> 
> pointPlot
> 
> tv <- fitTimeVaryingTransformToPointSets( pointSets, timePoints = NULL,
>   numberOfTimeSteps = 3, domainImage = domainImage,
>   numberOfFittingLevels = 4, meshSize = 2, splineOrder = 3, displacementWeights = NULL,
>   numberOfCompositions = 100, compositionStepSize = 0.5, numberOfIntegrationSteps = 10,
>   sigma = 0.0, convergenceThreshold = 0.0, rasterizePoints = FALSE, verbose = TRUE )
> 
> # Warp between the endpoints
> 
> warpedRectangleToCircle <- applyAntsrTransformToPoint( tv$forwardTransform, rectangle )
> warpedCircleToRectangle <- applyAntsrTransformToPoint( tv$inverseTransform, circle )
> 
> allPointsX <- c( warpedRectangleToCircle[,1], circle[,1], warpedCircleToRectangle[,1], rectangle[,1] )
> allPointsY <- c( warpedRectangleToCircle[,2], circle[,2], warpedCircleToRectangle[,2], rectangle[,2] )
> pointsType <- c( rep( "RectangleToCircle", numberOfPoints ), rep( "Circle", numberOfPoints ),
>                  rep( "CircleToRectangle", numberOfPoints ), rep( "Rectangle", numberOfPoints ) )
> 
> pointPlotDf <- data.frame( X = allPointsX, Y = allPointsY, Type = pointsType )
> pointPlot <- ggplot( data = pointPlotDf ) +
>              geom_point( aes( x = X, y = Y, colour = Type, shape = Type ) ) +
>              geom_path( aes( x = X, y = Y, colour = Type ), alpha = 0.25 ) +
>              ggtitle( "Warping between endpoints" )
> 
> pointPlot
> 
> # Warp the endpoints (i.e., rectangle and circle) to the middle (i.e., square)
> 
> forwardXfrm <- createAntsrTransform( type = "DisplacementFieldTransform",
>         displacement.field = integrateVelocityField( tv$velocityField, 0, 0.5 ) )
> warpedRectangleToSquare <- applyAntsrTransformToPoint( forwardXfrm, rectangle )
> inverseXfrm <- createAntsrTransform( type = "DisplacementFieldTransform",
>         displacement.field = integrateVelocityField( tv$velocityField, 1.0, 0.5 ) )
> warpedCircleToSquare <- applyAntsrTransformToPoint( inverseXfrm, circle )
> 
> allPointsX <- c( rectangle[,1], square[,1], circle[,1], warpedRectangleToSquare[,1], warpedCircleToSquare[,1]  )
> allPointsY <- c( rectangle[,2], square[,2], circle[,2], warpedRectangleToSquare[,2], warpedCircleToSquare[,2] )
> pointsType <- c( rep( "Rectangle", numberOfPoints ), rep( "Square", numberOfPoints ), rep( "Circle", numberOfPoints ),
>                  rep( "RectangleToSquare", numberOfPoints ), rep( "CircleToSquare", numberOfPoints ) )
> 
> pointPlotDf <- data.frame( X = allPointsX, Y = allPointsY, Type = pointsType )
> pointPlot <- ggplot( data = pointPlotDf ) +
>              geom_point( aes( x = X, y = Y, colour = Type, shape = Type ) ) +
>              geom_path( aes( x = X, y = Y, colour = Type ), alpha = 0.25 ) +
>              ggtitle( "Warping to middle" )
> 
> pointPlot

Raw

b_4_ants_interfaces_misc_utilities.md

ANTsR / ANTsPy

Miscellaneous utilities

Bias correction and denoising

Python

>>> import ants
>>>
>>> r16 = ants.image_read(ants.get_ants_data("r16"))
>>> r16_mask = ants.get_mask(r16)
>>> ants.plot(r16, overlay=r16_mask, overlay_alpha = 0.5)
>>> 
>>> # Bias correction
>>> r16_n3 = ants.n3_bias_field_correction2(r16, r16_mask, verbose=True)
>>> r16_n4 = ants.n4_bias_field_correction(r16, r16_mask, verbose=True)
>>> 
>>> # Denoising
>>> r16_denoised = ants.denoise_image(r16_n3, r16_mask, v=1)

R

> library( ANTsR )
>
> r16 <- antsImageRead( getANTsRData( "r16" ) )
> r16Mask <- getMask( r16 )
> plot( r16, r16Mask, alpha = 0.5 )
> 
> # Bias correction
> r16N3 <- n3BiasFieldCorrection2( r16, r16Mask, verbose = TRUE )
> r16N4 <- n4BiasFieldCorrection( r16, r16Mask, verbose = TRUE )
> 
> # Denoising
> r16Denoised <- denoiseImage( r16N3, r16Mask, verbose = TRUE )

Atropos segmentation and quantification

Python

>>> import ants
>>>
>>> r16 = ants.image_read(ants.get_ants_data("r16"))
>>> r16_mask = ants.get_mask(r16)
>>> 
>>> # Atropos
>>> r16_atropos = ants.atropos(r16, r16_mask)
>>> ants.plot(r16, overlay=r16_atropos['segmentation'], overlay_alpha=0.5)
>>> 
>>> # Label measures
>>> ants.label_geometry_measures(r16_atropos['segmentation'])
   Label  VolumeInMillimeters  SurfaceAreaInMillimetersSquared  Eccentricity  ...  BoundingBoxUpper_y  IntegratedIntensity  WeightedCentroid_x  WeightedCentroid_y
0      1               2328.0                      2018.412135      0.579791  ...                 212                 2328          132.506014          131.398196
1      2               7319.0                      3645.703919      0.540221  ...                 211                14638          124.857084          125.767045
2      3               8370.0                      2155.330184      0.646684  ...                 209                25110          124.823775          131.760096

[3 rows x 17 columns]

R

> library( ANTsR )
>
> r16 <- antsImageRead( getANTsRData( "r16" ) )
> r16Mask <- getMask( r16 )
> 
> # Atropos
> r16Atropos <- atropos( r16, r16Mask )
> plot( r16, r16Mask, alpha = 0.5 )
> 
> # Label measures
> labelGeometryMeasures( r16Atropos$segmentation )
  Label VolumeInMillimeters SurfaceAreaInMillimetersSquared Eccentricity
1     1                2331                        2020.078    0.5797638
2     2                7311                        3646.259    0.5395833
3     3                8375                        2155.886    0.6470552
  Elongation Orientation Centroid_x Centroid_y AxesLength_x AxesLength_y
1   1.227318    1.565385   132.5324   131.4033     137.0855     168.2476
2   1.187744    1.497297   124.8014   125.7078     148.7605     176.6894
3   1.311573    1.447841   124.8623   131.8048     123.3535     161.7871
  BoundingBoxLower_x BoundingBoxUpper_x BoundingBoxLower_y BoundingBoxUpper_y
1                 57                195                 44                212
2                 58                194                 45                211
3                 60                192                 49                209
  IntegratedIntensity WeightedCentroid_x WeightedCentroid_y
1                2331           132.5324           131.4033
2               14622           124.8014           125.7078
3               25125           124.8623           131.8048

Histogram matching

The following examples demonstrate a variant of the well-known Nyul histogram standardization algorithm that permits masking and uses a B-spline approximation technique for histogram fitting. For the ITK implementation, use histogram_matching(...) in ANTsPy or histogramMatching(...) in ANTsR.

Python

>>> import ants
>>> from matplotlib import pyplot
>>> 
>>> source_image = ants.image_read(ants.get_ants_data("r16"))
>>> source_mask = ants.get_mask(source_image)
>>> 
>>> reference_image = ants.image_read(ants.get_ants_data("r64"))
>>> reference_mask = ants.get_mask(reference_image)
>>> 
>>> tx_source_image = ants.histogram_match_image2(source_image, reference_image, 
...                                               source_mask, reference_mask,
...                                               match_points=64,
...                                               transform_domain_size=255)
>>>
>>> source = source_image[source_mask != 0].flatten()
>>> reference = reference_image[reference_mask != 0].flatten()
>>> tx_source = tx_source_image[source_mask != 0].flatten()
>>>
>>> pyplot.hist([source, reference, tx_source], 64, alpha=0.5, label=["Source", "Ref", "Tx"])
>>> pyplot.legend(loc="upper right")
>>> pyplot.show()

R

> library( ANTsR )
> library( ggplot2 )
> 
> sourceImage <- antsImageRead( getANTsRData( "r16" ) )
> referenceImage <- antsImageRead( getANTsRData( "r64" ) )
> sourceMask <- getMask( sourceImage )
> referenceMask <- getMask( referenceImage )
> 
> txSourceImage <- histogramMatchImage2( txSourceImage, referenceImage, 
>                                        sourceMask, referenceMask,
>                                        matchPoints = 64, 
>                                        transformDomainSize = 255 )
> 
> sourceHist <- as.vector( sourceImage[sourceMask != 0] )
> referenceHist <- as.vector( referenceImage[referenceMask != 0] )
> txSourceHist <- as.vector( txSourceImage[sourceMask != 0] )
>
> histData <- data.frame( values = c( sourceHist, referenceHist, txSourceHist ),
>                         imageType = c( rep( "Src", length( sourceHist ) ),  
>                                        rep( "Ref", length( referenceHist ) ),   
>                                        rep( "Tx", length( txSourceHist ) ) ) 
>                       )
>
> histPlot <- ggplot( histData, aes( x = values, color = imageType ) ) +
>                     geom_density( alpha = 0.5, position = "identity")
> histPlot

Raw

b_5_ants_interfaces_deformation_based_morphometry.md

ANTsR / ANTsPy

Deformation-based morphometry

Python

>>> import ants
>>> import numpy as np
>>> import pandas as pd
>>> import random
>>> import statsmodels
>>> 
>>> # Build template
>>> r16 = ants.image_read(ants.get_ants_data('r16'))
>>> r27 = ants.image_read(ants.get_ants_data('r27'))
>>> r30 = ants.image_read(ants.get_ants_data('r30'))
>>> r62 = ants.image_read(ants.get_ants_data('r62'))
>>> r64 = ants.image_read(ants.get_ants_data('r64'))
>>> r85 = ants.image_read(ants.get_ants_data('r85'))
>>> 
>>> rimage_list = [r16, r27, r30, r62, r64, r85]
>>> number_of_subjects = len(rimage_list)
>>> rtemplate = ants.build_template(initial_template=None, image_list=rimage_list, type_of_transform="antsRegistrationSyNQuick[s]")
>>> 
>>> rtemplate_seg = ants.kmeans_segmentation(rtemplate, 3)['segmentation']
>>> rtemplate_gm = ants.threshold_image(rtemplate_seg, 2, 2, 1, 0)
>>> 
>>> # Generate log jacobian images and get brain volumes for later use
>>> log_jacobian_image_list = list()
>>> brain_volume = list()
>>> for i in range(number_of_subjects):
>>>     print("Registration subject", i, "to template.")
>>>     reg = ants.registration(fixed=rtemplate, moving=rimage_list[i], type_of_transform="antsRegistrationSyNQuick[b]")
>>>     log_jacobian_image_list.append(ants.create_jacobian_determinant_image(rtemplate, reg['fwdtransforms'][0], do_log=True))     
>>>     geom = ants.label_geometry_measures(ants.get_mask(rimage_list[i]))
>>>     brain_volume.append(geom['VolumeInMillimeters'][0])
>>> 
>>> # Create matrix to facilitate statistical analysis
>>> log_jacobian = ants.image_list_to_matrix(log_jacobian_image_list, rtemplate_gm)
>>> 
>>> # Create simulated outcomes and covariates
>>> neuro_condition = np.zeros((number_of_subjects,))
>>> age = np.zeros((number_of_subjects,))
>>> sex = np.zeros((number_of_subjects,))
>>> for i in range(number_of_subjects):
>>>     neuro_condition[i] = random.random()
>>>     age[i] = random.randint(40, 70)
>>>     sex[i] = random.randint(0, 1)
>>> data = {'neuro_condition' : neuro_condition, 'age' : age, 'sex' : sex, 'brain_volume' : brain_volume}
>>> covariates = pd.DataFrame(data)
>>> rformula = "neuro_condition ~ age + sex + brain_volume + log_jacobian"
>>> 
>>> # Run image-based regression
>>> dbm = ants.ilr(covariates, {"log_jacobian" : log_jacobian}, rformula)
>>> 
>>> # FDR correction
>>> log_jacobian_p_values = dbm['pValues']['pval_log_jacobian']
>>> log_jacobian_q_values = statsmodels.stats.multitest.fdrcorrection(log_jacobian_p_values, alpha=0.05, method='poscorr', is_sorted=False)[1]
>>> 
>>> # Create overlay q-value image
>>> log_jacobian_q_values_image = ants.matrix_to_images(np.reshape(log_jacobian_q_values, (1, len(log_jacobian_q_values))), rtemplate_gm)

R

> library( ANTsR )
> 
> # Build template
> r16 <- antsImageRead( getANTsRData( 'r16' ) )
> r27 <- antsImageRead( getANTsRData( 'r27' ) )
> r30 <- antsImageRead( getANTsRData( 'r30' ) )
> r62 <- antsImageRead( getANTsRData( 'r62' ) )
> r64 <- antsImageRead( getANTsRData( 'r64' ) )
> r85 <- antsImageRead( getANTsRData( 'r85' ) )
> 
> rimageList <- list( r16, r27, r30, r62, r64, r85 )
> numberOfSubjects <- length( rimageList )
> initialAverage <- antsAverageImages( rimageList )
> rtemplate <- buildTemplate( initialTemplate = initialAverage, imgList = rimageList, typeofTransform = "antsRegistrationSyNQuick[s]", verbose = TRUE )
> 
> rtemplateSeg <- thresholdImage( rtemplate, "Kmeans", 3 )
> rtemplateGM <- thresholdImage( rtemplateSeg, 2, 2, 1, 0 )
> 
> # Generate log jacobian images and get brain volumes for later use
> logJacobianImageList <- list()
> brainVolume <- c()
> for( i in seq.int( numberOfSubjects ) )
>   {
>   cat( "Registration subject", i, "to template.\n")
>   reg <- antsRegistration( fixed = rtemplate, moving = rimageList[[i]], typeofTransform = "antsRegistrationSyNQuick[b]" )
>   logJacobianImageList[[i]] <- createJacobianDeterminantImage( rtemplate, reg$fwdtransforms[[1]], doLog = TRUE )
> 
>   geom <- labelGeometryMeasures( getMask( rimageList[[i]] ) )
>   brainVolume[i] <- geom$VolumeInMillimeters[1] 
>   }
> 
> # Create matrix to facilitate statistical analysis
> logJacobian <- imageListToMatrix( logJacobianImageList, rtemplateGM )
> 
> # Create simulated outcomes and covariates
> neuroCondition <- rnorm( numberOfSubjects )
> age <- sample( 40:70, numberOfSubjects, replace = TRUE )
> sex <- sample( c( 0, 1 ), numberOfSubjects, replace = TRUE )
> covariates <- data.frame( neuroCondition = neuroCondition, age = age, sex = sex, brainVolume = brainVolume)
> rformula <- "neuroCondition ~ age + sex + brainVolume + logJacobian"
> 
> # Run image-based regression
> dbm <- ilr( covariates, list( logJacobian = logJacobian ), rformula )
>
# FDR correction
> logJacobianPvalues <- dbm$pValue[5,]
> logJacobianQvalues <- p.adjust( logJacobianPvalues, method = "fdr" )
> 
> # Create overlay q-value image
> logJacobianQvaluesImage <- matrixToImages( matrix( data = logJacobianQvalues, nrow = 1 ), rtemplateGM )

Raw

c_0_antsxnet_basic.md

ANTsXNet

Installation

Python: git clone https://github.com/ANTsX/ANTsPyNet; cd ANTsPyNet; python3 setup.py install
R: devtools::install_github( "ANTsX/ANTsRNet" )

A couple important notes

Primarily geared towards applications
No need for a GPU

Raw

c_1_antsxnet_applications_1_brain.md

ANTsXNet

Brain applications

Brain extraction

T1 (multiple "flavors")
T2
T2 star
T1/T2 infant
mean bold
FA

Python

>>> import ants
>>> import antspynet
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> 
>>> # ANTs-flavored
>>> seg = antspynet.brain_extraction(t1, modality="t1", verbose=True)
>>> ants.plot(t1, overlay=seg, overlay_alpha=0.5)
>>>
>>> # FreeSurfer-flavored
>>> seg = antspynet.brain_extraction(t1, modality="t1nobrainer", verbose=True)
>>> ants.plot(t1, overlay=seg, overlay_alpha=0.5)
>>>
>>> # Combined
>>> seg = antspynet.brain_extraction(t1, modality="t1combined", verbose=True)
>>> ants.plot(t1, overlay=seg, overlay_alpha=0.5)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1 <- antsImageRead( getANTsXNetData( 'mprage_hippmapp3r' ) )
> 
> # ANTs-flavored
> seg <- brainExtraction( t1, modality = "t1", verbose = TRUE )
> plot( t1, seg, alpha = 0.5 )
>
> # FreeSurfer-flavored
> seg <- brainExtraction( t1, modality = "t1nobrainer", verbose = TRUE )
> plot( t1, seg, alpha = 0.5 )
>
> # Combined
> seg <- brainExtraction( t1, modality = "t1combined", verbose = TRUE )
> plot( t1, seg, alpha = 0.5 )

Deep Atropos

Python

>>> import ants
>>> import antspynet
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> seg = antspynet.deep_atropos(t1, verbose=True)
>>> ants.plot(t1, overlay=seg['segmentation_image'], overlay_alpha=0.75)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1 <- antsImageRead( getANTsXNetData( 'mprage_hippmapp3r' ) )
> seg <- deepAtropos( t1, verbose = TRUE )
> plot( t1, seg$segmentationImage, alpha = 0.75 )

Cortical thickness

Python

>>> import ants
>>> import antspynet
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> kk = antspynet.cortical_thickness(t1, verbose=True)
>>> ants.plot(t1, overlay=kk['thickness_image'], overlay_alpha=0.75)
>>>
>>> # Also see antspynet.longitudinal_cortical_thickness(...)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1 <- antsImageRead( getANTsXNetData( 'mprage_hippmapp3r' ) )
> kk <- corticalThickness( t1, verbose = TRUE )
> plot( t1, kk$thicknessImage, alpha = 0.75 )
>
> # Also see longitudinalCorticalThickness(...)

Desikan-Killiany-Tourville parcellation

Python

>>> import ants
>>> import antspynet
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> dkt = antspynet.desikan_killiany_tourville_labeling(t1, do_lobar_parcellation=True, verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1 <- antsImageRead( getANTsXNetData( 'mprage_hippmapp3r' ) )
> dkt <- desikanKillianyTourvilleLabeling( t1, doLobarParcellation = TRUE, verbose = TRUE )

Harvard-Oxford Atlas subcortical parcellation

Python

>>> import ants
>>> import antspynet
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('adni'))
>>> hoa = antspynet.harvard_oxford_atlas_labeling(t1, verbose=True)

DeepFLASH (medial temporal lobe parcellation)

T1
T1/T2

Python

>>> import ants
>>> import antspynet
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> df = antspynet.deep_flash(t1, verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1 <- antsImageRead( getANTsXNetData( 'mprage_hippmapp3r' ) )
> df <- deepFlash( t1, verbose = TRUE )

Hippmapp3r (hippocampal segmentation)

Python

>>> import ants
>>> import antspynet
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> hipp = antspynet.hippmapp3r_segmentation(t1, verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1 <- antsImageRead( getANTsXNetData( 'mprage_hippmapp3r' ) )
> hipp <- hippMapp3rSegmentation( t1, verbose = TRUE )

Brain Age

Python

>>> import ants
>>> import antspynet
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> age = antspynet.brain_age(t1, number_of_simulations=3, sd_affine=0.01, verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1 <- antsImageRead( getANTsXNetData( 'mprage_hippmapp3r' ) )
> age <- brainAge( t1, numberOfSimulations = 3, sdAffine = 0.01, verbose = TRUE )

Claustrum segmentation

Python

>>> import ants
>>> import antspynet
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> seg = antspynet.claustrum_segmentation(t1, verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1 <- antsImageRead( getANTsXNetData( 'mprage_hippmapp3r' ) )
> seg <- claustrumSegmentation( t1, verbose = TRUE )

Hypothalamus segmentation

Python

>>> import ants
>>> import antspynet
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('kirby'))
>>> seg = antspynet.hypothalamus_segmentation(t1, verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1 <- antsImageRead( getANTsXNetData( 'kirby' ) )
> seg <- hypothalamusSegmentation( t1, verbose = TRUE )

White matter hyperintensities segmentation (SYSU)

Python

>>> import ants
>>> import antspynet
>>> import tensorflow as tf
>>>
>>> t1_file = tf.keras.utils.get_file(fname="t1.nii.gz", origin="https://figshare.com/ndownloader/files/40251796", force_download=True)
>>> t1 = ants.image_read(t1_file)
>>> flair_file = tf.keras.utils.get_file(fname="flair.nii.gz", origin="https://figshare.com/ndownloader/files/40251793", force_download=True)
>>> flair = ants.image_read(flair_file)
>>>
>>> wmh = antspynet.sysu_media_wmh_segmentation(flair, t1, verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1File <- tensorflow::tf$keras$utils$get_file( fname = "t1.nii.gz", origin = "https://figshare.com/ndownloader/files/40251796", force_download = TRUE )
> t1 <- antsImageRead( t1File )
> flairFile <- tensorflow::tf$keras$utils$get_file( fname = "flair.nii.gz", origin = "https://figshare.com/ndownloader/files/40251793", force_download = TRUE )
> flair <- antsImageRead( flairFile )
> 
> wmh <- sysuMediaWmhSegmentation( flair, t1, verbose = TRUE )

White matter hyperintensities segmentation (Hypermapp3r)

Python

>>> import ants
>>> import antspynet
>>> import tensorflow as tf
>>>
>>> t1_file = tf.keras.utils.get_file(fname="t1.nii.gz", origin="https://figshare.com/ndownloader/files/40251796", force_download=True)
>>> t1 = ants.image_read(t1_file)
>>> flair_file = tf.keras.utils.get_file(fname="flair.nii.gz", origin="https://figshare.com/ndownloader/files/40251793", force_download=True)
>>> flair = ants.image_read(flair_file)
>>>
>>> wmh = antspynet.hypermapp3r_segmentation(t1, flair, verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1File <- tensorflow::tf$keras$utils$get_file( fname = "t1.nii.gz", origin = "https://figshare.com/ndownloader/files/40251796", force_download = TRUE )
> t1 <- antsImageRead( t1File )
> flairFile <- tensorflow::tf$keras$utils$get_file( fname = "flair.nii.gz", origin = "https://figshare.com/ndownloader/files/40251793", force_download = TRUE )
> flair <- antsImageRead( flairFile )
> 
> wmh <- hyperMapp3rSegmentation( t1, flair, verbose = TRUE )

White matter hyperintensities segmentation (SHIVA)

Python

>>> import ants
>>> import antspynet
>>> import tensorflow as tf
>>>
>>> t1_file = tf.keras.utils.get_file(fname="t1.nii.gz", origin="https://figshare.com/ndownloader/files/40251796", force_download=True)
>>> t1 = ants.image_read(t1_file)
>>> flair_file = tf.keras.utils.get_file(fname="flair.nii.gz", origin="https://figshare.com/ndownloader/files/40251793", force_download=True)
>>> flair = ants.image_read(flair_file)
>>>
>>> wmh = antspynet.shiva_wmh_segmentation(flair, t1, which_model = "all", verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1File <- tensorflow::tf$keras$utils$get_file( fname = "t1.nii.gz", origin = "https://figshare.com/ndownloader/files/40251796", force_download = TRUE )
> t1 <- antsImageRead( t1File )
> flairFile <- tensorflow::tf$keras$utils$get_file( fname = "flair.nii.gz", origin = "https://figshare.com/ndownloader/files/40251793", force_download = TRUE )
> flair <- antsImageRead( flairFile )
> 
> wmh <- shivaWmhSegmentation( flair, t1, whichModel = "all", verbose = TRUE )

White matter hyperintensities segmentation (ANTsXNet)

Python

>>> import ants
>>> import antspynet
>>> import tensorflow as tf
>>>
>>> t1_file = tf.keras.utils.get_file(fname="t1.nii.gz", origin="https://figshare.com/ndownloader/files/40251796", force_download=True)
>>> t1 = ants.image_read(t1_file)
>>> t1 = ants.resample_image(t1, (240, 240, 64), use_voxels=True)
>>> flair_file = tf.keras.utils.get_file(fname="flair.nii.gz", origin="https://figshare.com/ndownloader/files/40251793", force_download=True)
>>> flair = ants.image_read(flair_file)
>>> flair = ants.resample_image(flair, (240, 240, 64), use_voxels=True)
>>>
>>> wmh = antspynet.wmh_segmentation(flair, t1, use_combined_model=True, verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1File <- tensorflow::tf$keras$utils$get_file( fname = "t1.nii.gz", origin = "https://figshare.com/ndownloader/files/40251796", force_download = TRUE )
> t1 <- antsImageRead( t1File )
> t1 <- resampleImage( t1, c( 240, 240, 64 ), useVoxels = TRUE )
> flairFile <- tensorflow::tf$keras$utils$get_file( fname = "flair.nii.gz", origin = "https://figshare.com/ndownloader/files/40251793", force_download = TRUE )
> flair <- antsImageRead( flairFile )
> flair <- resampleImage( flair, c( 240, 240, 64 ), useVoxels = TRUE )
> 
> wmh <- wmhSegmentation( t1, flair, useCombinedModel = TRUE, verbose = TRUE )

Perivascular spaces segmentation (SHIVA)

Python

>>> import ants
>>> import antspynet
>>> import tensorflow as tf
>>>
>>> t1_file = tf.keras.utils.get_file(fname="pvs_t1.nii.gz", origin="https://figshare.com/ndownloader/files/48675367", force_download=True)
>>> t1 = ants.image_read(t1_file)
>>> flair_file = tf.keras.utils.get_file(fname="pvs_flair.nii.gz", origin="https://figshare.com/ndownloader/files/48675352", force_download=True)
>>> flair = ants.image_read(flair_file)
>>>
>>> pvs = antspynet.shiva_pvs_segmentation(t1, flair, which_model = "all", verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1File <- tensorflow::tf$keras$utils$get_file( fname = "pvs_t1.nii.gz", origin = "https://figshare.com/ndownloader/files/48675367", force_download = TRUE )
> t1 <- antsImageRead( t1File )
> flairFile <- tensorflow::tf$keras$utils$get_file( fname = "pvs_flair.nii.gz", origin = "https://figshare.com/ndownloader/files/48675352", force_download = TRUE )
> flair <- antsImageRead( flairFile )
> 
> pvs <- shivaPvsSegmentation( t1, flair, whichModel = "all", verbose = TRUE )

Cerebellum morphology

Python

>>> import ants
>>> import antspynet
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data("mprage_hippmapp3r"))
>>>
>>> # Computing the thickness image will take relatively more time than the other	
>>> # outputs so the recommendation would be to run it initially without computing the	
>>> # thickness image to just get a sense of the application.	
>>>
>>> cereb = antspynet.cerebellum_morphology(t1, compute_thickness_image=True, verbose=True)	
>>> 	
>>> # possible refinement with cerebellum estimate (comment out since it's not needed for this image).	
>>> # mask = ants.threshold_image(cereb['cerebellum_probability_image'], 0.5, 1, 1, 0)	
>>> # cereb = antspynet.cerebellum_morphology(t1, cerebellum_mask=mask, verbose=True)	
>>> 
>>> ants.image_write(cereb['cerebellum_probability_image'], "cerebellum_probability_mask.nii.gz")	
>>> ants.image_write(cereb['thickness_image'], "kk.nii.gz")	
>>> ants.image_write(cereb['parcellation_segmentation_image'], "parcelation.nii.gz")

R

> library( ANTsR )	
> library( ANTsRNet )	
>	
> t1 <- antsImageRead( getANTsXNetData( "mprage_hippmapp3r" ) )	
>	
> # Computing the thickness image will take relatively more time than the other	
> # outputs so the recommendation would be to run it initially without computing the	
> # thickness image to just get a sense of the application.	
> cereb <- cerebellumMorphology( t1, computeThicknessImage = TRUE, verbose = TRUE )	
>	
> # possible refinement with cerebellum estimate (comment out since it's not needed for this image).	
> # mask <- thresholdImage( cereb$cerebellumProbabilityImage, 0.5, 1, 1, 0 )	
> # cereb <- cerebellumMorphology( t1, cerebellumMask = mask, verbose = TRUE )	
> 	
> # Write output to disk	
> antsImageWrite( cereb$cerebellumProbabilityImage, "cerebellumProbabilityMask.nii.gz" )	
> antsImageWrite( cereb$tissueSegmentationImage, "tissue.nii.gz" )	
> antsImageWrite( cereb$thicknessImage, "kk.nii.gz" )	
> antsImageWrite( cereb$parcellationSegmentationImage, "parcelation.nii.gz" )

Brain tumor segmentation

Python

>>> import ants
>>> import antspynet
>>> import tensorflow as tf
>>>
>>> flair_file = tf.keras.utils.get_file(fname="flair.nii.gz", origin="https://figshare.com/ndownloader/files/42385077", force_download=True)
>>> flair = ants.image_read(flair_file)
>>> t1_file = tf.keras.utils.get_file(fname="t1.nii.gz", origin="https://figshare.com/ndownloader/files/42385071", force_download=True)
>>> t1 = ants.image_read(t1_file)
>>> t1_contrast_file = tf.keras.utils.get_file(fname="t1_contrast.nii.gz", origin="https://figshare.com/ndownloader/files/42385068", force_download=True)
>>> t1_contrast = ants.image_read(t1_contrast_file)
>>> t2_file = tf.keras.utils.get_file(fname="t2.nii.gz", origin="https://figshare.com/ndownloader/files/42385074", force_download=True)
>>> t2 = ants.image_read(t2_file)
>>>
>>> bt = antspynet.brain_tumor_segmentation(flair, t1, t1_contrast, t2, patch_stride_length=32, verbose=True)
>>> # ants.image_write(bt['segmentation_image'], "brain_tumor_segmentation.nii.gz")

R

> library( ANTsR )
> library( ANTsRNet )
>
> flairFile <- tensorflow::tf$keras$utils$get_file( fname = "flair.nii.gz", origin = "https://figshare.com/ndownloader/files/42385077", force_download = TRUE )
> flair <- antsImageRead( flairFile )
> t1File <- tensorflow::tf$keras$utils$get_file( fname = "t1.nii.gz", origin = "https://figshare.com/ndownloader/files/42385071", force_download = TRUE )
> t1 <- antsImageRead( t1File )
> t1ContrastFile <- tensorflow::tf$keras$utils$get_file( fname = "t1_contrast.nii.gz", origin = "https://figshare.com/ndownloader/files/42385068", force_download = TRUE )
> t1Contrast <- antsImageRead( t1ContrastFile )
> t2File <- tensorflow::tf$keras$utils$get_file( fname = "t2.nii.gz", origin = "https://figshare.com/ndownloader/files/42385074", force_download = TRUE )
> t2 <- antsImageRead( t2File )
> 
> bt <- brainTumorSegmentation( flair, t1, t1Contrast, t2, patchStrideLength = 32, verbose = TRUE )
> #antsImageWrite( bt$segmentationImage, "brainTumorSegmentation.nii.gz" )

MRA-TOF vessel segmentation

Python

>>> import ants
>>> import antspynet
>>> import tensorflow as tf
>>>
>>> mmbop_file = tf.keras.utils.get_file(fname="mra.nii.gz", origin="https://figshare.com/ndownloader/files/46406755", force_download=True)
>>> mra = ants.image_read(mra_file)
>>> vessels = antspynet.mra_brain_vessel_segmentation(mra, verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> mmbopFile <- tensorflow::tf$keras$utils$get_file( fname = "mra.nii.gz", origin = "https://figshare.com/ndownloader/files/46406755", force_download = TRUE )
> mra <- antsImageRead( mmbopFile )
> vessels <- mraBrainVesselSegmentation( mra, verbose = TRUE )

Lesion segmentation (WIP)

Python

>>> import ants
>>> import antspynet
>>> import tensorflow as tf
>>>
>>> t1_file = tf.keras.utils.get_file(fname="t1w_with_lesion.nii.gz", origin="https://figshare.com/ndownloader/files/44053868", force_download=True)
>>> t1 = ants.image_read(t1_file)
>>>
>>> probability_mask = antspynet.lesion_segmentation(t1, do_preprocessing=True, verbose=True)
>>> ants.image_write(probability_mask, "lesion_probability_mask.nii.gz")

R

> library( ANTsR )
> library( ANTsRNet )
> library( tensorflow )
>
> t1File <- tensorflow::tf$keras$utils$get_file( fname = "t1w_with_lesion.nii.gz", origin = "https://figshare.com/ndownloader/files/44053868", force_download = TRUE )
> t1 <- antsImageRead( t1File )
>
> probabilityMask <- lesionSegmentation( t1, doPreprocessing = TRUE, verbose = TRUE )
> antsImageWrite( probabilityMask, "lesion_probability_mask.nii.gz" )

Whole head inpainting

Python

>>> import ants
>>> import antspynet
>>> import tensorflow as tf
>>>
>>> t1_file = tf.keras.utils.get_file(fname="t1w_with_lesion.nii.gz", origin="https://figshare.com/ndownloader/files/44053868", force_download=True)
>>> t1 = ants.image_read(t1_file)
>>>
>>> probability_mask = antspynet.lesion_segmentation(t1, do_preprocessing=True, verbose=True)
>>> lesion_mask = ants.threshold_image(probability_mask, 0.5, 1.1, 1, 0)
>>> t1_inpainted = antspynet.whole_head_inpainting(t1, roi_mask=lesion_mask, modality="t1", mode="axial", verbose=True)
>>> ants.image_write(t1_inpainted, "t1_repaired.nii.gz")

R

> library( ANTsR )
> library( ANTsRNet )
> library( tensorflow )
>
> t1File <- tensorflow::tf$keras$utils$get_file( fname = "t1w_with_lesion.nii.gz", origin = "https://figshare.com/ndownloader/files/44053868", force_download = TRUE )
> t1 <- antsImageRead( t1File )
>
> probabilityMask <- lesionSegmentation( t1, doPreprocessing = TRUE, verbose = TRUE )
> lesionMask <- thresholdImage( probabilityMask, 0.5, 1.1, 1, 0 )
> t1Inpainted <- wholeHeadInpainting( t1, roiMask = lesionMask, modality = "t1", mode = "axial", verbose = TRUE )
> antsImageWrite( t1Inpainted, "t1_repaired.nii.gz" )

Raw

c_1_antsxnet_applications_2_lung.md

Lung applications

Lung extraction

CT (left, right, and airways)
Proton (left, right)
Ventilation-based whole lung segmentation
X-ray (left, right)

Python

>>> import ants
>>> import antspynet
>>> import tensorflow as tf
>>>
>>> # CT
>>> ct_file = tf.keras.utils.get_file(fname="ctLung.nii.gz", origin="https://figshare.com/ndownloader/files/42934234", force_download=True)
>>> ct = ants.image_read(ct_file)
>>> lung_ex = antspynet.lung_extraction(ct, modality="ct", verbose=True)
>>> ants.plot(ants.iMath_normalize(ct), lung_ex['segmentation_image'], alpha=0.5)
>>>
>>> # Proton (whole lung)
>>> proton_file = tf.keras.utils.get_file(fname="protonLung.nii.gz", origin="https://figshare.com/ndownloader/files/42934228", force_download=True)
>>> proton = ants.image_read(proton_file)
>>> lung_ex = antspynet.lung_extraction(proton, modality="proton", verbose=True)
>>> ants.plot(ants.iMath_normalize(proton), lung_ex['segmentation_image'], alpha=0.5)
>>>
>>> # Proton (lobes)
>>> lung_ex = antspynet.lung_extraction(proton, modality="protonLobes", verbose=True)
>>> ants.plot(ants.iMath_normalize(proton), lung_ex['segmentation_image'], alpha=0.5)
>>>
>>> # Lobes from lung mask
>>> lung_mask = ants.threshold_image(lung_ex['segmentation_image'], 0, 0, 0, 1 )
>>> lung_ex = antspynet.lung_extraction(lung_mask, modality="maskLobes", verbose=True)
>>> ants.plot(ants.iMath_normalize(proton), lung_ex['segmentation_image'], alpha=0.5)
>>>
>>> # Chest x-ray
>>> cxr_file = tf.keras.utils.get_file(fname="cxr.nii.gz", origin="https://figshare.com/ndownloader/files/42934237", force_download=True)
>>> cxr = ants.image_read(cxr_file)
>>> lung_ex = antspynet.lung_extraction(cxr, modality="xray", verbose=True)
>>> ants.plot(ants.iMath_normalize(cxr), lung_ex['segmentation_image'], alpha=0.5)
>>>
>>> # Ventilation
>>> ventilation_file = tf.keras.utils.get_file(fname="ventilationLung.nii.gz", origin="https://figshare.com/ndownloader/files/42934231", force_download=True)
>>> ventilation = ants.image_read(ventilation_file)
>>> lung_ex = antspynet.lung_extraction(ventilation, modality="ventilation", verbose=True)
>>> ants.plot(ants.iMath_normalize(ventilation), ants.threshold_image(lung_ex, 0.5, 1, 1, 0), alpha=0.5)

R

> library( ANTsR )
> library( ANTsRNet )
>
> # CT
> ctFile <- tensorflow::tf$keras$utils$get_file( fname = "ctLung.nii.gz", origin = "https://figshare.com/ndownloader/files/42934234", force_download = TRUE )
> ct <- antsImageRead( ctFile )
> lungEx <- lungExtraction( ct, modality = "ct", verbose = TRUE )
> plot( ct, lungEx$segmentationImage, alpha = 0.5 )
>
> # Proton (whole lung)
> protonFile <- tensorflow::tf$keras$utils$get_file( fname = "protonLung.nii.gz", origin = "https://figshare.com/ndownloader/files/42934228", force_download = TRUE )
> proton <- antsImageRead( protonFile )
> lungEx <- lungExtraction( proton, modality = "proton", verbose = TRUE )
> plot( proton, lungEx$segmentationImage, alpha = 0.5 )
>
> # Proton (lobes)
> lungEx <- lungExtraction( proton, modality = "protonLobes", verbose = TRUE )
> plot( proton, lungEx$segmentationImage, alpha = 0.5 )
>
> # Lobes from lung mask
> lungMask <- thresholdImage( lungEx$segmentationImage, 0, 0, 0, 1 )
> lungEx <- lungExtraction( lungMask, modality = "maskLobes", verbose = TRUE )
> plot( proton, lungEx$segmentationImage, alpha = 0.5 )
>
> # Chest x-ray
> cxrFile <- tensorflow::tf$keras$utils$get_file( fname = "cxr.nii.gz", origin = "https://figshare.com/ndownloader/files/42934237", force_download = TRUE )
> cxr <- antsImageRead( cxrFile )
> lungEx <- lungExtraction( cxr, modality = "xray", verbose = TRUE )
> plot( cxr, lungEx$segmentationImage, alpha = 0.5 )
>
> # Ventilation
> ventilationFile <- tensorflow::tf$keras$utils$get_file( fname = "ventilationLung.nii.gz", origin = "https://figshare.com/ndownloader/files/42934231", force_download = TRUE )
> ventilation <- antsImageRead( ventilationFile )
> lungEx <- lungExtraction( ventilation, modality = "ventilation", verbose = TRUE )
> plot( ventilation, thresholdImage( lungEx, 0.5, 1, 1, 0 ), alpha = 0.5 )

Functional lung segmentation

Python

>>> import ants
>>> import antspynet
>>> import tensorflow as tf
>>>
>>> # Process the proton image to get the lung mask
>>> proton_file = tf.keras.utils.get_file(fname="protonLung.nii.gz", origin="https://figshare.com/ndownloader/files/42934228", force_download=True)
>>> proton = ants.image_read(proton_file)
>>> lung_ex = antspynet.lung_extraction(proton, modality="proton", verbose=True)
>>> lung_mask = ants.threshold_image(lung_ex['segmentation_image'], 0, 0, 0, 1 )
>>>
>>> # Use 'el bicho' to label the four ventilation levels
>>> ventilation_file = tf.keras.utils.get_file(fname="ventilationLung.nii.gz", origin="https://figshare.com/ndownloader/files/42934231", force_download=True)
>>> ventilation = ants.image_read(ventilation_file)
>>> eb = antspynet.el_bicho(ventilation, lung_mask, verbose=True)
>>> ants.plot(ants.iMath_normalize(ventilation), eb['segmentation_image'], alpha=0.5)

R

> library( ANTsR )
> library( ANTsRNet )
>
> # Process the proton image to get the lung mask
> protonFile <- tensorflow::tf$keras$utils$get_file( fname = "protonLung.nii.gz", origin = "https://figshare.com/ndownloader/files/42934228", force_download = TRUE )
> proton <- antsImageRead( protonFile )
> lungEx <- lungExtraction( proton, modality = "proton", verbose = TRUE )
> lungMask <- thresholdImage( lungEx$segmentationImage, 0, 0, 0, 1 )
>
> # Use 'el bicho' to label the four ventilation levels
> ventilationFile <- tensorflow::tf$keras$utils$get_file( fname = "ventilationLung.nii.gz", origin = "https://figshare.com/ndownloader/files/42934231", force_download = TRUE )
> ventilation <- antsImageRead( ventilationFile )
> eb <- elBicho( ventilation, lungMask, verbose = TRUE )
> plot( ventilation, eb$segmentationImage, alpha = 0.5 )

Pulmonary artery segmentation

Python

>>> import ants
>>> import antspynet
>>> import tensorflow as tf
>>>
>>> ct_file = tf.keras.utils.get_file(fname="ctLung.nii.gz", origin="https://figshare.com/ndownloader/files/42934234", force_download=True)
>>> ct = ants.image_read(ct_file)
>>> arteries = antspynet.lung_pulmonary_artery_segmentation(ct, verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> ctFile <- tensorflow::tf$keras$utils$get_file( fname = "ctLung.nii.gz", origin = "https://figshare.com/ndownloader/files/42934234", force_download = TRUE )
> ct <- antsImageRead( ctFile )
> arteries <- lungPulmonaryArterySegmentation( ct, verbose = TRUE )

Pulmonary airway segmentation

Python

>>> import ants
>>> import antspynet
>>> import tensorflow as tf
>>>
>>> ct_file = tf.keras.utils.get_file(fname="ctLung.nii.gz", origin="https://figshare.com/ndownloader/files/42934234", force_download=True)
>>> ct = ants.image_read(ct_file)
>>> airways = antspynet.lung_airway_segmentation(ct, verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> ctFile <- tensorflow::tf$keras$utils$get_file( fname = "ctLung.nii.gz", origin = "https://figshare.com/ndownloader/files/42934234", force_download = TRUE )
> ct <- antsImageRead( ctFile )
> airways <- lungAirwaySegmentation( ct, verbose = TRUE )

CheXNet

Python

>>> import ants
>>> import antspynet
>>> import deepsimlr
>>>
>>> # This particular subject is diagnosed as having cardiomegaly 
>>> cxr_file = tf.keras.utils.get_file(fname="cxr.nii.gz", origin="https://figshare.com/ndownloader/files/42429156", force_download=True)
>>> cxr = ants.image_read(cxr_file)
>>> 
>>> # Pytorch version adapted from https://github.com/jrzech/reproduce-chexnet
>>> deepsimlr.chexnet(cxr)
   Atelectasis  Cardiomegaly  Effusion  Infiltration      Mass    Nodule  Pneumonia  Pneumothorax  Consolidation     Edema  Emphysema  Fibrosis  Pleural_Thickening    Hernia
0     0.027259      0.745355  0.042421      0.049225  0.003613  0.026083   0.005946      0.003007       0.005099  0.002063   0.017297  0.059956             0.04354  0.004266
>>>
>>> # Keras implementation in ANTsXNet (original)
>>> antspynet.chexnet(cxr, use_antsxnet_variant=False)
   Atelectasis  Cardiomegaly  Effusion  Infiltration      Mass    Nodule  Pneumonia  Pneumothorax  Consolidation     Edema  Emphysema  Fibrosis  Pleural_Thickening   Hernia
0     0.133284       0.64578  0.085462       0.09486  0.019405  0.018975   0.013067      0.021075       0.028001  0.005851   0.018955    0.0555            0.034219  0.00972
>>>
>>> # Keras implementation in ANTsXNet (variation)
>>> antspynet.chexnet(cxr, use_antsxnet_variant=True, include_tuberculosis_diagnosis=True)
   Atelectasis  Cardiomegaly  Effusion  Infiltration      Mass    Nodule  Pneumonia  Pneumothorax  Consolidation     Edema  Emphysema  Fibrosis  Pleural_Thickening    Hernia  Tuberculosis
0     0.069623      0.898534  0.039589      0.055751  0.005931  0.030309   0.006876      0.013999       0.010621  0.005677   0.048858  0.021668            0.013867  0.008546      0.000285

R

> library( ANTsR )
> library( ANTsRNet )
>
> cxrFile <- tensorflow::tf$keras$utils$get_file( fname = "cxr.nii.gz", origin = "https://figshare.com/ndownloader/files/42429156", force_download = TRUE )
> cxr <- antsImageRead( cxrFile )
> 
> # Keras implementation in ANTsXNet (original)
> chexnet( cxr, useANTsXNetVariant = FALSE )
  Atelectasis Cardiomegaly   Effusion Infiltration       Mass     Nodule
1   0.1065275    0.7038919 0.09204669    0.1006764 0.02943837 0.02722431
   Pneumonia Pneumothorax Consolidation       Edema  Emphysema   Fibrosis
1 0.01293034   0.01917528    0.02944043 0.006727949 0.01533213 0.06293546
  Pleural_Thickening     Hernia
1         0.03341614 0.01022097
> 
> # Keras implementation in ANTsXNet (variation)
> chexnet( cxr, useANTsXNetVariant = TRUE, includeTuberculosisDiagnosis = TRUE )
  Atelectasis Cardiomegaly  Effusion Infiltration       Mass     Nodule
1  0.07443799    0.9660769 0.1168589   0.08752631 0.01058979 0.03241358
   Pneumonia Pneumothorax Consolidation      Edema Emphysema  Fibrosis
1 0.01211402   0.02042338    0.01615181 0.01346725  0.052293 0.0357863
  Pleural_Thickening      Hernia Tuberculosis
1          0.0235217 0.008928899 0.0002852035
>

Evaluation on CheXNet testing data

Compare with original implementation

AUC values

                     Pytorch     Keras  ANTsXNet
Atelectasis         0.832871  0.822537  0.863367
Cardiomegaly        0.946268  0.938183  0.962786
Consolidation       0.816950  0.815728  0.855411
Edema               0.912380  0.915611  0.938205
Effusion            0.893323  0.891513  0.916355
Emphysema           0.906865  0.887093  0.920824
Fibrosis            0.834109  0.827323  0.868360
Hernia              0.882469  0.824022  0.919963
Infiltration        0.770469  0.769772  0.805857
Mass                0.875950  0.834911  0.893959
Nodule              0.788035  0.764380  0.823646
Pleural_Thickening  0.828046  0.814205  0.854000
Pneumonia           0.758576  0.751609  0.800324
Pneumothorax        0.880895  0.870271  0.910933

Raw

c_1_antsxnet_applications_3_mouse.md

Mouse applications

Mouse brain extraction

Python

>>> import ants
>>> import antspynet
>>>
>>> mouse_t2_file = tf.keras.utils.get_file(fname="mouse.nii.gz", origin="https://figshare.com/ndownloader/files/45289309", force_download=True)
>>> mouse_t2 = ants.image_read(mouse_t2_file)
>>> mouse_t2_n4 = ants.n4_bias_field_correction(mouse_t2, 
                                                rescale_intensities=True,
                                                shrink_factor=2, 
                                                convergence={'iters': [50, 50, 50, 50], 'tol': 0.0}, 
                                                spline_param=20, verbose=True)
>>> mask = antspynet.mouse_brain_extraction(mouse_t2_n4, modality='t2', verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> mouseT2File <- tensorflow::tf$keras$utils$get_file( fname="mouse.nii.gz",  origin = "https://figshare.com/ndownloader/files/45289309", force_download = TRUE )
> mouseT2 <- antsImageRead( mouseT2File )
> mouseT2N4 <- n4BiasFieldCorrection( mouseT2, 
                                      rescaleIntensities = TRUE,
                                      shrinkFactor = 2, 
                                      convergence = list( iters = c( 50, 50, 50, 50 ), tol = 0.0 ), 
                                      splineParam = 20, verbose = TRUE )
>>> mask <- mouseBrainExtraction( mouseT2N4, modality = 't2', verbose = TRUE )

Mouse brain parcellation

Python

>>> import ants
>>> import antspynet
>>>
>>> mouse_t2_file = tf.keras.utils.get_file(fname="mouse.nii.gz", origin="https://figshare.com/ndownloader/files/45289309", force_download=True)
>>> mouse_t2 = ants.image_read(mouse_t2_file)
>>> mouse_t2_n4 = ants.n4_bias_field_correction(mouse_t2, 
                                                rescale_intensities=True,
                                                shrink_factor=2, 
                                                convergence={'iters': [50, 50, 50, 50], 'tol': 0.0}, 
                                                spline_param=20, verbose=True)
>>> parc_nick = antspynet.mouse_brain_parcellation(mouse_t2_n4, 
                                                   mask=None, 
                                                   which_parcellation="nick",      
                                                   return_isotropic_output=True,  
                                                   verbose=True)
>>> parc_tct = antspynet.mouse_brain_parcellation(mouse_t2_n4, 
                                                  mask=None, 
                                                  which_parcellation="tct",      
                                                  return_isotropic_output=True,  
                                                  verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> mouseT2File <- tensorflow::tf$keras$utils$get_file( fname="mouse.nii.gz",  origin = "https://figshare.com/ndownloader/files/45289309", force_download = TRUE )
> mouseT2 <- antsImageRead( mouseT2File )
> mouseT2N4 <- n4BiasFieldCorrection( mouseT2, 
                                      rescaleIntensities = TRUE,
                                      shrinkFactor = 2, 
                                      convergence = list( iters = c( 50, 50, 50, 50 ), tol = 0.0 ), 
                                      splineParam = 20, verbose = TRUE )
> parcNick <- mouseBrainParcellation( mouseT2N4, 
                                      mask = NULL,
                                      whichParcellation = 'nick', 
                                      returnIsotropicOutput = TRUE,
                                      verbose = TRUE )
> parcTct <- mouseBrainParcellation( mouseT2N4, 
                                     mask = NULL,
                                     whichParcellation = 'tct', 
                                     returnIsotropicOutput = TRUE,
                                     verbose = TRUE )

Mouse cortical thickness

Python

>>> import ants
>>> import antspynet
>>>
>>> mouse_t2_file = tf.keras.utils.get_file(fname="mouse.nii.gz", origin="https://figshare.com/ndownloader/files/45289309", force_download=True)
>>> mouse_t2 = ants.image_read(mouse_t2_file)
>>> mouse_t2_n4 = ants.n4_bias_field_correction(mouse_t2, 
                                                rescale_intensities=True,
                                                shrink_factor=2, 
                                                convergence={'iters': [50, 50, 50, 50], 'tol': 0.0}, 
                                                spline_param=20, verbose=True)
>>> kk = antspynet.mouse_cortical_thickness(mouse_t2_n4, 
                                            mask=None, 
                                            return_isotropic_output=True,                                    
                                            verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> mouseT2File <- tensorflow::tf$keras$utils$get_file( fname="mouse.nii.gz",  origin = "https://figshare.com/ndownloader/files/45289309", force_download = TRUE )
> mouseT2 <- antsImageRead( mouseT2File )
> mouseT2N4 <- n4BiasFieldCorrection( mouseT2, 
                                      rescaleIntensities = TRUE,
                                      shrinkFactor = 2, 
                                      convergence = list( iters = c( 50, 50, 50, 50 ), tol = 0.0 ), 
                                      splineParam = 20, verbose = TRUE )
>>> kk <- mouseCorticalThickness( mouseT2N4, 
                                  mask = NULL,
                                  returnIsotropicOutput = TRUE,
                                  verbose = TRUE )

Raw

c_1_antsxnet_applications_4_general.md

General applications

MRI super resolution

Python

>>> import ants
>>> import antspynet
>>>
>>> # This method is deprecated.  Please see https://github.com/stnava/siq.
>>> 
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> t1_lr = ants.resample_image(t1, (4, 4, 4), use_voxels=False)
>>> ants.plot(t1_lr)
>>> t1_sr = antspynet.mri_super_resolution(t1_lr, verbose=True)
>>> ants.plot(t1_sr)

R

> library( ANTsR )
> library( ANTsRNet )
>
> # This method is deprecated.  Please see https://github.com/stnava/siq.
> 
> t1 <- antsImageRead( getANTsXNetData( 'mprage_hippmapp3r' ) )
> t1Lr <- resampleImage( t1, c( 4, 4, 4 ), useVoxels = FALSE )
> plot( t1Lr )
> t1Sr <- mriSuperResolution( t1Lr, verbose = TRUE )
> plot( t1Sr )

No reference image quality assesment using TID

Python

>>> import ants
>>> import antspynet
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> qa = antspynet.tid_neural_image_assessment(t1, verbose=True)

R

> library( ANTsR )
> library( ANTsRNet )
>
> t1 <- antsImageRead( getANTsXNetData( 'mprage_hippmapp3r' ) )
> qa <- tidNeuralImageAssessment( t1, verbose = TRUE )

Full reference image quality assessment

Python

>>> import ants
>>> import antspynet 
>>> 
>>> r16 = ants.image_read(ants.get_data("r16"))
>>> r64 = ants.image_read(ants.get_data("r64"))
>>> psnr_value = antspynet.psnr(r16, r64)
>>> ssim_value = antspynet.ssim(r16, r64)

R

> library( ANTsR )
> library( ANTsRNet ) 
> 
> r16 <- antsImageRead( getANTsRData( "r16" ) )
> r64 <- antsImageRead( getANTsRData( "r64" ) )
> psnrValue <- PSNR( as.array( r16 ), as.array( r64 ) )
> ssimValue <- SSIM( as.array( r16 ), as.array( r64 ) )

Raw

c_2_antsxnet_data_augmentation.md

Data augmentation

Noise

Python

>>> import ants
>>> import antspynet
>>> import random
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> 
>>> t1 = ants.iMath( t1, "Normalize" )
>>> noise_parameters = (0.0, random.uniform(0, 0.05))
>>> t1 = ants.add_noise_to_image(t1, noise_model="additivegaussian", noise_parameters=noise_parameters)
>>>

R

> library( ANTsR )
> library( ANTsRNet )
> 
> t1 <- antsImageRead( getANTsXNetData( 'mprage_hippmapp3r' ) )
> 
> t1 <- t1 %>% iMath( "Normalize" )
> noiseParameters <- c( 0.0, runif( 1, 0.0, 0.05 ) )
> t1 <- addNoiseToImage( t1, noiseModel = "additivegaussian", noiseParameters = noiseParameters )

Histogram intensity warping

Python

>>> import ants
>>> import antspynet
>>> import random
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> 
>>> break_points = [0.2, 0.4, 0.6, 0.8]
>>> displacements = list()
>>> for b in range(len(break_points)):
>>>     displacements.append(abs(random.gauss(0, 0.05)))
>>>     if random.sample((True, False), 1)[0]:
>>>         displacements[b] *= -1
>>>
>>> t1 = antspynet.histogram_warp_image_intensities(t1, 
                                                    break_points=break_points, 
                                                    clamp_end_points=(True, False), 
                                                    displacements=displacements)

R

> library( ANTsR )
> library( ANTsRNet )
> 
> t1 <- antsImageRead( getANTsXNetData( 'mprage_hippmapp3r' ) )
> 
> breakpoints <- c( 0.2, 0.4, 0.6, 0.8 )
> displacements <- c()
> for( b in seq.int( length( breakpoints ) ) )
>   {
>   displacements[b] <- abs( runif( 1, 0, 0.05 ) )
>   if( sample( c( TRUE, FALSE ), 1 ) )
>     {
>     displacements[b] = displacements[b] * -1
>     }
>   } 
> 
> t1 <- histogramWarpImageIntensities( t1, 
                                       breakPoints = breakpoints, 
                                       clampEndPoints = c( TRUE, FALSE ), 
                                       displacements = displacements )

Simulate bias field

Python

>>> import ants
>>> import numpy as np
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> log_field = antspynet.simulate_bias_field(t1, 
                                              number_of_points=10, 
                                              sd_bias_field=1.0, 
                                              number_of_fitting_levels=2, 
                                              mesh_size=10)
>>> log_field = log_field.iMath("Normalize")
>>> log_field_array = np.power(np.exp(log_field.numpy()), 4)
>>> t1 = t1 * ants.from_numpy(log_field_array, origin=t1.origin, spacing=t1.spacing, direction=t1.direction)

R

> library( ANTsR )
> 
> t1 <- antsImageRead( getANTsXnetData( "mprage_hippmapp3r" ) )
> logField <- simulateBiasField(image, 
                                numberOfPoints = 10, 
                                sdBiasField = 1.0, 
                                numberOfFittingLevels = 2, 
                                meshSize = 10 ) %>% iMath( "Normalize" )
> logField <- ( exp( logField ) )^4
> t1 <- t1 * logField

Random spatial transformations

Python

>>> import ants
>>> import antspynet
>>>
>>> t1 = ants.image_read(antspynet.get_antsxnet_data('mprage_hippmapp3r'))
>>> segmentation = ants.threshold_image(t1, "Otsu", 3)
>>>
>>> data_augmentation = antspynet.randomly_transform_image_data(t1,
             [[t1]],
             [segmentation],
             number_of_simulations=2,
             transform_type='affineAndDeformation',
             sd_affine=0.01,
             deformation_transform_type="bspline",
             number_of_random_points=1000,
             sd_noise=2.0,
             number_of_fitting_levels=4,
             mesh_size=1,
             sd_smoothing=4.0,
             input_image_interpolator='linear',
             segmentation_image_interpolator='nearestNeighbor')
>>>
>>> simulated_image1 = data_augmentation['simulated_images'][0][0]
>>> simulated_segmentation1 = data_augmentation['simulated_segmentation_images'][0]
>>> simulated_image2 = data_augmentation['simulated_images'][1][0]
>>> simulated_segmentation2 = data_augmentation['simulated_segmentation_images'][1]

R

> library( ANTsR )
> library( ANTsRNet )
> 
> t1 <- antsImageRead( getANTsXNetData( 'mprage_hippmapp3r' ) )
> segmentation <- thresholdImage( t1, "Otsu", 3 )
 Otsu Thresh with 3 thresholds
> 
> dataAugmentation <- randomlyTransformImageData( t1,
              list( list( t1 ) ),
              list( segmentation ),
              numberOfSimulations = 2,
              transformType = 'affineAndDeformation',
              sdAffine = 0.01,
              deformationTransformType = "bspline",
              numberOfRandomPoints = 1000,
              sdNoise = 2.0,
              numberOfFittingLevels = 4,
              meshSize = 1,
              sdSmoothing=4.0,
              inputImageInterpolator = 'linear',
              segmentationImageInterpolator = 'nearestNeighbor' )
> 
> simulatedImage1 <- dataAugmentation$simulatedImages[0][0]
> simulatedSegmentation1 <- dataAugmentation$simulatedSegmentationImages[0]
> simulatedImage2 <- dataAugmentation$simulatedImages[1][0]
> simulatedSegmentation2 <- dataAugmentation$simulatedSegmentationImages[1]

Combined

Python

>>> import antspynet
>>> import ants
>>> 
>>> t1s = [[ants.image_read(antspynet.get_antsxnet_data("kirby"))]]
>>> aug = antspynet.data_augmentation(t1s, segmentation_image_list=None, pointset_list=None, 
>>>                                     number_of_simulations=10, 
>>>                                     reference_image=None, 
>>>                                     transform_type='affineAndDeformation', 
>>>                                     noise_model='additivegaussian', 
>>>                                     noise_parameters=(0.0, 0.05), 
>>>                                     sd_simulated_bias_field=0.05, 
>>>                                     sd_histogram_warping=0.05, 
>>>                                     sd_affine=0.05, 
>>>                                     output_numpy_file_prefix=None, 
>>>                                     verbose=True)
>>> for i in range(len(aug['simulated_images'])):
>>>     print("Writing simulated image ", str(i))
>>>     ants.image_write(aug['simulated_images'][i][0], "simulated_image" + str(i) + ".nii.gz")

R

> library( ANTsR )
> library( ANTsRNet )
> 
> t1s <- list( list( antsImageRead( getANTsXNetData( 'kirby' ) ) ) )
> aug <- dataAugmentation( t1s, segmentationImageList = NULL,
              pointsetList = NULL,
              numberOfSimulations = 10, 
              referenceImage = NULL,
              transformType = 'affineAndDeformation',
              sdAffine = 0.05,
              noiseModel = 'additivegaussian',
              noiseParameters = c( 0.0, 0.05 ),
              sdSimulatedBiasField = 1.0,
              sdHistogramWarping = 0.05,
              outputNumpyFilePrefix = NULL,
              verbose = TRUE )
> for( i in seq.int( length( aug$simulatedImages[[0]] ) )
>   {
>   cat( "Writing simulated image", i, "\n" )
>   antsImageWrite( aug$simulatedImages[[i][[0]]], paste0( "simulatedImages", i, ".nii.gz" ) )
>   }