H&E to Xenium SIFT-based Alignment (HEX-SIFT)

About

Xenium is a partially-destructive, highly multiplexed in-situ hybridization method in which the finished slide can still be used for certain staining (H&E and some antibodies (not CODEX)) protocols. In order to best interpret the output of these slides with each other it is best if we can align the post-Xenium stained image (such as H&E) with the dapi image generated by Xenium. This is done using the SIFT feature detection algorithm published by David G Lowe (https://doi.org/10.1023/B:VISI.0000029664.99615.94). The US-patent on SIFT expired in 2020 (https://patents.google.com/patent/US6711293B1/en). Accurate image alignment will then allow you to evaluate the accuracy of cell segmentation as well as better interpret localization of transcripts in the tissue morphology. Depending on how you run it this tool can take anywhere from 5 minutes to 20 minutes to run for a section that takes a quarter to half a Xenium slide. Most of that time is just writing the image output files.

Preparing input files

To run this pipeline you need only 2 input files:

1. the dapi image from the xenium run. This is found in the Xenium output run folder and is called either morphology_mip.ome.tif or morphology_focus/morphology_focus_000X.ome.tif (replace X with an integer 1 through 4).
2. the image of the stained slide. The input image can be a .qptiff, a .tiff, a .tif, or a .ome.tif file. I haven't tried other image formats. For consistency I am currently using .ome.tif files generated in the following manner.
   • In some instances if an image of the entire slide is taken and you may want to first crop out the region of interest with the desired piece of tissue on it. To do this follow the below instructions sourced from https://forum.image.sc/t/how-to-crop-an-image/73550:
     • Open the H&E image with QuPath . If given the option to "generate a pyramid dynamically" from the image when opening it in QuPath do so.
     • Select the rectangular annotation drawing tool in the top left (alternatively: from the top menu drop-down Tools > rectangle)
     • Draw a rectangular box around the tissue that matches the desired Xenium sample. Do not click anywhere else in the image. The box needs to stay selected (yellow). If it is not selected (yellow), then you need to click the edge of the box so it is selected (yellow). Alternatively you can go to the Annotations tab on the left and click on the annotation that outlines the region you are trying to export.
     • While the box is still selected go to File > Export Images > OME TIFF.
     • Change the settings to the following:
       • Compression type: ZLIB (lossless)
       • Pyramidal downsample: 2
       • Tile size: 1024 px
       • Parallelize export (if this is checked then it exports faster but takes more memory)
   • if you initially align the one image and then retake the image with the a microsope you will need to re-align the second image to the Xenium dapi image and the results will differ.
   • If you want to re-take the image of a specific region of interest instead of the full slide, then you should still be able to align it with this tool by saving the image as a raw tif file, but I have not tested it. You will need to take a image large enough that sufficient keypoint matches are found and I do not know what that size might be. I typically use this script to align 20X magnification H&E images. I have not tested it on 40X images. The main reason for using 20X images is to save on memory and storage space.

Example commands

To run this tool you need a python environment with numpy, OpenCV (cv2), tifffile, and optparse, (future plan to also use ome_types package). I have tested this in python 3.9 and python 3.10. You can clone the environment I have previously used to run the tools as follows:

conda create --name hexsift --file py3.10_spec_file_2024-05-06.txt
conda env create -f py3.10_environment.yml --name hexsift

If you do not yet have Anaconda installed on your lab cluster then please see the Conda Installation on Linux guide for how to do so. This tools requires one of the following two script to run. The current version now enables alignment to multimodal images. It allows alignment using either just the dapi channel from the multimodal image (strongly recommended) or the the full multi-modal image. For non-multimodal images the only major change is in how the dapi and H&E images are read into the script.

image_alignment_HE_to_xen_dapi_post_xoa_v3.0.py

For alignment of H&E images to Xenium that are of non-multimodal samples prior to the Xenium Onboard Analysis (XOA) update to v3.0 you can still use this script:

image_alignment_HE_to_xen_dapi_pre_xoa_v3.0.py

Options

There are currently 2 flags that you will need to use every time you run the command and one additional flag you may need to use depending on your image:

• -i or --image is the path to the image file that you will be aligning to the Xenium dapi image. For now this should be an H&E image. If it is an image of a fluorescent antibody stain, then channel that contains dapi needs to be the first image (channel) stored in the input image file (i.e the first channel in axis 3 he_image[:,:,0]). For now if it is not then you will have to re-order the image channels and re-save the image. One way of doing this is by using the provided writing_ome_tif_v9.py script. python writing_ome_tif_v9.py -i <your_image_file.tif> -s <name_you_want_the_output_file_to_be>.
• -d or --dapi is used to specify the path to the Xenium dapi image from the xenium output files. This is the image that output image is aligned to. This image is not warped or shifted during the alignment process. For xenium output from XOA (Xenium Onboard Analysis) vrom v1.9 or older, it should be a uint16 (unsigned 16-bit integer) scaled image file called morphology_mip.ome.tif. You can open the analysis_summary.html file if you are unsure which XOA version you are working with. For Xenium output from XOA v2 or newer, it should be a uint16 (unsigned 16-bit integer) scaled image file located at morphology__focus/morphology_focus_0000.ome.tif. If you are working with a multimodal image and you are soft-linking the image file to a new folder and providing the path to the link instead of the path to the original folder then you may need to also softlink all of the other morphology_focus_000X.ome.tif files to the same folder as well. Do not rename the files from multimodal images. If you provide the -d flag at anypoint in the script call then it will only use the dapi channel of a multimodal image.
• --multimodal_path is used to specify the path to the folder containing the multimodal images from xenium platform. This is the folder that contains the 4 files: morphology_focus_0000.ome.tif, morphology_focus_0001.ome.tif, morphology_focus_0002.ome.tif, and morphology_focus_0003.ome.tif. Do not rename the files if using this option. To align to the full multimodal image provide just the path to the multimodal image folder with the --multimodal_path flag. DO NOT USE THIS FLAG AT THE SAME TIME AS THE -d FLAG. Currently alignment with the full multimodal image does not work well because I probably need to play with the weights of the different channels when grayscaling to create. If anyone else wants to modify how the grayscale is done at this line xen_grey = np.mean(xen_color, axis=0) and finds weights that work consistently please let me know.
• -f or --flip is the flag you will need to use if your image also needs to be flipped. SIFT only allows for your image to be translated (moved up/down/left/right), rotated, or scaled (shrink/stretch). It does not allow for images to be flipped or transposed (transposed = rotated and then flipped). If your image needs to be flipped or transposed then you need to use this flag when running your pipeline. If you do not use this flag then your output will look extremely warped or like you are falling down in a cliff while looking up at the rising edges of said cliff. This is the first thing to check if the result of running this pipeline looks like the above. However, you can still get the above result when you do not need to use this flag and you don't use it. More on that below. There are 3 optional flags that can improve your results.
• -b is the flag you will use when you want to specify for the script to only use the blue channel of the image. This results in the alignment only relying on the first image channel (blue) of the BGR (Blue, Green, Red) channels that all bright-field images include. In otherwords the output will only rely on he_image[:,:,0]. Generally I find using this flag improves the overall alignment for all samples I have tried it with so far. In addition though for some samples it is needed to ensure a high-quality alignment. This is because the morphology_mip.ome.tif only contains a dapi stain of the nuclei on a Xenium slide. The H&E stain of that same slide stains the slide pink and purple with most of that purple being present in the nuclei. Therefore the To align the slide we need to create a grayscale image and if one part of that image is a intense red it is because there are few to no nuclei there and that signal will not be present anywhere in the Xenium image you are aligning to. To get around this we can filter the image to only the blue pixel values since the most intense blue pixel values are typically found in the cell nuclei of the H&E slide. If parts of your aligned output do not properly overlap the Xenium dapi image then use this flag if you are not already. If your input H&E image is already very blue in regions that are not nuclei then using this flag may or may not help.
• -s or --scaling is the flag you use to specify how much the image is down-scaled when the image is aligned. Specify a whole number from 0-7. Default is 4. It is not recommended that you change this value for an initial alignment beyond those values suggested below in the example commands. This is the level of dapi and H&E image downscaling that will be used to calculate the Homograpphy matrix for image alignment. A level of 4 means that the Y and X axes of the H&E image will be divided by 2^4 (i.e. 16) during initial homography matrix calculation. The Homography matrix will then be scaled back up by the resulting scale matrix for use in the full size H&E alignment. Large amounts of tissue detachment (1/5th or more of the total area of the image) will decrease image alignment quality. If there is tissue detachment, it is still sometimes possible to align the H&E with a scaling power of 4 or less, however you may need to increase to 5,6 or 7. Images alignments generating with scaling of 5,6, or 7 are often of lower quality.
• -m or --matches_limit is used to specify any whole number less than the total number of keypoint matches found. Default is 20. It is not recommended that you change this value for an initial alignment beyond those values suggested below in the example commands. Keypoints and their descriptors of an image are identified using the SIFT algorithm and they are then matched using a FLANN matcher. The top n keypoint matches are then used for image alignment. If there are a large number of potentially good keypoint matches then increasing this number can improve the alignment. Otherwise decreasing it can hurt the alignment. When the script runs normally part of the stdout states something along the lines Number of good matches: 248. So long as the value of -m is less than the value printed to stdout in that line your alignment should be of decent quality and the script will not crash. If you decrease the downscaling of the image prior to alignment (use a smaller value for -s then you should generally provide a larger value with -m). If large portions of the tissue have detached from the slide during or prior to the Xenium image generation, or during H&E generation, then it is likely that the number of matches will be lower. Tissue detachment will decrease successful image alignment

Command

The following are some example commands:

NOTE: To date I have found the highest quality alignments to be generated when using the flags: -b -s 0 -m 60. However the total size of all the output files combined is much larger, it can use a lot more memory (200GB+ due to the number of keypoints being identified on the imge) and the script takes much longer to run (most of this additional time is just from writing all of the output files), also during SIFT keypoint detection the memory can easily take up half of a 500GB node if this option is used. Most of the output files are just QC and used to keep track of which parts of the image are used to generate the alignment. The equivalent to these commands that the script uses by default is -s 4 -m 20 (uses closer to 30-70GB). A good middle ground is -b -s 2 -m 40. Increasing the number of matches used while keeping the scaling power (e.g. going from -b -s 2 -m 40 to -b -s 2 -m 60) the same does not increase memory requirements significantly. Including more keypoint matches and using the full resolution image is not always a good thing as incorrect matches can still be found which is why we don't just use all keypoint matches. We only want to take the highest quality keypoint matches. The large memory problem is something that I am planning on addressing in the future. Note that the extreme memory usage can probably be addressed by tiling over the image when identifying keypoints (would then need to loop over all keypoints and check to make sure they are all unique) and I just haven't had time to do that yet, but plan to do so in the future.

Commands post XOA v3

It is recommended that you use this version for all future alignments.

Dapi only alignment example (recommended for most use cases):

cd /folder/where/your/want/the/output/to/be/
python3 /path/to/HEX-SIFT/image_alignment_HE_to_xen_dapi_post_xoa_v3.0.py -i /path/to/image_RAW.ome.tif -d /path/to/Xenium/output/folder/morphology_focus/morphology_focus_0000.ome.tif -b -s 4 -m 20

# this gives the same result as if you leave off the `-s 4 -m 20`
cd /folder/where/your/want/the/output/to/be/
python3 /path/to/HEX-SIFT/image_alignment_HE_to_xen_dapi_post_xoa_v3.0.py -i /path/to/image_RAW.ome.tif -d /path/to/Xenium/output/folder/morphology_focus/morphology_focus_0000.ome.tif -b

cd /folder/where/your/want/the/output/to/be/
python3 /path/to/HEX-SIFT/image_alignment_HE_to_xen_dapi_post_xoa_v3.0.py -i /path/to/image_RAW.ome.tif -d /path/to/Xenium/output/folder/morphology_focus/morphology_focus_0000.ome.tif -b -s 0 -m 60

Multi-modal alignment example: (less likely to be successful at lower -s values. While it is still possible to use the full multimodal image and only the blue channel of the H&E image as input it is not recommended that you do so (i.e. do not use -b)

cd /folder/where/your/want/the/output/to/be/
python3 /path/to/HEX-SIFT/image_alignment_HE_to_xen_dapi_post_xoa_v3.0.py -i /path/to/image_RAW.ome.tif -multimodal_path /path/to/Xenium/output/folder/morphology_focus/ -s 4 -m 20

Example when the image needs to be flipped

cd /folder/where/your/want/the/output/to/be
python3 /path/to/HEX-SIFT/image_alignment_HE_to_xen_dapi_post_xoa_v3.0.py -i /path/to/image_RAW.tif -d /path/to/Xenium/folder/morphology_focus/morphology_focus_0000.ome.tif -f 

Commands pre XOA v3

cd /folder/where/your/want/the/output/to/be/
python3 /diskmnt/Projects/Users/austins2/tools/xenium_image_alignment/image_alignment_HE_to_xen_dapi_pre_xoa_v3.0.py -i /path/to/image_RAW.ome.tif -d /path/to/Xenium/folder/morphology_mip.ome.tif -f 

cd /folder/where/your/want/the/output/to/be/
python3 /diskmnt/Projects/Users/austins2/tools/xenium_image_alignment/image_alignment_HE_to_xen_dapi_pre_xoa_v3.0.py -i /path/to/image_RAW.ome.tif -d /path/to/Xenium/folder/morphology_mip.ome.tif -f -b

cd /folder/where/your/want/the/output/to/be/
python3 /diskmnt/Projects/Users/austins2/tools/xenium_image_alignment/image_alignment_HE_to_xen_dapi_pre_xoa_v3.0.py -i /path/to/image_RAW.ome.tif -d /path/to/Xenium/folder/morphology_mip.ome.tif -b -s 2 -m 40

cd /folder/where/your/want/the/output/to/be/
python3 /diskmnt/Projects/Users/austins2/tools/xenium_image_alignment/image_alignment_HE_to_xen_dapi_pre_xoa_v3.0.py -i /path/to/image_RAW.ome.tif -d /path/to/Xenium/folder/morphology_mip.ome.tif -b -s 0 -m 60

Evaluating the output quality

There are 17 files that are generated by this tool. The vast majority are for the purposes of QC and debugging when necessary. The 3 files that you will need as output are as follows:

• Dapi_stacked_with_HE_aligned.ome.tif - This file contains two tiled, pyramidal, zlib-compressed images in it. Each image is stored as a separate channel. The first channel is the dapi image that is input when the script is run using the -d flag. The second image is the aligned greyscale of the H&E image (or whatever image is input with the -i flag. The dapi image is not the original dapi image however. The original dapi image is scaled from 0 to 65535 (the unsigned 16-bit integer limit = (2^16)-1). The dapi image in this file has been rescaled from 0 to 255 (the unsigned 8-bit integer (uint8) limit). This allows the image to be viewed in imageJ and QuPath and at the same time compared to the other image in the file which is scaled from 0 to 255 (uint8 limit). To view this image open Fiji. The go to Plugins > Bio-Formats > Bio-Formats Importer. Then use the browser to navigate to and select the Dapi_stacked_with_HE_aligned.ome.tif image. You will then see a menu of options. On the top left change the drop down to "Hyperstack" (the channels field should then automatically change to "XYCZT"). If you want to see image metadata then check the box that says "Open Image Metadata". Then click okay and when you see a menu of image options appear uncheck the box next to "Series 1" If you do not unckeck the box next to Series 1 there is a large chance that the image will fail to open if you system does not have enough memory. Scroll until you find an image size less than ~14,000 x ~14,000 pixels and check the box next to that Series. You can then open the image using the button at the bottom. From here you can use Fiji to navigate around the image using the hand tool to move left/right/up/down, the up/down arrow keys to zoom in and out, the left and right arrow keys to swap back and forth from the H&E image to the dapi image.
• he_aligned.ome.tif - This file is an ome tif that contains the full color aligned H&E image. You can inspect the aligned image in QuPath. If you want to use the aligned H&E image for other analysis in Xenium Explorer, Python or R this is also the file that you should use. See below for how to view this image in 10X Genomics Xenium Explorer.
• he_aligned.imageJ.ome.tif - This file contains the same aligned H&E image as before but it is saved in a slightly different format so that the image can be read by imageJ. You will not need this image unless you want to look at your image in imageJ. If you do then follow the instructions under the Dapi_stacked_with_HE_aligned.ome.tif file to open the image in Fiji.
• channels.pickle - This is a pickle dump of the following list: [dapi_image_scale_0_255, found_large]. dapi_image_scale_0_255 is the dapi image rescaled from np.uint16 (0-65535) to np.uint8 (0-255) in a numpy array. found_large is the aligned full resolution aligned H&E image in a numpy array. This file contains the image data for the aligned H&E and xenium dapi image. If there is a problem with either of the he_aligned.ome.tif file or the he_aligned.imageJ.ome.tif, then the channels.pickle file can be used to re-write a new version of the images. A brief description of the remaining files:
• found.tif - This file contains the downscaled version of the grayscale H&E image that is the output of the shrunk image that is used to find the homography matrix that is used in the final scaled up alignment. If you use -s 0 when you run the script then this file is just the full size grayscaled H&E image.
• found_full_res.tif - the full resolution version of the aligned grayscaled H&E image. This is a grayscale of the final he_aligned.ome.tif that is used for subsequent analysis.
• matched.tif - This image file shows the matching keypoints that were used to calculate the homography between the dapi image and the H&E image. All keypoints are circled. Matching keypoints are connected by a line.
• dapi_gray_scale_target.tif - This is the version of the dapi image that the H&E image is aligned to. It is re-scaled to the unsigned 8-bit integer limit (values go from 0 to 255).
• he_blue.tif - This is the blue pixel values of the H&E image. If the -b flag is used this is the image that is used identify keypoints and generate the homography matrix.
• he_red.tif - This is the red pixel values of the H&E image. If the -b flag is used this channel is not used to identify keypoints and generate the homography matrix. If any large region of this image is near black then it is a good idea to use the -b flag.
• he_green.tif - This is the green pixel values of the H&E image. If the -b flag is used this channel is not used to identify keypoints and generate the homography matrix. If any large region of this image is near black then it is a good idea to use the -b flag.
• he_small.tif - This is the grayscaled H&E image pre-image thresholding and alignment. If the -b flag is used this is just the blue pixel/channel values of the H&E image.
• he_thresh.tif - This is the blurred ((Gaussian blur (3,3)), thresholded (cv2.THRESH_BINARY_INV + cv2.THRESH_OTSU) image generated from the grayscaled H&E image.
• he_close.tif - This is morphological close (disk 9x9) of the blurred, thresholded, grayscale H&E image.
• he_erode_2.tif - This is the mask that is applied to the grayscaled H&E image. It is generated by eroding the H&E he_close.tif image 2 times using a disk ellipse (2,2).
• he_masked.tif - This is the masked grayscale H&E image.
• he_aligned.pickle - This is a pickle dump of the aligned full-resolution, RGB-colored H&E image. If something is wrong with the he_aligned.ome.tif file. Then this pickle file can be used to re-write the image.
• scaled_M.csv - This is the homography matrix as it calculated from OpenCV2 that has been scaled up for use on the full size image via the scale matrix. It has then been inversed for compatibility with skimage since skimage.transform.warp() is what is now used for the full size H&E image warp. It is used to calculate the position of the dapi image corner points so the keypoints of the dapi image match the H&E image.
• perspectiveM_large.csv - This is the perspective transformation matrix that is then used to warp the H&E image to the dapi image. It is calculated based on the transformation of the dapi destination points. If you need the matrix for warping the H&E image to the dapi for any downstream stuff I think this is it.

Importing aligned images into Xenium Browser

If you plan on importing the image into Xenium Explorer then you should do so using the he_aligned.ome.tif image. To do so open the experiment.xenium file from the Xenium output files. Then select the "Image" drop down. Click the "Add image" button. Click "Import Image". Navigate to your downloaded he_aligned.ome.tif file and select it. Wait for a the browser to show a window with a small version of your H&E image (this may take a minute or two). Give the H&E image a name and click "Continue". At the next pop-up select "Yes, skip alignment" and click "Done".

Troubleshooting

• I usually run the alignment 3 times with the following settings and then select the one that gives the best alignment in the region of interest: -b -s 4 -m 20, -b -s 2 -m 40, and -b -s 0 -m 60. I then review the he_aligned.ome.tif file in Xenium explorer to select the best one.
• If you are working with a large tissue section (takes up close to 1/2 of the usable area of the Xenium slide) and or is nuclei dense it can be helpful to use -m higher than 60. In such cases I find -m 100 to be helpful.
• Sidenote: When the tissue is dehydrated during H&E staining after the Xenium data generation it is not uncommon for tissue to shrink or shift slightly. The result is that some region of a image are often aligned slightly better than others. The result on the data analysis side is that regions in the center of the image are often better aligned than regions on the edges of the tissue.
• If you need to work with a region specifically on the edge of the tissue in a section near where a large amount of the tissue has detached, HEX-SIFT on the full image may not give a high quality image alignment. In this case you can try croping the H&E image to an area larger than you ROI that excludes the detached tissue (do the same for the dapi image - make sure to save a copy and not overwrite the original file). The cropped images may then be aligned with HEX-SIFT.
• Yes the scale bars shown when viewing the image in QuPath are broken. I have tried to fix this by writing the image metadata multiple different ways but have yet to be successful. I aim to fix this in the future.