Ylaine/2009-7-15: Difference between revisions

From ZhangLabWiki
Jump to navigation Jump to search
>Ylaine
No edit summary
>Ylaine
 
(6 intermediate revisions by the same user not shown)
Line 1: Line 1:
==Comparing Maq/Maq and Bwa/Sam==
==Maq/Sam==
Outputs from the two aligning/consensus calling methods often have different read numbers associated with SNPs at the same location, suggesting that the mapping algorithms produce different results. This is confirmed by direct examination of the output 'pileup' files.
Outputs from the two aligning/consensus calling methods (Maq/Maq and Bwa/Sam) often have different read numbers associated with SNPs at the same location, suggesting that the mapping algorithms produce different results. This is confirmed by direct examination of the output 'pileup' files.
* Run Sam on Maq output. Created new directory "MaqSamTarget"
* Run Sam on Maq output. Created new directory "MaqSamTarget"
  /Users/kunzhang/Downloads/samtools/misc/maq2sam-long ../MaqTarget/all.map > all.sam
  /Users/kunzhang/Downloads/samtools/misc/maq2sam-long ../MaqTarget/all.map > all.sam
Line 7: Line 7:
  /Users/kunzhang/Downloads/samtools/samtools rmdupse all.sorted.bam all.unique.bam
  /Users/kunzhang/Downloads/samtools/samtools rmdupse all.sorted.bam all.unique.bam
  /Users/kunzhang/Downloads/samtools/samtools pileup -c -f /Users/kunzhang/WorkSpace/Exome/targets/ccdsExon_unique_miRNA_200bp.fa.combined all.unique.bam > all.pileup
  /Users/kunzhang/Downloads/samtools/samtools pileup -c -f /Users/kunzhang/WorkSpace/Exome/targets/ccdsExon_unique_miRNA_200bp.fa.combined all.unique.bam > all.pileup
  /Users/kunzhang/Downloads/samtools/misc/samtools.pl varFilter all.pileup > snp.txt
  /Users/kunzhang/Downloads/samtools/misc/samtools.pl varFilter –D 50000 -d2 all.pileup > snp.txt
  /Users/kunzhang/WorkSpace/Exome/Solexa/scripts/parseExomeCnsSNP.pl /Users/kunzhang/WorkSpace/Exome/targets/ccdsExon_unique_miRNA_200bp.fa.seqStartTable.txt snp.txt > snp.chr.txt
  /Users/kunzhang/WorkSpace/Exome/Solexa/scripts/parseExomeCnsSNP.pl /Users/kunzhang/WorkSpace/Exome/targets/ccdsExon_unique_miRNA_200bp.fa.seqStartTable.txt snp.txt > snp.chr.txt
* Output is still different from MAQ (lower read numbers)
* Try without removing duplicates:
/Users/kunzhang/Downloads/samtools/samtools pileup -c -f /Users/kunzhang/WorkSpace/Exome/targets/ccdsExon_unique_miRNA_200bp.fa.combined all.sorted.bam > nonunique.pileup
/Users/kunzhang/Downloads/samtools/misc/samtools.pl varFilter –D 50000 -d2 nonunique.pileup > snp.nonunique.txt
/Users/kunzhang/WorkSpace/Exome/Solexa/scripts/parseExomeCnsSNP.pl /Users/kunzhang/WorkSpace/Exome/targets/ccdsExon_unique_miRNA_200bp.fa.seqStartTable.txt snp.nonunique.txt > snp.nonunique.chr.txt
* Read numbers tend to be +1 higher than for Maq/Maq
* Found 24532 SNPs (vs 8376 for Maq easyrun; 79100 for varFilter file from Dr. Zhang)
* SNPs for which Maq reported a depth of 255 reads are missing altogether from Maq/Sam
===SNP Analysis===
* Put dbSNP matching method inside separate script
./matchCnsSnp2Ref_dbSNP.pl MaqSamTarget/snp.nonunique.chr.txt > snp.nonunique.dbSNP
* 10428/24532 in dbSNP (42.51%)
* Compare to 1000 Genome and HapMap calls (modified matchCnsSnp2Ref_AG.pl)
grep 'rs' snp.nonunique.dbSNP | ../matchCnsSnp2Ref_AG.pl > snp.nonunique.compare.in_dbSNP
grep -v 'rs' snp.nonunique.dbSNP | ../matchCnsSnp2Ref_AG.pl > snp.nonunique.compare.no_dbSNP
* Result in table below. FP rate is 49.51%
{| {{table}}
| align="center" style="background:#f0f0f0;"|'''1K Genome'''
| align="center" style="background:#f0f0f0;"|'''In dbSNP'''
| align="center" style="background:#f0f0f0;"|'''Percent'''
| align="center" style="background:#f0f0f0;"|'''Not in dbSNP'''
| align="center" style="background:#f0f0f0;"|'''Percent'''
|-
| Match||8038||77.08%||228||1.62%
|-
| Miss||1154||11.07%||13706||97.18%
|-
| Mismatch||1236||11.85%||80||0.57%
|-
| Total||10428||||14104||
|-
|
|}
* Discovered problem with extracting columns using 'awk'
** If column entry is (tab) in a tab-delimited file, then the entry is merged with the delimiters on either side
** Should re-do statistics from other methods :(
[[Image:Maqsam_nodbSNP.jpg‎]]

Latest revision as of 23:19, 15 July 2009

Maq/Sam[edit]

Outputs from the two aligning/consensus calling methods (Maq/Maq and Bwa/Sam) often have different read numbers associated with SNPs at the same location, suggesting that the mapping algorithms produce different results. This is confirmed by direct examination of the output 'pileup' files.

  • Run Sam on Maq output. Created new directory "MaqSamTarget"
/Users/kunzhang/Downloads/samtools/misc/maq2sam-long ../MaqTarget/all.map > all.sam
/Users/kunzhang/Downloads/samtools/samtools import /Users/kunzhang/WorkSpace/Exome/targets/ccdsExon_unique_miRNA_200bp.fa.combined.fai all.sam all.bam
/Users/kunzhang/Downloads/samtools/samtools sort all.bam all.sorted
/Users/kunzhang/Downloads/samtools/samtools rmdupse all.sorted.bam all.unique.bam
/Users/kunzhang/Downloads/samtools/samtools pileup -c -f /Users/kunzhang/WorkSpace/Exome/targets/ccdsExon_unique_miRNA_200bp.fa.combined all.unique.bam > all.pileup
/Users/kunzhang/Downloads/samtools/misc/samtools.pl varFilter –D 50000 -d2 all.pileup > snp.txt
/Users/kunzhang/WorkSpace/Exome/Solexa/scripts/parseExomeCnsSNP.pl /Users/kunzhang/WorkSpace/Exome/targets/ccdsExon_unique_miRNA_200bp.fa.seqStartTable.txt snp.txt > snp.chr.txt
  • Output is still different from MAQ (lower read numbers)
  • Try without removing duplicates:
/Users/kunzhang/Downloads/samtools/samtools pileup -c -f /Users/kunzhang/WorkSpace/Exome/targets/ccdsExon_unique_miRNA_200bp.fa.combined all.sorted.bam > nonunique.pileup
/Users/kunzhang/Downloads/samtools/misc/samtools.pl varFilter –D 50000 -d2 nonunique.pileup > snp.nonunique.txt
/Users/kunzhang/WorkSpace/Exome/Solexa/scripts/parseExomeCnsSNP.pl /Users/kunzhang/WorkSpace/Exome/targets/ccdsExon_unique_miRNA_200bp.fa.seqStartTable.txt snp.nonunique.txt > snp.nonunique.chr.txt
  • Read numbers tend to be +1 higher than for Maq/Maq
  • Found 24532 SNPs (vs 8376 for Maq easyrun; 79100 for varFilter file from Dr. Zhang)
  • SNPs for which Maq reported a depth of 255 reads are missing altogether from Maq/Sam

SNP Analysis[edit]

  • Put dbSNP matching method inside separate script
./matchCnsSnp2Ref_dbSNP.pl MaqSamTarget/snp.nonunique.chr.txt > snp.nonunique.dbSNP
  • 10428/24532 in dbSNP (42.51%)
  • Compare to 1000 Genome and HapMap calls (modified matchCnsSnp2Ref_AG.pl)
grep 'rs' snp.nonunique.dbSNP | ../matchCnsSnp2Ref_AG.pl > snp.nonunique.compare.in_dbSNP
grep -v 'rs' snp.nonunique.dbSNP | ../matchCnsSnp2Ref_AG.pl > snp.nonunique.compare.no_dbSNP
  • Result in table below. FP rate is 49.51%
1K Genome In dbSNP Percent Not in dbSNP Percent
Match 8038 77.08% 228 1.62%
Miss 1154 11.07% 13706 97.18%
Mismatch 1236 11.85% 80 0.57%
Total 10428 14104
  • Discovered problem with extracting columns using 'awk'
    • If column entry is (tab) in a tab-delimited file, then the entry is merged with the delimiters on either side
    • Should re-do statistics from other methods :(

File:Maqsam nodbSNP.jpg