Tweet
Dear all,
I am mapping one end mRNA sequencing on TOPHAT with human reference genome hg19.
It looks like script worked partially and It has generated .BAM file of 7 million MB. But geneexpression.txt are empty.
Here is the error file:
[Fri Oct 28 18:13:55 2011] Beginning TopHat run (v1.3.0)
-----------------------------------------------
[Fri Oct 28 18:13:55 2011] Preparing output location /home/ah644/tophat_output/intestinal/mRNA/C1055/
[Fri Oct 28 18:13:55 2011] Checking for Bowtie index files
[Fri Oct 28 18:13:55 2011] Checking for reference FASTA file
[Fri Oct 28 18:13:55 2011] Checking for Bowtie
Bowtie version: 0.12.7.0
[Fri Oct 28 18:13:55 2011] Checking for Samtools
Samtools Version: 0.1.17
[Fri Oct 28 18:13:55 2011] Generating SAM header for /home/ah644/hg19andhs37/hg19
[Fri Oct 28 18:13:56 2011] Preparing reads
format: fastq
quality scale: phred33 (default)
[Fri Oct 28 18:13:56 2011] Reading known junctions from GTF file
Left reads: min. length=76, count=151037396
[Fri Oct 28 19:33:26 2011] Mapping left_kept_reads against hg19 with Bowtie
[Fri Oct 28 21:46:41 2011] Processing bowtie hits
[Fri Oct 28 23:43:18 2011] Mapping left_kept_reads_seg1 against hg19 with Bowtie (1/3)
[Sat Oct 29 00:38:57 2011] Mapping left_kept_reads_seg2 against hg19 with Bowtie (2/3)
[Sat Oct 29 01:34:48 2011] Mapping left_kept_reads_seg3 against hg19 with Bowtie (3/3)
[Sat Oct 29 02:21:31 2011] Searching for junctions via segment mapping
[Sat Oct 29 05:17:18 2011] Retrieving sequences for splices
[Sat Oct 29 05:22:36 2011] Indexing splices
[Sat Oct 29 05:23:14 2011] Mapping left_kept_reads_seg1 against segment_juncs with Bowtie (1/3)
[Sat Oct 29 05:48:31 2011] Mapping left_kept_reads_seg2 against segment_juncs with Bowtie (2/3)
[Sat Oct 29 06:11:28 2011] Mapping left_kept_reads_seg3 against segment_juncs with Bowtie (3/3)
[Sat Oct 29 06:28:37 2011] Joining segment hits
[Sat Oct 29 08:53:02 2011] Reporting output tracks
-----------------------------------------------
Run complete [16:42:32 elapsed]
/opt/torque/mom_priv/jobs/80964.bulldogn-rocks.wss.yale.edu.SC: line 15: coverageBed: command not found
Pipeline:
#!/bin/bash
# usage: qsub -v SAMPLE="name" intestinaltophat1lane.pbs
#PBS -l walltime=100:00:00
#PBS -m ae
#PBS -M [email protected]
### Set the queue name, given to you when you get a reservation.
#PBS -l nodes=1: ppn=8
ROOT=/home/ah644/projects
cd $ROOT/$SAMPLE
ROOTB=/home/ah644/tophat_output/
TOPHATOUTPUT=$ROOTB$SAMPLE
mkdir $TOPHATOUTPUT
SEQFILES=`ls *_sequence.txt`
tophat -p 8 -r 50 --microexon-search -G /home/ah644/hg19andhs37/Homo_sapiens.GRCh37.62.gtf -o $TOPHATOUTPUT /home/ah644/hg19andhs37/hg19 $SEQFILES
coverageBed -abam $TOPHATOUTPUT/accepted_hits.bam -b /home/ah644/hg19andhs37/hg19ensembl.BED > $TOPHATOUTPUT/geneexpression.txt
cut -f 1,2,3,4,13,14,15,16 $TOPHATOUTPUT/geneexpression.txt > $TOPHATOUTPUT/temp
cut -f 5 /home/ah644/hg19andhs37/hg19ensemblgeneinfo.txt > $TOPHATOUTPUT/temp2
paste $TOPHATOUTPUT/temp $TOPHATOUTPUT/temp2 > $TOPHATOUTPUT/reducedgeneexpression.txt
rm $TOPHATOUTPUT/temp*
please help me in correcting the script.
many thanks,
Adam
I am mapping one end mRNA sequencing on TOPHAT with human reference genome hg19.
It looks like script worked partially and It has generated .BAM file of 7 million MB. But geneexpression.txt are empty.
Here is the error file:
[Fri Oct 28 18:13:55 2011] Beginning TopHat run (v1.3.0)
-----------------------------------------------
[Fri Oct 28 18:13:55 2011] Preparing output location /home/ah644/tophat_output/intestinal/mRNA/C1055/
[Fri Oct 28 18:13:55 2011] Checking for Bowtie index files
[Fri Oct 28 18:13:55 2011] Checking for reference FASTA file
[Fri Oct 28 18:13:55 2011] Checking for Bowtie
Bowtie version: 0.12.7.0
[Fri Oct 28 18:13:55 2011] Checking for Samtools
Samtools Version: 0.1.17
[Fri Oct 28 18:13:55 2011] Generating SAM header for /home/ah644/hg19andhs37/hg19
[Fri Oct 28 18:13:56 2011] Preparing reads
format: fastq
quality scale: phred33 (default)
[Fri Oct 28 18:13:56 2011] Reading known junctions from GTF file
Left reads: min. length=76, count=151037396
[Fri Oct 28 19:33:26 2011] Mapping left_kept_reads against hg19 with Bowtie
[Fri Oct 28 21:46:41 2011] Processing bowtie hits
[Fri Oct 28 23:43:18 2011] Mapping left_kept_reads_seg1 against hg19 with Bowtie (1/3)
[Sat Oct 29 00:38:57 2011] Mapping left_kept_reads_seg2 against hg19 with Bowtie (2/3)
[Sat Oct 29 01:34:48 2011] Mapping left_kept_reads_seg3 against hg19 with Bowtie (3/3)
[Sat Oct 29 02:21:31 2011] Searching for junctions via segment mapping
[Sat Oct 29 05:17:18 2011] Retrieving sequences for splices
[Sat Oct 29 05:22:36 2011] Indexing splices
[Sat Oct 29 05:23:14 2011] Mapping left_kept_reads_seg1 against segment_juncs with Bowtie (1/3)
[Sat Oct 29 05:48:31 2011] Mapping left_kept_reads_seg2 against segment_juncs with Bowtie (2/3)
[Sat Oct 29 06:11:28 2011] Mapping left_kept_reads_seg3 against segment_juncs with Bowtie (3/3)
[Sat Oct 29 06:28:37 2011] Joining segment hits
[Sat Oct 29 08:53:02 2011] Reporting output tracks
-----------------------------------------------
Run complete [16:42:32 elapsed]
/opt/torque/mom_priv/jobs/80964.bulldogn-rocks.wss.yale.edu.SC: line 15: coverageBed: command not found
Pipeline:
#!/bin/bash
# usage: qsub -v SAMPLE="name" intestinaltophat1lane.pbs
#PBS -l walltime=100:00:00
#PBS -m ae
#PBS -M [email protected]
### Set the queue name, given to you when you get a reservation.
#PBS -l nodes=1: ppn=8
ROOT=/home/ah644/projects
cd $ROOT/$SAMPLE
ROOTB=/home/ah644/tophat_output/
TOPHATOUTPUT=$ROOTB$SAMPLE
mkdir $TOPHATOUTPUT
SEQFILES=`ls *_sequence.txt`
tophat -p 8 -r 50 --microexon-search -G /home/ah644/hg19andhs37/Homo_sapiens.GRCh37.62.gtf -o $TOPHATOUTPUT /home/ah644/hg19andhs37/hg19 $SEQFILES
coverageBed -abam $TOPHATOUTPUT/accepted_hits.bam -b /home/ah644/hg19andhs37/hg19ensembl.BED > $TOPHATOUTPUT/geneexpression.txt
cut -f 1,2,3,4,13,14,15,16 $TOPHATOUTPUT/geneexpression.txt > $TOPHATOUTPUT/temp
cut -f 5 /home/ah644/hg19andhs37/hg19ensemblgeneinfo.txt > $TOPHATOUTPUT/temp2
paste $TOPHATOUTPUT/temp $TOPHATOUTPUT/temp2 > $TOPHATOUTPUT/reducedgeneexpression.txt
rm $TOPHATOUTPUT/temp*
please help me in correcting the script.
many thanks,
Adam