Reference data preparation - FunGen-xQTL Consortium

Prepare the genome, annotation, alignment indices and optional genomic reference products used by downstream xQTL workflows.

Miniprotocol Timing¶

This is the total duration for all phases; module-specific timings appear on their respective pages.

Timing: TBD

Overview¶

This mini-protocol walks through how the reference data used throughout the xQTL Protocol are downloaded, formatted, and indexed. It is an introductory, top-to-bottom guide: each step calls a single workflow from reference_data_preparation.ipynb, with additional workflows from VCF_QC.ipynb, generalized_TADB.ipynb, and rss_ld_sketch.ipynb.

The commands are organized as selectable routes rather than one mandatory 13-step chain. Steps 1–10 prepare the core genome and the RNA-seq, RSEM, or Picard resources needed by the selected route. Steps 11–13 are independent optional workflows for dbSNP annotation, TAD-based association windows, and RSS LD sketches.

Steps¶

Choose a route before running commands; the 13 commands are not one mandatory chain.

Analysis goal	Commands to run, in order	Inputs
Basic genome and gene annotation	1 → 2 → 5 → 6	-
Bulk RNA-seq alignment with STAR	1 → 2 → 3 → 5 → 6 → 7 → 8	-
Transcript quantification with RSEM	1 → 2 → 3 → 5 → 6 → 7 → 8 → 9	-
Picard RNA-seq quality control	1 → 2 → 3 → 5 → 6 → 7 → 8 → 10	-
Add dbSNP rsIDs to a study VCF	4 → 11	Internet access and a study VCF
TAD-based association windows	12	Tissue-specific TAD calls and gene coordinates
RSS fine-mapping or TWAS	13	LD-block definitions and genotype VCFs

Run only the row matching the intended analysis goal, following its commands in numerical order.

1. Download the human genome¶

What it does: Download the GRCh38 reference FASTA used by sequence-aware tools.

sos run pipeline/reference_data_preparation.ipynb download_hg_reference --cwd output/reference_data

2. Download the gene annotation¶

What it does: Download the Ensembl gene annotation used to define genes and transcripts.

sos run pipeline/reference_data_preparation.ipynb download_gene_annotation --cwd output/reference_data

3. Download the ERCC reference¶

What it does: Download ERCC spike-in sequences and annotations for RNA-seq reference construction.

sos run pipeline/reference_data_preparation.ipynb download_ercc_reference --cwd output/reference_data

4. Download dbSNP¶

What it does: Download the dbSNP variant resource used when rsID annotation is required.

sos run pipeline/reference_data_preparation.ipynb download_dbsnp --cwd output/reference_data

5. Format and index the genome¶

What it does: Remove unsupported alternate sequences, append ERCC sequences and create FASTA indices.

sos run pipeline/reference_data_preparation.ipynb hg_reference --cwd output/reference_data --ercc-reference output/reference_data/ERCC92.fa --hg-reference output/reference_data/GRCh38_full_analysis_set_plus_decoy_hla.fa

6. Format the gene annotation¶

What it does: Add chromosome prefixes and construct the protocol-compatible gene models.

sos run pipeline/reference_data_preparation.ipynb hg_gtf --cwd output/reference_data --hg-gtf output/reference_data/Homo_sapiens.GRCh38.103.chr.gtf --hg-reference output/reference_data/GRCh38_full_analysis_set_plus_decoy_hla.noALT_noHLA_noDecoy.fasta --stranded

7. Combine gene and ERCC annotations¶

What it does: Combine the processed human and ERCC annotations into the GTF used downstream.

sos run pipeline/reference_data_preparation.ipynb gene_annotation --cwd output/reference_data --ercc-gtf tests/fixtures/reference_data_preparation/ERCC92.gtf --hg-gtf output/reference_data/Homo_sapiens.GRCh38.103.chr.gtf --hg-reference output/reference_data/GRCh38_full_analysis_set_plus_decoy_hla.noALT_noHLA_noDecoy.fasta --stranded

8. Build the STAR index¶

What it does: Build the STAR genome index used for RNA-seq alignment.

sos run pipeline/reference_data_preparation.ipynb STAR_index --cwd output/reference_data --hg-reference output/reference_data/GRCh38_full_analysis_set_plus_decoy_hla.noALT_noHLA_noDecoy_ERCC.fasta --numThreads 10 --mem 40G

9. Build the RSEM index¶

What it does: Build the RSEM reference used for transcript-level expression quantification.

sos run pipeline/reference_data_preparation.ipynb RSEM_index --cwd output/reference_data --hg-reference output/reference_data/GRCh38_full_analysis_set_plus_decoy_hla.noALT_noHLA_noDecoy_ERCC.fasta --hg-gtf output/reference_data/Homo_sapiens.GRCh38.103.chr.reformatted.ERCC.gtf

10. Generate RefFlat annotation¶

What it does: Convert the processed GTF into the RefFlat annotation used by Picard RNA-seq QC.

sos run pipeline/reference_data_preparation.ipynb RefFlat_generation --cwd output/reference_data --hg-gtf output/reference_data/Homo_sapiens.GRCh38.103.chr.reformatted.ERCC.gtf

11. Annotate a study VCF with dbSNP identifiers (optional)¶

What it does: Add rsIDs to a study VCF when downstream tools require named variants.

sos run pipeline/VCF_QC.ipynb dbsnp_annotate \
    --genoFile tests/fixtures/vcf_qc/protocol_example.genotype.chr22.vcf.gz \
    --cwd output/vcf_qc

12. Construct generalized TAD boundaries (optional)¶

What it does: Combine tissue-specific TAD calls into customized cis-association windows.

sos run pipeline/generalized_TADB.ipynb default \
    --tad-input tests/fixtures/generalized_TADB/protocol_example.brain_TADs.txt \
    --gene-coords tests/fixtures/generalized_TADB/protocol_example.gene_start_end.tsv \
    --cwd output/tadb

13. Construct an LD sketch for RSS analysis (optional)¶

What it does: Generate, process and merge LD-sketch products for summary-statistics regression.

sos run pipeline/rss_ld_sketch.ipynb generate_W --n-samples 60 --output-dir output/rss_ld_sketch --B 50 --seed 123 --cwd output/rss_ld_sketch

sos run pipeline/rss_ld_sketch.ipynb process_block --ld-block-file tests/fixtures/rss_ld_sketch/protocol_example.ld_blocks.bed --chrom 22 --vcf-base output/rss_ld_sketch/rss_ld_sketch --vcf-prefix protocol_example.genotype. --output-dir output/rss_ld_sketch --W-matrix output/rss_ld_sketch/W_B50.rds --B 50 --cohort-id protocol_example --cwd output/rss_ld_sketch

sos run pipeline/rss_ld_sketch.ipynb merge_chrom --output-dir output/rss_ld_sketch --cohort-id protocol_example --chrom 22 --cwd output/rss_ld_sketch

Output¶

Route	Output products and relative paths
Downloads	`output/reference_data/GRCh38_full_analysis_set_plus_decoy_hla.fa`; `output/reference_data/Homo_sapiens.GRCh38.103.chr.gtf`; `output/reference_data/ERCC92.fa`; `tests/fixtures/reference_data_preparation/ERCC92.gtf`; `output/reference_data/00-All.vcf.gz`; `output/reference_data/00-All.vcf.gz.tbi`
Basic genome	`output/reference_data/GRCh38_full_analysis_set_plus_decoy_hla.noALT_noHLA_noDecoy.fasta`; `output/reference_data/GRCh38_full_analysis_set_plus_decoy_hla.noALT_noHLA_noDecoy_ERCC.fasta`; `output/reference_data/GRCh38_full_analysis_set_plus_decoy_hla.noALT_noHLA_noDecoy_ERCC.dict` and FASTA index
Basic annotation	`output/reference_data/Homo_sapiens.GRCh38.103.chr.reformatted.gtf`; `output/reference_data/Homo_sapiens.GRCh38.103.chr.reformatted.collapse_only.gene.gtf`; `output/reference_data/Homo_sapiens.GRCh38.103.chr.reformatted.ERCC.gtf`
STAR	`output/reference_data/STAR_Index/`
RSEM	`output/reference_data/RSEM_Index/`
Picard QC	`output/reference_data/Homo_sapiens.GRCh38.103.chr.reformatted.ERCC.ref.flat`
dbSNP study-VCF annotation	`output/reference_data/00-All.variants.gz` lookup plus the annotated VCF under `output/vcf_qc/`
TAD windows	`tests/fixtures/generalized_TADB/expected/generalized_TAD.tsv`; `tests/fixtures/generalized_TADB/expected/generalized_TADB.tsv`; `tests/fixtures/generalized_TADB/expected/TADB_enhanced_cis.bed`; `tests/fixtures/generalized_TADB/expected/extended_TADB.bed`
RSS LD sketch	`output/rss_ld_sketch/W_B50.rds`; block-level `output/rss_ld_sketch/protocol_example.<block>.dosage.gz`; chromosome-level `output/rss_ld_sketch/protocol_example.22.pgen` products

Anticipated Results¶

After the selected route among steps 1–10, the reference directory contains the formatted genome and the annotations required by that analysis. The basic route produces the processed genome and gene annotation; the RNA-seq routes additionally produce the combined ERCC annotation and, as selected, the STAR index, RSEM index, or RefFlat annotation.

Optional step 11 produces a study VCF annotated with dbSNP rsIDs. Optional step 12 produces generalized TAD and TAD-boundary files such as generalized_TAD.tsv, generalized_TADB.tsv, TADB_enhanced_cis.bed, and extended_TADB.bed. Optional step 13 produces the random-projection matrix and the block- and chromosome-level LD-sketch products consumed by RSS workflows.

Proceed to the molecular-phenotype or genotype-preprocessing mini-protocol required by the selected analysis route.

Command Interface¶

List the workflows and parameters available in the reference data preparation module.

sos run pipeline/reference_data_preparation.ipynb -h