## 2 Materials and Methods
### Registration
\(R\) defines a minimization process such that Equation (2) is unfolded as  
\[T_{i\to j}^{*}=\underset{T_{i\to j}}{argmin}\mu\big{(}I_{i},\ I_{j}\ \circ\ T_{i\to j}\big{)}+\lambda\rho\big{(}T_{i\to j}\big{)} \tag{5}\]  
where \(\mu\) is the similarity metric (SSD, NMI, or NCC, the latter with the kernel size, for instance, NCC[\(3\times 3\)]), \(\lambda\) is a scalar parameter, and \(\rho\) is an optional regularization term. Initially, affine registration was performed between the image pairs, using an optimization of the dissimilarity metric based on a limited-memory Broyden-Fletcher-Goldfarb-Shanno (L-BFGS) algorithm [41], denoted by  
\[A_{i\to j}\ =\ R_{aff}\big{(}I_{i}\ \to\ I_{j};\ \mu,\ A_{0}\big{)} \tag{6}\]  
where \(A_{0}\) is the initial rigid transformation between the images. The initial transformation was obtained using a brute force search, where 4500 pairs of rigid transformations (which captured all possibly combinations of random rotations and translations for the specific dataset) were applied to the moving image and the combination, and the best NCC metric value was saved as \(A_{0}\). Specifically, a standard deviation of 180\({}^{\circ}\) for the angle (ensuring all rotations are sampled) and the standard deviation of the random displacement in each coordinate was equal to 10% of the input image width used, which was large enough to showcase deformation but at the same time small enough to mitigate folding in the dataset due to extreme deformations. Figure 6 illustrates the error in the landmarks as a function of the number of random iterations for the given dataset. This brute force search, performed at the highest pyramid level and not requiring computation of metric gradients, had significant impact on robustness and was relatively fast, i.e., contributing only a few seconds to the total registration time. Figure 7 illustrates an example result from the application of these steps before the actual affine registration. Following the affine registration, the diffeomorphic registration of slice \(j\) to \(i\) was applied:  
\[\varphi_{i\to j}\ =\ R_{diff}\big{(}I_{i}\ \to\ I_{j};\ \mu,\ \sigma_{s}, \sigma_{t},\ N\big{)} \tag{7}\]  
where \(\sigma_{s}\) and \(\sigma_{t}\) are the regularization parameters for the registration and \(N\) is the number of iterations required at each multi-resolution pyramid, e.g., \(N\) = {100,50,10} refers to 100 iterations at 4\(\times\), 50 at 2\(\times\), and 10 at full resolution. Larger values of \(\sigma_{s}\) result in more smoothing, and larger values of \(\sigma_{t}\) amount to less overall deformation. Furthermore, _Greedy_ uses an optimized smoothing of the deformation fields on the basis of the ITK recursive Gaussian smoothing classes [42]. The actual registration was computed in an iterative manner using the update equations [43]:  
\[\psi^{\gamma}=Id+\varepsilon^{\gamma}\cdot\left[ \begin{array}{c}G_{\sigma_{5}} \ast\ D_{\sigma_{i}^{\gamma}\to j}\mu\left(I_{i},\ I_{j}\ o\ \varphi_{i\to j}^{\gamma}\right)\end{array} \right] \tag{8}\]  
\[\varphi_{i\to j}^{\gamma+1}=\ G_{\sigma_{i}}\ast\left(\varphi_{i\to j}^{ \gamma}\ o\ \psi^{\gamma}\right) \tag{9}\]  
Figure 6: Landmark error over the number of random iterations for the initial transformation. Figure 7: Example results of our affine registration step. (**A**) The affine registration estimated and applied on the resampled and padded images; (**B**) application of the registration applied to the original full-scale images; (**C**) the NCC response map between source and target before and after affine registration. \[\varphi^{0}_{i\to j}=Id \tag{10}\]  
where \(\gamma\) is the current iteration, \(D_{\varphi_{i}^{T}\to j}\mu\) is the gradient of the metric with respect to \(\varphi\), \(\varepsilon^{\gamma}\) is the step size, \(G_{\sigma}*\varphi\) denotes the convolution of \(\varphi\) with an isotropic Gaussian kernel with a standard deviation of \(\sigma\), and \(Id\) is the identity transformation.