Subscribed unsubscribe Subscribe Subscribe

CPSAM.org

computer, programming, statistics and more!

マイナスから始めるExome解析

0から始めるエクソームデータ解析
を参考にさせていただきました

とりあえず現行では

sudo apt-get install ant
sudo apt-get install gradle
git clone https://github.com/broadinstitute/picard.git
cd picard
git clone https://github.com/samtools/htsjdk.git
./gradle -jar
cd ../
ant

でpicardのbuildまでは成功する


Broad instituteからbudle(reference genome, interval fileなどがセットになったもの)の入手
ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/b2.8/b37より適宜入手する

Reference genome: human_g1k_v37.fasta, human_g1k_v37.fasta.fai
Variant file: dbsnp_138.b37.vcf, dbsnp_138.b37.vcf.idx
Insert-deletion file: 1000G_phase1.indels.b37.vcf, 1000G_phase1.indels.b37.vcf.idx


サンプルデータの入手(http://www.ncbi.nlm.nih.gov/pubmed/21936905)

wget ftp://ftp.ddbj.nig.ac.jp/ddbj_database/dra/fastq/SRA044/SRA044780/SRX091464/SRR330441_1.fastq.bz2
wget ftp://ftp.ddbj.nig.ac.jp/ddbj_database/dra/fastq/SRA044/SRA044780/SRX091464/SRR330441_2.fastq.bz2
bunzip2 *.fastq.bz2

Reference genomeにindexをつける

bwa index -p human_g1k_v37.fasta -a bwtsw human_g1k_v37.fasta

picardでreferenceのdictionaryを作る

java -jar ../tools/picard.jar CreateSequenceDictionary \
 R=human_g1k_v37.fasta \
 O=human_g1k_v37.dict

アラインメント

bwa aln -t 2 human_g1k_v37.fasta SRR330441_1.fastq > SRR330441_1.sai
bwa aln -t 2 human_g1k_v37.fasta SRR330441_2.fastq > SRR330441_2.sai
bwa sampe -P human_g1k_v37.fasta \
 -r '@RG\tID:01\tSM:s6\tPL:Illumina' \ 
 -f SRR330441.sam \
 SRR330441_1.sai SRR330441_2.sai \
 SRR330441_1.fastq SRR330441_2.fastq

bamへの変換とソート、インデックス付与

samtools view -bS SRR330441.sam > SRR330441.bam
samtools sort SRR330441.bam -o SRR330441.sorted.bam
samtools index SRR330441.sorted.bam

Picardによる重複リードの除去とこれに対するインデックス付与

java -Xmx2G -jar ../tools/picard.jar MarkDuplicates \
 ASSUME_SORTED=true \
 REMOVE_DUPLICATES=true \
 INPUT=SRR330441.sorted.bam \
 OUTPUT=SRR330441.removed.sorted.bam \
 METRICS_FILE=duplicate \
 VALIDATION_STRINGENCY=LENIENT

samtools index SRR330441.removed.sorted.bam


bamファイル名を記載した"files.list"を作る

echo 'SRR330441.removed.sorted.bam' > files.list

GATKのRealignerTargetCreatorを使って既知のinsert/deletion付近の再マッピング用のファイルを作る

java -jar ../tools/GenomeAnalysisTK.jar \
 -T RealignerTargetCreator \
 -R human_g1k_v37.fasta \
 -I files.list \
 --known ../bundle/1000G_phase1.indels.b37.vcf \
 -log intervals.log \
 -L ../sureselect/S07604514_AllTracks.bed
 -o SRR330441.intervals
Remove all ads