Assay performance
The nCounter analysis data system was selected with possible use by hatchery managers in mind. Nanostring Technologies designs the custom CodeSets for nCounter® analysis based on submitted target sequences and conducts the genomic analyses required to avoid non-specific hybridization, and provides free software, nSolver, to quality check, normalize, and analyze the data. In addition, the system can include over 800 genes and new probes can easily be exchanged in the CodeSet. In our assay the mean raw reads for a given gene were consistent across the three Groups or populations with an average CV of 18.6% (Additional file 1: Table S4D) comparing among the three studies. The high abundance of reads for the five mitochondrial genes limited read values for many nuclear genes with 14 below acceptable limits including two of the reference genes ef1a and ppia (Additional file 1: Table S4D). The CodeSet can be designed to attenuate high abundance genes to alleviate this problem. The cause for apparent instability of b-actin in A1 is not known but it was also found to be unstable among the 20 samples from A1 previously analyzed by RNA-Seq (Ma et al. 2019). Nevertheless, the reliability of reference genes can be questionable with egg quality in which the proportion of mitochondrial and nuclear transcripts can vary with the quality of the eggs (Ma et al. 2019) and the efficiency of methods to capture polyadenylated transcripts can vary for shorter poly(A) tail lengths as with stored maternal mRNAs. Whereas the majority of cytosolic nuclear transcripts in most cells are polyadenylated with a poly(A) tail greater than 80 nucleotides (Villalba et al. 2011; Curanovic et al. 2013), stored maternal nuclear transcripts possess a short poly(A) tail of around 15–40 nucleotides that are elongated to over 80 nucleotides through cytoplasmic polyadenylation during activation (Cabada et al. 1977; Bachvarova 1992; Richter 1999; Villalba et al. 2011; Gohin et al. 2014). Oligo(dT) capture approaches in general are not very efficient at capturing mRNAs with shorter poly(A) tails (Cabada et al. 1977; Meijer et al. 2007; Blower et al. 2013). Assay performance and normalization may be improved with the use of an RNA spike-in as described by Bettwgowda et al. (2006) for use in bovine oocytes. An RNA with a poly(A) tail with a minimum of 25 nucleotides that would likely be efficiently capture by most oligo(dt) methods, and can be used starting at homogenization, should be considered.
Eyeing Rate And Early Embryonic Viability
Previous studies have characterized early mortality in line A. We characterized mortality in the 20 families from population A1 when used in our previous analysis of transcript abundance (Ma et al. 2019). In these 20 families almost all the mortality observed by eyeing, took place between the 8- and 32-cell stages. Stoddard et al. (2005) studying the same line determined most of the mortality in subfertile embryos took place by the second cleavage interval. Although the timing was slightly different, in both studies most of the mortality took place before the MBT and therefore before ZGA. Although the timing of the ZGA has not been determined for rainbow trout, the major wave of ZGA in most fish investigated takes place after the 64-cell stage (Kane and Kimmel 1993; Zamir et al. 1997; Hall et al. 2004; Mathavan et al. 2005; Kleppe et al. 2012). Furthermore, the fish species investigated develop more rapidly than rainbow trout and in general the number of cleavage divisions completed before ZGA is greater in animals that develop more slowly (Marlow 2010). We normally see a more prolonged period of mortality in embryos from our NCCCWA broodstock as observed with the families selected for the present study. Unfortunately, we sampled the embryos in population B earlier than when we sampled A1, and therefore were not able to determine what percentage of mortality took place before the 32-cell stage which was the end of most pre-eyeing mortality in A1. Nevertheless, there was considerable mortality between fertilization and the 8-cell stage; the 8-cell stage and the 16-32-cell stages at about 20 h; ~20 h and streak; and between streak and eyeing (Table 1). As mentioned, the streak rate is a very rough estimate, but the values support at least about 28% of mortality by eyeing took place after the ZGA, suggesting the causes for mortality are likely, at least in part, different for the A and B populations.
Gene Expression
Mitochondrial genes
The rainbow trout mitochondrial genome encodes 13 polypeptides, two rRNAs, 22 tRNAs and a non-coding region (Zardoya et al. 1995). All 13 polypeptide genes and the mt-dlp region transcripts were found to be reduced to a similar degree in eggs with low eyeing rates by Ma et al. (2019). The present assay included five mitochondrial genes, one for a polypeptide from each of the four complexes of the electron transport chain, and the non-coding mt-dlp region. In the zebrafish, mitochondria are required for oxidative phosphorylation to generate ATP for early cleavage events as the maternal ATP pool is insufficient to sustain the proteasomal pathway required for protein degradation needed to advance beyond the 32-cell stage (Dutta and Sinha 2017). Only two of the mitochondrial genes were found to be differentially expressed among the three studies, with mt-cyb being increased with eyeing rate in A1 and A2 and mt-dlp transcripts being negatively correlated with eyeing rate in B (Table 2–5). The P-values reported in Ma et al. (2019) for the differences in expression for the mitochondrial genes were generally much higher than for the nuclear genes selected to be in the present assay, which may also contribute to why a greater proportion of the selected nuclear genes were identified as DEGs in the present groups. Nevertheless, all the mitochondrial genes in all groups trended towards increasing with eyeing rate except for mt-dlp transcripts in group B, consistent with decreased mitochondrial gene expression being a common feature of eggs with reduced developmental competence. Monitoring a suite of mitochondrial genes for small but similar differences in expression may be a reliable approach to identifying eggs with compromised mitochondrial function.
Nuclear Genes
A total of 48 nuclear genes were at detectable levels for all three groups. Transcript abundance was correlated with eyeing rate in 30 of these nuclear genes in at least one of the three groups (Table 5). All the DEGs except ctsz, cycB and mr-1 were among the 49 nuclear genes included in the assay because they were identified as DEGs in our previous RNA-seq study on Group A1. Of the 17 genes included in the assay that had previously been reported as DEGs in studies of rainbow trout other than just Ma et al. (2019) (Aegerter et al. 2004; 2005; Bonnet et al. 2007b), eight were identified as a DEG in at least one of the present groups, seven including igf-2 were not significantly altered, and two (igf-1 and npm2) were below the detection limit in all groups. Therefore, among the genes that were represented at detectable levels, about half of the transcripts selected from our previous study on Group A1, and half of the transcripts identified from studies by different investigators, served as markers for egg quality in at least one of our study groups supporting a fair degree of overlap in genes altered with egg quality among rainbow trout populations.
The highest number of nuclear DEGs, 25, was in group A1 (Table 2) as would be expected since most of the genes in the assay were identified as DEGs using samples from this group (Ma et al. 2019). About half as many genes were DEGs in A2 (Table 3), with 13 of the 14 being DEGs in both A1 and A2 and only fbxo5 significant in A2 and not A1 (Table 5). All 26 genes identified as DEGs for either A1 or A2 in the present study, including fbxo5, were DEGs in Ma et al. (2019). Transcript levels of all 14 genes that were DEGs in both A1 and A2 were positively correlated with eyeing rate and only myo1b and galnt3 were negatively correlated in A1. The mitochondrial gene mt-cyb was also increased with eyeing rate in both A1 and A2. This uniformity in the profiles of maternal transcript dysregulation with poor egg quality for both year classes for population A supports a consistent cause for reduced egg quality in this line. A consistent pattern of dysregulation increases the chances that further investigation can elucidate physiological reasons and an underlying condition for variation in egg quality among females of a broodstock. Again, these transcripts derive from a global transcriptome analysis of egg quality in this line (Ma et al. 2019).
Despite being from a different population than the one from which most of the genes in the assay were selected, Group B (Table 4) had a similar number of nuclear genes with expression levels that correlated with eyeing rate as A2, 11; of which seven were shared with A1 and A2 (Table 5). The four genes that were unique DEGs to group B were also the only four that decreased in abundance with increased eyeing rate. One of the four genes that were unique DEGs to Group B, rpl30, was identified as a DEG in A1 following RNA-Seq analysis (Ma et al. 2019). Therefore, whereas all the genes identified as DEGs in A1 and A2 were previously reported in Ma et al. 2019, Group B had three genes that were not. These genes included the two highest expressed nuclear genes in the assay, cycB and ctsz. Transcript abundance of cycB increased with post-ovulatory ageing and malformation rate at yolk sac resorption, but not survival at eyeing (Aegerter et al. 2004); and ctsz increased with post-ovulatory ageing and decreased with eyeing rate (Aegerter et al. 2005). The abundance of the remaining unique DEG, mr-1, decreased in eggs from fish induced to ovulate with hormone injections; a treatment that also results in more deformed embryos at yolk sac resorption (Bonnet et al. 2007b). Whereas ctsz and cycB generally increased with measures or treatments associated with reduced egg quality among the studies, mr-1 which increased with eyeing rate in our studies, decreased in response to hormone induced spawning.
Two genes, tubb which increased with eyeing rate in A1 and A2, and krt18, which increased with eyeing rate in A1, were also identified as DEGs with eyeing rate by Aegerter et al. (2005). Bonnet et al. (2007b) found rpl24, myo1b, and apoc1, which were DEGs in A1, to have altered expression in response to treatments that changed spawning time, which were in turn shown to increase malformation rates at yolk sac resorption. However, the direction of the effects did not agree among the studies. Transcripts for tubb were increased with eyeing rate in both studies. On the other hand, krt18 decreased with increasing eyeing rate in Aegerter et al. (2005) but was increased with increased egg eyeing rate in Group A1. Nevertheless, the R2 was very low for A1, 0.0453 (Table 2), and the families with medium fertility had the highest mean transcript levels in Groups A1, A2, and B (Table 2–4). Transcript abundance for myo1b was negatively correlated with eyeing rates in the present study and increased in eggs from fish induced to spawn with photoperiod shifting (Bonnet et al. 2007b), whereas apoc1 and rpl24 were positively correlated with eyeing rates in the present study but were increased in eggs of fish injected to spawn with a gonadotropin-releasing hormone analog (Bonnet et al. 2007b).
Many genes that were not significantly correlated with eyeing rate or were below detection in the present study were identified as DEGs in the previous studies of egg quality in rainbow trout. Genes for which no correlation between transcript abundance and eyeing rates were found in the present study included genes of the IGF system, igf-2 and igfr1b. Aegerter et al. (2004) found igf-1, which was below detection in the present assay, and igf-2, to be associated with eyeing rate, and although igfr1b was not associated with eyeing rate, it was decreased with increased malformation rates at yolk sac resorption. Transcript abundance of krt8 and ptgs2 were not altered in our groups although they were negatively correlated with eyeing rate by Aegerter et al. (2005). Transcript abundance of npm2 and igf-1 were increased with eyeing rate in the same study but were below detection in our assay. Detection limits for several of these genes in the present study impaired a comparison among studies. Transcript abundance of pyc and ntan1 were not altered in our groups but were increased in fish exposed to a shifted photoperiod and hormone injection to affect time of spawning, respectively (Bonnet et al. 2007b).
In all, differences in transcript profiles can be found among all studies and populations in which a relation between egg quality and maternal transcriptome have been investigated in rainbow trout. The bases for these differences among studies are not known but there were differences in when mortality took place and treatments associated with decreased developmental competence among the studies that likely influenced or were influenced by the transcript profiles. In the present study almost all the mortality that would take place by eyeing likely took place by 24 h in line A based on studies by Ma et al. (2019) on Group A1 and Stoddard et al. (2015), looking at an earlier year class of the same line. Later mortalities in Group B suggests at least some different factors associated with mortality before eyeing. Further analyses including time of mortality within population B families may suggest if there are transcript profiles associated with early and later embryonic mortality. Aegerter and colleagues (2004; 2005) used post ovulatory ageing to decrease egg quality whereas it was purposely avoided as in the present study. Genes identified as DEGs in the studies on post ovulatory ageing that were not DEGs in the present study or had opposite correlations might suggest those genes that are more specific or responsive to post ovulatory ageing. This includes krt8, krt18, ptgs2, and possibly igf-2 or the IGF system in general. Different expression profiles between the present study and the study by Bonnet et al. 2007b may not be surprising considering the transcript levels in Bonnet et al. 2007b were in response to treatments that were in turn associated with increased malformations at yolk sac resorption as a measure of egg quality. Thus, the transcript profiles in this previous study were not only not in association with the same egg quality metric as in the present study, they were not even in direct association to that egg quality metric, malformations at yolk sac resorption. Nevertheless, the differences in responses of many genes among the studies shows these transcripts cannot be used as general markers of egg quality but does not negate important roles for these transcripts in embryo development.
The identification of six DEGs shared among our three groups, dcaf11, impa2, mrpl39_like, senp7, tfip11 and uchl1, supports some consistency among populations in transcripts that are altered with the egg quality metric eyeing rate. Considering the six genes were identified as DEGs in our previous study that included in-depth genomic analyses within the context of all the identified DEGs in Group A1 (Ma et al. 2019), we will not again discuss possible functions of these genes in detail. The genes do however belong to both disparate and overlapping functional pathways. Several are involved in proteolysis with dcaf11 and uchl1 functioning in ubiquitin pathways, and senp7 and impa2 having hydrolase activity. On the other hand, tfip11 is involved in RNA processing and mrpl39 at least, is involved in mitochondrial function. Unfortunately, they were not included as target genes in previous studies on rainbow trout in which eyeing rate was a metric (Aegerter et al. 2004; 2005). Nevertheless, they were dysregulated in all three of our groups and therefore are the most promising as markers for poor egg quality. Unfortunately, even for these genes the regression coefficients and R2 were low and therefore the expression of any one gene would not be effective at predicting quality of a single batch of eggs. As previously suggested, a suite of genes would likely be required to identify eggs of different developmental competence (Chapman et al. 2014; Sullivan et al. 2015). Still more data would be required if the bases or causes of gene dysregulation, or more specific quality outcomes are to be predicted.