RSEM-PR...FERENCE(1)

πŸ“¦
version 1.3.3 Β· 5 versions
, in ubuntu noble, ubuntu resolute

πŸ’– RSEM-PR...FERENCE(1) πŸ’–

okay bestie, welcome to the User Contributed Perl Documentation glow-up of rsem-prepare-reference. yes she's long. yes she's worth it. grab a lil snack πŸ“


✨ NAME ✨

rsem-prepare-reference β€” get your transcript references all cute and ready for RSEM, and optionally build BOWTIE / BOWTIE2 / STAR / HISAT2 (transcriptome) indices while you're at it. no cap, she does a lot 🫢


✨ SYNOPSIS ✨

rsem-prepare-reference [options] reference_fasta_file(s) reference_name

✨ ARGUMENTS ✨

  • reference_fasta_file(s) β€” okay so this is EITHER a comma-separated list of Multi-FASTA formatted files OR just a directory name. If you hand her a directory, RSEM reads every file in it that ends in .fa or .fasta. Now β€” what should those files contain?? plot twist: transcript sequences, or an entire genome. which one depends on whether you're using the --gtf option. πŸ’…
  • reference_name β€” the name of your reference. RSEM generates a bunch of reference-related files and prefixes them all with this name, so pick something you'll recognize. she can absolutely contain path info too, like '/ref/mm9'. ✨

✨ OPTIONS ✨

  • --gtf <file> β€” if this is on, RSEM assumes reference_fasta_file(s) is holding the sequence of a genome, and she'll extract transcript reference sequences using the gene annotations in <file>, which needs to be in GTF format. 🧬

If both this and --gff3 are off, RSEM assumes reference_fasta_file(s) holds the reference transcripts. in that case, she assumes the name of each sequence in the Multi-FASTA files is its transcript_id. periodt.

(Default: off)

  • --gff3 <file> β€” same energy, but the annotation file is in GFF3 format instead of GTF. RSEM first converts it to GTF format under the filename reference_name.gtf. ⚠️ please please please make sure that reference_name.gtf does not already exist. (Default: off)
  • --gff3-RNA-patterns <pattern> β€” <pattern> is a comma-separated list of transcript categories, like "mRNA,rRNA". only transcripts matching the <pattern> get extracted. βœ‚οΈ (Default: "mRNA")
  • --gff3-genes-as-transcripts β€” this one's for the untypical girlies, like viruses, whose GFF3 files only contain genes. RSEM will assume each gene is a unique transcript when she converts the GFF3 file into GTF format. 🦠
  • --trusted-sources <sources> β€” <sources> is a comma-separated list of trusted sources, e.g. "ENSEMBL,HAVANA". only transcripts coming from these sources get extracted. if this option is off, all sources are accepted. (Default: off)
  • --transcript-to-gene-map <file> β€” use info from <file> to map transcript (isoform) ids β†’ gene ids. each line of <file> should look like this:
gene_id transcript_id

with the two fields separated by a tab character.

if you're using a GTF file for the "UCSC Genes" gene set from the UCSC Genome Browser, then the "knownIsoforms.txt" file (grabbed from the "Downloads" section of the UCSC Genome Browser site) is exactly this format. iykyk 🫢

if this option is off, then the isoform→gene mapping depends on whether --gtf is specified. if --gtf is specified, RSEM uses the "gene_id" and "transcript_id" attributes in the GTF file. otherwise, RSEM assumes every sequence in the reference sequence files is a separate gene.

(Default: off)

  • --allele-to-gene-map <file> β€” use info from <file> to provide gene_id and transcript_id info for each allele-specific transcript. each line of <file> should look like:
gene_id transcript_id allele_id

with the fields separated by a tab character.

this option is designed for quantifying allele-specific expression. 🧬✨ it is only valid if the --gtf option is not specified. allele_id should be the sequence names presented in the Multi-FASTA-formatted files.

(Default: off)

  • --polyA β€” add poly(A) tails to the end of all reference isoforms. the length of the poly(A) tail added is specified by the --polyA-length option. STAR aligner users may not want this one, jsyk. (Default: do not add poly(A) tail to any of the isoforms)
  • --polyA-length <int> β€” the length of the poly(A) tails to be added. (Default: 125)
  • --no-polyA-subset <file> β€” only meaningful if --polyA is specified. do not add poly(A) tails to those transcripts listed in <file>. <file> is a file containing a list of transcript_ids. (Default: off)
  • --bowtie β€” build Bowtie indices. (Default: off)
  • --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)
  • --bowtie2 β€” build Bowtie 2 indices. (Default: off)
  • --bowtie2-path <path> β€” the path to the Bowtie 2 executables. (Default: the path to Bowtie 2 executables is assumed to be in the user's PATH environment variable)
  • --star β€” build STAR indices. (Default: off)
  • --star-path <path> β€” the path to STAR's executable. (Default: the path to STAR executable is assumed to be in user's PATH environment variable)
  • --star-sjdboverhang <int> β€” length of the genomic sequence around annotated junction. it's only used for STAR, to build the splice junctions database, and it's not needed for Bowtie or Bowtie2. it gets passed as the --sjdbOverhang option to STAR. according to STAR's manual, its ideal value is max(ReadLength)-1 β€” e.g. for 2x101 paired-end reads, the ideal value is 101-1=100. in most cases, the default value of 100 will work as well as the ideal value. 🌈 (Default: 100)
  • --hisat2-hca β€” build HISAT2 indices on the transcriptome according to the Human Cell Atlas (HCA) SMART-Seq2 pipeline. (Default: off)
  • --hisat2-path <path> β€” the path to the HISAT2 executables. (Default: the path to HISAT2 executables is assumed to be in the user's PATH environment variable)
  • -p/--num-threads <int> β€” number of threads to use for building STAR's genome indices. πŸ”₯ (Default: 1)
  • -q/--quiet β€” suppress the output of logging information. 🀫 (Default: off)
  • -h/--help β€” show help information.

✨ PRIOR-ENHANCED RSEM OPTIONS ✨

  • --prep-pRSEM β€” a Boolean indicating whether to prepare reference files for pRSEM, including building Bowtie indices for a genome and selecting training set isoforms. the index files will be used for aligning ChIP-seq reads in prior-enhanced RSEM, and the training set isoforms will be used for learning prior. πŸŽ€ a path to Bowtie executables and a mappability file in bigWig format are required when this option is on. currently, Bowtie2 is not supported for prior-enhanced RSEM. (Default: off)
  • --mappability-bigwig-file <string> β€” full path to a whole-genome mappability file in bigWig format. this file is required for running prior-enhanced RSEM. it's used for selecting a training set of isoforms for prior-learning. this file can either be downloaded from the UCSC Genome Browser or generated by GEM (Derrien et al., 2012, PLoS One). (Default: "")

✨ DESCRIPTION ✨

this program extracts/preprocesses the reference sequences for RSEM and prior-enhanced RSEM. πŸ’… she can optionally build Bowtie indices (with the --bowtie option) and/or Bowtie 2 indices (with the --bowtie2 option) using their default parameters. she can also optionally build STAR indices (with the --star option) using parameters from ENCODE3's STAR-RSEM pipeline. 🌟

for prior-enhanced RSEM, she can build Bowtie genomic indices and select training set isoforms (with the options --prep-pRSEM and --mappability-bigwig-file <string>).

if an alternative aligner is to be used, indices for that particular aligner can be built from either reference_name.idx.fa or reference_name.n2g.idx.fa (see OUTPUT for details πŸ‘€). this program is used in conjunction with the rsem-calculate-expression program. bestie duo fr fr. 🫢


✨ OUTPUT ✨

this program will generate reference_name.grp, reference_name.ti, reference_name.transcripts.fa, reference_name.seq, reference_name.chrlist (if --gtf is on), reference_name.idx.fa, reference_name.n2g.idx.fa, optional Bowtie/Bowtie 2 index files, and optional STAR index files. so many files. she's an overachiever. πŸ‘‘

reference_name.grp, reference_name.ti, reference_name.seq, and reference_name.chrlist are used by RSEM internally. 🀫

reference_name.transcripts.fa contains the extracted reference transcripts in Multi-FASTA format. poly(A) tails are not added, and it may contain lower case bases in its sequences if the corresponding genomic regions are soft-masked. ✨

reference_name.idx.fa and reference_name.n2g.idx.fa are used by aligners to build their own indices. in these two files, all sequence bases are converted into upper case. in addition, poly(A) tails are added if the --polyA option is set. the only difference between reference_name.idx.fa and reference_name.n2g.idx.fa is that reference_name.n2g.idx.fa additionally converts all 'N' characters to 'G' characters. this conversion is in particular desired for aligners (e.g. Bowtie) that do not allow reads to overlap with 'N' characters in the reference sequences. otherwise, reference_name.idx.fa should be used to build the aligner's index files. RSEM uses reference_name.idx.fa to build Bowtie 2 indices and reference_name.n2g.idx.fa to build Bowtie indices. 🎯 for visualizing the transcript-coordinate-based BAM files generated by RSEM in IGV, reference_name.idx.fa should be imported as a "genome" (see Visualization section in README.md for details).

if the whole genome is indexed for prior-enhanced RSEM, all the index files will be generated with prefix as reference_name_prsem. selected isoforms for the training set are listed in the file reference_name_prsem.training_tr_crd 🌸


🌸 EXAMPLES 🌸

okay bestie, this is the fun part β€” real-life scenarios, real commands, copy-paste energy. let's go πŸ’…


✨ Example 1: mouse RNA-Seq + UCSC mm9 + Bowtie ✨

So picture this: we got mouse RNA-Seq data and we wanna use the UCSC mm9 version of the mouse genome. We've already downloaded the UCSC Genes transcript annotations in GTF format (saved as mm9.gtf) using the Table Browser, plus the knownIsoforms.txt file for mm9 from the UCSC Downloads. We also have all the chromosome files for mm9 sitting in the directory /data/mm9. 🐭

We want the generated reference files to land under /ref with the name mouse_0. And we're not adding any poly(A) tails.

⚠️ heads up bestie: GTF files generated from UCSC's Table Browser do not contain isoform-gene relationship information. For the UCSC Genes annotation, that info can be obtained from the knownIsoforms.txt file. iykyk.

Also we wanna build Bowtie indices, and the Bowtie executables live in /sw/bowtie.

There are two ways to write this command, so pick your fighter:

 rsem-prepare-reference --gtf mm9.gtf \
                        --transcript-to-gene-map knownIsoforms.txt \
                        --bowtie \
                        --bowtie-path /sw/bowtie \
                        /data/mm9/chr1.fa,/data/mm9/chr2.fa,...,/data/mm9/chrM.fa \
                        /ref/mouse_0

OR

 rsem-prepare-reference --gtf mm9.gtf \
                        --transcript-to-gene-map knownIsoforms.txt \
                        --bowtie \
                        --bowtie-path /sw/bowtie \
                        /data/mm9 \
                        /ref/mouse_0

✨ Example 2: same vibe, but now with Bowtie 2 too ✨

Suppose we also wanna build Bowtie 2 indices in the above example, and Bowtie 2 executables live in /sw/bowtie2. The command becomes:

 rsem-prepare-reference --gtf mm9.gtf \
                        --transcript-to-gene-map knownIsoforms.txt \
                        --bowtie \
                        --bowtie-path /sw/bowtie \
                        --bowtie2 \
                        --bowtie2-path /sw/bowtie2 \
                        /data/mm9 \
                        /ref/mouse_0

two aligners, one command, we love to see it 🫢


✨ Example 3: STAR era ✨

Now we wanna build STAR indices in the above example and save the index files under /ref with the name mouse_0. Assuming the STAR executable is /sw/STAR, the command will be:

 rsem-prepare-reference --gtf mm9.gtf \
                        --transcript-to-gene-map knownIsoforms.txt \
                        --star \
                        --star-path /sw/STAR \
                        -p 8 \
                        /data/mm9/chr1.fa,/data/mm9/chr2.fa,...,/data/mm9/chrM.fa \
                        /ref/mouse_0

OR

 rsem-prepare-reference --gtf mm9.gtf \
                        --transcript-to-gene-map knownIsoforms.txt \
                        --star \
                        --star-path /sw/STAR \
                        -p 8 \
                        /data/mm9
                        /ref/mouse_0

STAR genome index files will be saved under /ref/. slay ⭐


✨ Example 4: prior-enhanced RSEM ✨

Suppose we wanna prepare references for prior-enhanced RSEM in the above example. In this scenario, both STAR and Bowtie are required to build genomic indices β€” STAR for RNA-seq reads and Bowtie for ChIP-seq reads. Assuming their executables are under /sw/STAR and /sw/Bowtie, respectively. Also assuming the mappability file for the mouse genome is /data/mm9.bigWig. The command will be:

 rsem-prepare-reference --gtf mm9.gtf \
                        --transcript-to-gene-map knownIsoforms.txt \
                        --star \
                        --star-path /sw/STAR \
                        -p 8 \
                        --prep-pRSEM \
                        --bowtie-path /sw/Bowtie \
                        --mappability-bigwig-file /data/mm9.bigWig \
                        /data/mm9/chr1.fa,/data/mm9/chr2.fa,...,/data/mm9/chrM.fa \
                        /ref/mouse_0

OR

 rsem-prepare-reference --gtf mm9.gtf \
                        --transcript-to-gene-map knownIsoforms.txt \
                        --star \
                        --star-path /sw/STAR \
                        -p 8 \
                        --prep-pRSEM \
                        --bowtie-path /sw/Bowtie \
                        --mappability-bigwig-file /data/mm9.bigWig \
                        /data/mm9
                        /ref/mouse_0

Both STAR and Bowtie's index files will be saved under /ref/. Bowtie files will have the name prefix mouse_0_prsem. πŸŽ€


✨ Example 5: EST tags, poly(A) tails, no Bowtie ✨

Okay so last one: suppose we only have transcripts from EST tags stored in mm9.fasta, and the isoform-gene information stored in mapping.txt. We wanna add 125bp long poly(A) tails to all transcripts. The reference_name is set as mouse_125.

On top of that, we do not wanna build Bowtie/Bowtie 2 indices, and we'll use an alternative aligner to align reads against either mouse_125.idx.fa or mouse_125.idx.n2g.fa:

 rsem-prepare-reference --transcript-to-gene-map mapping.txt \
                        --polyA \
                        mm9.fasta \
                        mouse_125

no bowtie, no problem πŸ’–


and that's a wrap on this whole page bestie!! this doc came to you courtesy of perl v5.38.2, dated 2024-04-14, and it's part of the User Contributed Perl Documentation ✨ now go build those references, you absolute legend πŸ‘‘πŸ¦‹πŸŒˆ

πŸ“¦ where this page comes from