I'm trying the tuxedo protocol as well as analysis using DESeq2, I wanted to see if I can see in the difference of DE genes. DESeq2 uses something called count matrix ,which I generated using both ht-seq and featurecounts .Now I want to know how do I transform those counts to make it usable for DESeq2 analysis
I see the count table from ht-seq and featurecount whats the difference between both the output.? for ht-seq I see a long list of ensemble ID with the counts where as for featurecount there is lot more information .
I have two 4 samples Control vs Test with replicate of each of them ,so I have the all together 4 count table.
So how can i use those counts into my DESeq2 ?
Any help and suggestion would be highly appreciated.
•
written
by
1769mkc
2 answers
The output of both HT-seq and featurecounts can almost be fed directly in DESeq2.
- for FeatureCounts :
Output looks like :
Geneid Chr Start End Strand Length cond1 cond2
SPBC460.05 I 16470 18062 + 1593 1 12 24
SPBC460.02c II;II 8856;9651 9365;9803 -;- 663 329
SPAC212.11 I 1 5662 - 5662 0 0 0
First you need to remove 5 firsts columns and name the rows as gene names
cond1 cond2
SPBC460.05 12 24
SPBC460.02c 663 329
SPAC212.11 0 0
then just use the DESeqDataSetFromMatrix function of DESeq2 to convert it.
- for HTseqCounts : see comments below
u mean HTseq right ?
okay so if I remove the first 5 columns then I can use that as input to DESeq2 , I will do that But i have 4 featurecount data sets so can i just remove first 5 columns and simply join them?
i would like to do the same with HTseqCounts data which gives just enslebleID and count value so how can I put all my samples in a single count matrix ?
I would use cut and paste in the terminal. In this case, I read the 16 counts files in the counts directory (I have 16 conditions), retrieve the useful columns and save the final matrix that I can load into R>DESeq.
paste counts/* | cut -f 1,2,4,6,8,10,12,14,16,18,20,22,24,26,28,30,32 > counts/merged_counts.mat
Note that the paste process the files in alphanumerical order, which may not be the desired order. Reorering the columns in R afterward might be needed.
I just pasted code I used in a project with 16 samples. U can adapt it for 4 conditions.
paste cond1 cond2 cond3 cond4 | cut -f 1,2,4,6,8 > merged.txt
That can be a good advice except if your bam file is made of paired-end reads and is not name-sorted. FeatureCounts can sort the bam files automatically but it can be very time-consuming.
There are many different ways to do such basic data manipulation. I like to do it in the terminal when I can because it can be faster, requires less memory and can be syntaxically shorter. One alternative using only R would require to read all the files, save them into memory, merge them, then cut the unwanted columns and clear the original tables from memory. For me it is less efficient (especially with 16 samples), but it also works.
You can also use featureCounts directly from R, which will give you an R object which you can use directly as input to DESeq2 and similar tools, see the Rsubread package. This avoids tampering around with inputfiles.
Log in to answer this question.
-
Batch correction for Nanopore RNAseq
written by Nico80 •
-
DESeq2 for different design and normalized counts
written by t.ru •
-
remove rRNA reads before DE analysis?
written by anubhav9 •
-
DESeq2 input file from featureCount
written by Bioinf_Questions •
-
Dual RNA-seq differential gene expression analysis
written by Karla •
-
How can I get a count of mRNA + lncRNAs from Gencode + lncipedia?
written by njk639 •
-
deseq2 read counts vs featurecount read counts
written by Adeler001 •
-
STAR-HT-Seq/featureCount got much more gene expression count than RSEM did
written by CY
-
cufflink issues in Tuxedo protocol
written by 1769mkc
-
Rna -Seq cummeRbund plot issues
written by 1769mkc
