Matt:LabNotes/2014-7-30: Difference between revisions
Jump to navigation
Jump to search
>Mzcai mNo edit summary |
>Mzcai mNo edit summary |
||
Line 1: | Line 1: | ||
==RNA-SeQC== | ==RNA-SeQC== | ||
*Tophat aligned to hg19 BAM input | |||
===Tial 1=== | |||
====Add Read Groups and '''coordinate sort'''==== | |||
*Picard.AddOrReplaceReadGroups | |||
java -Xmx2g -jar /home/kunzhang/softwares/picard-tools-1.79/AddOrReplaceReadGroups.jar INPUT='accepted_hits.bam' OUTPUT='accepted_hits_RanHex.coorsortRG.bam' SORT_ORDER=coordinate RGID='1' RGLB='RTprimer' RGPL='ILLUMINA' RGPU='flowcell-barcode.lane' RGSM='RanHex' | java -Xmx2g -jar /home/kunzhang/softwares/picard-tools-1.79/AddOrReplaceReadGroups.jar INPUT='accepted_hits.bam' OUTPUT='accepted_hits_RanHex.coorsortRG.bam' SORT_ORDER=coordinate RGID='1' RGLB='RTprimer' RGPL='ILLUMINA' RGPU='flowcell-barcode.lane' RGSM='RanHex' | ||
====Index reference fasta and create sequence dictionary==== | |||
*reference fasta index already exists | |||
/home/kunzhang/HsGenome/hg19/HsGenome19.fa.fai | /home/kunzhang/HsGenome/hg19/HsGenome19.fa.fai | ||
*Picard.CreateSequenceDictionary | |||
java -Xmx2g -jar /home/kunzhang/softwares/picard-tools-1.79/CreateSequenceDictionary.jar REFERENCE=/home/kunzhang/HsGenome/hg19/HsGenome19.fa OUTPUT=/home/mzcai/Genomes/HsGenome19.dict | java -Xmx2g -jar /home/kunzhang/softwares/picard-tools-1.79/CreateSequenceDictionary.jar REFERENCE=/home/kunzhang/HsGenome/hg19/HsGenome19.fa OUTPUT=/home/mzcai/Genomes/HsGenome19.dict | ||
====Match order of contigs in Bam file to reference file==== | |||
*Picard.ReorderSam | |||
java -Xmx2g -jar /home/kunzhang/softwares/picard-tools-1.79/ReorderSam.jar INPUT=accepted_hits_RanHex.coorsortRG.bam OUTPUT=accepted_hits_RanHex.coorsortRGreorder.bam REFERENCE=/home/mzcai/Genomes/HsGenome19.fa | java -Xmx2g -jar /home/kunzhang/softwares/picard-tools-1.79/ReorderSam.jar INPUT=accepted_hits_RanHex.coorsortRG.bam OUTPUT=accepted_hits_RanHex.coorsortRGreorder.bam REFERENCE=/home/mzcai/Genomes/HsGenome19.fa | ||
====Index bam file==== | |||
samtools index accepted_hits_RanHex.coorsortRGreorder.bam | samtools index accepted_hits_RanHex.coorsortRGreorder.bam | ||
====Create sample file==== | |||
*~/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_samples.txt | |||
Sample ID Bam File Notes | Sample ID Bam File Notes | ||
RanHex /home/mzcai/scratch/RanHex_EnrichmRNA/inVitroRTSeq/tophat_hg19unmask_Indx26_RanHex/accepted_hits_RanHex.coorsortRGreorder.bam Indx26RanHex | RanHex /home/mzcai/scratch/RanHex_EnrichmRNA/inVitroRTSeq/tophat_hg19unmask_Indx26_RanHex/accepted_hits_RanHex.coorsortRGreorder.bam Indx26RanHex | ||
====Run RNA-SeQC v1.1.7.jar==== | |||
java -jar /home/kunzhang/softwares/RNA-SeQC_v1.1.7.jar -bwa /home/kunzhang/softwares/bwa-0.7.5a/bwa -BWArRNA /home/mzcai/scratch/RanHex_EnrichmRNA/mart_export_exons_totalrRNA.fa -o /home/mzcai/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_output -r /home/mzcai/Genomes/HsGenome19.fa -s /home/mzcai/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_samples.txt -singleEnd -t /home/mzcai/LTS/Genomes/Homo_sapiens.GRCh37.75.chr.gtf -ttype 2 | java -jar /home/kunzhang/softwares/RNA-SeQC_v1.1.7.jar -bwa /home/kunzhang/softwares/bwa-0.7.5a/bwa -BWArRNA /home/mzcai/scratch/RanHex_EnrichmRNA/mart_export_exons_totalrRNA.fa -o /home/mzcai/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_output -r /home/mzcai/Genomes/HsGenome19.fa -s /home/mzcai/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_samples.txt -singleEnd -t /home/mzcai/LTS/Genomes/Homo_sapiens.GRCh37.75.chr.gtf -ttype 2 | ||
The required transcript_id attribute was not found on line chr1 pseudogene gene 11869 14412 . + . gene_id "ENSG00000223972"; gene_name "DDX11L1"; gene_source "ensembl_havana"; gene_biotype "pseudogene"; | The required transcript_id attribute was not found on line chr1 pseudogene gene 11869 14412 . + . gene_id "ENSG00000223972"; gene_name "DDX11L1"; gene_source "ensembl_havana"; gene_biotype "pseudogene"; | ||
===Download resources from [http://www.broadinstitute.org/cancer/cga/rnaseqc_download broadinstitute]=== | ===Trial 2=== | ||
*Use provided gencode gtf annotation file | |||
**Also use provided rRNA ref file and GC definitions file to be consistent with example | |||
====Download resources from [http://www.broadinstitute.org/cancer/cga/rnaseqc_download broadinstitute]==== | |||
*~/Genomes/RNA-SeQCResources | *~/Genomes/RNA-SeQCResources | ||
*Download human rRNA reference file | *Download human rRNA reference file | ||
Line 32: | Line 43: | ||
gencode.v7.gc.txt | gencode.v7.gc.txt | ||
====Run RNA-SeQC v1.1.7.jar==== | |||
java -jar /home/kunzhang/softwares/RNA-SeQC_v1.1.7.jar -bwa /home/kunzhang/softwares/bwa-0.7.5a/ -BWArRNA ~/Genomes/RNA-SeqCResources/human_all_rRNA.fasta -s ~/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_samples.txt -t ~/Genomes/RNA-SeqCResources/gencode.v7.annotation.gtf -r ~/Genomes/HsGenome19.fa -o ~/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_output/ -strat gc -gc ~/Genomes/RNA-SeqCResources/gencode.v7.gc.txt -singleEnd | java -jar /home/kunzhang/softwares/RNA-SeQC_v1.1.7.jar -bwa /home/kunzhang/softwares/bwa-0.7.5a/ -BWArRNA ~/Genomes/RNA-SeqCResources/human_all_rRNA.fasta -s ~/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_samples.txt -t ~/Genomes/RNA-SeqCResources/gencode.v7.annotation.gtf -r ~/Genomes/HsGenome19.fa -o ~/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_output/ -strat gc -gc ~/Genomes/RNA-SeqCResources/gencode.v7.gc.txt -singleEnd | ||
org.broadinstitute.sting.utils.exceptions.UserException$LexicographicallySortedSequenceDictionary: Lexicographically sorted human genome sequence detected in reads. | org.broadinstitute.sting.utils.exceptions.UserException$LexicographicallySortedSequenceDictionary: Lexicographically sorted human genome sequence detected in reads. | ||
For safety's sake the GATK requires human contigs in karyotypic order: 1, 2, ..., 10, 11, ..., 20, 21, 22, X, Y with M either leading or trailing these contigs. | For safety's sake the GATK requires human contigs in karyotypic order: 1, 2, ..., 10, 11, ..., 20, 21, 22, X, Y with M either leading or trailing these contigs. | ||
This is because all distributed GATK resources are sorted in karyotypic order, and your processing will fail when you need to use these files. | This is because all distributed GATK resources are sorted in karyotypic order, and your processing will fail when you need to use these files. | ||
You can use the ReorderSam utility to fix this problem: http://www.broadinstitute.org/gsa/wiki/index.php/ReorderSam | You can use the ReorderSam utility to fix this problem: http://www.broadinstitute.org/gsa/wiki/index.php/ReorderSam | ||
===Trial 3=== | |||
===Reorder HsGenome19.fa to karyotypic order=== | ====Reorder HsGenome19.fa to karyotypic order==== | ||
csplit -f c HsGenome19.fa /'>chr9'/ /'>chr5'/ /'>chr2'/ /'>chr17'/ /'>chr14'/ /'>chr4'/ /'>chr16'/ /'>chr21'/ /'>chrM'/ /'>chr7'/ /'>chr3'/ /'>chr18'/ /'>chr12'/ /'>chrX'/ /'>chr13'/ /'>chr15'/ /'>chr8'/ /'>chr22'/ /'>chr11'/ /'>chr10'/ /'>chr20'/ /'>chr19'/ /'>chr6'/ /'>chr1'/ /'>chrY'/ | csplit -f c HsGenome19.fa /'>chr9'/ /'>chr5'/ /'>chr2'/ /'>chr17'/ /'>chr14'/ /'>chr4'/ /'>chr16'/ /'>chr21'/ /'>chrM'/ /'>chr7'/ /'>chr3'/ /'>chr18'/ /'>chr12'/ /'>chrX'/ /'>chr13'/ /'>chr15'/ /'>chr8'/ /'>chr22'/ /'>chr11'/ /'>chr10'/ /'>chr20'/ /'>chr19'/ /'>chr6'/ /'>chr1'/ /'>chrY'/ | ||
cat c24 c03 c11 c06 c02 c23 c10 c17 c01 c20 c19 c13 c15 c05 c16 c07 c04 c12 c22 c21 c08 c18 c14 c25 c09 > HsGenome19_karyorder.fa | cat c24 c03 c11 c06 c02 c23 c10 c17 c01 c20 c19 c13 c15 c05 c16 c07 c04 c12 c22 c21 c08 c18 c14 c25 c09 > HsGenome19_karyorder.fa | ||
rm HsGenome19.* | |||
mv HsGenome19_karyorder.fa ./HsGenome19.fa | |||
{| {{table}} | {| {{table}} | ||
Line 101: | Line 115: | ||
| >chrY||/'>chrY'/||c25||c09||chrM | | >chrY||/'>chrY'/||c25||c09||chrM | ||
|} | |} | ||
==== |
Revision as of 17:58, 1 August 2014
RNA-SeQC
- Tophat aligned to hg19 BAM input
Tial 1
Add Read Groups and coordinate sort
- Picard.AddOrReplaceReadGroups
java -Xmx2g -jar /home/kunzhang/softwares/picard-tools-1.79/AddOrReplaceReadGroups.jar INPUT='accepted_hits.bam' OUTPUT='accepted_hits_RanHex.coorsortRG.bam' SORT_ORDER=coordinate RGID='1' RGLB='RTprimer' RGPL='ILLUMINA' RGPU='flowcell-barcode.lane' RGSM='RanHex'
Index reference fasta and create sequence dictionary
- reference fasta index already exists
/home/kunzhang/HsGenome/hg19/HsGenome19.fa.fai
- Picard.CreateSequenceDictionary
java -Xmx2g -jar /home/kunzhang/softwares/picard-tools-1.79/CreateSequenceDictionary.jar REFERENCE=/home/kunzhang/HsGenome/hg19/HsGenome19.fa OUTPUT=/home/mzcai/Genomes/HsGenome19.dict
Match order of contigs in Bam file to reference file
- Picard.ReorderSam
java -Xmx2g -jar /home/kunzhang/softwares/picard-tools-1.79/ReorderSam.jar INPUT=accepted_hits_RanHex.coorsortRG.bam OUTPUT=accepted_hits_RanHex.coorsortRGreorder.bam REFERENCE=/home/mzcai/Genomes/HsGenome19.fa
Index bam file
samtools index accepted_hits_RanHex.coorsortRGreorder.bam
Create sample file
- ~/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_samples.txt
Sample ID Bam File Notes RanHex /home/mzcai/scratch/RanHex_EnrichmRNA/inVitroRTSeq/tophat_hg19unmask_Indx26_RanHex/accepted_hits_RanHex.coorsortRGreorder.bam Indx26RanHex
Run RNA-SeQC v1.1.7.jar
java -jar /home/kunzhang/softwares/RNA-SeQC_v1.1.7.jar -bwa /home/kunzhang/softwares/bwa-0.7.5a/bwa -BWArRNA /home/mzcai/scratch/RanHex_EnrichmRNA/mart_export_exons_totalrRNA.fa -o /home/mzcai/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_output -r /home/mzcai/Genomes/HsGenome19.fa -s /home/mzcai/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_samples.txt -singleEnd -t /home/mzcai/LTS/Genomes/Homo_sapiens.GRCh37.75.chr.gtf -ttype 2
The required transcript_id attribute was not found on line chr1 pseudogene gene 11869 14412 . + . gene_id "ENSG00000223972"; gene_name "DDX11L1"; gene_source "ensembl_havana"; gene_biotype "pseudogene";
Trial 2
- Use provided gencode gtf annotation file
- Also use provided rRNA ref file and GC definitions file to be consistent with example
Download resources from broadinstitute
- ~/Genomes/RNA-SeQCResources
- Download human rRNA reference file
human_all_rRNA.fasta
- Download GENCODE GTF file
gencode.v7.annotation.gtf
- Download GENCODE GC definitions file
gencode.v7.gc.txt
Run RNA-SeQC v1.1.7.jar
java -jar /home/kunzhang/softwares/RNA-SeQC_v1.1.7.jar -bwa /home/kunzhang/softwares/bwa-0.7.5a/ -BWArRNA ~/Genomes/RNA-SeqCResources/human_all_rRNA.fasta -s ~/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_samples.txt -t ~/Genomes/RNA-SeqCResources/gencode.v7.annotation.gtf -r ~/Genomes/HsGenome19.fa -o ~/scratch/RanHex_EnrichmRNA/inVitroRTSeq/RNA-SeQC_output/ -strat gc -gc ~/Genomes/RNA-SeqCResources/gencode.v7.gc.txt -singleEnd org.broadinstitute.sting.utils.exceptions.UserException$LexicographicallySortedSequenceDictionary: Lexicographically sorted human genome sequence detected in reads. For safety's sake the GATK requires human contigs in karyotypic order: 1, 2, ..., 10, 11, ..., 20, 21, 22, X, Y with M either leading or trailing these contigs. This is because all distributed GATK resources are sorted in karyotypic order, and your processing will fail when you need to use these files. You can use the ReorderSam utility to fix this problem: http://www.broadinstitute.org/gsa/wiki/index.php/ReorderSam
Trial 3
Reorder HsGenome19.fa to karyotypic order
csplit -f c HsGenome19.fa /'>chr9'/ /'>chr5'/ /'>chr2'/ /'>chr17'/ /'>chr14'/ /'>chr4'/ /'>chr16'/ /'>chr21'/ /'>chrM'/ /'>chr7'/ /'>chr3'/ /'>chr18'/ /'>chr12'/ /'>chrX'/ /'>chr13'/ /'>chr15'/ /'>chr8'/ /'>chr22'/ /'>chr11'/ /'>chr10'/ /'>chr20'/ /'>chr19'/ /'>chr6'/ /'>chr1'/ /'>chrY'/ cat c24 c03 c11 c06 c02 c23 c10 c17 c01 c20 c19 c13 c15 c05 c16 c07 c04 c12 c22 c21 c08 c18 c14 c25 c09 > HsGenome19_karyorder.fa rm HsGenome19.* mv HsGenome19_karyorder.fa ./HsGenome19.fa
Old chr order | Add // | csplit output files | cat order | New chr order |
>chr9 | /'>chr9'/ | c01 | c24 | chr1 |
>chr5 | /'>chr5'/ | c02 | c03 | chr2 |
>chr2 | /'>chr2'/ | c03 | c11 | chr3 |
>chr17 | /'>chr17'/ | c04 | c06 | chr4 |
>chr14 | /'>chr14'/ | c05 | c02 | chr5 |
>chr4 | /'>chr4'/ | c06 | c23 | chr6 |
>chr16 | /'>chr16'/ | c07 | c10 | chr7 |
>chr21 | /'>chr21'/ | c08 | c17 | chr8 |
>chrM | /'>chrM'/ | c09 | c01 | chr9 |
>chr7 | /'>chr7'/ | c10 | c20 | chr10 |
>chr3 | /'>chr3'/ | c11 | c19 | chr11 |
>chr18 | /'>chr18'/ | c12 | c13 | chr12 |
>chr12 | /'>chr12'/ | c13 | c15 | chr13 |
>chrX | /'>chrX'/ | c14 | c05 | chr14 |
>chr13 | /'>chr13'/ | c15 | c16 | chr15 |
>chr15 | /'>chr15'/ | c16 | c07 | chr16 |
>chr8 | /'>chr8'/ | c17 | c04 | chr17 |
>chr22 | /'>chr22'/ | c18 | c12 | chr18 |
>chr11 | /'>chr11'/ | c19 | c22 | chr19 |
>chr10 | /'>chr10'/ | c20 | c21 | chr20 |
>chr20 | /'>chr20'/ | c21 | c08 | chr21 |
>chr19 | /'>chr19'/ | c22 | c18 | chr22 |
>chr6 | /'>chr6'/ | c23 | c14 | chrX |
>chr1 | /'>chr1'/ | c24 | c25 | chrY |
>chrY | /'>chrY'/ | c25 | c09 | chrM |