There are two main approaches to flatten a warped specimen.
- Either, during reconstruction, by making use of the markers present on the specimen boundaries; this is what we refer as flattening on the fly.
- Or by flattening an already reconstructed volume.
While the first approach is very attractive because it is a pseudo-automated method (no additional tracking or segmenting is necessary), the results remain strongly dependent on the initial marker distribution. The second approach can be more precise, but requires extra manual steps since a model for the specimen boundaries needs to be provided.
TxBR flattening on the fly
A specific choice for the bundle adjustment gauge ambiguity
As we mentioned in the how to run TxBR section, the bundle adjustment solution to the alignment problem is only known within an invertible 3D transformation T. Given a set of solution for the projection map Pθ and the 3D marker positions Xτ, another possible solution is:
θ denotes for a tilt and τ for a marker.
We can search for T as a polynomial deformation applied along the normal to the specimen (Z-axis):
so that the markers are brought on two parallel planes. The coefficients of h are calculated by minimizing the following quantity:
with respect to hij, Z1 and Z2. This is equivalent to searching for the best-fit polynomial surfaces going through the set of gold markers.
This method is suited to cases where markers remain on the specimen surfaces, such as for plastic sections. We nevertheless expect it to hold as well when the markers are found within the specimen.
Running TxBR flattening on-the-fly
In order to enable the on-the-fly flattening capability when running TxBR, the --flatten-order flag should be added to the runtxbr.py command line script (or txbr_align.py when just running the alignment step), followed by the order n of h, for e.g.:
runtxbr.py -b Naoko3A9 --flatten-order 6
We applied the on-the-fly flattening on a dendrite specimen from a mouse hippocampus (see figure 1). A series of 61 tilts was taken by increment of 2 degrees. From the little animation below, it is slightly apparent that the specimen is warped.
Fig. 1: Tilt series of a dendrite specimen.
Indeed, by looking at the reconstruction without on-the-fly flattening (third image on figure 2), the sample appears quite distorted, especially from the top view (XZ view). It is not even clear if the specimen boundaries will be well describe with low order polynomial function. We used a flattening-order of 6 (this is a very big value for this method), results are displayed on the other snapshots of figure 2.
Fig. 2: On-the-fly flattening of a warped dendrite specimen. Click an image to get an enlarged view.
In the first two views of figure 2, the 3D positions of the fiducial markers (80 of them) are displayed with white beads, before and after applying the flattening transformation (at alignment). The surfaces in the first view represent the best-fit 6th-order polynomial function crossing the set of markers.In the last two views of figure 2, a XYZ snapshot for both the warped and flattened reconstruction volume is shown. The final result is fairly flat, even though some slight bumps on the surface still remain — this particular example is especially difficult to flatten. For better results, the method outlined in the section Flattening a warped volume in 3D should be used.
High order projection maps and on-the-fly flattening
The on-the-fly flattening can be applied in conjunction with a non-linear projection map. We have applied this feature to the test case of a mitochondria study, by typing the following command line:
runtxbr.py -b wtped3 --n2 2 --flatten-order 2
Here, the --n2 and --flatten-order flags are respectively set to 2 and 2. It is important to note that while the order of the projection map is quadratic, its effective order in term of computational operations is 4, the product n2 ‚òì n between the order of the projection map (n2=2) and the flattening order (here n=2). This means in practice, one will always to use low order polynomial function h to describe the specimen boundaries, so the effective projection map order remains also low.
In this particular example, the tilt series is made of 3 by 3 mosaics obtained using the image shift capability of an electron microscope. Micrographs were stitched together using the blendmont utility from IMOD. As very obvious on the leftmost plot of figure 3, the specimen is warped in a saddle shape manner. Nevertheless, spatial distribution of the markers allowed for a very good description of the specimen boundaries: 138 gold markers (white beads) were initially tracked on this series of 116 tilts (1 degree increment). The on-the-fly method is very effective as it is clear on the second plot of figure 3.
Fig. 3: On-the-fly flattening of a warped specimen featuring a mitochondria. Click an image to get an enlarged view.
Finally the last two snapshots in figure 3 are XYZ view of the specimen reconstruction when the on-the-fly_ option if OFF (warped case) and ON (flattened case). Let's mention that the original size of an image from the tilt series is approximately 12k ‚òì 12k pixels.
Flattening a warped volume in 3D
As an alternative approach to the previous section, flattening is performed directly on a 3D reconstruction [1].
Local Expansion and Compression along the sample normal axis
The volume density u is then readjusted along the Z axis (normal to the specimen or direction of flattening) according to the following rule:
where
Equations Z=f1(X,Y) and Z=f2(X,Y) describe the boundary surfaces of the specimen. They are evaluated after a manual step to segment the boundary surfaces has been done – see next paragraph. Z1 and Z2 represent the locations of the specimen boundaries on the Z axis after flattening has occurred. In the default setup, f1 and f2 are approximated with polynomial functions.
How to run the flattening routine on a 3D volume
The first step in this flattening procedure is to create an IMOD model file that fully characterizes the specimen boundaries. The 3dmod software can be used for this purpose. To complete this step, a TxBR reconstruction of the warped sample is necessary. Two objects, one for each of the specimen boundaries, need to be created in the model. Then, points that belong to the boundary surfaces should be manually picked on a roughly regular XY grid, and added to their corresponding object. It does not matter if all the points of an object belong to only one contour or to several different contours.
Distribution of the points characterizing a boundary surface does not need to be strictly regular; less points are needed on the slowly varying part of the surface, but overall there should not be any region (in X and Y) free of "markers". The optimal average point density directly depends on the order of the polynomial function describing the surfaces. In the default configuration (order 10,10), it is at least a 10 by 10 grid, so at least 100 points; the more the better!
Creating an accurate boundary model file is critical to obtain the flattest result. Do not hesitate to work from the three view (XY, XZ and YZ) at the same time when you build the model!
In Figure 5, a boundary model for the dendrite specimen used in a previous section is displayed. On the left side, the model stands by itself while on the right side it is shown on top of the warped reconstruction volume. The reconstruction file (Naoko3A9_z_-60.0.out) and the model file (Naoko3A9.volmod) can be downloaded from here.
Fig 5. Model for the specimen boundary surfaces.
To implement the flattening, just type the following instruction as a command line:
txbr_flatten.py -B Naoko3A9.volmod -V Naoko3A9_z_-60.0.out
The program will then output a flattened version of the reconstruction named Naoko3A9_z_-60.0.out.uwrpd; the extension uwrpd is appended to the name of the volume file to process. As is visible from figure 6, the flattening result is fairly good.
Fig. 6: Flattened reconstruction of the dendrite specimen. Click on the left to enlarge.
A few options are available for the txbr_flatten.py script. Order of the polynomial map describing the surfaces can be changed by using the --order flag. Its default value of 10,10 can be replace for instance to 5,5 with typing:
txbr_flatten.py -B boundary_model -V volume_to_flatten --order 5,5
In that case, fewer points are necessary in the boundary model.
Polynomials are sometimes not the best choice to describe accurately the boundary surfaces, in particular when high frequency variations occur in their topography. An --interpolation flag is available to handle this type of situation. It means that on top of the polynomial description, an interpolation scheme will ensure the calculated surface to cross every point of the boundary model.
txbr_flatten.py -B Naoko3A9.volmod -V Naoko3A9_z_-60.0.out --interpolation
It is also possible to specify the direction where to apply the flattening procedure with the --D flag. Default direction is Z.
Finally, a choice regarding the width evaluation (that is the quantity Z2-Z1) of the flattened volume can be done by setting the flag --width-approximation. Three options are available: if set to min (respectively max), the final width will correspond to the minimum (respectively maximum) width of the warped specimen in the flattening direction; if set to isov, volumes of the flattened and warped specimen are equal. Default value for the --width-approximation option is isov.
Light Microscopy example
The 3D flattening routine described here is not just restricted to the Electron Tomography field. In this paragraph, we display results obtained for volumes created with light microscopy techniques [2]. In this large scale mosaic brain study, 6 consecutive sections of 500 mosaics had to be stacked together. However as it is visible in the left side of figure 7, the sections do not exhibit a perfect slab-like shape, making the serial-stacking step difficult to implement.
Fig. 7: on the left, the boundary model used to flatten one section of this mouse brain study;
on the right, the shape of the brain section before (in yellow) and after (in blue) being flattened.
Prior to stacking, sections were flattened to ease this process. On the right side of figure 7, the shape of a brain section is displayed before and after flattening.
Note that in this study, the two boundary surfaces do not exhibit the same spatial variations and are not parallel; this underlines the importance of describing the two boundary surfaces independently of each other in the model.
Finally, in figure 8, the original and flattened volume of a brain section are displayed side by side, while slices in the normal direction to the specimen are scrolling back and forth.
Fig. 8: Scrolling in an original (left) and flattened (right) mouse brain section slice by slice.
References
- S. Phan, M. Terada and A. Lawrence, “Serial Reconstruction and Montaging from Large-Field Electron Microscope Tomograms”, accepted in 31st Annual International Conference of the IEEE Engineering in Medicine and Biology Society (2009)
- M. L. Berlanga, E. A. Bushong, S. Lamont, S. Wu, B. S. Phung, M. E. Martone, M. H. Ellisman, Three-dimensional reconstruction of serial mouse brain sections using high-resolution large-scale mosaic images, Program No. 484.11. 2009 Neuroscience Meeting Planner. Chicago, IL: Society for Neuroscience, 2009. Online