ORACLE RNA-seq pipeline

ORACLE is an end-to-end workflow for paired-end bulk RNA-seq. It starts with FASTQ files and produces quality-control reports, aligned reads, gene-level counts, differential-expression results, pathway enrichment, genome-browser tracks, and exploratory transcript assemblies.

The workflow supports human, mouse, and rat references and runs through Snakemake with isolated Conda environments. New users should begin with the complete tutorial; experienced users can use the short run guide.

Scope

This pipeline is designed for:

• paired-end bulk RNA-seq;
• replicated, metadata-driven differential-expression experiments;
• human, mouse, or rat reference genomes;
• gene-level DESeq2 analysis from featureCounts or Salmon/tximport; and
• GO, KEGG, and preranked GSEA interpretation.

It does not currently implement single-end processing, single-cell RNA-seq, transcript-level differential expression, differential transcript usage, or within-species sample-identity testing. Salmon produces transcript estimates, but the implemented DE path summarizes them to genes. StringTie assemblies are exploratory annotations, not a transcript-level statistical test.

Workflow

paired-end FASTQ
  |
  |-- 01  FastQC + fastp                    raw-read QC and trimming
  |-- 01b FastQ Screen, optional/default-off genome-panel screening
  v
trimmed FASTQ
  |
  |-- 02  STAR two-pass alignment
  v
coordinate-sorted BAM
  |
  |-- 03  duplicate marking, RSeQC, Picard, Qualimap, preseq, flagstat
  |        + confidence-aware strandedness inference
  |
  |-- 04  featureCounts                     gene-level counts
  |-- 05  Salmon                            transcript estimates; optional DE source
  |-- 06  StringTie + gffcompare            strand-aware exploratory assembly
  |-- 07  export_for_r                      counts, metadata, tx2gene, resolved config
  |-- 08  MultiQC                           aggregate QC report
  `-- 09  DESeq2 + shrinkage + GO/KEGG/GSEA

Important safeguards include:

• strict validation of sample IDs, FASTQ pairs, design variables, reference levels, comparisons, replicate structure, and design-matrix rank;
• a fully resolved configuration saved with the results and passed unchanged to the R analysis;
• confidence-aware strand inference: ambiguous automatic calls stop the run instead of silently becoming unstranded;
• paired-fragment featureCounts settings that require both mates and exclude chimeric pairs (-B -C);
• fail-closed Salmon selection: requesting Salmon counts cannot silently fall back to featureCounts;
• shrinkage that must succeed rather than silently returning unshrunk effects;
• GO/KEGG over-representation tests using the tested-gene universe; and
• GSEA ranked by the signed Wald statistic by default.

Quick start

git clone https://github.com/challagandla/oracle-rna-seq.git
cd oracle-rna-seq

bash setup.sh
bash prepare_references.sh human

# Replace the template rows and review the full analysis design.
${EDITOR:-nano} config/samples.tsv
${EDITOR:-nano} config/config.yaml

bash run.sh --dry-run
bash run.sh --cores 8

Do not treat the quick start as permission to skip experimental-design review. Before the dry run, confirm that biological replicates are independent, batches are not confounded with the comparison, and every modeled variable is recorded correctly in the sample sheet. The tutorial provides a full preflight.

Installation

setup.sh uses an existing Conda installation or installs a checksum-verified Miniforge release without editing shell startup files. It creates the small oracle-rna-runner environment and smoke-tests the 17 rule environments under workflow/envs/.

bash setup.sh                 # install runner and rule environments
bash setup.sh --check         # verify an existing installation
bash setup.sh --runner-only   # install the runner; create rule envs on demand

No shell activation is required. run.sh locates the runner, selects the appropriate Snakemake Conda-deployment option, and prevents personal R/Python library overrides from leaking into the analysis. The environment YAML files are the declared software specifications; they are not exact platform-specific lockfiles.

Reference indexing is resource intensive. A mammalian STAR index commonly needs tens of gigabytes of RAM and the complete reference plus indexes needs substantial disk space. Actual requirements depend on the assembly and local Conda package cache.

Required inputs

config/samples.tsv must contain one row per biological library and these columns:

sample_id    fastq_r1    fastq_r2

Add every variable used by the statistical model, such as condition, batch, donor, sex, or age. Categorical and continuous variables are declared explicitly in config/config.yaml; numeric-looking category codes must not be allowed to change meaning accidentally.

Set de.primary_factor to the main categorical factor of interest. It anchors the reference level and is the default factor for comparisons that omit one.

The checked-in configuration defaults to human. Build matching references with:

bash prepare_references.sh human

Use mouse or rat only after selecting the same species in config/config.yaml. Reference paths are defined in config/species.yaml.

Optional FastQ Screen

FastQ Screen is disabled by default because its reference panel is large and is not shipped with the repository. To use it:

bash run.sh --cores 8 fastq_screen_genomes
# Then set fastq_screen.enabled: true in config/config.yaml.
bash run.sh --cores 8 contamination

It can reveal reads matching other species, rRNA, vectors, PhiX, or other panel references. It cannot establish that a human sample belongs to a particular human donor and therefore cannot rule out within-species sample swaps.

Common commands

bash run.sh --dry-run
bash run.sh --cores 16
bash run.sh --cores 16 --rerun-incomplete

bash run.sh --cores 8 qc_raw
bash run.sh --cores 8 contamination
bash run.sh --cores 16 align_star
bash run.sh --cores 16 quant
bash run.sh --cores 16 novel_transcripts
bash run.sh --cores 8 export_for_r
bash run.sh --cores 4 run_r
bash run.sh --cores 8 multiqc

bash run.sh --cores 8 --configfile path/to/project-config.yaml

export_for_r materializes the gene-count matrix, metadata, transcript-to-gene map, and resolved configuration. run_r runs the statistical stage and does not require bigWig generation. When de.use_salmon: true, it also requires each configured sample's completed quant.sf.

Main outputs

• results/multiqc/multiqc_report.html — aggregate technical QC.
• results/qc/strandedness.tsv — per-sample strand fractions, calls, and confidence information.
• results/r_inputs/ — the exact inputs handed to the statistical analysis, including resolved_config.yaml.
• results/r_outputs/DE_<comparison>.csv — complete per-gene DESeq2 result for each configured comparison.
• results/r_outputs/analysis_summary.csv — comparison-level summary.
• results/r_outputs/design_matrix.tsv and metadata_used.tsv — the validated model matrix and ordered metadata actually analyzed.
• results/r_outputs/results_names.txt and count_source.txt — resolved DESeq2 result names and the featureCounts or Salmon/tximport source.
• results/r_outputs/sessionInfo.txt — R and package session record.
• results/stringtie/gffcompare.* — exploratory assembled-transcript comparison.
• results/tracks/<sample>.bw — normalized genome-browser coverage.

The R rule log captures both standard output and errors. A completed command is only the start of interpretation: review MultiQC and sample-level diagnostics before trusting differential-expression or enrichment results.

Interpreting DE results

padj controls the false-discovery rate across the tested genes. The configured absolute log2-fold-change cutoff is a separate biological effect-size screen; it does not change the null hypothesis tested by DESeq2. Report both the adjusted p-value and effect size, and inspect expression/QC context rather than reducing the result to a binary label.

GO/KEGG over-representation analysis uses all tested, mappable genes as the background. GSEA uses the signed Wald statistic by default and asks whether a gene set is concentrated toward either end of the full ranked list. Neither method proves pathway activation or causality.

Supported species

Species	Assembly/annotation	OrgDb	KEGG code
Human	GRCh38 / GENCODE v46	org.Hs.eg.db	hsa
Mouse	GRCm39 / GENCODE vM34	org.Mm.eg.db	mmu
Rat	mRatBN7.2 / Ensembl 112	org.Rn.eg.db	rno

Adding another species requires more than a path block: add reference metadata, download/index logic, the correct OrgDb package to the DE environment, KEGG and identifier mappings, then validate reference preparation, DAG resolution, and a small end-to-end run.

Documentation

• Complete beginner-to-advanced tutorial
• Concise run guide
• Best-practice checklist

License and usage

The pipeline's own code is MIT — see LICENSE. It bundles no third-party source code or reference data; tools are installed through Conda and references are downloaded at run time. See THIRD_PARTY_LICENSES.md for the dependency and data license breakdown.

KEGG pathway data and UCSC command-line tools have commercial-use restrictions. Academic and non-commercial users should cite the relevant sources. Commercial users must obtain the applicable licenses or replace/disable those components as described in the third-party license notice.

Contact

Anil Challagandla — [email protected]