TSS profiles ChIP seq

Question

[Lila M]

Hi everyone! I'm performing some ChIP seq experiments and I want to read depth coverage over the region -500 to 3000pb of the TSS of all Ensemble protein coding transcript.

I have different approach, but as I am totally naive, I appreciate very much your advises.

a)First of all I've normalized the bam files using bamCoverage (deepTools), the output is a bedgraph, then I will do the anntotaion peak over TSS using ChIPseeker (R)

b)Also, I can do the peakannot over the bed file (after bam file conversion) and I will normalize the peakannot when I will plot it.

Other questions is, Is necessary make the peakcalling to calculate the read depth over TSS?

I'm looking forward to hear from you

Than you!!!

Answer 1

[Carlo Yague]

a) looks good to me.

b) I don't really understand what you want to achieve here. What is 'doing a peakannot' ? Do you refer to the annotatePeak function of the ChIPseeker R package ?

It is not necessary to make a peakcalling to calculate read depth over TSS. You only need genome coverage and TSS annotation.

[Lila M]

Hi again, As ChIPseeker gives me some error, you know any other way in which I can get the read depth coverage over TSS onces I have the normalized file ? Thanks!!

[Carlo Yague]

There are many many tools to do that. Since you are already familiar with deeptools, you could try their computeMatrix tool.

[Lila M]

Thank you, Do you know what file can I use as TSS.bed in computeMatrix? I downloaded one from USC but the program gives to me an error message

Warning: chrY:9748406-9748406 is an invalid BED interval! Ignoring it.

Thanks!

[Carlo Yague]

Could it be that your ChIP-seq data is from a female and that there is no Y ? This is a warning, not an error, and it could be ok to proceed if the non-Y coordinates are considered correct by the script.

More generally, check that the chromosome names in the bed file are the same as in the bedgraph file.

Carlo

Answer 2

[Lila M]

Thank you very much, Other question that I have, for bamCoverage I have two options, --nomralizeto1x --normalizeUsingRPKN

Is there any way in which I could specify that I want a normalization equal to 25 milliom reads?

Thanks!

[Carlo Yague]

a normalization equal to 25 million reads

Do you mean that you want to normalize as if you had 25 million reads mapped ? If so, you need to use a third option : --scaleFactor set to 25000000/number_of_reads_mapped. To calculate the number of reads mapped, you can just use samtools flagstat on the .bam file.

[Lila M]

One more question, as I can't introduce in scale factor 2500000/31567890, the argument should be: --sacaleFactor 25000000 or --scaleFactor 0.64?

Thanks!

[Carlo Yague]

25000000/31567890 = 0.7919 : this is your scaleFactor for this particular bam file :)

[Lila M]

last question, when I've tried to run reads <- import.bed(con="NGS974Norm")

NGSNorm is the bedgraph file, 4 cols and size ~ 400MB

cov <- coverage(reads) peakAnno <- annotatePeak(reads, tssRegion=c(-500, 2000), TxDb=TxDb.Hsapiens.UCSC.hg19.knownGene, annoDb="org.Hs.eg.db")

I have the next error Error in rsqlite_send_query(conn@ptr, statement) : out of memory

any ideas?

[Lila M]

Thank you! I'm working on that :)