use Getopt::Long;
use Pod::Usage;
-use File::Basename;
+use FindBin;
+use lib $FindBin::Bin;
use strict;
+use rsem_perl_utils;
+
my $status;
my $gtfF = "";
$size = scalar(@list);
}
+pod2usage(-msg => "reference_fasta_file(s) is empty! Please check if you provide the correct folder name or file suffixes!", -exitval => 2, -verbose => 2) if ($size <= 0);
+
if ($no_polyA) { $polyAChoice = 1 }
elsif ($subsetFile ne "") { $polyAChoice = 2; }
if ($bowtie_path ne "") { $bowtie_path .= "/"; }
-my ($fn, $dir, $suf) = fileparse($0);
+my $dir = "$FindBin::Bin/";
my $command = "";
if ($type == 0) {
if ($mappingF ne "") { $command .= " 1 $mappingF"; }
else { $command .= " 0"; }
$command .= " @list";
- print "$command\n";
- $status = system($command);
- if ($status != 0) {
- print "rsem-extract-reference-transcripts failed! Please check if you provide correct parameters/options for the pipeline!\n";
- exit(-1);
- }
- print "\n";
+ &runCommand($command);
}
else {
$"=" ";
if ($mappingF ne "") { $command .= " 1 $mappingF"; }
else { $command .= " 0"; }
$command .= " @list";
- print "$command\n";
- $status = system($command);
- if ($status != 0) {
- print "rsem-synthesis-reference-transcripts failed! Please check if you provide correct parameters/options for the pipeline!\n";
- exit(-1);
- }
- print "\n";
+ &runCommand($command);
}
$command = $dir."rsem-preref $ARGV[1].transcripts.fa $polyAChoice $ARGV[1] -l $polyALen";
if ($polyAChoice == 2) { $command .= " -f $subsetFile"; }
if ($no_ntog) { $command .= " --no-ntog"; }
if ($quiet) { $command .= " -q"; }
-
-print "$command\n";
-$status = system($command);
-if ($status != 0) {
- print "rsem-preref failed! Please check if you provide correct parameters/options for the pipeline!\n";
- exit(-1);
-}
-print "\n";
+
+&runCommand($command);
if (!$no_bowtie) {
$command = $bowtie_path."bowtie-build -f";
if ($quiet) { $command .= " -q"; }
$command .= " $ARGV[1].idx.fa $ARGV[1]";
-
- print "$command\n";
- $status = system($command);
- if ($status != 0) {
- print "bowtie-build failed! Please check if you provide correct parameters/options for the pipeline!\n";
- exit(-1);
- }
- print "\n";
+
+ &runCommand($command);
}
__END__
=head1 SYNOPSIS
-=over
-
- rsem-prepare-reference [options] reference_fasta_file(s) reference_name
-
-=back
+rsem-prepare-reference [options] reference_fasta_file(s) reference_name
=head1 ARGUMENTS
gene_id transcript_id
with the two fields separated by a tab character.
-
+
If you are using a GTF file for the "UCSC Genes" gene set from the UCSC Genome Browser, then the "knownIsoforms.txt" file (obtained from the "Downloads" section of the UCSC Genome Browser site) is of this format.
If this option is off, then the mapping of isoforms to genes depends on whether the --gtf option is specified. If --gtf is specified, then RSEM uses the "gene_id" and "transcript_id" attributes in the GTF file. Otherwise, RSEM assumes that each sequence in the reference sequence files is a separate gene.
=item B<--no-polyA>
-Do not add poly(A) tails to the end of reference isoforms. (Default: add poly(A) tails to all transcripts)
+Do not add poly(A) tails to the end of reference isoforms. (Default: adding poly(A) tails to all transcripts)
=item B<--no-polyA-subset> <file>
=item B<--bowtie-path> <path>
-The path to the bowtie executables. (Default: the path to bowtie executables is assumed to be in the user's PATH environment variable)
+The path to the Bowtie executables. (Default: the path to Bowtie executables is assumed to be in the user's PATH environment variable)
=item B<--no-bowtie>
'reference_name.grp', 'reference_name.ti', 'reference_name.seq', 'reference_name.idx.fa', and 'reference_name.chrlist' are used by RSEM internally.
-B<'reference_name.transcripts.fa'> contains the extracted reference transcripts in FASTA format. Poly(A) tails are not added.
+B<'reference_name.transcripts.fa'> contains the extracted reference transcripts in FASTA format. Poly(A) tails are added unless '--no-polyA' is set.
=head1 EXAMPLES
-1) Suppose we have mouse RNA-Seq data and want to use the UCSC mm9 version of the mouse genome. We have downloaded the UCSC Genes transcript annotations in GTF format (as mm9.gtf) using the Table Browser and the knownIsoforms.txt file for mm9 from the UCSC Downloads. We also have all chromosome files for mm9 in the directory '/data/mm9'. We want to put the generated reference files under '/ref' with name 'mm9'. We'll add poly(A) tails with length 125. Please note that GTF files generated from UCSC's Table Browser do not contain isoform-gene relationship information. For the UCSC Genes annotation, this information can be obtained from the knownIsoforms.txt file. Suppose we want to build Bowtie indices and Bowtie executables are found in '/sw/bowtie'.
+1) Suppose we have mouse RNA-Seq data and want to use the UCSC mm9 version of the mouse genome. We have downloaded the UCSC Genes transcript annotations in GTF format (as mm9.gtf) using the Table Browser and the knownIsoforms.txt file for mm9 from the UCSC Downloads. We also have all chromosome files for mm9 in the directory '/data/mm9'. We want to put the generated reference files under '/ref' with name 'mouse_125'. We'll add poly(A) tails with length 125. Please note that GTF files generated from UCSC's Table Browser do not contain isoform-gene relationship information. For the UCSC Genes annotation, this information can be obtained from the knownIsoforms.txt file. Suppose we want to build Bowtie indices and Bowtie executables are found in '/sw/bowtie'.
There are two ways to write the command:
--transcript-to-gene-map knownIsoforms.txt \
--bowtie-path /sw/bowtie \
/data/mm9/chr1.fa,/data/mm9/chr2.fa,...,/data/mm9/chrM.fa \
- /ref/mm9
+ /ref/mouse_125
OR
--transcript-to-gene-map knownIsoforms.txt \
--bowtie-path /sw/bowtie \
/data/mm9 \
- /ref/mm9
+ /ref/mouse_125
-2) Suppose we only have transcripts from EST tags in 'mm9.fasta'. In addition, we also have isoform-gene information in 'mapping.txt'. We do not want to add any poly(A) tails. The reference_name will be set to 'mm9'. In addition, we do not want to build Bowtie indices, and will use an alternative aligner to align reads against the 'mm9.idx.fa' output file:
+2) Suppose we only have transcripts from EST tags in 'mm9.fasta'. In addition, we also have isoform-gene information in 'mapping.txt'. We do not want to add any poly(A) tails. The reference_name will be set to 'mouse_0'. In addition, we do not want to build Bowtie indices, and will use an alternative aligner to align reads against the 'mouse_0.idx.fa' output file:
rsem-prepare-reference --transcript-to-gene-map mapping.txt \
--no-polyA \
--no-bowtie \
mm9.fasta \
- mm9
+ mouse_0
=cut
projects / rsem.git / blobdiff