Dinh/Dinh 2012/NOTES/2012-2-13: Difference between revisions
Jump to navigation
Jump to search
>Dinh No edit summary |
>Dinh No edit summary |
||
Line 36: | Line 36: | ||
java -Xmx2G -jar /usr/local/fp4/FindPeaks.jar -input tmp -aligner bed -output . -wig_step_size 1 - dist_type 0 350 -subpeaks 0.5 -landerwaterman 0.001 -name $f -one_per | java -Xmx2G -jar /usr/local/fp4/FindPeaks.jar -input tmp -aligner bed -output . -wig_step_size 1 - dist_type 0 350 -subpeaks 0.5 -landerwaterman 0.001 -name $f -one_per | ||
done; | done; | ||
* Number of peaks called for each dataset: | |||
41109 normalized.R.cac.BED_fixed_350_subpeaks.peaks | |||
820 normalized.R.dnmt3a.BED_fixed_350_subpeaks.peaks | |||
12382 normalized.R.fc.BED_fixed_350_subpeaks.peaks | |||
14825 normalized.R.hmc.BED_fixed_350_subpeaks.peaks | |||
28826 normalized.R.mc.BED_fixed_350_subpeaks.peaks | |||
1301 normalized.R.tdg.BED_fixed_350_subpeaks.peaks | |||
2177 normalized.R.tet2.BED_fixed_350_subpeaks.peaks | |||
61289 normalized.S.cac.BED_fixed_350_subpeaks.peaks | |||
1431 normalized.S.dnmt3a.BED_fixed_350_subpeaks.peaks | |||
9107 normalized.S.fc.BED_fixed_350_subpeaks.peaks | |||
51385 normalized.S.hmc.BED_fixed_350_subpeaks.peaks | |||
52640 normalized.S.mc.BED_fixed_350_subpeaks.peaks | |||
1222 normalized.S.tdg.BED_fixed_350_subpeaks.peaks | |||
1090 normalized.S.tet2.BED_fixed_350_subpeaks.peaks | |||
* Things to note: For all properly paired reads there are fewer peaks for ChIP data than for dIP data. This is expected since the ChIP are peaks for DNA binding proteins (they are specific). Whereas for the dIP data, it is possible to have more signals since there are many location for cytosines to be oxydized. | |||
* I also used the same normalization and peak calling method (except with dist_type 0 110) for the improperly paired reads. | |||
16759 normalized.R.cac.BED_fixed_110_subpeaks.peaks | |||
95 normalized.R.dnmt3a.BED_fixed_110_subpeaks.peaks | |||
4250 normalized.R.fc.BED_fixed_110_subpeaks.peaks | |||
1720 normalized.R.hmc.BED_fixed_110_subpeaks.peaks | |||
2047 normalized.R.mc.BED_fixed_110_subpeaks.peaks | |||
94 normalized.R.tdg.BED_fixed_110_subpeaks.peaks | |||
130 normalized.R.tet2.BED_fixed_110_subpeaks.peaks | |||
16148 normalized.S.cac.BED_fixed_110_subpeaks.peaks | |||
112 normalized.S.dnmt3a.BED_fixed_110_subpeaks.peaks | |||
1688 normalized.S.fc.BED_fixed_110_subpeaks.peaks | |||
1097 normalized.S.hmc.BED_fixed_110_subpeaks.peaks | |||
1212 normalized.S.mc.BED_fixed_110_subpeaks.peaks | |||
109 normalized.S.tdg.BED_fixed_110_subpeaks.peaks | |||
104 normalized.S.tet2.BED_fixed_110_subpeaks.peaks | |||
* Things to note: much fewer peaks found for ChIP data, which is expected, since DNMT3a, TDG, and Tet2 are not known to bind in repeats. Up to 10x greater signals are observed for caC marks in both resting and stimulated neurons. | |||
* I need to validate the peaks called by going back to the raw data. Using FindPeak 4.0, I converted un-normalized reads to wig files. | |||
java -Xmx2G -jar /usr/local/fp4/FindPeaks.jar -input ../sorted.Indx1.BED -aligner bed -output R.cac/ -duplicatefilter -wig_step_size 10 -max_pet_size 1000 -dist_type 3 -name R.cac -one_per |
Revision as of 00:38, 15 February 2012
24 ChIP-seq - Stimulated and Resting Neurons
- Analysis at repeats.
Downloaded RepeatMasker track fro mm9. Then filtered and get Non-overlapping windows: 2072 awk '{if($12 ~ /DNA/) print $0}' RepeatMasker_mm9.txt | ./getNonOverlapping.pl > DNA_repeats.nonoverlapping & 2073 awk '{if($12 ~ /LINE/) print $0}' RepeatMasker_mm9.txt | ./getNonOverlapping.pl > LINE_repeats.nonoverlapping & 2074 awk '{if($12 ~ /SINE/) print $0}' RepeatMasker_mm9.txt | ./getNonOverlapping.pl > SINE_repeats.nonoverlapping & 2075 awk '{if($12 ~ /LTR/) print $0}' RepeatMasker_mm9.txt | ./getNonOverlapping.pl > LTR_repeats.nonoverlapping & 2076 awk '{if($12 ~ /Satellite/) print $0}' RepeatMasker_mm9.txt | ./getNonOverlapping.pl > Satellite_repeats.nonoverlapping & File:GetNonOverlapping.txt
- Note that some repeat segments overlaps between these difference repeat class.
cat *_repeats.nonoverlapping | sort -k1,1 -k2,3n > sorted.Repeats (3,188,316 entries) ./removeOverlapping < sorted.Repeats > tmp mv tmp sorted.Repeats (3,187,848 entries) File:RemoveOverlapping.txt
- Remove chr*_random entries
awk '{if($1 !~ /_random/) print $0}' sorted.Repeats > sorted.Repeats_norandom
- Final repeat regions file: 3,127,206 entries
Started with 3,188,316 non-overlapping entries for DNA, LINE, SINE, LTR, and Satellite repeats Then removed entries/segments which have over-lapping labels -> 3,187,848 entries left Then removed entries on chr*_random -> 3,127,206 entries -> Non-overlapping repeat regions
- Strand specific normalization with both IgG and Input controls. Normalization of read depth is taken care of by regression model.
for chr in chr1 chr2 chr3 chr4 chr5 chr6 chr7 chr8 chr9 chr10 chr11 chr12 chr13 chr14 chr15 chr16 chr17 chr18 chr19 chrX chrY chrM do ./getCHR.pl $chr < Indx1.BED | sort -k2,3n > treat.BED ./getCHR.pl $chr < Indx7.BED | sort -k2,3n > input.BED ./getCHR.pl $chr < Indx4.BED | sort -k2,3n > iggcl.BED ftumch2 -i treat.BED -c input.BED iggcl.BED -o normalized.BED -rl 350 -cs 100000 -ws 5 -iv 5 cat normalized.BED >> normalized.R.cac.BED done >> log.Rcac * repeat this for all the samples.
- Peak calling using FindPeak 4.0. The reason why MACS1.4 wouldn't work because MACS14 "shifts" read toward the center when reads are pileup on the forward strand and reverse strand. This shifting model wouldn't work for finding enrichment of certain regions in the genome (unlike motif finding, the size of the region is unknown and may be very large).
for f in normalized.*BED do awk '{print $1"\t"$2"\t"$3"\t"$4"\t0\t"$6}' $f > tmp java -Xmx2G -jar /usr/local/fp4/FindPeaks.jar -input tmp -aligner bed -output . -wig_step_size 1 - dist_type 0 350 -subpeaks 0.5 -landerwaterman 0.001 -name $f -one_per done;
- Number of peaks called for each dataset:
41109 normalized.R.cac.BED_fixed_350_subpeaks.peaks 820 normalized.R.dnmt3a.BED_fixed_350_subpeaks.peaks 12382 normalized.R.fc.BED_fixed_350_subpeaks.peaks 14825 normalized.R.hmc.BED_fixed_350_subpeaks.peaks 28826 normalized.R.mc.BED_fixed_350_subpeaks.peaks 1301 normalized.R.tdg.BED_fixed_350_subpeaks.peaks 2177 normalized.R.tet2.BED_fixed_350_subpeaks.peaks 61289 normalized.S.cac.BED_fixed_350_subpeaks.peaks 1431 normalized.S.dnmt3a.BED_fixed_350_subpeaks.peaks 9107 normalized.S.fc.BED_fixed_350_subpeaks.peaks 51385 normalized.S.hmc.BED_fixed_350_subpeaks.peaks 52640 normalized.S.mc.BED_fixed_350_subpeaks.peaks 1222 normalized.S.tdg.BED_fixed_350_subpeaks.peaks 1090 normalized.S.tet2.BED_fixed_350_subpeaks.peaks
- Things to note: For all properly paired reads there are fewer peaks for ChIP data than for dIP data. This is expected since the ChIP are peaks for DNA binding proteins (they are specific). Whereas for the dIP data, it is possible to have more signals since there are many location for cytosines to be oxydized.
- I also used the same normalization and peak calling method (except with dist_type 0 110) for the improperly paired reads.
16759 normalized.R.cac.BED_fixed_110_subpeaks.peaks 95 normalized.R.dnmt3a.BED_fixed_110_subpeaks.peaks 4250 normalized.R.fc.BED_fixed_110_subpeaks.peaks 1720 normalized.R.hmc.BED_fixed_110_subpeaks.peaks 2047 normalized.R.mc.BED_fixed_110_subpeaks.peaks 94 normalized.R.tdg.BED_fixed_110_subpeaks.peaks 130 normalized.R.tet2.BED_fixed_110_subpeaks.peaks 16148 normalized.S.cac.BED_fixed_110_subpeaks.peaks 112 normalized.S.dnmt3a.BED_fixed_110_subpeaks.peaks 1688 normalized.S.fc.BED_fixed_110_subpeaks.peaks 1097 normalized.S.hmc.BED_fixed_110_subpeaks.peaks 1212 normalized.S.mc.BED_fixed_110_subpeaks.peaks 109 normalized.S.tdg.BED_fixed_110_subpeaks.peaks 104 normalized.S.tet2.BED_fixed_110_subpeaks.peaks
- Things to note: much fewer peaks found for ChIP data, which is expected, since DNMT3a, TDG, and Tet2 are not known to bind in repeats. Up to 10x greater signals are observed for caC marks in both resting and stimulated neurons.
- I need to validate the peaks called by going back to the raw data. Using FindPeak 4.0, I converted un-normalized reads to wig files.
java -Xmx2G -jar /usr/local/fp4/FindPeaks.jar -input ../sorted.Indx1.BED -aligner bed -output R.cac/ -duplicatefilter -wig_step_size 10 -max_pet_size 1000 -dist_type 3 -name R.cac -one_per