I need to count the expression of a particular gene from rna-seq data. I've mapped the reads using tophat and then I use cuffdiff for analyse the differential gene expression.
I want to see if also I use a raw count I obtain the same results. I mean with cufflinks I see the sample A have low quantity of a particula gene, but using the raw estimation. I see the oppostive results using my script.
I merged the hg19 to obatain a one dimension for each gene and the I use intersectbed.
Then I use a script in perl with I make this procedure:
(Count reads for a each gene/kb of that gene) /number of fragment.
example( reads count on TP53/kb dimension TP53/number of fragment mapped)
What it is wrong?
Thanks for any help!
thanks for any suggestion?
I want to see if also I use a raw count I obtain the same results. I mean with cufflinks I see the sample A have low quantity of a particula gene, but using the raw estimation. I see the oppostive results using my script.
I merged the hg19 to obatain a one dimension for each gene and the I use intersectbed.
Code:
intersectBed -bed -abam example.bam -b ../RefSeq/hg19_merged.bed -wa -wb | cut -f 16 | sort | uniq -c > example.counts.txt
(Count reads for a each gene/kb of that gene) /number of fragment.
example( reads count on TP53/kb dimension TP53/number of fragment mapped)
What it is wrong?
Thanks for any help!
thanks for any suggestion?
Comment