For this experiment, we selected a female Pengze crucian carp, a male Pengze crucian carp, and a male sea carp (Cyprinus acutidorsalis), all of which were artificially selected from Guangzhou Jianbo Fish Farm Co., Ltd., as our breeding stock. Artificial insemination was performed on 13 March 2023. Artificial insemination was conducted on 13 March 2023. The female’s eggs were evenly divided into two portions, each fertilized separately with the sperm from the male Pengze crucian carp and the male sea carp. The fertilized eggs from the homologous and heterologous groups were evenly distributed onto hatching net plates and placed in two concrete hatching pools, each measuring 6 m × 4 m, within the same room. These pools were equipped with a circulating water system, and the water temperature was maintained between 24.5°C and 25.8°C with continuous aeration for 24 h. After hatching, the fry were initially fed boiled egg yolk for 3 days. Once the fry reached a length of 3-4 cm, we randomly selected 100 individuals from each group and transferred them to two separate circular tanks, each with a diameter of 3 m and a height of 1.5 m, for breeding. The daily diet for these fish consisted of feed containing 38% protein.

Parental muscle tissue samples from two Pengze crucian carps and one male sea carp were collected on 13 March 2023. The samples were immediately placed into microcentrifuge tubes pre-cooled with liquid nitrogen and then stored at −80°C for preservation. The parental samples were labeled as follows: PF (female Pengze crucian carp), PM (male Pengze crucian carp), and CM (male sea carp). On 18 September 2023, nine homologous offspring and nine heterologous offspring were randomly selected for the same sampling procedure. The samples from the homologous offspring were designated as PP1, PP2, PP3, PP4, PP5, PP6, PP7, PP8, and PP9, while those from the heterologous offspring were labeled as PC1, PC2, PC3, PC4, PC5, PC6, PC7, PC8, and PC9.

On 19 September 2023, 21 frozen samples were sent to Wuhan Grandomics Biosciences Co., Ltd. (https://www.grandomics.com/) for third-generation Nanopore resequencing. Samples were collected, and high molecular weight genomic DNA was prepared by the CTAB method and followed by purification with QIAGEN^® Genomic kit (Cat#13343, QIAGEN) for regular sequencing, according to the standard operating procedure provided by the manufacturer. After the sample was qualified, size-select of long DNA fragments were performed using the PippinHT system (Sage Science, United States). Next, the ends of DNA fragments were repaired, and A-ligation reaction were conducted with NEBNext Ultra II End Repair/dA-tailing Kit (Cat# E7546). The adapter in the SQK-LSK109 (Oxford Nanopore Technologies, UK) was used for further ligation reaction and DNA liabrary was measured by Qubit^® 4.0 Fluorometer (Invitrogen, United States). About 700 ng DNA library was constructed and performed on a Nanopore PromethION sequencer instrument (Oxford Nanopore Technologies, UK) at the Genome Center of Grandomics (Wuhan, China).

Raw data from Nanopore sequencing were base-called and converted into the fastq format using Drado v0.3.4 (https://github.com/nanoporetech/dorado) for subsequent analysis. The fastq data of the 21 samples underwent adapter and low-quality sequence removal using Porechop v0.2.4 (Bonenfant et al., 2023). Subsequently, with Nanofilt v2.8.0 (De Coster et al., 2018), we trimmed 90 bp from the beginning and 80 bp from the end of each read, filtering out reads shorter than 500 bp or with a quality score below 7. The cleaned data from each sample were aligned to the published triploid Carassius gibelio genome (NCBI RefSeq assembly: GCA_023724105.1) using Minimap2 v2.2.6 (Li, 2018) with default parameters. Alignment information for each sample was assessed using samtools v1.19 (Danecek et al., 2021) and mosdepth v0.3.4 (Pedersen and Quinlan, 2018).

To ensure ensure subsequent differential analysis between females and males as well as the offspring genomes, a specific strategy for SV detection was employed (Supplementary Figure S1):

1. SVs in the PF data were detected using default parameters of Sniffles2 v2.2 (Smolka et al., 2024), cuteSV v2.0.3 (Jiang et al., 2022), and SVIM v2.0.0 (Heller and Vingron, 2019). The results from these three tools were compiled to create a composite dataset of PF’s SVs.

Only DEL, INS, INV, and DUP were analyzed.

To detect DSVs between PF and offsprings, as well as PF and males, the following approach was adopted. Firstly, SVs from males and offspring, filtered earlier, were merged with SV results of PF using SURVIVOR v1.0.6 (Jeffares et al., 2017). The parameters set for this were: “SURVIVOR merge sample.files 1,000 2 1 1 0 50 output.vcf.” If males and offspring SVs could be merged with PF’s SVs, these were eliminated using a Python script. This process was repeated until no consistent SVs with PF were detectable. The remaining DSVs were then screened to exclude those with read coverage in PF less than 4 in the DSV region. This resulted in the identification of genuine and significant DSVs between PF and offspring, and between PF and males.

For accurate and effective detection of MSSVs in offspring, the approach first involved excluding DSVs larger than 10 Kb. Then, SURVIVOR was used to detect MSSVs by merging homologous offsprings with the the homologous male and heterologous offsprings with the heterologous male, using the parameters: “SURVIVOR merge sample.files 100 2 1 1 0 50 output.vcf.” Subsequently, a separate analysis was conducted on both sets of nine homologous and heterologous offspring to identify common MSSVs in each group, with the parameters: “SURVIVOR merge sample.files 100 9 1 1 0 50 output.vcf.”

Gene annotation for MSSVs detected in each offspring sample was carried out using ANNOVAR v2020.06.07 (Wang et al., 2010). This process was aimed at understanding the potential impacts of these variations. A local database was established using the annotation database corresponding to the crucian carp reference genome, and the MSSVs in offspring were annotated using the default parameters.

To validate MSSVs in offsprings, a focus was placed on a specific MSSV within the exonic regions of protein-coding genes that could cause frameshift mutations. Using samtools v1.19, reads from the bam alignment file that encompass the specified region were extracted and saved as a new bam file. This was followed by visualization using Geneious Prime v2022.2.2 (Deng et al., 2020), for further analysis. Primers for PCR validation of the MSSV were designed using Primer3 v2.3.7 (Untergasser et al., 2012). The DNA template with primers underwent a series of thermal cycles: initial denaturation at 95°C for 5 min, followed by 30 cycles of amplification(denaturation at 94°C for 30 s, annealing at 56°C for 45 s, and extension at 72°C for 30 s), and a final extension phase at 72°C for 10 min.

In this investigation, we employed third-generation Nanopore sequencing technology to resequence the whole genomes of three parental and eighteen offspring samples, obtaining approximately 348.69 Gb of raw data, which was reduced to about 338.80 Gb after data cleaning. The filtered dataset for all samples was roughly 15 Gb, with an N50 length varying from 17 Kb to 32 Kb. The average read quality surpassed a score of 11, achieving a Q7 score of 100% (Supplementary Table S1). Alignment with the reference genome demonstrated a match rate exceeding 99%, with an average depth of approximately 9x (Supplementary Table S2). To maximize the replication of the female SVs and the accuracy of males and offspring SVs, we utilized multiple algorithms for PF and only Sniffles2 for males and offspring (Cuenca-Guardiola et al., 2023). Detection of PF’s SVs through tools such as Sniffles2, cuteSV, and SVIM resulted in 474,628, 1,012,709, and 87,059 SVs, respectively. These data were aggregated to represent SVs of PF (Table 1). Detection of SVs in males and offspring retained between 48.2% and 64.9% of the SVs after filtered (Supplementary Figure S2; Supplementary Table S3). CM showed a total of 492,833 SVs spanning 227.11 Mb, while PM and offspring ranged between 217,203 and 277,826 SVs spanning 246.25 Mb to 367.57 Mb (Figures 2, 3). DEL and INS were predominant in number, while INV and DUP, though fewer, were significant in total length, representing larger segment variations. The larger genomic sequence differences between the sea carp and the Pengze crucian carp are determined through the calculation of DSVs (Supplementary Table S4). CM with higher SV numbers but lower total length, indicating smaller fragment SVs. Subsequent DSVs detection showed CM had 44–67 times more SVs than PM and offspring (Figure 4; Supplementary Figure S3), reflecting genomic sequence evolution driven by its unique growth environment. PM and offspring exhibited consistency in the number and type of SVs, indicating genetic stability and limited changes in environmental adaptability.

To detect DSVs within the genomes of male and offspring samples, we employed a combination of SURVIVOR and custom Python scripts. This involved iterative processes of merging the SVs detected in the PF with those of the offspring, followed by rigorous filtering to ensure that identified DSVs were unique and not present in PF. In the offspring, DSVs were detected in a range of 4,783–7,665, constituting approximately 2.19%–2.66% of their filtered SVs (Figure 4). The observed variation indicates the emergence of novel SVs, which may result from a variety of mechanisms. These include sporadic homologous recombination and high-frequency gene conversion during the maturation of oocytes (Wang et al., 2022), the suppression of DNA double-strand break formation and homologous recombination in oocytes (Lu et al., 2022), or the incorporation of male genetic fragments during the gynogenetic process (Zhang et al., 1999; Zhao et al., 2004a; Mao et al., 2020; Ding et al., 2021; Zhao et al., 2021). The process of recombination, depending on its accuracy, plays a dual role: inaccurate recombination can lead to the introduction of new genetic variations, whereas precise recombination has the potential to correct existing genetic errors (Burssed et al., 2022). Consequently, it is hypothesized that the DSVs observed in offspring may originate from both standard and aberrant recombination repair pathways, as well as directly from male genetic contributions.

In further analysis, we distinguished MSSVs in offspring by merging DSVs from homologous pairs with those from the PM and DSVs from heterologous pairs with the CM. Findings revealed that 21.4% to 23.4% of the DSVs in homologous offspring were MSSVs, with counts between 1,195 and 1709, while heterologous offspring had 4.3% to 5.0% MSSVs, with counts between 236 and 350 (Figure 5). This finding suggests that a significant proportion of the DSVs in the offspring can be traced back to male genetic contributions. Remarkably, even though the CM exhibited a DSV count that was 52 times higher than that of the PM (Supplementary Figure S3), the offspring derived from homologous sperm stimulation displayed a substantially higher prevalence of MSSVs—ranging from 3 to 7 times greater—compared to those conceived through heterologous sperm stimulation.

To identify common MSSVs across groups, we conducted a merged analysis for the nine homologous and nine heterologous offspring separately. In the gynogenetic triploid crucian carp, which possesses 150 basic chromosomes forming two sets of haplotypes (Wang et al., 2022), this analysis unveiled 14 common MSSVs across 11 chromosomes in homologous offspring (A4, A11, A12, A18, A22, A25, B1, B5, B14, B22, and B25) and only 4 common MSSVs across four chromosomes in heterologous offspring (A3, B14, B20, B21) (Supplementary Figure S4). This outcome suggests that the distribution of MSSVs is essentially random across chromosomes, without significant commonality between the two groups of offspring.

It is generally believed that MSSVs in offspring may have potential effects, we utilized the ANNOVAR tool for comprehensive gene annotation. Analysis showed that homologous offspring exhibited between 483 and 342 intergenic MSSVs, accounting for roughly 26.20% to 31.84% of their total MSSVs, and between 1,167 and 840 MSSVs in gene-associated regions, representing about 68.16%–73.80%. However, heterologous offspring presented with 58–90 intergenic MSSVs, making up approximately 21.53%–26.87%, and 180 to 277 MSSVs within gene-associated regions, equating to about 73.13%–78.47% (Supplementary Table S5).

Generally, intergenic variations, located more than 2 Kb from known functional gene domains, are thought to exert minimal direct influence on gene functionality (Tropp, 2012). In contrast, MSSVs in proximate gene-associated regions (including upstream, downstream, exonic, UTR3, UTR5, ncRNA_intronic, ncRNA_exonic, and ncRNA_splicing areas) are positioned to significantly affect gene operation. The majority of MSSVs were identified in or near gene sequences (Figure 6), suggesting they could have considerable impacts, irrespective of the sperm stimulation method. Predominantly, impactful MSSVs were discovered within introns—non-coding segments known to contain sequences that might enhance gene expression (Rose, 2018). Following the exclusion of intergenic MSSVs, a detailed statistical analysis was conducted, revealing a significant number of genes affected in both homologous (777–1,064) and heterologous (185–284) offspring, primarily within protein-coding and lncRNA gene categories (Figure 7). This suggests these MSSVs predominantly affect biological types of gene like protein coding and lncRNA. This indicates a substantial influence of these MSSVs on crucial gene types, such as protein-coding and lncRNA genes, which are essential for protein synthesis and in regulating gene expression, epigenetics, and cell differentiation, respectively (Alberts et al., 2014; Quinn and Chang, 2016).

Distinctly, the annotation of common MSSVs varied between groups. In homologous offspring, among the 14 common MSSVs, distribution included 2 intergenic, 8 intronic, 2 upstream, 1 across both upstream and downstream, and 1 in UTR3 areas. For heterologous offspring, 4 common MSSVs were identified: 1 intergenic and 3 intronic (Supplementary Figure S5). Targeted analysis on exon region variations, crucial for encoding protein sequences and thus directly influencing protein structure and function, identified exon MSSVs in 6 genes among homologous offspring, all classified as INS. Notably, the sequence and position of these MSSVs in genes such as polyunsaturated fatty acid 5-lipoxygenase, transmembrane protease serine 13-like, and UNC93-like protein MFSD11 were consistent across several offspring (Table 2). While no universally conserved MSSVs were detected across all homologous offspring, specific exon MSSVs recurrently appeared in 3-4 offspring, hinting at a selective presence in distinct regions. In contrast, the identification of exon MSSVs in heterologous offspring was limited due to stringent filtering challenges.

Confirming the reliability of detected MSSVs in offspring. We focused on validating an MSSV within the polyunsaturated fatty acid 5-lipoxygenase gene, chosen due to its significant mutation length and prevalence across multiple offspring samples. This mutation, situated on chromosome A7, constitutes a frameshift mutation that potentially alters the gene’s coding sequence, thereby modifying the amino acid sequence during protein synthesis (Table 2). Polyunsaturated fatty acid 5-lipoxygenase, a dioxygenase enzyme, plays a crucial role in catalyzing the conversion of polyunsaturated fatty acids, such as linoleic acid and arachidonic acid, into peroxides. This process is pivotal for initiating inflammatory responses and activating various cell death pathways, including apoptosis, ferroptosis, and necroptosis (Sun et al., 2019). While research on this gene’s function in fish remains sparse, its mutations are believed to significantly impact its enzymatic activity, potentially altering cell death mechanisms in offspring. The study of this gene in fish is currently limited, but its mutation can affect its function and consequently influence cell death patterns in offspring.

To accurately pinpoint the mutation’s insertion site, we utilized Geneious Prime for sequence alignment and correction. The insertion was identified at position 28,699,675, spanning a sequence length of approximately 5174 bp to 5210 bp (Figure 8; Table 2). The observed sequence length variability is likely due to the inherent error rates associated with sequencing technology and alignment processes. Notably, within PF sample, this mutation was absent, suggesting a homozygous state. Conversely, in PM and offspring samples PP1, PP3, PP4, and PP5, the mutation was supported by at least four reads, indicating heterozygosity. A striking match of four consecutive bases “AGATCC” at the end of the insertion sequence with the original sequence at the mutation site further substantiates this finding (Figure 9).

For further validation, primers were designed using Primer3, targeting both the region preceding the insertion point and the sequence immediately before the variant. Subsequent DNA electrophoresis revealed distinct bands in PM and offspring samples PP1, PP3, PP4, and PP5, but not in PF, unequivocally confirming the inheritance of this MSSV in a subset of the offspring (Figure 9; Supplementary Figure S6).