eRegulon difference between two groups in specific cell type

[ianhsu99]

Hi,

First of all, thanks for this excellent software.

I have a processed dataset with 13 cell types and two groups (disease and control). I want to see the regulon difference between the two groups in Fibroblasts.

I wondered if it makes sense to run the SCENIC+ on Fibroblasts and use groups (e.g. disease and control) as the variable for the run_scenicplus() function.

Cheers,

Ian

[SeppeDeWinter]

Hi @ianhsu99

Thank you for your kind words.

I converted the issue to a discussion.

Yes that does make sense. This variable does not affect the SCENIC+ inference itself but it does influence the downstream analysis that is done (i.e. RSS calculation and calculating the correlation between TF expression and eRegulon AUC values).

Best,

Seppe

[mary77]

Hi Seppe,

I have the same situation. My question is that at which step should we subset our data? Right before running scenicplus or at the beginning of the pipeline?

[SeppeDeWinter]

It depends a bit on your research question.

SCENIC+ depends on having a good contrast, i.e. you need a population of cells with enough heterogeneity to be able to find biologically representative motifs for each cell state. However subsetting your cells before hand might give you more detailed differences, with the cost of missing the regulons which are nog shared by all cells.

I hope this brings some clarity.

Best,

Seppe

[zfliu01]

Hi Seppe,

Follow up questions: so in practice if we're asking the same question that Ian was asking, how do we look at the cell type specific difference in eRegulons (i.e. differences in fibroblast between control and diseased samples)? If we set diseased vs control in run_scenicplus(), are we going to get eRegulon across all cell types that is different between control and diseased samples?

We're trying to look at difference between normal and diseased sample but we only care about one specific cell type, and we're unsure how to proceed.

Thank you so much for all the work and help!

Best,
Jeff

[SeppeDeWinter]

Hi @zfliu01

Yes, I think it is reasonable to try.
Already with the pycisTopic analysis you should be able to see wether you are picking up any differences (are there any differential topics? DARs? Are those enriched for different motifs?).

Good luck with the analysis.

All the best,

Seppe

p.s. thanks for joining the webinar, looking forward!

[zfliu01]

Hi Seppe,

Thank you so much for the developmental branch, which makes it a lot more streamlined! I was hoping to give a faster response but unfortunately it seems like now after pycistopic analysis it only outputs bed files.

I planned to do three comparisons to see which one makes sense:

1. using all samples and all cell types, contrast diseased vs normal state (which is probably going to pick up all cells, but hepatocytes are a major portion of this dataset).
2. separate diseased and normal samples, contrast hepatocytes vs other cells, and compare the eRegulon afterward.
3. subset the initial scRNA-seq dataset for only hepatocytes, and then contrast diseased vs normal state.

I was able to go through the first one (after 2 weeks of battling with our cluster admin; the peak usage is 1.1T ram since this is a 70K+ cell dataset). I saw our expected TF showing up multiple times in the cistrome HTML result, but not in the DEM HTML result, albeit we do see its paralog and cofactor). But both HTMLs contain a huge list. Do you think it is a good sign? Also, looking at the final mudata and comparing diseased vs normal does give us a high direct_++ (0.57) and direct+- score (0.55) in the diseased sample and high direct-+ (0.63) and direct-_- (0.60) in the normal sample for the TF that we're looking at.

I will also try the other two in the coming days and see how they differ. But after running prerequisites before scenic+ it seems like some of the bed files are empty. I think scenic+ wouldn't take empty bed files? Would changing to some less stringent parameters for pycistopic analysis work? If so, what are the parameters I should look to change?

I also tried to analyze the mudata object at the end, but I'm not sure how to repeat a lot of the analysis that was in the original tutorial since the data structure is now different. Specifically, I was trying to create graphs of interaction network. Let me know if you think this should be asked in a new question/issue.

Best,
Jeff

[SeppeDeWinter]

Hi @zfliu01

Thank you for the nice feedback.

it seems like some of the bed files are empty.
For what step is this?

I also tried to analyze the mudata object at the end, but I'm not sure how to repeat a lot of the analysis that was in the original tutorial since the data structure is now different.

I added a function last week to convert the mudata object to the old style scenicplus_obj: f216792, this you can use with the old functions.

All the best,

Seppe

[zfliu01]

Hi Seppe,

Thank you so much for your help along the way! We're putting the analysis together and preparing a manuscript now.

In the end, I did my analysis by subsetting only hepatocytes and compared diseased vs control cells. As we're looking back at the eRegulon that our TF is in, we noticed that the gene list of the eRegulon doesn't overlap previously experimentally validated target gene list. Obviously there are multiple possible reasons such as cell type differences. One of our colleagues is concerned that we might be looking at false positives. I think what we're seeing in the eRegulon is just differences coming out from comparing disease vs normal, since canonical targets are going to be activated by the TF in both diseased and normal cells, subtracted as the background.

I'm wondering if you have any input on this question. Or is there anything that we should check to make sure our analysis is correct before submitting it to the journal? Any insight will be greatly appreciated.

Best,
Jeff

[uqnsarke]

following