Compare gene expression across two experiments

[adilysyn]

I have two STOmics experiments - control and mutant.
I would like to compare gene expression across the two different chips.
Can you please instruct me how to integrate them first using Harmony batch-correction method, and how to perform clustering/differential expression on the integrated level for a region of interest within clusters?

Thank you in advance for your help,
Adily

[Jinghao-Chen]

Hi, it is recommended that you try using Stereopy for batch effect correction and clustering analysis, as Harmony has been pre-integrated into Stereopy.
https://stereopy.readthedocs.io/en/latest/Tutorials%28Multi-sample%29/Batch_Effect.html#Integrating

[adilysyn]

Thank you very much!
So after the integration step (ms_data.tl.batches_integrate) and clustering (ms_data.tl.leiden), how I can compare the gene expression between the control and mutant?

Best regards,
Adily

[Jinghao-Chen]

Thank you very much! So after the integration step (ms_data.tl.batches_integrate) and clustering (ms_data.tl.leiden), how I can compare the gene expression between the control and mutant?

Best regards, Adily

Hi, if you want to compare the gene expression between two chips, you could refer to the following code. For example, I set cluster_res_key='batch' in my code, which means performing differential analysis based on different batches (chips). The case_groups parameter specifies the mutant samples labeled as '1', and control_groups specifies the control samples labeled as '0'.

ms_data.tl.find_marker_genes(cluster_res_key='batch',case_groups='1',control_groups='0',
                             res_key='batch_diff_1-vs-0',
                             use_highly_genes=False, use_raw =False)

If you want to compare the differences of the same cluster between two chips, you could refer to the following code. Here, I use batch-leiden to represent the same cluster across different batches (chips). In my code, I compared the differences between cluster 10 of chip '0' and cluster 10 of chip '1'.

new_label='batch-leiden'

ms_data.obs[new_label] = ms_data.obs['batch'].astype(str) + '-' + ms_data.obs['leiden'].astype(str)
ms_data.merged_data.tl.result[new_label] = ms_data.merged_data.tl.result['batch'].astype(str) + '-' + ms_data.merged_data.tl.result['leiden'].astype(str)

ms_data.tl.find_marker_genes(cluster_res_key='batch-leiden',case_groups='0-10',control_groups='1-10',
                             use_highly_genes=False, use_raw =False)

Finally, the following code would be helpful for obtaining the results of differential expression analysis.

ms_data.plt.marker_genes_scatter(res_key = 'batch_diff_1-vs-0')
ms_data.tl.result['scope_[0,1]']['batch_diff_1-vs-0']

[adilysyn]

Thank you for the detailed response!
My problem is that the regions of interest cluster differently in the mutant and in the control. For example: 15, 9, 6, 2, 1 in the control vs 14, 9, 16 in the mutant. How I can compare gene expression in my ROI?

Best regards,
Adily

[Jinghao-Chen]

Thank you for the detailed response! My problem is that the regions of interest cluster differently in the mutant and in the control. For example: 15, 9, 6, 2, 1 in the control vs 14, 9, 16 in the mutant. How I can compare gene expression in my ROI?

Best regards, Adily

Hi, you may try the following code: compile the clusters of interest from the same sample into a list as the input. The default value for res_key is marker_genes; alternatively, you may assign it to other value.

ms_data.tl.find_marker_genes(cluster_res_key='batch-leiden',case_groups=['1-14','1-9','1-16'],control_groups=['0-15','0-9','0-6','0-2','0-1'],use_highly_genes=False, use_raw =False)

[adilysyn]

Thank you!

Based on my understanding, there are two possible approaches:

Option 1:
Perform batch correction on the GEM files from each experiment, and then use your last command with a list of clusters for each batch.
At this point, I would like to clarify one thing: does ms_data retain the original clusters, or are they recalculated by the following command?

ms_data.tl.leiden(
    scope=slice_generator[:],
    mode='integrate',
    neighbors_res_key='neighbors_integrated',
    res_key='leiden'
)

Option 2:
Perform batch correction on the GEM files derived from the StereoMap lasso ROIs, and then run the command you shared initially:

ms_data.tl.find_marker_genes(
    cluster_res_key='batch',
    case_groups='1',
    control_groups='0',
    res_key='batch_diff_1-vs-0',
    use_highly_genes=False,
    use_raw=False
)

Please let me know which option you would recommend, or if there is a preferred workflow between the two.

Best regards,
Adily

[Jinghao-Chen]

@adilysyn Hi, for Option 1, clustering will be recalculated in the code. Since batch effects have been removed, the two chips will be integrated for reclustering; in addition, GEM files themselves do not contain clustering result information and thus require manual supplementation of relevant information.

If the goal is to compare differences after integrating the GEM files of the two ROIs, Option 2 is more appropriate.

Best regards,
Jinghao