statOmics / tradeSeq

TRAjectory-based Differential Expression analysis for SEQuencing data
Other
237 stars 27 forks source link

Unexpectedly high number of genes associated with pseudotime #213

Open nickhir opened 1 year ago

nickhir commented 1 year ago

Hello,

I am using tradeSeq downstream of slingshot. slingshot has identified 6 different trajectories, that upon closer inspection seem biologically meaningful. I wanted to use tradeSeq to identify genes whose expression is associated with the pseudotime of the respective trajectory. Afterwards I used the plotSmoothers function to visually inspect some of the results. Here I noticed for a couple of genes, that the associationTest detected a lot of genes as being associated with pseudotime, while visually it does not look like it. I have posted an extreme example below:

image

Here, associationTest indicates, that all lineages have a p.value of ~ 0, except the lineage 5. Visually it looks like only for lineage 4 there should be an association.

This is the code I used to run the analysis:

slingshot_sce  <- readRDS("SLINGSHOT_RESULT")

geneFilter <- apply(as.matrix(counts(slingshot_sce)),1,function(x){
  sum(x >= 10) >= 5
})

slingshot_sce <- slingshot_sce[geneFilter,]

counts <- as.matrix(BiocGenerics::counts(slingshot_no_end_isc_start))
pseudotime <- slingPseudotime(slingshot_no_end_isc_start, na = FALSE)
cellWeights <- slingCurveWeights(slingshot_no_end_isc_start)
sce <- fitGAM(counts = counts, pseudotime = pseudotime, cellWeights = cellWeights,
                 nknots = 7, verbose = FALSE)

association_result <- tradeSeq::associationTest(sce, lineages = TRUE)

I would greatly appreciate, if somebody has an idea why I observe such phenomena, or if I misunderstand the test that tradeSeq is doing.

Cheers!

koenvandenberge commented 1 year ago

Hi @nickhir

Thanks for reporting, this is quite peculiar indeed. Do you mind sharing the sce object after running fitGAM, and also let us know what version of tradeSeq you are using?

Thanks, Koen

nickhir commented 1 year ago

Sorry for getting back so late. I talked to the person who generated the data and he is fine with me sharing it. What do you think is the best way to transfer the data? The object is 60MB in size.

koenvandenberge commented 1 year ago

Hi @nickhir

Apologies for coming back to this so late. Feel free to send to koen.vdberge at gmail dot com if still relevant, with a description of which genes may be acting funny. You may also want to check out #209.