Dear Joshua,

I am not sure why you need to separate both components? What’s wrong with simply running EasyReg iteratively? It won’t be any slower, since the affine transform is computed almost instantaneously from the segmentations?

Cheers,

/E

 

--

Juan Eugenio Iglesias

http://www.jeiglesias.com

 

From: freesurfer-bounces@nmr.mgh.harvard.edu <freesurfer-bounces@nmr.mgh.harvard.edu> on behalf of Joshua K Lee <jkilee@ucdavis.edu>
Date: Saturday, September 7, 2024 at 12:59 PM
To: Freesurfer support list <freesurfer@nmr.mgh.harvard.edu>
Subject: [Freesurfer] using easy_reg to build unbiased template: separating linear and non linear components

        External Email - Use Caution        

apologies if this is a duplicate, I seem to be having some trouble with the listserv.


I've really impressed with mri_easy_reg and I was hoping to use it in a similar way as ANTs built_template scripts, which requires separation of affine from nonlinear transforms for their unbiased template updating algorithm. However, mri_easy_reg only returns a single transform object, combining affine and nonlinear components.

I added an argument as an option for producing separate transforms

```
def main():

    parser.add_argument("--separate_transforms", default=None, help="(optional) Produce Affine and NonLinear Fields Separately, if 'a' produce combined transform too. If 'b', do not produce combined transform")

```

then include
```
        if args.separate_transforms is not None:

            sepAII3 = II2

            sepAJJ3 = JJ2

            sepAKK3 = KK2

            sepAaffine = torch.tensor(np.matmul(Mflo, atlas_aff), device='cpu')

            sepARAS_X = affine[0, 0] * sepAII3 + affine[0, 1] * sepAJJ3 + affine[0, 2] * sepAKK3 + affine[0, 3]

            sepARAS_Y = affine[1, 0] * sepAII3 + affine[1, 1] * sepAJJ3 + affine[1, 2] * sepAKK3 + affine[1, 3]

            sepARAS_Z = affine[2, 0] * sepAII3 + affine[2, 1] * sepAJJ3 + affine[2, 2] * sepAKK3 + affine[2, 3]
```
and
```
            if args.separate_transforms is not None:

                #create separate nonlinear naked field

                ndepth, nheight, nwidth, _ = f2r_field.shape

                # Create meshgrid indices for each dimension

                nII, nJJ, nKK = np.meshgrid(np.arange(ndepth), np.arange(nheight), np.arange(nwidth), indexing='ij')

                # Convert numpy arrays to PyTorch tensors

                nII2 = torch.tensor(nII, dtype=torch.float32, device='cpu')  # or device='cuda' if using GPU

                nJJ2 = torch.tensor(nJJ, dtype=torch.float32, device='cpu')

                nKK2 = torch.tensor(nKK, dtype=torch.float32, device='cpu')

                # Now use these tensors as inputs to your interpolation function

                nFIELD = fast_3D_interp_field_torch(f2r_field, nII2, nJJ2, nKK2)

                nII3 = nII2 + nFIELD[:, :, :, 0]

                nJJ3 = nJJ2 + nFIELD[:, :, :, 1]

                nKK3 = nKK2 + nFIELD[:, :, :, 2]
```
```
and saves them
        if args.separate_transforms is not None:

            print("args.separate_transforms is not None")

            # Correctly combine the prefix with the filename

            fwd_field_dir = os.path.dirname(args.fwd_field)

            fwd_field_basename = os.path.basename(args.fwd_field)

            nonlinear_fwd_field = os.path.join(fwd_field_dir, nonlinear_prefix + fwd_field_basename)

            print(f"  Saving nonlinear field to: {nonlinear_fwd_field}")

            try:

                nonlinear_field = torch.stack([nII3, nJJ3, nKK3], axis=-1)

                print(f"  Nonlinear field shape: {nonlinear_field.shape}")

                save_volume(nonlinear_field, Raff, Rh, nonlinear_fwd_field, n_dims=3)

                print("  Nonlinear field saved successfully")

            except Exception as e:

                print(f"  Error saving nonlinear field: {str(e)}")

 

            # Do the same for the affine field

            affine_fwd_field = os.path.join(fwd_field_dir, affine_prefix + fwd_field_basename)

            print(f"  Saving affine field to: {affine_fwd_field}")

            try:

                save_volume(torch.stack([sepARAS_X, sepARAS_Y, sepARAS_Z], axis=-1), Raff, Rh, affine_fwd_field, n_dims=3)

                print("  Affine field saved successfully")

            except Exception as e:

                print(f"  Error saving affine field: {str(e)}")

```

full python file attached.

but it doesn't work. When I try applying them to floating image one at a time, the images go haywire.
On the left is after applying the affine transform my code generated, on the right is applying the nonlinear transform onto the image on the left.

 

Josh Lee

Assistant Professional Researcher

Department of Psychiatry and Behavioral Sciences

University of California, MIND Institute

jkilee@ucdavis.edu