Dinh 2011/NOTES/2011-11-14
Jump to navigation
Jump to search
Yi Zhang UNC-5caC/5fC[edit]
8,506,009 reads for 5fC 39,354,845 reads for 5caC 12,518,460 reads for IgG 25,206,212 reads for Input
Mapped with bowtie: (# alignments)[edit]
bowtie -n 2 --best --strata -l 28 -k 1 -m 1 -p 16 --phred64-quals --un $f.unmappable.txt /media/1TB_store1/Rui_Scratch/iGenome/Mus_musculus/UCSC/mm9/Sequence/BowtieIndex/genome $f $f.bowtie.map 2> $f.log > YZ_UNC-5caC-DIP.fastq.log <== # reads with at least one reported alignment: 20985742 (53.32%) # reads that failed to align: 2218974 (5.64%) # reads with alignments suppressed due to -m: 16150129 (41.04%) Reported 20985742 alignments to 1 output stream(s)
> YZ_UNC-5fC-DIP.fastq.log <== # reads with at least one reported alignment: 3740642 (43.98%) # reads that failed to align: 970771 (11.41%) # reads with alignments suppressed due to -m: 3794596 (44.61%) Reported 3740642 alignments to 1 output stream(s)
> YZ_UNC-IgG-DIP.fastq.log <== # reads with at least one reported alignment: 8023681 (64.09%) # reads that failed to align: 1072660 (8.57%) # reads with alignments suppressed due to -m: 3422119 (27.34%) Reported 8023681 alignments to 1 output stream(s)
> YZ_UNC-Input.fastq.log <== # reads with at least one reported alignment: 16508323 (65.49%) # reads that failed to align: 3733811 (14.81%) # reads with alignments suppressed due to -m: 4964078 (19.69%) Reported 16508323 alignments to 1 output stream(s)
Balance hits[edit]
2191 head -3740642 YZ_UNC-5caC-DIP.fastq.bowtie.map > sampled_YZ_UNC-5caC-DIP.fastq.bowtie.map & 2192 head -3740642 YZ_UNC-IgG-DIP.fastq.bowtie.map > sampled_YZ_UNC-IgG-DIP.fastq.bowtie.map & 2193 head -3740642 YZ_UNC-Input.fastq.bowtie.map > sampled_YZ_UNC-Input-DIP.fastq.bowtie.map &
Peak calling using MACS[edit]
macs14 -t sampled_YZ_UNC-5caC-DIP.fastq.bowtie.map -c sampled_YZ_UNC-IgG-DIP.fastq.bowtie.map -f BOWTIE -g mm -n 5caC-vsIgG -w --call-subpeaks > macs14_5caCvsIgG_balanced & macs14 -t sampled_YZ_UNC-5caC-DIP.fastq.bowtie.map -c sampled_YZ_UNC-Input-DIP.fastq.bowtie.map -f BOWTIE -g mm -n 5caC-vsInput -w --call-subpeaks > macs14_5caCvsInput_balanced & macs14 -t YZ_UNC-5fC-DIP.fastq.bowtie.map -c sampled_YZ_UNC-IgG-DIP.fastq.bowtie.map -f BOWTIE -g mm -n 5fC-vsIgG -w --call-subpeaks > macs14_5fCvsIgG_balanced & macs14 -t YZ_UNC-5fC-DIP.fastq.bowtie.map -c sampled_YZ_UNC-Input-DIP.fastq.bowtie.map -f BOWTIE -g mm -n 5fC-vsInput -w --call-subpeaks > macs14_5fCvsInput_balanced &
Remapping using BWA to check enrichment in repeats[edit]
- For ChIP-seq experiments, it is okay to use reads mapped to multiple locations, since a multi-mapped read can be split to its multiple locations, and we have data from the control experiment to normalize the enrichment results.
- I decided to re-map using BWA since it have better sensitivity and I think better sensitivity is good for mapping to repetitive regions of the genome.
- Use reads which map uniquely in one experiment -> unique genome test
- Use reads which map up to many different regions in the genome -> repetitive genome test
- Building mm9 bwa index: I copied genome.fa from Rui's iGenome folder:
cp /media/1TB_store1/Rui_Scratch/iGenome/Mus_musculus/UCSC/mm9/Sequence/WholeGenomeFasta/genome.fa BWAIndex/ /usr/local/bwa-0.5.9/bwa index -p mm9 -a bwtsw genome.fa &
- Understanding BWA tags:
- By default, aln finds all best hits and all suboptimal hits (within the given constraints)
- Using samse, each read will be aligned and a best hit will be reported at random with the XA tag which will list up to a maximal number of alternative hits (default 3). Reads which map to more than (3) optimal+suboptimal locations will not have the XA tag but a best hit will be reported at random.
- XT tag will indicate whether a read mapped uniquely or in repeats.
- X0 tag will indicate how many best hits were found.
- Understanding BWA flags:
awk '{print $2;}' *.sam | sort -u
- 0 = mapped, hit given on forward strand
- 4 = unmapped read
- 16 = mapped, hit given on reverse strand
- 20 = 4+16, read is unmapped, hit given on reverse strand but the read falls off the end of the chromosome.
- Piping script:
for f in *.fastq do /usr/local/bwa-0.5.9/bwa aln -l 28 -k 2 -t 16 -M 15 -I /media/disk-1/Dinh_Scratch/mm9Annotations-iGenomeUCSC/BWAIndex/mm9 $f > $f.sai /usr/local/bwa-0.5.9/bwa samse -n 1000 /media/disk-1/Dinh_Scratch/mm9Annotations-iGenomeUCSC/BWAIndex/mm9 $f.sai $f | ./parseIPtags-bwaSAM.pl $f.sam done;
for f in *.sam do /usr/local/samtools-0.1.17/samtools view -uSt /media/disk-1/Dinh_Scratch/mm9Annotations-iGenomeUCSC/BWAIndex/genome.fai $f | /usr/local/samtools-0.1.17/samtools sort - sorted.$f /usr/local/samtools-0.1.17/samtools rmdup -S sorted.$f.bam rmdup.sorted.$f.bam /usr/local/samtools-0.1.17/samtools flagstat rmdup.sorted.$f.bam > rmdup.$f.LOG /usr/local/samtools-0.1.17/samtools flagstat sorted.$f.bam > sorted.$f.LOG done;
File:ParseIPtags-bwaSAM.txt
- RepeatTags (after rmdup)
5caC 171,839 (198,681) 5fC 177,485 (73,294) IgG 12,904 (12,113) Input 43,686 (36,396)
- My original repeatTags files were really large, I have to figure out why samtools sort removed many lines of repeatTags:
The number of lines in each file ( I checked that there are no empty lines printout by the perl script. 19139320 YZ_UNC-5caC-DIP.fastq.repeatTags.sam 4854616 YZ_UNC-5fC-DIP.fastq.repeatTags.sam 4079042 YZ_UNC-IgG-DIP.fastq.repeatTags.sam 6198485 YZ_UNC-Input.fastq.repeatTags.sam
- UniqueTags (after rmdup)
5caC 22,864,806 (20,399,989) 5fC 4,147,394 (3,546,988) IgG 8,010,605 (7,419,742) Input 16,825,683 (16,295,873)
Balance unique hits[edit]
- Use sort -R to get random hits
2611 samtools view rmdup.sorted.5fC.uniqueTags.bam > rmdup.sorted.5fC.uniqueTags.sam & 2615 samtools view rmdup.sorted.5caC.uniqueTags.bam | sort -R | head -3546988 > rmdup.sorted.5caC.uniqueTags.sam & 2616 samtools view rmdup.sorted.IgG.uniqueTags.bam | sort -R | head -3546988 > rmdup.sorted.IgG.uniqueTags.sam & 2617 samtools view rmdup.sorted.Input.uniqueTags.bam | sort -R | head -3546988 > rmdup.sorted.Input.uniqueTags.sam &
Peak calling of unique hits using MACS[edit]
macs14 -t rmdup.sorted.5caC.uniqueTags.sam -c rmdup.sorted.IgG.uniqueTags.sam -f SAM -g mm -n 5caC-vsIgG -w --call-subpeaks > macs14_5caCvsIgG_balanced & macs14 -t rmdup.sorted.5caC.uniqueTags.sam -c rmdup.sorted.Input.uniqueTags.sam -f SAM -g mm -n 5caC-vsInput -w --call-subpeaks > macs14_5caCvsInput_balanced & macs14 -t rmdup.sorted.5fC.uniqueTags.sam -c rmdup.sorted.IgG.uniqueTags.sam -f SAM -g mm -n 5fC-vsIgG -w --call-subpeaks > macs14_5fCvsIgG_balanced & macs14 -t rmdup.sorted.5fC.uniqueTags.sam -c rmdup.sorted.Input.uniqueTags.sam -f SAM -g mm -n 5fC-vsInput -w --call-subpeaks > macs14_5fCvsInput_balanced &