starsolo-spliced-unspliced
Quantifies spliced, unspliced, and ambiguous count matrices from 10x FASTQs using STARsolo's Velocyto mode — a faster drop-in replacement for velocyto.py plus CellRanger. Use when preparing upstream i
SKILL.mdskills/starsolo-spliced-unspliced/SKILL.md
GitHubSTARsolo: Spliced/Unspliced Count Quantification
Purpose: Generate spliced, unspliced, and ambiguous count matrices from FASTQ files using STARsolo's Velocyto mode — a drop-in replacement for velocyto.py that is ~10x faster than CellRanger.
When to Use
- You need spliced/unspliced/ambiguous layers for RNA velocity (scVelo, CellRank, Chronocell)
- Starting from FASTQ files (not BAMs) — STARsolo maps + quantifies in one pass
- Replacing a velocyto.py + CellRanger pipeline with a single STAR run
- Need CellRanger-compatible output format with velocity layers added
Not for: If you only have BAMs (no FASTQs), use velocyto.py directly or STARsolo's BAM-input mode.
Key Concept: soloFeatures Velocyto
The critical parameter is:
--soloFeatures Gene Velocyto
This calculates Spliced, Unspliced, and Ambiguous counts per cell per gene, similar to velocyto.py (La Manno et al.). Gene is required alongside Velocyto.
For maximum information, request all features in one run:
--soloFeatures Gene GeneFull SJ Velocyto
Quick Start: 10X Chromium V3
1. Build Genome Index (one-time)
STAR --runMode genomeGenerate \
--runThreadN 32 \
--genomeDir /path/to/star_index \
--genomeFastaFiles /path/to/genome.fa \
--sjdbGTFfile /path/to/genes.gtf
Use CellRanger-filtered GTF for best concordance with CellRanger counts. Add --genomeSAsparseD 3 to reduce RAM to ~16GB (at 30-50% speed cost).
2. Run STARsolo with Velocyto
STAR \
--genomeDir /path/to/star_index \
--readFilesIn Read2.fastq.gz Read1.fastq.gz \
--readFilesCommand zcat \
--runThreadN 32 \
--soloType CB_UMI_Simple \
--soloCBwhitelist /path/to/3M-february-2018.txt \
--soloCBstart 1 --soloCBlen 16 \
--soloUMIstart 17 --soloUMIlen 12 \
--soloFeatures Gene Velocyto \
--soloCBmatchWLtype 1MM_multi_Nbase_pseudocounts \
--soloUMIfiltering MultiGeneUMI_CR \
--soloUMIdedup 1MM_CR \
--clipAdapterType CellRanger4 \
--outFilterScoreMin 30 \
--soloCellFilter EmptyDrops_CR \
--outSAMtype BAM SortedByCoordinate \
--outSAMattributes NH HI nM AS CR UR CY UY CB UB \
--outFileNamePrefix /path/to/output/sample1_
Critical read order: --readFilesIn cDNA_Read2 Barcode_Read1 (Read2 first for 10X).
For multiple lanes:
--readFilesIn R2_L1.fq.gz,R2_L2.fq.gz R1_L1.fq.gz,R1_L2.fq.gz
3. Output Structure
sample1_Solo.out/
Gene/
raw/ # Unfiltered: barcodes.tsv, features.tsv, matrix.mtx
filtered/ # Cell-called: barcodes.tsv, features.tsv, matrix.mtx
Velocyto/
raw/
spliced.mtx # Exonic reads
unspliced.mtx # Intronic reads
ambiguous.mtx # Reads overlapping exon-intron boundaries
barcodes.tsv
features.tsv
filtered/
spliced.mtx
unspliced.mtx
ambiguous.mtx
barcodes.tsv
features.tsv
Load into Python (scVelo / Scanpy)
import scanpy as sc
import scvelo as scv
from scipy.io import mmread
import pandas as pd
import numpy as np
from scipy.sparse import csr_matrix
def load_starsolo_velocity(solo_out_dir, use_filtered=True):
"""Load STARsolo Velocyto output into AnnData with spliced/unspliced layers."""
subdir = "filtered" if use_filtered else "raw"
# Load standard gene counts
gene_dir = f"{solo_out_dir}/Gene/{subdir}"
adata = sc.read_mtx(f"{gene_dir}/matrix.mtx").T
adata.obs_names = pd.read_csv(f"{gene_dir}/barcodes.tsv", header=None)[0].values
features = pd.read_csv(f"{gene_dir}/features.tsv", sep="\t", header=None)
adata.var_names = features[0].values
if features.shape[1] > 1:
adata.var["gene_name"] = features[1].values
# Load velocity layers
velo_dir = f"{solo_out_dir}/Velocyto/{subdir}"
barcodes_velo = pd.read_csv(f"{velo_dir}/barcodes.tsv", header=None)[0].values
# Ensure barcode order matches
assert np.array_equal(adata.obs_names, barcodes_velo), \
"Barcode mismatch between Gene and Velocyto outputs"
adata.layers["spliced"] = csr_matrix(mmread(f"{velo_dir}/spliced.mtx").T)
adata.layers["unspliced"] = csr_matrix(mmread(f"{velo_dir}/unspliced.mtx").T)
adata.layers["ambiguous"] = csr_matrix(mmread(f"{velo_dir}/ambiguous.mtx").T)
return adata
# Usage
adata = load_starsolo_velocity("/path/to/sample1_Solo.out")
# Verify layers
print(f"Spliced: {adata.layers['spliced'].sum():.0f} counts")
print(f"Unspliced: {adata.layers['unspliced'].sum():.0f} counts")
print(f"Ambiguous: {adata.layers['ambiguous'].sum():.0f} counts")
# Proceed with scVelo
scv.pp.filter_and_normalize(adata, min_shared_counts=20, n_top_genes=2000)
scv.pp.moments(adata, n_pcs=30, n_neighbors=30)
scv.tl.velocity(adata, mode="dynamical")
scv.tl.velocity_graph(adata)
Assay-Specific Presets
10X Chromium V3 (default)
--soloCBstart 1 --soloCBlen 16 --soloUMIstart 17 --soloUMIlen 12
--soloCBwhitelist 3M-february-2018.txt
10X Chromium V2
--soloCBstart 1 --soloCBlen 16 --soloUMIstart 17 --soloUMIlen 10
--soloCBwhitelist 737K-august-2016.txt
10X 5' (R2 only)
--soloCBstart 1 --soloCBlen 16 --soloUMIstart 17 --soloUMIlen 10
--soloStrand Reverse
--soloCBwhitelist 737K-august-2016.txt
10X 5' Paired-End (R1 longer than 39nt)
--soloBarcodeMate 1 --clip5pNbases 39 0
--soloCBstart 1 --soloCBlen 16 --soloUMIstart 17 --soloUMIlen 10
10X Multiome (Gene Expression)
Same as V3 but use the ARC Multiome whitelist.
Drop-seq / SeqWell
--soloCBstart 1 --soloCBlen 12 --soloUMIstart 13 --soloUMIlen 8
SHARE-seq
--soloCBstart 1 --soloCBlen 24 --soloUMIstart 25 --soloUMIlen 10
Matching CellRanger Versions
CellRanger >= 4.x (recommended)
--clipAdapterType CellRanger4 \
--outFilterScoreMin 30 \
--soloCBmatchWLtype 1MM_multi_Nbase_pseudocounts \
--soloUMIfiltering MultiGeneUMI_CR \
--soloUMIdedup 1MM_CR
CellRanger 3.x
--soloCBmatchWLtype 1MM_multi_Nbase_pseudocounts \
--soloUMIfiltering MultiGeneUMI_CR \
--soloUMIdedup 1MM_CR
CellRanger 2.x (STAR defaults)
No extra flags needed.
Cell Filtering Options
EmptyDrops (CellRanger 3+ style, recommended)
--soloCellFilter EmptyDrops_CR
Parameters: nExpectedCells maxPercentile maxMinRatio indMin indMax umiMin umiMinFracMedian candMaxN FDR simN
Defaults: 3000 0.99 10 45000 90000 500 0.01 20000 0.01 10000
Knee filtering (CellRanger 2 style)
--soloCellFilter CellRanger2.2 3000 0.99 10
No filtering (get raw, filter yourself)
--soloCellFilter None
Re-filter existing raw matrix (no re-mapping)
STAR --runMode soloCellFiltering \
/path/to/Solo.out/Gene/raw/ \
/path/to/output_prefix \
--soloCellFilter EmptyDrops_CR
Multi-Gene Read Handling
For genes with overlapping loci or paralog families:
--soloMultiMappers EM # Maximum Likelihood (best accuracy)
--soloMultiMappers Uniform # Simple fractional assignment
--soloMultiMappers Rescue # First EM iteration
--soloMultiMappers PropUnique # Proportional to unique mappers
Output: UniqueAndMult-EM.mtx etc. alongside standard matrix.mtx (unique only).
BAM Input Mode
If starting from existing BAM files (e.g., CellRanger output):
STAR \
--genomeDir /path/to/star_index \
--readFilesIn input.bam \
--readFilesType SAM SE \
--readFilesCommand samtools view -F 0x100 \
--soloType CB_UMI_Simple \
--soloInputSAMattrBarcodeSeq CR UR \
--soloInputSAMattrBarcodeQual CY UY \
--soloCBwhitelist /path/to/whitelist.txt \
--soloFeatures Gene Velocyto \
--outSAMtype BAM SortedByCoordinate \
--outSAMattributes NH HI nM AS CR UR CY UY CB UB
Use --readFilesSAMattrKeep None to avoid duplicated SAM attributes.
Batch Processing Script
#!/usr/bin/env bash
# run_starsolo_batch.sh - Process multiple 10X V3 samples
set -euo pipefail
GENOME_DIR="/path/to/star_index"
WHITELIST="/path/to/3M-february-2018.txt"
OUTPUT_BASE="/path/to/results"
THREADS=32
declare -A SAMPLES=(
["sample1"]="/path/to/sample1_fastqs"
["sample2"]="/path/to/sample2_fastqs"
)
for sample in "${!SAMPLES[@]}"; do
fastq_dir="${SAMPLES[$sample]}"
out_dir="${OUTPUT_BASE}/${sample}"
mkdir -p "${out_dir}"
# Collect FASTQ files
R1=$(ls "${fastq_dir}"/*_R1_001.fastq.gz | sort | paste -sd,)
R2=$(ls "${fastq_dir}"/*_R2_001.fastq.gz | sort | paste -sd,)
echo "Processing ${sample}..."
STAR \
--genomeDir "${GENOME_DIR}" \
--readFilesIn "${R2}" "${R1}" \
--readFilesCommand zcat \
--runThreadN "${THREADS}" \
--soloType CB_UMI_Simple \
--soloCBwhitelist "${WHITELIST}" \
--soloCBstart 1 --soloCBlen 16 \
--soloUMIstart 17 --soloUMIlen 12 \
--soloFeatures Gene GeneFull Velocyto \
--soloCBmatchWLtype 1MM_multi_Nbase_pseudocounts \
--soloUMIfiltering MultiGeneUMI_CR \
--soloUMIdedup 1MM_CR \
--clipAdapterType CellRanger4 \
--outFilterScoreMin 30 \
--soloCellFilter EmptyDrops_CR \
--outSAMtype BAM SortedByCoordinate \
--outSAMattributes NH HI nM AS CR UR CY UY CB UB \
--outFileNamePrefix "${out_dir}/"
echo "Done: ${sample} -> ${out_dir}/Solo.out/"
done
Resource Requirements
| Component | Recommendation |
|---|---|
| RAM | 32-40GB (genome loading); add --genomeSAsparseD 3 at index build to reduce to ~16GB |
| CPUs | 16-32 threads (--runThreadN) |
| Disk | ~500GB per sample (BAM + intermediate); output matrices are small |
| STAR version | >= 2.7.9a (Velocyto improvements); recommended 2.7.11b |
| Time | ~20-40 min per sample (vs. hours for CellRanger) |
Whitelist Files
| Chemistry | Whitelist | UMI length |
|---|---|---|
| 10X V3 / V3.1 | 3M-february-2018.txt |
12 |
| 10X V2 | 737K-august-2016.txt |
10 |
| 10X 5' | 737K-august-2016.txt |
10 |
| 10X Multiome | ARC whitelist | 12 |
Download from: https://github.com/10XGenomics/cellranger/tree/master/lib/python/cellranger/barcodes/
Prebuilt References (Cloud/Terra)
| Keyword | Description |
|---|---|
GRCh38-2020-A |
Human (GENCODE v32/Ensembl 98) |
mm10-2020-A |
Mouse (GENCODE vM23/Ensembl 98) |
GRCh38-and-mm10-2020-A |
Human + Mouse (barnyard) |
Troubleshooting
| Issue | Solution |
|---|---|
| Empty Velocyto output | Ensure Gene is included: --soloFeatures Gene Velocyto (Velocyto requires Gene) |
| Barcode mismatch | Check read order: R2 (cDNA) first, R1 (barcode) second for 10X |
| Low unspliced fraction | Normal for mature mRNA-enriched protocols; check with adata.layers['unspliced'].sum() / adata.layers['spliced'].sum() — expect 5-20% for 3' protocols |
| RAM exceeded | Use sparse suffix array: --genomeSAsparseD 3 at genome generation step |
| Slow BAM sorting | Increase --limitBAMsortRAM or use --outSAMtype BAM Unsorted if BAM not needed |
| Version mismatch counts | Match GTF exactly between STAR index and CellRanger reference |
soloBarcodeReadLength error |
Set --soloBarcodeReadLength 0 to skip length check |
When not to use
- Do not use if you only have BAM files — use velocyto.py or STARsolo BAM-input mode instead.
- Do not use for standard gene-only counting where velocity layers are unneeded — plain CellRanger or STARsolo Gene is lighter.
See also
rna-velocity-trajectory
Upstream docs: https://github.com/alexdobin/STAR/blob/master/docs/STARsolo.md