Lake Trout Genomics: Comparative Analysis of Lean and Siscowet Ecotypes

🌐 Project website: robertslab.github.io/project-lake-trout — the annotation & phenotype landing page, top candidate genes, GO enrichment, interpretation guardrails, and links to both interactive genome browsers.

Overview

Comprehensive genomic analysis of Salvelinus namaycush (lake trout) comparing two distinct ecotypes: lean and siscowet. This repository contains multiple integrated analyses including:

• RNAseq differential expression analysis using parasitized/non-parasitized liver tissue
• PacBio HiFi DNA methylation profiling and differential methylation analysis
• Presence-Absence Variation (PAV) analysis to identify structural genomic variations
• Functional annotation & phenotype interpretation linking DMRs and PAVs to candidate genes
• Ecotype genome assemblies (PacBio HiFi, hifiasm) for lean and siscowet, with Liftoff annotation and synteny
• Two-genome comparison: structural-variant map, reciprocal PAV, gene-content sets, and a refined evidence model that retires the single-reference caveat
• Interactive genome browser for visualizing genomic features

Reference Genome

• Assembly: GCF_016432855.1 (SaNama_1.0)
• Species: Salvelinus namaycush (Lake Trout)
• Source: NCBI Genome

---

Key Analyses

1. RNAseq Differential Expression

Analysis of liver RNAseq data from parasitized and non-parasitized samples (NCBI BioProject PRJNA316738) to identify:

• Differentially expressed genes (DEGs) between subspecies
• Differentially expressed transcripts (DETs) and alternative isoforms
• Expression differences related to parasite status

Key Results:

• 202 differentially expressed transcripts (p < 0.05)
• Analysis performed using Ballgown
• See analyses/README.md for detailed results

Analysis Files:

• code/01-ballgown-analysis.Rmd - Primary differential expression analysis
• code/02-gene-explore.qmd - Gene-specific exploration

2. PacBio HiFi DNA Methylation Analysis

Whole-genome DNA methylation profiling using PacBio HiFi sequencing with 5mC modification calling:

• Sample-level methylation profiles for both ecotypes
• Differential methylation analysis between lean and siscowet
• Identification of Differentially Methylated Regions (DMRs)

Key Results:

• 540,040 CpG sites tested
• 4,440 significant differentially methylated cytosines (DMCs, p < 0.05)
• 302 Differentially Methylated Regions (DMRs)
  • 20 hypermethylated in siscowet
  • 282 hypomethylated in siscowet

Analysis Files:

• code/04-pacbio/ - PacBio workflow (alignment, QC, methylation calling)
• code/10-mCG-call.Rmd - Methylation calling
• code/13.3-hifiasm-differential-methylation-plan.md - Plan for extending differential methylation analysis to ecotype-specific hifiasm assemblies
• code/14-diff-meth.Rmd - Differential methylation analysis
• code/14-diff-meth.py - Python implementation for DMR identification

3. Presence-Absence Variation (PAV)

Genome-wide structural variation analysis identifying insertions and deletions between ecotypes:

• Coverage-based detection of absent regions (deletions)
• CIGAR-based detection of novel insertions
• Ecotype-specific and shared structural variants

Key Results:

• Lean-specific: 996,228 variants (770,891 insertions + 225,337 deletions)
• Siscowet-specific: 1,332,705 variants (1,086,799 insertions + 245,906 deletions)
• Shared: 878,372 variants common to both ecotypes

Analysis Files:

• code/11-pav.Rmd - PAV identification analysis
• code/12-pav.py - Python implementation
• code/15-diff-pav.py - Differential PAV analysis

4. Interactive Genome Browsers

Web-based genome browsers for exploring PAV and methylation data across the SaNama_1.0 assembly.

Features:

• Visualize ecotype-specific insertions and deletions
• View differential methylation tracks (CpG methylation, 8 samples; DMRs)
• Gene annotations with interactive navigation
• Mobile-responsive design

Live Demos:

• IGV.js (quick exploration): robertslab.github.io/project-lake-trout/genome-browser/
• JBrowse 2 (advanced analysis): robertslab.github.io/project-lake-trout/jbrowse/

Documentation: genome-browser/README.md

5. Functional Annotation & Phenotype Interpretation

The interpretive layer that turns DMR and PAV coordinates into genes and plausible ecotype phenotypes. A genome-wide RefSeq annotation backbone is built, differential features are assigned to genes with positional context, candidates are ranked, GO over-representation is tested, and the results are synthesized into hypothesized phenotype axes.

Key Results:

• 46,359 genes annotated (46,231 with a product, 34,367 with ≥1 GO term)
• 2,036 candidate genes within 5 kb of a DMR/DMC/stringent-PAV; 4 convergent (DMR and stringent siscowet deletion), led by znf883-like (LOC120032414)
• Exonic siscowet deletions in lipid-metabolism genes (angptl5, mogat2, epoxide hydrolase 1)
• Most defensible GO enrichment (deletion set): calcium ion transport (FDR 3×10⁻³), with gene-length and lean-reference caveats carried throughout
• 0 DMCs survive q < 0.1 — interpretation leads with DMR-level and stringent-PAV sets

Caveat (2025, now partly retired): the reference is a lean-background doubled-haploid genome, so siscowet-specific deletions were divergence-inflated. Section 6 re-tests these results on each ecotype's own genome: the calcium signal survives, the lipid-axis deletions and 3 of the 4 convergent genes do not. All links remain associations, not validated mechanisms.

Analysis Files:

• code/18-diff-annotation-phenotype-plan.md - Plan of work
• code/18-build-gene-function-table.py - Annotation backbone
• code/18.1-assign-features-to-genes.py - DMR/PAV → gene assignment
• code/18.2-integrate-candidates.py - Ranked candidate integration
• code/18.3-go-enrichment.py - GO over-representation
• code/18-diff-annotation-phenotype.Rmd - Phenotype synthesis report
• analyses/18-annotation/README.md - Outputs & provenance

---

6. Ecotype Genome Assemblies & Two-Genome Comparison

Each ecotype was assembled from its own PacBio HiFi reads (hifiasm, purge_dups), annotated by Liftoff from the reference, and compared to the reference and to each other. This layer replaces single-reference association with bidirectional evidence.

Assemblies: lean 2.58 Gb / 12,700 contigs / 55,437 lifted genes; siscowet 2.84 Gb / 8,159 contigs / 59,968 lifted genes; 50,379 genes shared.

Key Results:

• SyRI structural-variant map for lean↔ref, siscowet↔ref and lean↔siscowet; ~18–20 % of reference-based PAVs fall inside SV-confirmed regions
• Reciprocal PAV (each ecotype's reads → the other's genome): 340 siscowet-present and 343 lean-present regions (1.57 Mb), 115 genes, free of lean-reference inflation
• Gene content: 1,474 lean-only and 3,040 siscowet-only genes corroborated by both assemblies; 2,359 copy-number-divergent genes (dominated by tRNA / snRNA / histone arrays)
• Refined evidence model (5 independent lines, tiers A–C): 77 tier-A genes with read-level native support (66 flag-free), 4,011 tier-B. The calcium-transport GO signal survives after tandem-cluster collapsing (voltage-gated calcium channel activity, FDR 0.03); no lipid term does. Three of the four 2025 "convergent" candidates and the angptl5 / mogat2 lipid deletions were reference-distance artifacts.

Analysis Files:

• code/13.1-hifiasm-genome-assembly-*.Rmd - hifiasm assemblies
• code/20-ecotype-genomes.md, code/20.1-liftoff-annotation.Rmd - annotation plan and Liftoff
• code/21.1-gene-anchored-synteny.Rmd - purge_dups + MCScanX synteny
• code/22-synteny-browser-integration.qmd - synteny tracks and JBrowse synteny view
• code/23-next-phase-research-plan.md - two-genome plan and progress log
• code/23.1-genome-sv-map.Rmd, code/23.2-sv-pav-crossvalidate.py, code/23.3-sisco-sv-repair.py - SV map
• code/23.4-reciprocal-pav.py, code/23.5-gene-set-ops.py, code/23.6-integrate-phase2.py - reciprocal PAV, gene sets, integration
• code/23.7-refine-evidence-model.py - refined evidence lines, tiers, cluster-collapsed GO
• analyses/23-integration/refined/README.md - refined results and tables

Open phases: native differential methylation (Phase 3), BRAKER3 de novo annotation (Phase 4), browser tracks and manuscript (Phase 5).

---

Repository Structure

project-lake-trout/
├── index.html / index.qmd  # Project website landing page (GitHub Pages)
├── code/                    # Analysis scripts and notebooks
│   ├── 01-ballgown-analysis.Rmd       # RNAseq differential expression
│   ├── 02-gene-explore.qmd            # Gene exploration
│   ├── 04-pacbio/                     # PacBio HiFi analysis workflow
│   ├── 05-pacbio-align.Rmd            # PacBio alignment
│   ├── 07-pacbio-QC.Rmd               # PacBio quality control
│   ├── 10-mCG-call.Rmd                # Methylation calling
│   ├── 11-pav.Rmd                     # PAV analysis
│   ├── 13.3-hifiasm-differential-methylation-plan.md # Plan for hifiasm-based differential methylation
│   ├── 14-diff-meth.Rmd/py            # Differential methylation
│   ├── 15-diff-pav.py                 # Differential PAV
│   ├── 18-*                           # Annotation, candidate integration, GO, phenotype
│   ├── 13.1-*, 20-*, 21.1-*, 22-*     # Ecotype assemblies, Liftoff, synteny, browser integration
│   └── 23-*                           # Two-genome SV map, reciprocal PAV, gene sets, evidence model
├── data/                    # Raw data and metadata
│   ├── SraRunTable.csv                # RNAseq sample information
│   ├── ballgown-metadata.csv          # Ballgown metadata
│   └── *.bed                          # Gene annotations
├── output/                  # Assembly, Liftoff, synteny outputs (large files on gannet)
├── analyses/                # Analysis outputs and results
│   ├── DEG-*.csv                      # Differentially expressed genes
│   ├── DET-*.csv                      # Differentially expressed transcripts
│   ├── 04-pacbio/                     # PacBio analysis outputs
│   ├── 14-diff-meth/                  # Methylation results
│   ├── 15-diff-pav/                   # PAV results
│   ├── 18-annotation/                 # Functional annotation & candidates
│   ├── 23-reciprocal-pav/, 23-gene-sets/ # Two-genome results
│   └── 23-integration/refined/        # Five-line evidence model, tiers, collapsed GO
├── genome-browser/          # Interactive IGV.js genome browser
│   ├── index.html                     # Browser interface
│   ├── prepare_data.py                # Data preparation script
│   └── data/                          # Browser data files
├── jbrowse/                 # JBrowse 2 genome browser
└── figures/                 # Generated figures and plots

See README files in each subdirectory for detailed information about specific analyses.

---

Sample Information

RNAseq Samples

• Lean Nonparasitized: NPLL32, NPLL34, NPLL44, NPLL46, NPLL56, NPLL61
• Lean Parasitized: PLL20, PLL31, PLL43, PLL55, PLL59, PLL62
• Siscowet Nonparasitized: NPSL15, NPSL24, NPSL29, NPSL36, NPSL50, NPSL58
• Siscowet Parasitized: PSL13, PSL16, PSL35, PSL49, PSL53, PSL63

See data/SraRunTable.csv for complete RNAseq sample metadata.

PacBio HiFi Samples (for methylation and PAV analysis)

• Lean: bc2041, bc2068, bc2069, bc2070
• Siscowet: bc2071, bc2072, bc2073, bc2096

---

Pre-Analysis Data Processing

The following data processing steps were performed prior to analyses in this repository:

• SRA RNAseq data retrieval, trimming, and QC

• RNAseq alignment and splice site identification using Hisat2 and StringTie

• Convert NCBI genome GFF to BED

---

Technologies Used

• R/RStudio: Statistical analysis and visualization
  • Ballgown: Differential expression analysis
  • tidyverse: Data manipulation
• Python: Data processing and analysis pipelines
  • pysam, pandas, numpy: Data manipulation
  • modbampy: Modified base parsing
• PacBio Tools: HiFi sequencing analysis
  • pbmm2: Read alignment
  • pb-CpG-tools: Methylation calling
• IGV.js: Interactive genome visualization
• Quarto/RMarkdown: Reproducible analysis notebooks

---

Citation

If you use data or methods from this repository, please cite:

• Project website: robertslab.github.io/project-lake-trout
• Lake Trout RNAseq data: NCBI BioProject PRJNA316738
• Reference genome: GCF_016432855.1 (SaNama_1.0)

---

Contact

Roberts Lab
School of Aquatic and Fishery Sciences
University of Washington

For questions or issues, please open a GitHub issue in this repository.