Pipelines – omics.co.in

DNA-Seq • Germline Authority: Broad Institute / nf-core/sarek

1. GATK4 Whole-Genome / Exome Germline Variant Discovery

The international reference standard for mapping high-throughput Illumina reads, performing Base Quality Score Recalibration (BQSR), and discovering high-confidence germline SNPs and indels in single samples or joint multi-sample cohorts.

Flow: Raw FASTQ → fastp (QC/Trimming) → BWA-MEM2 (Alignment) → MarkDuplicatesSpark → BaseRecalibrator → HaplotypeCaller (GVCF) → GenotypeGVCFs → DeepVariant / VQSR

# 1. Quality Control & Adapter Trimming
fastp -i sample_R1.fastq.gz -I sample_R2.fastq.gz \
      -o trimmed_R1.fastq.gz -O trimmed_R2.fastq.gz \
      --detect_adapter_for_pe --thread 8 --html fastp_qc.html

# 2. Ultra-Fast Read Mapping with BWA-MEM2 (GRCh38)
bwa-mem2 mem -t 16 -R "@RG\tID:SM01\tSM:Sample01\tPL:ILLUMINA" \
    GRCh38_noalt.fasta trimmed_R1.fastq.gz trimmed_R2.fastq.gz | \
    samtools sort -@ 8 -m 4G -o sample.sorted.bam -

# 3. Mark Duplicate Reads
gatk MarkDuplicatesSpark \
    -I sample.sorted.bam \
    -O sample.dedup.bam \
    -M metrics.txt --conf 'spark.executor.cores=8'

# 4. Base Quality Score Recalibration (BQSR)
gatk BaseRecalibrator \
    -R GRCh38_noalt.fasta -I sample.dedup.bam \
    --known-sites dbsnp_146.hg38.vcf.gz \
    --known-sites Mills_and_1000G_gold_standard.indels.hg38.vcf.gz \
    -O recal_data.table

gatk ApplyBQSR \
    -R GRCh38_noalt.fasta -I sample.dedup.bam \
    --bqsr-recal-file recal_data.table -O sample.recal.bam

# 5. Germline Variant Calling with HaplotypeCaller (GVCF mode)
gatk HaplotypeCaller \
    -R GRCh38_noalt.fasta -I sample.recal.bam \
    -O sample.g.vcf.gz -ERC GVCF --native-pair-hmm-threads 8

DNA-Seq • Oncology Authority: TCGA / GATK Mutect2 / ICGC

2. Somatic Cancer Mutation Calling & Neoantigen Profiling

Paired Tumor-Normal somatic mutation detection designed to identify sub-clonal somatic single nucleotide variants (SNVs), small indels, cross-sample contamination, and functional variant effect annotations.

Flow: Tumor/Normal BAMs → GATK Mutect2 → GetPileupSummaries → CalculateContamination → FilterMutectCalls → Ensembl VEP → vcf2maf

# 1. Somatic Variant Calling with GATK Mutect2
gatk Mutect2 \
    -R GRCh38.fasta \
    -I tumor.bam -tumor Tumor_Sample \
    -I normal.bam -normal Normal_Sample \
    --germline-resource af-only-gnomad.hg38.vcf.gz \
    --panel-of-normals pon.hg38.vcf.gz \
    -O somatic_unfiltered.vcf.gz \
    --f1r2-tar-gz f1r2.tar.gz

# 2. Learn Read Orientation Model (Artifact filtering)
gatk LearnReadOrientationModel -I f1r2.tar.gz -O read-orientation-model.tar.gz

# 3. Estimate Cross-Sample Contamination
gatk GetPileupSummaries -I tumor.bam -V common_variants.vcf.gz -L common_variants.vcf.gz -O tumor_pileups.table
gatk GetPileupSummaries -I normal.bam -V common_variants.vcf.gz -L common_variants.vcf.gz -O normal_pileups.table
gatk CalculateContamination -I tumor_pileups.table -matched normal_pileups.table -O contamination.table

# 4. Filter Somatic Calls
gatk FilterMutectCalls \
    -R GRCh38.fasta -V somatic_unfiltered.vcf.gz \
    --contamination-table contamination.table \
    --ob-priors read-orientation-model.tar.gz \
    -O somatic_filtered.vcf.gz

# 5. Annotation with Ensembl Variant Effect Predictor (VEP)
vep --input_file somatic_filtered.vcf.gz --output_file annotated_somatic.vcf \
    --format vcf --vcf --symbol --terms SO --cache --assembly GRCh38 --fork 8

Long-Reads • SVs Authority: Human Pangenome Consortium / PacBio / ONT

3. Long-Read Structural Variant (SV) Discovery & Phasing

Uncovers large deletions, insertions, inversions, duplications, and chromosomal translocations (>50 bp) using PacBio HiFi or Oxford Nanopore reads with Sniffles2 and CuteSV.

Flow: HiFi/ONT FASTQ → minimap2 (-ax map-hifi) → samtools sort → Sniffles2 (SV calling) → CuteSV → SURVIVOR merge → sniffles genotype

# 1. Long-Read Alignment with Minimap2
minimap2 -ax map-hifi -t 16 --MD -Y GRCh38.fasta hifi_reads.fastq.gz | \
    samtools sort -@ 8 -m 4G -o hifi_aligned.bam -
samtools index hifi_aligned.bam

# 2. Structural Variant Calling with Sniffles2
sniffles --input hifi_aligned.bam \
         --vcf hifi_sniffles_sv.vcf.gz \
         --reference GRCh38.fasta \
         --threads 16 --minsvlen 50 --phase

# 3. Orthogonal SV Calling with CuteSV
cuteSV hifi_aligned.bam GRCh38.fasta hifi_cutesv.vcf work_dir/ \
       --threads 16 --min_size 50 --min_support 3 --genotype

# 4. Consensus Merging with SURVIVOR
ls hifi_sniffles_sv.vcf.gz hifi_cutesv.vcf > vcf_list.txt
SURVIVOR merge vcf_list.txt 1000 2 1 1 0 50 merged_high_confidence_sv.vcf

Transcriptomics • Bulk RNA Authority: nf-core/rnaseq / ENCODE

4. Bulk RNA-Seq Quantification & Differential Expression

Splice-aware read alignment, gene and transcript quantification with STAR and Salmon, followed by negative binomial generalized linear modeling with DESeq2.

Flow: RNA FASTQ → fastp → STAR (2-pass alignment) → featureCounts → DESeq2 (Wald test) → clusterProfiler (GO/KEGG)

# 1. Splice-Aware Alignment with STAR
STAR --runThreadN 16 \
     --genomeDir /ref/STAR_hg38_index \
     --readFilesIn sample_R1.fastq.gz sample_R2.fastq.gz \
     --readFilesCommand zcat \
     --outFileNamePrefix sample_star_ \
     --outSAMtype BAM SortedByCoordinate \
     --outSAMstrandField intronMotif \
     --quantMode GeneCounts

# 2. Exon-Level Counting with featureCounts
featureCounts -T 16 -p -B -C \
              -a gencode.v46.annotation.gtf \
              -o gene_expression_matrix.txt \
              *.bam

# 3. Differential Expression in R with DESeq2
# R script:
library(DESeq2)
counts <- read.table("gene_expression_matrix.txt", header=TRUE, row.names=1, skip=1)[, 6:11]
colData <- data.frame(condition = factor(c("ctrl","ctrl","ctrl","treat","treat","treat")))
dds <- DESeqDataSetFromMatrix(countData = counts, colData = colData, design = ~ condition)
dds <- DESeq(dds)
res <- results(dds, contrast=c("condition","treat","ctrl"), alpha=0.05)
write.csv(as.data.frame(res), "differential_expression_results.csv")

Single-Cell • scRNA-Seq Authority: 10x Genomics / Seurat v5 / Scanpy

5. Single-Cell RNA-Seq (scRNA-Seq) End-to-End Processing

De-multiplexing 10x droplet microfluidic libraries, cell-barcode filtering, Scrublet doublet identification, SCTransform normalization, Harmony batch integration, and CellTypist automated annotation.

Flow: Raw 10x FASTQs → Cell Ranger count / STARsolo → Quality Filter (Mito < 15%) → Scrublet → SCTransform / Harmony → UMAP → CellTypist

# 1. 10x Genomics Cell Ranger Count Pipeline
cellranger count --id=PBMC_10k_scRNA \
                 --transcriptome=/ref/refdata-gex-GRCh38-2020-A \
                 --fastqs=/data/pbmc_fastqs \
                 --sample=pbmc_10k \
                 --localcores=32 --localmem=128

# 2. Python Scanpy Single-Cell Analysis Script
import scanpy as sc
import scrublet as scr

# Load Cell Ranger filtered matrix
adata = sc.read_10x_mtx('PBMC_10k_scRNA/outs/filtered_feature_bc_matrix/')
sc.pp.filter_cells(adata, min_genes=200)
sc.pp.filter_genes(adata, min_cells=3)
adata.var['mt'] = adata.var_names.str.startswith('MT-')
sc.pp.calculate_qc_metrics(adata, qc_vars=['mt'], percent_top=None, log1p=False, inplace=True)
adata = adata[adata.obs.pct_counts_mt < 15, :]

# Doublet Scoring
scrub = scr.Scrublet(adata.X)
doublet_scores, predicted_doublets = scrub.scrub_doublets()
adata.obs['is_doublet'] = predicted_doublets
adata = adata[~adata.obs['is_doublet'], :]

# Normalization & Dimensionality Reduction
sc.pp.normalize_total(adata, target_sum=1e4)
sc.pp.log1p(adata)
sc.pp.highly_variable_genes(adata, min_mean=0.0125, max_mean=3, min_disp=0.5)
sc.pp.pca(adata, svd_solver='arpack')
sc.pp.neighbors(adata, n_neighbors=15, n_pcs=30)
sc.tl.umap(adata)
sc.tl.leiden(adata, resolution=0.8)

Spatial • Visium / Xenium Authority: 10x Genomics Space Ranger / Seurat Spatial

6. Spatial Transcriptomics & Cellular Spot Deconvolution

Mapping spatially resolved gene expression in intact histological tissue sections, identifying spatially variable genes via Moran’s I, and deconvolving multi-cell spots with scRNA-seq reference atlases using RCTD.

Flow: Visium FASTQs + H&E Image → Space Ranger → Seurat Load10X_Spatial → SCTransform → FindSpatiallyVariableFeatures → RCTD Deconvolution

# 1. 10x Space Ranger Alignment & Image Registration
spaceranger count --id=Cortex_Visium \
                  --transcriptome=/ref/refdata-gex-GRCh38-2020-A \
                  --fastqs=/data/visium_fastqs \
                  --sample=cortex_sample \
                  --image=cortex_he_image.tif \
                  --slide=V10M12-001 --area=A1

# 2. R Seurat Spatial Processing
library(Seurat)
spatial_obj <- Load10X_Spatial(data.dir = "Cortex_Visium/outs/")
spatial_obj <- SCTransform(spatial_obj, assay = "Spatial", verbose = FALSE)
spatial_obj <- RunPCA(spatial_obj, assay = "SCT", verbose = FALSE)
spatial_obj <- FindNeighbors(spatial_obj, reduction = "pca", dims = 1:30)
spatial_obj <- FindClusters(spatial_obj, verbose = FALSE)
spatial_obj <- RunUMAP(spatial_obj, reduction = "pca", dims = 1:30)

# Identify Spatially Variable Genes
spatial_obj <- FindSpatiallyVariableFeatures(spatial_obj, assay = "SCT", 
                                             features = VariableFeatures(spatial_obj)[1:1000], 
                                             selection.method = "moransi")

Epigenomics • ATAC-Seq Authority: ENCODE Project / nf-core/atacseq

7. ATAC-Seq Open Chromatin & Regulatory Footprinting

Assay for Transposase-Accessible Chromatin to profile open promoters, distant enhancers, and transcription factor footprints with Tn5 transposase shift corrections (+4 bp/-5 bp).

Flow: ATAC FASTQ → Bowtie2 (-X 2000) → Filter Mito Reads → alignmentSieve (--ATACshift) → MACS3 callpeak → HOMER motif analysis

# 1. Alignment with Bowtie2 (allowing fragment lengths up to 2kb)
bowtie2 -p 16 -X 2000 --very-sensitive \
        -x GRCh38_bowtie2_index \
        -1 atac_R1.fastq.gz -2 atac_R2.fastq.gz | \
        samtools view -bS -q 30 -F 1804 - | \
        samtools sort -@ 8 -o atac_aligned.bam -

# 2. Filter Mitochondrial Reads & Deduplicate
samtools idxstats atac_aligned.bam | cut -f 1 | grep -v chrM | xargs samtools view -b atac_aligned.bam > atac_noMito.bam
picard MarkDuplicates I=atac_noMito.bam O=atac_dedup.bam M=dup_metrics.txt REMOVE_DUPLICATES=true

# 3. Correct for Tn5 Transposase Insertion Offset (+4 / -5 bp)
alignmentSieve --numberOfProcessors 16 --ATACshift -b atac_dedup.bam -o atac_shifted.bam
samtools sort -@ 8 atac_shifted.bam -o atac_shifted.sorted.bam && samtools index atac_shifted.sorted.bam

# 4. Narrow Peak Calling with MACS3
macs3 callpeak -t atac_shifted.sorted.bam -f BAMPE -n sample_atac \
               -g hs --nomodel --shift -100 --extsize 200 -q 0.01 --keep-dup all

# 5. Transcription Factor Motif Discovery with HOMER
findMotifsGenome.pl sample_atac_peaks.narrowPeak hg38 homer_motifs_out/ -size 200

Epigenomics • ChIP-Seq Authority: ENCODE Consortium Uniform Pipeline

8. ChIP-Seq Transcription Factor & Histone Mark Analysis

Crosslinked chromatin immunoprecipitation sequencing for pinpointing transcription factor binding sites (narrow peaks) and histone modification domains (broad peaks, e.g. H3K27me3, H3K4me1).

Flow: ChIP/Input FASTQ → Bowtie2 → Cross-Correlation QC (phantompeakqualtools) → MACS3 callpeak → deepTools bamCoverage → plotHeatmap

# 1. Alignment to Reference
bowtie2 -p 16 -x hg38_index -U chip_sample.fastq.gz | samtools sort -@ 8 -o chip.sorted.bam -
bowtie2 -p 16 -x hg38_index -U input_control.fastq.gz | samtools sort -@ 8 -o input.sorted.bam -

# 2. Peak Calling (Narrow Peaks for Transcription Factors)
macs3 callpeak -t chip.sorted.bam -c input.sorted.bam \
               -f BAM -g hs -n TF_sample -q 0.01 --outdir macs3_results/

# 3. Peak Calling (Broad Peaks for Histone Modifications like H3K27me3)
macs3 callpeak -t chip.sorted.bam -c input.sorted.bam \
               -f BAM -g hs -n Histone_sample --broad --broad-cutoff 0.05

# 4. Generate Normalized Signal BigWig Tracks with deepTools
bamCoverage -b chip.sorted.bam -o chip_normalized.bw \
            --binSize 10 --normalizeUsing RPKM --effectiveGenomeSize 2913022398

Epigenomics • Methylation Authority: Babraham Bioinformatics / nf-core/methylseq

9. Whole-Genome Bisulfite Sequencing (WGBS) Methylation

Base-resolution cytosine methylation detection across CpG, CHG, and CHH contexts, distinguishing 5-methylcytosine from unmethylated cytosines using 3-letter in silico bisulfite alignments.

Flow: Bisulfite FASTQ → Trim Galore → Bismark / bwa-meth → deduplicate_bismark → bismark_methylation_extractor → methylKit DMR analysis

# 1. Bisulfite Read Alignment with Bismark
bismark --multicore 8 --genome /ref/bismark_hg38/ \
        -1 sample_val_1.fq.gz -2 sample_val_2.fq.gz

# 2. PCR Duplicate Removal
deduplicate_bismark --paired sample_val_1_bismark_bt2_pe.bam

# 3. Methylation Extraction Across Contexts (CpG, CHG, CHH)
bismark_methylation_extractor --paired-end --no_overlap \
                              --comprehensive --merge_non_CpG \
                              --bedGraph --counts --buffer_size 16G \
                              sample_val_1_bismark_bt2_pe.deduplicated.bam

Metagenomics • WGS MAGs Authority: Joint Genome Institute / nf-core/mag

10. Shotgun Metagenomic Assembly & MAG Reconstruction

De novo assembly of microbial communities from complex environmental or host samples, multi-algorithm binning (MetaBAT2, MaxBin2), CheckM2 quality validation, and GTDB-Tk taxonomic resolution.

Flow: Raw FASTQ → Bowtie2 (Host Depletion) → MEGAHIT Co-Assembly → MetaBAT2/MaxBin2 → DASTool → CheckM2 (>90% comp, <5% cont) → GTDB-Tk

# 1. Host Read Depletion with Bowtie2
bowtie2 -p 16 -x human_hg38_index \
        -1 raw_R1.fastq.gz -2 raw_R2.fastq.gz \
        --un-conc-gz clean_microbe_R%.fastq.gz -S /dev/null

# 2. De Novo Metagenome Co-Assembly with MEGAHIT
megahit -1 clean_microbe_R1.fastq.gz -2 clean_microbe_R2.fastq.gz \
        -t 32 --min-contig-len 1500 -o megahit_assembly_out/

# 3. Contig Coverage Depth Calculation
bwa-mem2 index megahit_assembly_out/final.contigs.fa
bwa-mem2 mem -t 16 megahit_assembly_out/final.contigs.fa clean_microbe_R1.fastq.gz clean_microbe_R2.fastq.gz | \
    samtools sort -@ 8 -o mapped.sorted.bam -
jgi_summarize_bam_contig_depths --outputDepth depth.txt mapped.sorted.bam

# 4. Binning with MetaBAT2
metabat2 -i megahit_assembly_out/final.contigs.fa -a depth.txt -o bins/bin -t 16

# 5. Bin Quality Benchmark with CheckM2
checkm2 predict --threads 16 --input bins/ -x fa --output-directory checkm2_results/

# 6. Taxonomic Placement with GTDB-Tk
gtdbtk classify_wf --genome_dir bins/ --out_dir gtdbtk_out/ --cpus 16 -x fa

Microbiome • 16S/18S Authority: QIIME 2 / DADA2 / SILVA

11. 16S / 18S / ITS Marker Gene ASV Inference & Diversity

Exact Amplicon Sequence Variant (ASV) resolution with DADA2, error model modeling, taxonomic classification against SILVA 138, and phylogenetic alpha/beta diversity quantification.

Flow: Demux FASTQs → cutadapt (primer trim) → qiime dada2 denoise-paired → qiime feature-classifier (SILVA 138) → qiime diversity core-metrics

# 1. Import Demultiplexed Paired-End FASTQs into QIIME 2 Artifact
qiime tools import \
  --type 'SampleData[PairedEndSequencesWithQuality]' \
  --input-path manifest.tsv \
  --output-path demux_paired.qza \
  --input-format PairedEndFastqManifestPhred33V2

# 2. Denoising with DADA2 to resolve exact Amplicon Sequence Variants (ASVs)
qiime dada2 denoise-paired \
  --i-demultiplexed-seqs demux_paired.qza \
  --p-trim-left-f 17 --p-trim-left-r 21 \
  --p-trunc-len-f 250 --p-trunc-len-r 220 \
  --o-table table.qza \
  --o-representative-sequences rep-seqs.qza \
  --o-denoising-stats denoising-stats.qza \
  --p-n-threads 16

# 3. Taxonomic Classification with SILVA 138 Classifier
qiime feature-classifier classify-sklearn \
  --i-classifier silva-138-99-nb-classifier.qza \
  --i-reads rep-seqs.qza \
  --o-classification taxonomy.qza \
  --p-n-jobs 8

# 4. Phylogenetic Tree & Diversity Metrics (Alpha & Beta Diversity)
qiime phylogeny align-to-tree-mafft-fasttree \
  --i-sequences rep-seqs.qza \
  --o-alignment aligned-rep-seqs.qza \
  --o-masked-alignment masked-aligned-rep-seqs.qza \
  --o-tree unrooted-tree.qza \
  --o-rooted-tree rooted-tree.qza

Assembly • Telomere-to-Telomere Authority: T2T Consortium / Earth BioGenome Project

12. High-C & HiFi Telomere-to-Telomere (T2T) De Novo Genome Assembly

Gapless chromosome-scale phased diploid genome assembly combining ultra-accurate PacBio HiFi reads with chromosome conformation capture (Hi-C) scaffolding using Hifiasm and YaHS.

Flow: PacBio HiFi + Hi-C → Hifiasm (Phased contigs) → YaHS (Hi-C scaffolding) → Inspector / Merqury QV → BUSCO v5

# 1. Phased Telomere-to-Telomere Contig Assembly with Hifiasm
hifiasm -o genome_asm -t 32 \
        --h1 hic_R1.fastq.gz --h2 hic_R2.fastq.gz \
        hifi_reads.fastq.gz

# Extract Primary and Phased Alternate Haplotypes to FASTA
awk '/^S/{print ">"$2"\n"$3}' genome_asm.hic.hap1.p_ctg.gfa > hap1_contigs.fa
awk '/^S/{print ">"$2"\n"$3}' genome_asm.hic.hap2.p_ctg.gfa > hap2_contigs.fa

# 2. Map Hi-C Reads to Contigs for Scaffolding
bwa-mem2 index hap1_contigs.fa
bwa-mem2 mem -5SP -t 16 hap1_contigs.fa hic_R1.fastq.gz hic_R2.fastq.gz | \
    samtools view -bS - | samtools sort -@ 8 -n -o hic_mapped_namesort.bam -

# 3. Chromosome Scaffolding with YaHS
yahs hap1_contigs.fa hic_mapped_namesort.bam -o yahs_scaffolds

# 4. Assembly Completeness Quality Audit with BUSCO v5
busco -i yahs_scaffolds_scaffolds_final.fa -l embryophyta_odb10 -o busco_qc -m genome --cpu 16

Pan-Genomics • Graph Genomes Authority: Human Pangenome Reference Consortium (HPRC)

13. Pan-Genome Graph Construction & Variation Graph Genotyping

Transitioning from reference bias toward non-linear multi-assembly variation graphs using Minigraph-Cactus, indexing with vg, and mapping short reads via vg giraffe.

Flow: Multiple FASTAs → minigraph-cactus → vg autoindex → vg giraffe (bias-free mapping) → vg call

# 1. Pangenome Graph Generation with Minigraph-Cactus
cactus-pangenome ./jobStore ./assemblies.seqfile \
                 --outDir ./pangenome_out \
                 --outName human_pangenome \
                 --reference GRCh38 \
                 --vcf --gbz --gfa

# 2. Graph Indexing with vg
vg autoindex --workflow giraffe -p pangenome_out/human_pangenome.gbz -g pangenome_index

# 3. Align Short Reads to the Variation Graph with vg giraffe (eliminates reference bias)
vg giraffe --gbz-name pangenome_out/human_pangenome.gbz \
           --dist-name pangenome_out/human_pangenome.dist \
           --min-name pangenome_out/human_pangenome.min \
           -f sample_R1.fq.gz -f sample_R2.fq.gz > sample_mapped.gam

# 4. Graph-Based Variant Calling
vg pack -x pangenome_out/human_pangenome.gbz -g sample_mapped.gam -Q 5 -o sample.pack
vg call pangenome_out/human_pangenome.gbz -k sample.pack > sample_graph_variants.vcf

Pathogen Genomics • Virology Authority: ARTIC Network / Nextstrain / CDC

14. Viral Pathogen Surveillance, Consensus & Clade Phylodynamics

Multiplexed tiled amplicon sequencing for respiratory, zoonotic, and emerging viruses, trimming primer binding sites with iVar, generating consensus FASTA sequences, and clade classification via Nextclade and Pangolin.

Flow: Amplicon FASTQ → BWA-MEM → iVar trim (primer removal) → iVar consensus → Nextclade CLI → Pangolin lineage

# 1. Align Reads to Viral Reference Genome
bwa mem -t 8 viral_ref.fasta viral_R1.fq.gz viral_R2.fq.gz | \
    samtools sort -@ 4 -o viral_mapped.bam -

# 2. Trim Tiled Primer Sequences with iVar
ivar trim -b primer_scheme_v4.bed -p viral_trimmed.bam -i viral_mapped.bam -q 15 -m 30 -e
samtools sort -@ 4 viral_trimmed.bam -o viral_trimmed.sorted.bam

# 3. Call SNVs & Generate Majority Consensus Sequence
samtools mpileup -aa -A -d 0 -B -Q 0 viral_trimmed.sorted.bam | \
    ivar consensus -p sample_consensus -q 20 -t 0.6 -m 10 -n N

# 4. Global Lineage & Clade Assignment with Nextclade
nextclade run --input-dataset viral_dataset/ \
              --output-csv nextclade_clade_report.csv \
              sample_consensus.fa

Annotation • Eukaryotic Authority: BRAKER3 / GeneMark-ETP / AUGUSTUS

15. Eukaryotic Structural & Functional Genome Annotation

Fully automated structural gene prediction combining RNA-seq spliced alignments and large-scale cross-species protein homology using GeneMark-ETP and AUGUSTUS, followed by InterProScan domain annotation.

Flow: Masked FASTA + RNA-seq BAM + Protein OrthoDB → BRAKER3 → TSEBRA → InterProScan → eggNOG-mapper

# 1. Softmask Repeats in Genome Assembly with RepeatMasker
RepeatModeler -database genome_db -threads 16
RepeatMasker -lib genome_db-families.fa -xsmall -pa 16 assembly.fa

# 2. Structural Gene Prediction with BRAKER3
braker.pl --genome=assembly.fa.masked \
          --bam=rnaseq_spliced_aligned.bam \
          --prot_seq=orthodb11_clade.fasta \
          --threads=16 --gff3 --useexisting

# 3. Functional Domain Annotation with InterProScan
interproscan.sh -i braker_output/braker.aa -f tsv,gff3 \
                -appl Pfam,PANTHER,Gene3D,SUPERFAMILY -goterms -pa -cpu 16

Structural AI • AlphaFold & Docking Authority: DeepMind / ColabFold / AutoDock Vina

16. AI-Driven 3D Biomolecular Structure Prediction & In Silico Docking

End-to-end atomic structure inference from sequence FASTA using ColabFold/AlphaFold, computing confidence metrics (pLDDT, PAE), and performing virtual screening for small-molecule chemical ligands via AutoDock Vina.

Flow: Protein FASTA → MMseqs2 MSA → ColabFold batch (GPU) → pLDDT filter (>80) → AutoDock Vina Grid → Binding Affinity Scoring

# 1. Fast Batch Structure Prediction with ColabFold (MMseqs2 + AlphaFold2)
colabfold_batch target_proteins.fasta colabfold_results/ \
                --amber --templates --num-recycle 3 --use-gpu-relax

# 2. Extract Structure with Highest pLDDT
python3 -c "
import json
with open('colabfold_results/predicted_scores.json') as f:
    scores = json.load(f)
print('Top Model Mean pLDDT:', max(scores['plddt']))
"

# 3. Molecular Docking with AutoDock Vina
# Prepare Receptor and Ligand (PDBQT format)
vina --receptor top_alphafold_model.pdbqt \
     --ligand drug_candidate.pdbqt \
     --center_x 15.2 --center_y 24.8 --center_z -8.4 \
     --size_x 20 --size_y 20 --size_z 20 \
     --out docked_conformations.pdbqt \
     --cpu 8 --exhaustiveness 32

Standard Genomics & Bioinformatics Analysis Pipelines