| \n", " | spot_id | \n", "x | \n", "y | \n", "
|---|---|---|---|
| 0 | \n", "pixel_14 | \n", "-14 | \n", "67 | \n", "
| 1 | \n", "pixel_15 | \n", "-15 | \n", "67 | \n", "
| 2 | \n", "pixel_16 | \n", "-16 | \n", "67 | \n", "
| 3 | \n", "pixel_17 | \n", "-17 | \n", "67 | \n", "
| 4 | \n", "pixel_18 | \n", "-18 | \n", "67 | \n", "
| \n", " | spot_id | \n", "mz-108.901945 | \n", "mz-108.902095 | \n", "mz-109.02966 | \n", "mz-109.07286 | \n", "
|---|---|---|---|---|---|
| 0 | \n", "pixel_14 | \n", "194565.2656 | \n", "147422.4375 | \n", "107558.41410 | \n", "0.00000 | \n", "
| 1 | \n", "pixel_15 | \n", "585243.8750 | \n", "386686.1875 | \n", "93422.78125 | \n", "0.00000 | \n", "
| 2 | \n", "pixel_16 | \n", "268511.0000 | \n", "199858.3125 | \n", "140381.84380 | \n", "79260.86719 | \n", "
| 3 | \n", "pixel_17 | \n", "346582.0938 | \n", "285972.9375 | \n", "117012.73440 | \n", "0.00000 | \n", "
| 4 | \n", "pixel_18 | \n", "487304.3750 | \n", "354217.7188 | \n", "29423.32617 | \n", "8445.78418 | \n", "
check_inputs -> perform_coreg -> make_spaceranger -> create_barcode_matrix \n",
"\n",
"### Parameter selection\n",
"\n",
"There are several user-selected parameters to be aware of in this file:\n",
"1) include_barcode_matrix = False specifies whether or not the final stage of the pipeline, where MSI observations are aggregated to form matching observations to Visium barcodes so that 1:1 mappings can be identified across both modalities. This stage is optional but can be very useful for downstream analysis. This part of the pipeline is the most time demanding so you may choose to run the pipeline without the step initially to check the success of the coregistration. The same step can also be performed within the semla package using the CreateMultiModalObject function.\n",
" \n",
"2) Within each stage you can set verbose = True/False to get more or less information about the different stages of the process as they are being executed.\n",
" \n",
"3) In the perform_coreg stage, a set of 3 parameters can be used no_HE_transform , MSI2HE_transform and HE2HE_transform , each of which can be set to 'affine' or 'TPS' depending on whether you want to use a linear or non-linear transformation. We discuss in the manuscript the different circumstances where each choice can be advantageous. In this example, only the no_HE_transform will be used as we do not have an intermediate MSI image and it is set here to affine. Because we have both measurements taken from the same tissue section, it may make sense to choose a linear transformation as we will only see small alterations between sections rather than large deformations.\n",
"\n",
"4) In the create_barcode_matrix step, the agg_fn parameter specifies how MSI observations (for each pixel) will be aggregated to form a single observation per Visium barcode. The options are 'mean' or 'sum' and in this case we will use the mean.\n",
"\n",
"5) Finally, in create_barcode_matrix you can choose whether only Visium barcodes within the tissue are returned or whether to not perform filtering at this stage. "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Once all of these parameters have been set, you can run the pipeline! "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Running the snakemake workflow"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"From within the same working directory, you can use the same following line to run the pipeline:\n",
"\n",
"```\n",
"snakemake -j [number of cores]\n",
"```"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The command line output will look something like this for this dataset:\n",
"\n",
"```\n",
"Building DAG of jobs...\n",
"Using shell: /bin/bash\n",
"Provided cores: 1 (use --cores to define parallelism)\n",
"Rules claiming more threads will be scaled down.\n",
"Conda environments: ignored\n",
"Job stats:\n",
"job count\n",
"--------------------- -------\n",
"all 1\n",
"check_inputs 1\n",
"create_barcode_matrix 1\n",
"make_spaceranger 1\n",
"perform_coreg 1\n",
"total 5\n",
"\n",
"Select jobs to execute...\n",
"\n",
"[Mon Feb 24 10:13:07 2025]\n",
"Job 4: Checking inputs.\n",
"Reason: Missing output files: output/summary.csv\n",
"\n",
"Working on sample V11L12-038_D1\n",
"Checking for intermediate MSI image...\n",
"Checking metadata and intensity tables exist and can be read...\n",
"Checking metadata contains spot_id, x and y columns...\n",
"Checking intensity table begins with metadata column...\n",
"Checking spot_id columns in metadata and intensity table match exactly\n",
"Storing dimensions and examples of data\n",
"Summary saved to output/summary.csv\n",
"[Mon Feb 24 10:13:17 2025]\n",
"Finished job 4.\n",
"1 of 5 steps (20%) done\n",
"Select jobs to execute...\n",
"\n",
"[Mon Feb 24 10:13:17 2025]\n",
"Job 3: Performing co-registration.\n",
"Reason: Missing output files: output/V11L12-038_D1/transformed.png, output/V11L12-038_D1/transformed.csv; Input files updated by another job: output/summary.csv\n",
"\n",
"No MSI H&E image identified.\n",
"Mapping MSI data to Visium H&E\n",
"Using affine transform to map MSI data to Visium H&E\n",
"Saving Visium H&E with new MSI coordinates overlaid...\n",
"Saving new MSI coordinates...\n",
"[Mon Feb 24 10:13:22 2025]\n",
"Finished job 3.\n",
"2 of 5 steps (40%) done\n",
"Select jobs to execute...\n",
"\n",
"[Mon Feb 24 10:13:22 2025]\n",
"Job 2: Creating spaceranger formatted output.\n",
"Reason: Missing output files: output/V11L12-038_D1/spaceranger/filtered_feature_bc_matrix.h5; Input files updated by another job: output/V11L12-038_D1/transformed.png, output/V11L12-038_D1/transformed.csv\n",
"\n",
"[Mon Feb 24 10:13:40 2025]\n",
"Finished job 2.\n",
"3 of 5 steps (60%) done\n",
"Select jobs to execute...\n",
"\n",
"[Mon Feb 24 10:13:40 2025]\n",
"Job 1: Generating aggregated data.\n",
"Reason: Missing output files: output/V11L12-038_D1/spaceranger_aggregated/filtered_feature_bc_matrix.h5; Input files updated by another job: output/V11L12-038_D1/spaceranger/filtered_feature_bc_matrix.h5\n",
"\n",
"Calculating distances between MSI pixels and Visium spots...\n",
"Grouping MSI pixels by Visium spot...\n",
"Number of MSI pixels matching to Visium spots: 2388\n",
"Number of Visium spots matching to MSI pixels: 2371\n",
"Number of MSI pixels not matching to Visium spots: 696\n",
"Number of Visium spots not matching to MSI pixels: 2621\n",
"Aggregating MSI pixels per Visium spot by mean...\n",
"Copying Visium coordinates which map to some MSI pixels...\n",
"The new MSI coordinate file has been saved, containing 2371 spots/pixels.\n",
"Creating JSON file...\n",
"Copying tissue image...\n",
"Spaceranger 'spatial' folder content ready!\n",
"Reading MSI peak data file...\n",
"Preparing MSI peak data...\n",
"Writing new feature-barcode matrix data...\n",
"Spaceranger 'filtered_feature_bc_matrix' folder content ready!\n",
"Writing filtered_feature_bc_matrix.h5 file...\n",
"Spaceranger 'filtered_feature_bc_matrix.h5' file ready!\n",
"Mock Space Ranger output completed for sample at path output/V11L12-038_D1/spaceranger_aggregated/\n",
"[Mon Feb 24 10:16:55 2025]\n",
"Finished job 1.\n",
"4 of 5 steps (80%) done\n",
"Select jobs to execute...\n",
"\n",
"[Mon Feb 24 10:16:55 2025]\n",
"localrule all:\n",
" input: output/V11L12-038_D1/spaceranger_aggregated/filtered_feature_bc_matrix.h5\n",
" jobid: 0\n",
" reason: Input files updated by another job: output/V11L12-038_D1/spaceranger_aggregated/filtered_feature_bc_matrix.h5\n",
" resources: tmpdir=/var/folders/vb/gvs9l8xj3q10l987dbjg2bkw0000gq/T\n",
"\n",
"[Mon Feb 24 10:16:55 2025]\n",
"Finished job 0.\n",
"5 of 5 steps (100%) done\n",
"Complete log: .snakemake/log/2025-02-24T101307.169550.snakemake.log\n",
"```"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We will now describe each of the outputs in turn corresponding to each stage of the pipeline."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"#### i) Input checking"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The output from the input checking stage ( check_inputs ) is a file output/summary.csv which looks like this:"
]
},
{
"cell_type": "code",
"execution_count": 21,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"| \n", " | Subfolder | \n", "MSI image | \n", "Num pixels | \n", "Name of first pixel | \n", "Num peaks | \n", "Name of first peak | \n", "
|---|---|---|---|---|---|---|
| 0 | \n", "V11L12-038_D1 | \n", "NaN | \n", "3084 | \n", "pixel_14 | \n", "3658 | \n", "mz-108.901945 | \n", "
perform_coreg ) is:\n",
"\n",
"- *transformed.csv* containing the newly mapped MSI coordinates which will now be in a comparable space to the Visium data\n",
"- *transformed.png* which is the resulting H&E image after mapping (in this case it will just be the original Visium image but if we had an intermediate MSI image it would be the adjusted version of this after coregistration)\n",
"- *transformed_withCoords_VisiumHE.png* overlays the transformed MSI coordinates onto the Visium H&E image so you can assess whether the coordinates overlap as you expect"
]
},
{
"cell_type": "code",
"execution_count": 24,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"| \n", " | spot_id | \n", "x | \n", "y | \n", "
|---|---|---|---|
| 0 | \n", "pixel_14 | \n", "1735.755610 | \n", "271.844118 | \n", "
| 1 | \n", "pixel_15 | \n", "1711.620147 | \n", "271.225193 | \n", "
| 2 | \n", "pixel_16 | \n", "1687.484684 | \n", "270.606269 | \n", "
| 3 | \n", "pixel_17 | \n", "1663.349221 | \n", "269.987345 | \n", "
| 4 | \n", "pixel_18 | \n", "1639.213758 | \n", "269.368421 | \n", "
| \n", " | visium_spot | \n", "MSI_spot | \n", "
|---|---|---|
| 0 | \n", "CGTCTCTCGCCGAGGC-1 | \n", "MSI_V11L12-038_D1_pixel_14 | \n", "
| 1 | \n", "AACCCGATAGGGCTTC-1 | \n", "MSI_V11L12-038_D1_pixel_15 | \n", "
| 2 | \n", "TGATGCTCACGTAGTC-1 | \n", "MSI_V11L12-038_D1_pixel_16 | \n", "
| 3 | \n", "CGGGCTACTTAAATTG-1 | \n", "MSI_V11L12-038_D1_pixel_17 | \n", "
| 4 | \n", "CTGAAACGGCCCTCAG-1 | \n", "MSI_V11L12-038_D1_pixel_18 | \n", "