Elementolab/ChIPseeqerGetReadDensityProfiles
From Icbwiki
(Difference between revisions)
| Revision as of 20:13, 20 October 2010 Eug2002 (Talk | contribs) ← Previous diff |
Revision as of 16:20, 8 November 2010 Eug2002 (Talk | contribs) Next diff → |
||
| Line 21: | Line 21: | ||
| -format eland -fraglen 0 > bcl6peaks.txt.centered2000.H3K4me3H3K9Ac | -format eland -fraglen 0 > bcl6peaks.txt.centered2000.H3K4me3H3K9Ac | ||
| - | Optiions are: | + | Options are: |
| -format STR [bed sam eland] | -format STR [bed sam eland] | ||
Revision as of 16:20, 8 November 2010
Back to Elementolab/ChIPseeqer_Tutorial
ChIPseeqerGetReadDensityProfiles
To run the tools directly from any folder, you need to add the $CHIPSEEQERDIR and $CHIPSEEQERDIR/SCRIPTS to your $PATH variable. Read How to set the CHIPSEEQERDIR variable.
- 1. If you haven't done it yet, run CS on ChIP-seq data to find peaks. Ignore this step if you have already located peaks.
ChIPseeqer.bin -chipdir BCL6CHIP -inputdir BCL6INPUT -t 15 -fold_t 2 -outfile bcl6peaks.txt
- 2. Run script to extract regions centered on peak summits (2kb windows used here, ie 1kb on each side)
extract_regions_around_peak_summits.pl --peakfile=bcl6peaks.txt --w=2000
This will create a file called bcl6peaks.txt.centered2000
- 3. Extract read densities from other ChIP-seq dataset (eg histone modification; input is ignored here only)
ChIPseeqerGetReadDensityProfiles.bin -intervals bcl6peaks.txt.centered2000 -chipdir H3K4ME3H3K9ACCHIP \ -format eland -fraglen 0 > bcl6peaks.txt.centered2000.H3K4me3H3K9Ac
Options are:
-format STR [bed sam eland] -fraglen INT [0 = no read extension otherwise extend to specified value] -nummtnorm INT [0/1, if 1, perform normalization by number of mapped nucleotides ] -uniquereads INT [0/1 if 1 collapose clonal reads; 1 = recommended for TF and histone modification, not for nucleosome positioning ] -ws INT [the window size, can be 10, 100 etc. Default is 10]
- 4. Load read densities in R, average them, and plot them
m <- read.table("bcl6peaks.txt.centered2000.H3K4me3H3K9Ac", row.names=1)
par(cex=1.3)
plot(10*(1:200)-1000, colMeans(m), type="l", lwd=5, col="red", ylab="Average H3K4me3+H3K9Ac read density",
xlab="Distance to maximum BCL6 peak height (bp)")
dev.copy2eps(file="H3K9Ac-H3K4me3-around-BCL6-peaks.eps")
