Phylogenomics and demographic history of O. formosanus.

(a) Geographic distribution and phenotypic variation of selected lineages within the Oncorhynchus masou species complex. (b) Species tree inferred from a concatenated alignment of 1,483 single-copy orthogroups. Blue horizontal bars representing the 95% confidence intervals. Numbers in red and blue denote the gained and lost orthogroups at each node, respectively. (c) Principal component analysis (PCA) showing genetic differentiation among masu salmon lineages. Numbers along dashed arrows represent mean pairwise homozygous SNP differences (in millions) between populations. (d) Time-calibrated diversification within Formosan (O. formosanus) and Japanese (Nagano Masu and Amago) salmon lineages visualized by DensiTree, indicating lineage splits occurring approximately 50,000–20,000 years ago. (e) Demographic histories inferred by PSMC, showing changes in effective population size for Japanese and Taiwanese salmon individuals.

Chromosome synteny and sex chromosome evolution in O. formosamous

(a) Genome-wide synteny comparison across Oncorhynchus species. Chromosomes highlighted in blue denote fusions unique to O. formosanus, and those in red indicate fusions in the common ancestor of masu salmon lineages. Colored lines represent pairwise single-copy orthologs. The scaffold containing the sex-determining gene (sdY) in O. m. masou follows data from Christensen et al. 28. (a) Normalized coverage of female (pink) and male (blue) sequencing reads mapped onto selected scaffolds from the male O. formosanus genome assembly. (c) Detailed annotation of repetitive elements and sequence alignments between putative X-linked (ptg00003l) and Y-specific scaffolds (ptg000741l and ptg002124l) from the male genome, aligned to chromosome X (scaffold13) from the female assembly. Numbers above shaded grey area indicate nucleotide identity between X chromosome sequences of male and female assemblies. Lines connecting Y-specific scaffolds (ptg000741l and ptg002124l) and scaffold13 represent homologous regions, with colors indicating amino acid sequence similarity. The exact location of the sdY gene on scaffold ptg000741l is indicated.

Population structure, genetic diversity and level of inbreeding of the Formosan salmon.

(a) PCA result of the population structure shows no significant differentiation among most populations, except for those from Luoyewei Creek and the downstream sites of Hehuan Creek (HEH_01, HEH_02, HEH_W1). Populations from Hehuan Creek occupy a broader range in PCA space (PC1 vs. PC2), indicating a higher genetic diversity. (b) Multilocus heterozygosity across individuals from the three creek systems. Hehuan Creek populations exhibit the highest genetic diversity, followed by Qijiawan Creek and lastly by Luoyewei Creek. (c) The level of inbreeding, estimated using runs of homozygosity (ROHs), is highest in Luoyewei Creek, intermediate in Qijiawan Creek, and lowest in Hehuan Creek.

Current effective population size Ne and demographic history of Formosan salmon populations.

(a) Estimated current Ne for Formosan salmon populations across multiple stream sites. The Qijiawan and Hehuan populations appear to maintain effective sizes above thresholds associated with inbreeding depression, while the Luoyewei Creek population, which likely comprises individuals from hatchery releases, has a very small effective population size. (b) Inferred demographic histories of populations from Qijiawan and Hehuan Creeks. The contrasting trajectories reflect differences in conservation management and habitat development across the two regions.

Population trajectories and extinction risk under varying typhoon frequencies.

(a-d) Simulated population dynamics over 20 years under different annual typhoon occurrence probabilities (p = 0.3, 0.5, 0.7, 0.9) over a 20-year simulation period. Panels a–d represent different disturbance scenarios, with each line indicating mean abundance (individuals aged 1-4, per 1,000 m²) for Qijiawan (QIJ; red) and Hehuan (HEH; blue) stream populations. Shaded bands represent 95% confidence intervals across 1,000 stochastic replicates per scenario. Simulations incorporated stage-specific survival rates, Leslie matrix projections, and probabilistic typhoon impact modeled using zero-inflated Gamma distributions derived from historical records.

Geographic distribution and sampling locations of Oncorhynchus formosanus populations in Taiwan’s Central Mountain Range.

The map shows the three creek systems harboring extant O. formosanus populations at elevations between 1,500-3,500 m. Qijiawan Creek (orange line) in the northern portion contains the historically protected population with sampling sites QIJ2, QIJ3, and QIJ4. Hehuan Creek (blue lines) in the southern portion includes both upstream sites (HEH01, HEHS1, HEHS2; dark blue circles) and downstream sites (HEH02, HEHW1; light blue circles) where reintroduction occurred in 2017. Luoyewei Creek (yellow line) in the eastern portion represents the hatchery-supplemented population (LUO). The grayscale gradient indicates elevation, with darker shades representing lower elevations. All populations are located within the upper tributaries of the Dajia River basin in central Taiwan.

Genome-wide Hi-C contact map of the O. formosanus genome assembly.

The Hi-C interaction heatmap illustrates the genome-wide three-dimensional organization and chromosomal-level scaffolding of the O. formosanus genome. Each blue box represents an individual chromosome-scale scaffold identified in this study. Stronger interaction signals are indicated in red, representing frequent physical interactions within chromosomes, while weaker or infrequent interactions appear in lighter colors.

Coalescent-based species tree inferred by ASTRAL for Oncorhynchus species.

Species phylogeny constructed using the coalescent-based approach implemented in ASTRAL from 1,483 single-copy orthologs. Branch support values represent local posterior probabilities, reflecting concordance among gene trees. This analysis provides the same topology and corroborates relationships among Oncorhynchus species compared to the species tree from concatenated alignments of Fig. 1b.

Pairwise genetic differentiation among masu salmon lineages.

Heatmap illustrating pairwise genetic differentiation (Fst, green shades) and absolute nucleotide divergence (Dxy, red shades) among five Oncorhynchus masou lineages.

Genealogical Divergence Index (GDI) analysis of Taiwanese and Japanese salmon lineages.

The genome-wide estimates of genealogical divergence index (GDI) reveals substantial genetic differentiation between Taiwanese (O. formosanus) and Japanese (O. masou) salmon lineages (the value was around 1.0; dark grey). The estimated GDI values between the two Japanese salmon lineages were highlighted in green and blue colors, where the values were around 0.7. High GDI values (>0.7) indicate strong genealogical divergence, implying long-term isolation and lineage differentiation. In contrast, the estimated GDI between the two Taiwanese population (highlighted in light grey) was smaller than 0.2, which suggests neglectable genetic divergence between the two populations.

Population phylogeny of Taiwanese and Japanese salmon lineages inferred from nuclear genome-wide and mitochondrial SNPs.

(a) Maximum likelihood phylogeny based on genome-wide SNP data from representative individuals of Taiwanese (O. formosanus: QIJ, HEH) and Japanese (Amago, Biwa, Nagono Masu and Hokkaido Masu) salmon populations, rooted with O. kisutch (coho salmon) as the outgroup. Branch lengths represent nucleotide substitutions per site. (b) Corresponding mitochondrial genome phylogeny inferred from mitogenome SNP data, also rooted with O. kisutch. Bootstrap support values (percent) for major branches are indicated. Together, these phylogenies illustrate consistent genetic differentiation patterns between nuclear and mitochondrial genomes across salmon lineages.

Relative proportions of loss-of-function to synonymous variants in homozygous and heterozygous states across salmon populations.

Scatter plot showing the ratio of loss-of-function (LOF) variants to synonymous variants in both homozygous (green) and heterozygous (orange) states among salmon populations: Amago, Biwa, QIJ (Qijiawan), HEH (Hehuan), Hokkaido Masu, and Nagano Masu. Horizontal bars indicate mean values for each population. Heterozygous LOF variants consistently exhibit higher ratios compared to homozygous LOF variants across all populations, aligning with expectations of increased accumulation of deleterious alleles and elevated homozygosity in populations with smaller effective population sizes. Asterisks denote statistical significance of the difference between homozygous and heterozygous states (*P<0.05, **P<0.01, ****P<0.0001).

Comparison of evolutionary rates (dN/dS ratios) across different chromosome categories in O. formosanus.

Boxplot showing the distribution of nonsynonymous-to-synonymous substitution rate ratios (dN/dS) for genes located on chromosomes classified as non-fused, fused in the last common ancestor of the O. masou species complex, or specifically fused in O. formosanus. "ns" denotes no significant difference, ** denote P < 0.01 in Wilcoxon rank sum test.

Alignment of the female-derived O. m. masou sdY-containing contig to chromosome 13.

Dot plot visualization illustrating the unique alignment of the female-derived O. m. masou contig containing the sex-determining gene (sdY) to the central genomic region (∼50 Mb) of chromosome 13. Alignment identity is indicated by a color scale, with higher identity in red. This confirms the precise chromosomal location of the sex-linked region at an ancestral chromosome fusion site.

Comparative alignment showing ancestral chromosome fusion between chromosome 13 of O. m. masou and chromosomes 14 and 25 of O. mykiss.

Visualization illustrating the genomic correspondence and alignment between chromosome 13 of O. m. masou and chromosomes 14 and 25 in rainbow trout (O. mykiss). The ancestral chromosome fusion site includes enriched regions of repetitive sequences, notably 281-mer satellite repeats, indicative of centromeric or neocentromeric activity and consistent with chromosome fusion events. Sequence similarity is depicted by color intensity.

Protein family (Pfam) domain expansions across salmonid genomes.

Heatmaps illustrate significant expansions in protein domains (Pfam) across selected salmonid genomes, highlighting enriched domains specifically abundant in (a) O. formosanus and (b) O. m. masou. Numbers in each cell indicate the total count of the corresponding Pfam domain within each salmonid genome. Colors represent the relative expansion (Z-score) of each domain across species, where positive Z-scores (green to blue shades) indicate higher-than-average domain expansions, and negative Z-scores (pink to red shades) represent lower-than-average copy numbers compared to other salmonids.

Phylogenetic relationships and divergence of expanded ice-binding lectin domain proteins (LbR_Ice_bind) in O. formosanus.

Maximum-likelihood phylogeny illustrating the evolutionary divergence of ice-binding lectin proteins (LbR_Ice_bind domain) expanded in O. formosanus, compared to characterized ice-structuring protein (ISP) types I, II, and IV (Bar Dolev et al. 2016a) identified in other Oncorhynchus species. The distinct cluster formed by O. formosanus sequences highlights potential lineage-specific functional specialization associated with adaptation to cold freshwater environments.

ddRAD-seq mapping results across the genome of O. formosanus.

SNPs were identified across the newly assembled O. formosanus genome using ddRAD-seq data with the window size of 1Mb. The mapping results show even genome-wide coverage, indicating the ddRAD-seq data effectively mapped variation across most regions of the reference genome.

Quasi-extinction probabilities and odds ratios under increasing typhoon disturbance.

(a) Quasi-extinction probability of Qijiawan (QIJ, light gray) and Hehuan (HEH, dark gray) across four typhoon frequencies. Bars represent means of 1,000 replicates; whiskers indicate exact 95% binomial confidence intervals. (b) Extinction odds ratios (QIJ / HEH) across typhoon scenarios. Gray points and error bars show Fisher’s exact test results with exact 95% confidence intervals; asterisks indicate significance levels (p < 0.05 *, p < 0.01 **, p < 0.001 ***, ns = not significant). The black line represents logistic GLM estimates of OR across scenarios, showing a monotonic decline with increasing typhoon frequency. The OR at p = 0.9 is omitted due to complete extinction in both populations.

Electrofishing vs. snorkel sampling frequency and detection ratio by length bin.

Stacked bars show the proportional frequency of fish in each corrected length bin (mm) as detected by electrofishing (purple) and snorkel surveys (orange) during the same sampling period in 2024 at Hehuan Creek. The green line with markers indicates the detection ratio (electrofishing frequency ÷ snorkel frequency) for each length bin, plotted against the secondary y-axis. Bins with a detection ratio >1 suggest underestimation by snorkel sampling, while ratios <1 suggest potential overestimation.

Total length distributions of electrofishing, raw snorkel, and corrected snorkel samples in Hehuan Creek.

Histograms show fish total length distributions (mm) obtained via electrofishing (blue), unadjusted snorkel observations (orange), and corrected snorkel estimates (green) during the 2024 sampling period in Hehuan Creek. Corrections for snorkel-derived lengths were applied using refractive index adjustment and detection ratio-based resampling to account for observational biases across length classes.

Monthly body-length distributions of Formosan land-locked salmon in two creeks (QIJ and HEH).

Relative-frequency histograms (grey bars) show the proportional length composition (mm) of snorkel-derived counts—after detection-ratio and refraction corrections—at seven study sites (Dam4, Dam3, Dam2, Salmon, S1, C1 and W1). Columns are ordered left-to-right by site within each tributary; the tributary names (QIJ, HEH) appear once above their respective site block. Rows correspond to individual monthly surveys (Jan, Jul, Oct, Nov, Dec 2024 and Jun 2025). Bars within each panel sum to 1, enabling direct comparison of length-class structure between site-months regardless of sample size. Blank panels indicate months or sites with no observations.

Monthly age-class composition of Formosan land-locked salmon in the QIJ and HEH tributaries of Hehuan Creek.

Stacked bars (shaded grey from light to dark for age-classes 0, 1, 2 and 3, respectively) show the posterior mean relative frequency of each age-class at seven study sites (Dam4, Dam3, Dam2, Salmon, S1, C1 and W1). Columns are ordered left-to-right by site within each tributary block (stream labels appear once above the block), while rows correspond to individual monthly surveys (January, July, October, November and December 2024, and June 2025). T-shaped lines denote the 95 % Bayesian credible interval around each bar. Within every panel the four bars sum to 1, allowing direct comparison of age structure across site-month combinations irrespective of sample size. Blank panels indicate months or sites with no observations.