Open qicaibiology opened 5 years ago
Hi ! can you paste the command you run for the DU analysis Thanks!
I used RNA-seq data from cells at different development stages and examine the different AS evetns:
My command for DU is this: du <- DUreportBinSplice( counts, targets, minGenReads = 10, minBinReads = 5, minRds = 0.05, contrast = NULL, forceGLM = FALSE, ignoreExternal = TRUE, ignoreIo = FALSE, ignoreI = FALSE, filterWithContrasted = FALSE)
my full command is follwoing: library(ASpli) library(GenomicFeatures) library(TxDb.Mmusculus.UCSC.mm10.knownGene) aTxDb <- TxDb.Mmusculus.UCSC.mm10.knownGene features <- binGenome(aTxDb) bamFiles <- c("DIVE8Rep1Aligned.sortedByCoord.out.bam", "DIVE8Rep2Aligned.sortedByCoord.out.bam", "DIVE8Rep3Aligned.sortedByCoord.out.bam", "DIVE8Rep4Aligned.sortedByCoord.out.bam", "DIV0Rep1Aligned.sortedByCoord.out.bam", "DIV0Rep2Aligned.sortedByCoord.out.bam", "DIV0Rep3Aligned.sortedByCoord.out.bam", "DIV7Rep1Aligned.sortedByCoord.out.bam", "DIV7Rep2Aligned.sortedByCoord.out.bam", "DIV7Rep3Aligned.sortedByCoord.out.bam", "DIV7Rep4Aligned.sortedByCoord.out.bam", "DIV7Rep5Aligned.sortedByCoord.out.bam") targets <- data.frame(row.names = c("DIVE8Rep1", "DIVE8Rep2", "DIVE8Rep3", "DIVE8Rep4", "DIV0Rep1", "DIV0Rep2", "DIV0Rep3", "DIV7Rep1", "DIV7Rep2", "DIV7Rep3", "DIV7Rep4", "DIV7Rep5"), bam = bamFiles, treat = c("DIVE8", "DIVE8", "DIVE8", "DIVE8", "DIV0", "DIV0", "DIV0", "DIV7", "DIV7", "DIV7", "DIV7", "DIV7"), stringAsFactors = FALSE) targets library(GenomicAlignments) bam <- loadBAM(targets) counts <- readCounts( features, bam, targets, cores = 8, readLength = 50L, maxISize = 500000) writeCounts(counts=counts, output.dir = "count_ESCDIV_4") writeRds(counts=counts, output.dir = "count_ESCDIV_4") as <- AsDiscover(counts, targets, features, bam, readLength=50L, threshold = 5) writeAS(as=as, output.dir="count_ESCDIV_4") du <- DUreportBinSplice( counts, targets, minGenReads = 10, minBinReads = 5, minRds = 0.05, contrast = NULL, forceGLM = FALSE, ignoreExternal = TRUE, ignoreIo = FALSE, ignoreI = FALSE, filterWithContrasted = FALSE) writeDU(du, output.dir="count_ESCDIV_4")
Thank you very much for quick response
Hi, Here you have the explanation in your case that you didnt define any specific contrast. Contrast is a vector representing the coefficients for each unique condition in the analysis. They are assigned in the order given by getConditions method The default value is NULL which correspond to a vector with value -1 in the first position, 1 in the second, and zero for the remaining positions. This is for the paired comparison of the second condition versus the first condition
So, in your case, logFC represent the log(mean coverage( reads )Condition 1 / mean coverage( reads)Ccondtion1)
thanks
by saying: "So, in your case, logFC represent the log(### mean coverage( reads )Condition 1 / mean coverage( reads)Ccondtion1)"
Do you mean: "mean coverage( reads )Condition 1 / mean coverage( reads)Ccondtion2" ?
Sorry, no. Cond2/Cond1
Got it!
Thanks,
So If I want to get the comparison between condition 3/condtion2 the Contrast=0?
Be careful with the order of your coefficients. Can you paste the output of: getConditions(targets) ?
I see. I ran it on server and did not write the output for the targets. I will do it step by step and see what is there.
Hi: Today I tried to compare AS events from 3 different developmental stages(progenitors, new_born neurons, Mature_neurons), I realized that I need to figure out what does logFC mean in the output file?
Thanks for any information.
Mike