SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
Contig vs contig or map against contig lib? JackieBadger Bioinformatics 1 05-30-2016 05:34 AM
High number of N in reads papori Illumina/Solexa 3 01-22-2013 05:21 AM
Number of Reads Rfriedman Bioinformatics 0 01-11-2012 01:43 PM
SRMA Problem SAMRecord contig does not match the current reference sequence contig gavin.oliver Bioinformatics 5 07-05-2011 05:28 AM
number of reads m_elena_bioinfo Bioinformatics 2 07-20-2010 08:47 AM

Reply
 
Thread Tools
Old 05-14-2013, 02:41 PM   #1
kmkocot
Member
 
Location: Alabama

Join Date: Jun 2009
Posts: 48
Default Easiest way to get number of reads per contig?

Hi all,

I'm trying to map 2 X 100 PE Illumina reads (transcriptome) back to a small number of contigs. I want to get the number of reads that align to each contig as quickly as possible. Here's the commands I've tried so far but it isn't working.

Code:
SCAFFOLD_INPUT_FILE=test_selected_sequences.fas
S1_FASTQ_INPUT_FILE=D13F4ACXX_s6_1_GSLv2-7_32_SL18114_test.fastq.gz
S2_FASTQ_INPUT_FILE=D13F4ACXX_s6_2_GSLv2-7_32_SL18114_test.fastq.gz
REF_MAP_BASENAME=test

BWA=`which bwa`
SAMTOOLS=`which samtools`
BCFTOOLS=`which bcftools`
VCFUTILS=`which vcfutils.pl`
BT2=`which bowtie2`
BT_BUILD=`which bowtie2-build`
TDSTAMP=`date | awk 'OFS="_" {print ($2,$3,$6,$4)}'`

${BT_BUILD} ${SCAFFOLD_INPUT_FILE} ${REF_MAP_BASENAME}
${BT2} -x ${REF_MAP_BASENAME} -1 ${S1_FASTQ_INPUT_FILE} -2 ${S2_FASTQ_INPUT_FILE} --qc-filter --no-unal -S ${REF_MAP_BASENAME}.bowtie2.mapped.sam
$SAMTOOLS view -bT ${SCAFFOLD_INPUT_FILE} ${REF_MAP_BASENAME}.bowtie2.mapped.sam > ${REF_MAP_BASENAME}.bowtie2.mapped.bam
$SAMTOOLS sort ${REF_MAP_BASENAME}.bowtie2.mapped.bam ${REF_MAP_BASENAME}.bowtie2.mapped.sorted
$SAMTOOLS idxstats ${REF_MAP_BASENAME}.bowtie2.mapped.sorted.bam
Here's the errors I get:
Code:
[samopen] SAM header is present: 28 sequences.
[bam_index_load] fail to load BAM index.
[bam_idxstats] fail to load the index.
The problem seems to be once the samtools apps start running. I'm using samtools version 0.1.19-44428cd. Any ideas on what I'm doing wrong? Is there an easier way?

Thanks!
Kevin
kmkocot is offline   Reply With Quote
Old 05-14-2013, 04:39 PM   #2
sdriscoll
I like code
 
Location: San Diego, CA, USA

Join Date: Sep 2009
Posts: 438
Default

revise:
Code:
${BT_BUILD} ${SCAFFOLD_INPUT_FILE} ${REF_MAP_BASENAME}
${BT2} -x ${REF_MAP_BASENAME} -1 ${S1_FASTQ_INPUT_FILE} -2 ${S2_FASTQ_INPUT_FILE} --qc-filter --no-unal -S ${REF_MAP_BASENAME}.bowtie2.mapped.sam
$SAMTOOLS view -bT ${SCAFFOLD_INPUT_FILE} ${REF_MAP_BASENAME}.bowtie2.mapped.sam > ${REF_MAP_BASENAME}.bowtie2.mapped.bam
$SAMTOOLS sort ${REF_MAP_BASENAME}.bowtie2.mapped.bam ${REF_MAP_BASENAME}.bowtie2.mapped.sorted
$SAMTOOLS idxstats ${REF_MAP_BASENAME}.bowtie2.mapped.sorted.bam
to

Code:
${BT_BUILD} ${SCAFFOLD_INPUT_FILE} ${REF_MAP_BASENAME}
${BT2} -x ${REF_MAP_BASENAME} -1 ${S1_FASTQ_INPUT_FILE} -2 ${S2_FASTQ_INPUT_FILE} --qc-filter --no-unal -S ${REF_MAP_BASENAME}.bowtie2.mapped.sam
$SAMTOOLS view -bT ${SCAFFOLD_INPUT_FILE} ${REF_MAP_BASENAME}.bowtie2.mapped.sam > ${REF_MAP_BASENAME}.bowtie2.mapped.bam
$SAMTOOLS sort ${REF_MAP_BASENAME}.bowtie2.mapped.bam ${REF_MAP_BASENAME}.bowtie2.mapped.sorted
$SAMTOOLS index ${REF_MAP_BASENAME}.bowtie2.mapped.sorted.bam
$SAMTOOLS idxstats ${REF_MAP_BASENAME}.bowtie2.mapped.sorted.bam
I added in the call to "samtools index" which is why the samtools idxstats call is failing.
sdriscoll is offline   Reply With Quote
Old 05-14-2013, 04:53 PM   #3
kmkocot
Member
 
Location: Alabama

Join Date: Jun 2009
Posts: 48
Default

Thanks! I'm testing it now.
kmkocot is offline   Reply With Quote
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 03:49 PM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO