Next Article in Journal
Matrix Metalloproteinases in Human Decidualized Endometrial Stromal Cells
Previous Article in Journal
Effect of Metformin and Simvastatin in Inhibiting Proadipogenic Transcription Factors
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Linking Pedigree Information to the Gene Expression Phenotype to Understand Differential Family Survival Mechanisms in Highly Fecund Fish: A Case Study in the Larviculture of Pacific Bluefin Tuna

1
Bioinformatics and Biosciences Division, Fisheries Stock Assessment Center, Fisheries Resources Institute, Japan Fisheries Research and Education Agency, 2-12-4 Fuku-ura, Kanazawa, Yokohama 236-8648, Kanagawa, Japan
2
Amami Field Station, Tuna Aquaculture Division, Aquaculture Research Department, Fisheries Technology Institute, Japan Fisheries Research and Education Agency, 955-5 Hyousakiyamahara, Setouchi 894-2414, Kagoshima, Japan
3
Highly Migratory Resources Division, Fisheries Stock Assessment Center, Fisheries Resources Institute, 5-7-1 Orido, Shimizu-Ku, Shizuoka-Shi 424-8633, Shizuoka, Japan
4
Shiogama Branch, Fisheries Resources Institute, Japan Fisheries Research and Education Agency, Shinhama 3-27-5, Shiogama 985-0001, Miyagi, Japan
5
Momoshima Field Station, Production Engineering Division, Fisheries Technology Institute, Japan Fisheries Research and Education Agency, 1760 Momoshima, Onomichi 722-0061, Hiroshima, Japan
*
Author to whom correspondence should be addressed.
Curr. Issues Mol. Biol. 2021, 43(3), 2098-2110; https://doi.org/10.3390/cimb43030145
Submission received: 27 October 2021 / Revised: 22 November 2021 / Accepted: 23 November 2021 / Published: 26 November 2021
(This article belongs to the Section Bioinformatics and Systems Biology)

Abstract

:
Mass spawning in fish culture often brings about a marked variance in family size, which can cause a reduction in effective population sizes in seed production for stock enhancement. This study reports an example of combined pedigree information and gene expression phenotypes to understand differential family survival mechanisms in early stages of Pacific bluefin tuna, Thunnus orientalis, in a mass culture tank. Initially, parentage was determined using the partial mitochondrial DNA control region sequence and 11 microsatellite loci at 1, 10, 15, and 40 days post-hatch (DPH). A dramatic proportional change in the families was observed at around 15 DPH; therefore, transcriptome analysis was conducted for the 15 DPH larvae using a previously developed oligonucleotide microarray. This analysis successfully addressed the family-specific gene expression phenotypes with 5739 differentially expressed genes and highlighted the importance of expression levels of gastric-function-related genes at the developmental stage for subsequent survival. This strategy demonstrated herein can be broadly applicable to species of interest in aquaculture to comprehend the molecular mechanism of parental effects on offspring survival, which will contribute to the optimization of breeding technologies.

1. Introduction

Early life stages in many fish species are characterized by high mortality. Fecundity of these fishes is accordingly big, and they typically broadcast their eggs through mass spawning [1]. When culturing these highly fecund fishes, mass spawning arbitrarily in a single tank or cage is a common practice [2]. This reproductive strategy often exhibits a highly imbalanced family structure in offspring, as revealed by genetic-marker-based parentage analyses in a number of aquaculture species, and thereby can lead to a high risk of rapid gene loss [3]. In an extreme case involving the Japanese flounder, Paralichthys olivaceus, more than 99% of larvae from a mass spawning event were sired by a single male in a cohort of six males [4]. Causes of highly variable family survival seem to vary across species, which are influenced by developmental stages, feed formulation, aquaculture facilities, and environmental factors. A better understanding of the mechanisms behind differential family survival will provide valuable information for improving breeding technology in aquaculture.
The genetic basis of natural variation in gene expression exists in organisms from yeast to humans [5], including fish [6,7]. This phenotypic trait, i.e., the gene expression phenotype, also exhibits familial aggregation and a heritable component [8,9,10,11,12]. Thus, linking pedigree information to the gene expression phenotype may provide clues for understanding the molecular mechanism underlying differential family survival in larviculture of aquaculture species.
To test this hypothesis, this study investigated the parental effects on survival and the gene expression phenotype in early stages of Pacific bluefin tuna (PBT), Thunnus orientalis, which were developed in a mass culture tank. Although in 2002, full-life-cycle aquaculture of PBT was developed at Kindai University in Japan, heavy mortality in the larviculture of PBT has been observed [13,14]. We initially examined parentage and family representation for progenies of PBT derived from the mass spawning of a broodstock group at 1, 10, 15, and 40 days post-hatch (DPH) in a mass culture tank. A remarkable degree of proportional change in the family representations occurred after 15 DPH; therefore, we performed transcriptome analysis of the 15 DPH larvae using a previously developed 44 K PBT oligonucleotide microarray (oligo-array) [15]. The analysis successfully addressed the causes of differential family survival by gene expression phenotypes.

2. Materials and Methods

2.1. Parentage Assessment

Fertilized eggs of PBT were obtained on August 31, 2013, from spontaneous spawns of broodstock (7 years old, n = 71) derived from wild-caught yearling tuna in a 40-m-diameter net pen at the Amami Field Station, Fisheries Research and Education Agency. Fertilized eggs were kept in a 500 L tank until 1 DPH, and the larvae at 1 DPH were transferred into a 50 kL mass culture tank (11,100 individuals per m3, water temperature 27–28 °C). Larvae and juveniles were fed the following: rotifer (Brachionus plicatilis sp. complex, 10 individuals per cc) from 2 to 18 DPH, Artemia nauplii from 12 to 24 DPH, yolk-sac larvae of the spangled emperor (Lethrinus nebulosus) from 15 to 39 DPH, and an artificial feed from 21 to 40 DPH (Figure 1). For parentage analysis, a random sample of whole larvae at 1, 10, 15, and 40 DPH (100–500 individuals each) was collected in 100% ethanol and stored at −20 °C until genomic DNA (gDNA) extraction. In addition, larvae of 15 DPH (200 individuals) were also preserved in RNAlater stabilized solution (Thermo Fisher Scientific, Waltham, MA, USA) for both parentage assessment (gDNA) and transcriptome analysis (RNA). The samples in RNAlater were stored at −20 °C until use. Next, 15 DPH larvae (96 individuals) were divided with a scalpel into the caudal portion for gDNA samples and the rest of the body, including the head and internal organs, for total RNA extraction. Genomic DNA samples were extracted from a random sample of 96 progenies at each sampling point using the QuickGene system (Kurabo, Neyagawa, Japan). In contrast, extraction of gDNA from broodstock was performed using a fin clip. Double-stranded sequencing of a region of mitochondrial DNA (mtDNA) encompassing the 3′ portion of cytochrome b to the left domain of the control region determined maternal lines (position 15,246 to 16,156 of a complete mtDNA sequence of Thunnus orientalis (AB185022)). The primers used for PCR and sequencing were CBCR-LT (CAC ATT AAA CCT GAA TGA TA) and MTCR-R1 (CAT TAT TGT ATT TGC ACT GTG A). The PCR profile has been described in another study [16]. PCR products were sequenced from both ends using the BigDye Terminator v.3.1 kit run on an ABI3730 sequencer (ABI, Foster City, CA, USA). Parentage analysis was performed based on the mtDNA control region and 11 microsatellite (MS) loci using the procedure previously described [17].

2.2. Transcriptome Analysis

Here, 4 individuals from each of the 15 DPH larvae from four full-sib families (16 samples) were used for transcriptome analysis. Total RNA was isolated using the RNeasy Plus Universal Mini Kit (Qiagen, Inc., Valencia, CA, USA), and samples were treated with 2 units of TURBO DNase from the TURBO DNase-free Kit (Life Technologies, Carlsbad, CA, USA). Transcriptome analysis was performed using a 44 K PBT oligo-array. Details of the oligo-array and experimental procedures have been described in another study [15]. Oligo-array experiments were performed using an Agilent One-Color platform following the One-Color Microarray-Based Gene Expression Analysis protocol version 6.5 (Agilent Technologies, Santa Clara, CA, USA). Oligo-array slides were scanned on an Agilent DNA Microarray Scanner with Surescan High-Resolution Technology (Agilent Technologies, Santa Clara, CA, USA), and features were extracted and quantified from scanned images using Feature Extraction software version 10.7.3.1 (Agilent Technologies, Santa Clara, CA, USA). Normalization and differential expression analyses were performed using GeneSpring version 12.6 (Agilent Technologies, Santa Clara, CA, USA). After 75th percentile shift normalization, the differentially expressed genes (DEGs) were considered significant according to one-way analysis of variance (ANOVA) followed by a Benjamini–Hochberg multiple testing correction and Tukey HSD post hoc test with a corrected p-value cut-off of 0.05 and a change of two-fold or greater in at least one full-sib family. Visualization of the DEGs was performed by constructing a hierarchical heatmap with Euclidean distance and complete linkage using GeneSpring version 12.6 (Agilent Technologies, Santa Clara, CA, USA).

2.3. RT-qPCR Assay

To validate oligo-array results, we performed an RT-qPCR assay for three pepsinogen genes (PG1, PG2, and PG3), two potassium-transporting ATPase subunit genes (ATP4A and ATP4B), and the β-actin gene (internal control) using the same RNA samples from the oligo-array analysis. PCR primers and FAM™ dye-labeled TaqMan® minor groove binder probes (Applied Biosystems, Foster City, CA, USA) are shown in Table 1.
cDNA was synthesized from 1 μg of total RNA using the PrimeScript™ RT reagent kit with gDNA Eraser (Takara, Otsu, Japan). RT-qPCR was performed using Premix Ex Taq™ (Probe qPCR; Takara, Otsu, Japan) in a Thermal Cycler Dice® Real-time System TP800 (Takara, Otsu, Japan) using the following program: 95 °C for 30 s, 40 cycles of 95 °C for 5 s, and 60 °C for 30 s. The relative expression level of each target gene was calculated using the 2−ΔΔCT method [18] and normalized to the expression of the β-actin gene. The qPCR data were analyzed using Welch’s t-test, and differences were considered statistically significant at p < 0.05.
An RT-qPCR assay for the PG2 gene was conducted for all genotyped 15 DPH progenies (96 individuals) by an absolute quantification method using standard curves derived from serial dilutions of the plasmid standard encoding the target sequence. cDNA synthesis and qPCR procedures were the same as above. Normalization of values was performed by dividing the number of copies of the PG2 gene by that of the β-actin gene. Statistical significance (p < 0.05) was determined by one-way ANOVA with Tukey post hoc analysis (among full-sib families) or by the Wilcoxon rank-sum test (between progenies of two females).

3. Results

3.1. The Influence of Parental Genetic Effects on Offspring Survival

A total of 44.4 × 104 larvae at 1 DPH were transferred into a 50 kL mass culture tank. Survival rates in the hatchery tank decreased sharply to 22.1% by 15 DPH (Figure 1). Finally, the survival rate was 0.5% at 40 DPH (Figure 1). A total of 376 progenies (96 each from 1, 10, and 15 DPH and 88 from 40 DPH) were successfully assigned to their parents. Among them, 14 full-sib families were generated by 7 females and 10 males. Although the family size was unequal at 1 DPH, a large variance in family size was observed dramatically after 15 DPH (Figure 2). In particular, ♀412 was the highest contributing dam at 1 DPH (44.8%) and 10 DPH (45.9%), but the contribution considerably decreased at 15 DPH (21.9%) and 40 DPH (13.6%) (Figure 2). In contrast, the contribution of ♀262 was much lower than that of ♀412 at 1 DPH (16.7%) and 10 DPH (28.2%), but the contribution dramatically increased at 15 DPH (58.4%) and 40 DPH (75.0%) (Figure 2). Families ♀262♂202 and ♀412♂202 shared the same sire (♂202), and the proportion of ♀262♂202 (31.8%) was 7-fold higher than that of ♀412♂202 (4.5%) at 40 DPH (Figure 2). Therefore, these results showed that the effect of dams on post-hatching survival is larger than that of sires, and dam ♀262 had a strong positive effect on offspring survivability during development. In addition, in ♀412-dammed families, a sire effect was also observed for ♂202, as the contribution was nearly unchanged from 5.2% (1 DPH) to 4.5% (40 DPH), while the other sibs in this maternal line markedly decreased from 19.8% (1 DPH) to 1.1% (40 DPH) for ♂297 and 19.8% (1 DPH) to 8.0% (40 DPH) for ♂387 (Figure 2).
Figure 1. Survival rate of hatchery-reared Pacific bluefin tuna larvae in a hatchery tank. DPH indicates days post-hatch. Food sources are also shown below the graph.
Figure 1. Survival rate of hatchery-reared Pacific bluefin tuna larvae in a hatchery tank. DPH indicates days post-hatch. Food sources are also shown below the graph.
Cimb 43 00145 g001

3.2. Parental Effect on Gene Expression Phenotype

Because a considerable change in family representation was observed after 15 DPH, transcriptome analysis was conducted for the PBT larvae at 15 DPH. Based on the parentage assessment (Figure 2), we chose four full-sib families (four individuals each), ♀262♂202, ♀262♂432 (increased family size), ♀412♂387 (decreased family size), and ♀412♂202 (nearly unchanged and small family size), for transcriptome analysis. The hierarchically clustered heat map clearly exhibited family-specific gene expression patterns with 5739 DEGs (Figure 3). In particular, the gene expression patterns were divided into two large clusters by dams: ♀262 and ♀412 (Figure 3). Thus, these results indicate the parental effect on the gene expression phenotype of offspring and a stronger maternal effect compared with the paternal effect at 15 DPH.
The gene expression profiles of three families (♀262♂202, ♀262♂432, and ♀412♂387) were compared with those of the nearly unchanged family, ♀412♂202. The top 20 most highly expressed genes in the three families are listed in Table 2. Notably, the acid-secretion-related genes, ATP4A and ATP4B, and pepsinogen genes, PG1, PG2, and PG3, were highly expressed in ♀262-dammed families and ♀412♂387, with more than a 10-fold difference (Table 2). These results were confirmed by reverse transcription–quantitative polymerase chain reaction (RT-qPCR; n = 4 for each family; Figure 4). ATP4A and ATP4B are subunits of H+/K+-ATPase, which is a gastric proton pump that facilitates the acidic environment necessary for acid hydrolysis of food in the stomach [19,20]. Pepsinogen is the major acidic protease in the stomach, which activates the active form of pepsin under acidic conditions [21]. In contrast, the expression of two potential gastric marker genes, mucin-5AC (MUC5AC) and the Fc fragment of IgG-binding protein (FCGBP), was high in the ♀262-dammed families and ♀412♂387 (Table 2). MUC5AC is a major component of the gastric mucus layer [22], while FCGBP is functionally related to gel-forming mucins [23], and is a highly expressed gene in the stomach of yellowtail (Seriola quinqueradiata) [24]. Furthermore, a trypsinogen (precursor of trypsin) and two solute carrier (SLC) superfamily genes (SLC22A31 and SLC26A9), which are feed-digestion- and nutrient-absorption-related genes, were highly expressed in ♀262-dammed families and ♀412♂387 (Table 2). Trypsin is a proteolytic digestive enzyme that plays an important role in the digestive ability of fish larvae [25], while the SLC superfamily genes encode a series of transporters that play important roles in the absorption of nutrients and ions [26]. It should be noted that murine SLC26A9 is abundantly expressed in the stomach and is an apical HCO3 transporter in gastric surface epithelial cells, playing an important role in protecting the gastric epithelium against acidic injury [27].
In contrast, ♀412-dammed families highly expressed immune-related genes, including GTPase IMAP family members [28], perforin [29], glucose-dependent insulinotropic receptor [30], Toll-like receptor 13 [31], and B cell receptor CD22 [32] (Table S1). It should also be noted that the ♀412-dammed families highly expressed transposable-element-related genes (Table S1). We also explored the highly expressed genes in ♀412♂202 (nearly unchanged and small family size) compared with ♀262♂202, ♀262♂432, and ♀412♂387 and identified eight candidate genes. Five of these eight genes also have putative functions in the immune system, including GIMP7, immune-associated nucleotide-binding protein 13 [33], heat shock protein 30 [34], Caspase-1 [35], and E3 ubiquitin-protein ligase RNF144A-A [36] (Table S2).

3.3. Correlation with the PG2 Expression Levels at 15 DPH and Larval Survival

Because pepsinogen genes showed highly variable expression among families (Table 2 and Figure 4), the gene expression levels of PG2 (the major pepsinogen in PBT [37]) were further analyzed for all 96 genotyped progenies. Comparisons of ♀262- and ♀412-dammed families clearly showed that the PG2 gene is expressed at significantly higher levels (p = 8 × 10−5) in ♀262-dammed families (increased family size) at 15 DPH (Figure 5a). Figure 5b demonstrates expression levels of PG2 among full-sib families. The strength of PG2 expression levels mostly corresponded to the family size at 15 DPH (Figure 2 and Figure 5b). The expression of the PG2 gene in the ♀412♂387 family was also higher than that in the ♀412♂202 family, as observed in the oligo-array results (Table 2), but the difference was not statistically significant (Figure 5b). Only two ♀262-dammed families, ♀262♂202 and ♀262♂432, had significantly higher PG2 expression (p < 0.05) than that of the ♀412♂202 family (Figure 5b).

4. Discussion

A dramatic proportional change in family representations with a strong effect of dams was observed during PBT larval development from 1 to 40 DPH, prominently from 10 to 15 DPH (Figure 2). These results suggested that critical events might be occurring that are responsible for subsequent larval survival at around 15 DPH, particularly in rising ♀262 and declining ♀412 progenies.
Transcriptome analysis of the 15 DPH larvae clearly exhibited that gene expression phenotypes are strongly correlated with the genetic architecture of their parents (Figure 3), with dams having a greater effect than sires, as expected by the parentage assessment (Figure 2). This transcriptome analysis highlighted the importance of the expression levels of genes related to the digestive function of the stomach of 15 DPH larvae for their subsequent survival (Table 2 and Figure 4). This finding was strongly supported by qRT-PCR analysis of PG2 using all genotyped larvae (n = 96). The gene expression strength was positively correlated with family survival and was significantly higher in ♀262-dammed families (increased families) than in ♀412-dammed families (declining families; Figure 2 and Figure 5). Although oligo-array results (n = 4) showed a high level of PG2 gene expression in ♀412♂387 (Table 2), qRT-PCR results (more quantitative methods and more sample number) did not highlight any statistically difference between the ♀412♂387 and ♀412♂202 families (Figure 5b). The ♀412♂387 family dramatically decreased from 1 to 40 DPH, and most of the larvae died before 40 DPH. Therefore, the results may have been biased by remaining survivors because of higher expression of these genes.
The 15 DPH larvae are at the postflexion stage, where gastric glands begin to form [38,39] and pepsin activity rapidly increases [39,40]. It should be noted that a dramatic increase in the expression of the gastric-function-related genes (Table 2), gastric proton pumps (ATP4A and ATP4B), pepsinogens (PG1, PG2, and PG3), and potential gastric markers (MUC5AC, FCGBP, and SLC26A9), has been observed from 13 to 15 DPH in a previous transcriptome profiling (same oligo-array used in this study) of PBT early developmental stages (1–25 DPH; unpublished data), suggesting the importance of these genes at 15 DPH. Thus, our transcriptome analysis suggested that the larval stomachs of low-mortality families are more functionally developed than those of high-mortality families at this developmental stage.
In larviculture of PBT, 15 DPH is approximately the timing of the initial feeding of fish larvae in the presence of rotifers and Artemia nauplii as prey in the culture tank (Figure 1). Subsequently, marked growth variations are observed, which are related to differences in the initiation of piscivory [14,41] and the use of yolk-sac larvae [42,43]. This growth variation leads to growth-selective mortality, and PBT larvae, which immediately use fish larvae, grow rapidly, and survive in the hatchery tank [14]. These findings suggest that piscivory is a key event for larval survival at around 15 DPH. This study strongly implied that the onset of gastric function in this stage is important for the digestion of large prey, such as fish larvae, which affects subsequent larval survival and is influenced by the genetic background of the parents. In contrast, Sabate et al. [44] reported that the chase behavior of PBT larvae, which is closely related to cannibalistic behavior, was first observed at 14 DPH at the flexion-to-postflexion stages and increased thereafter. Therefore, it is possible that PBT larvae with well-developed gastric functions (i.e., low-mortality families) are capable of cannibalistic attacks on relatively underdeveloped larvae (i.e., high-mortality families) and thereafter induce selective mortality.
Several immune- and transposable-element-related genes were more highly expressed in the ♀412-dammed families than in the ♀262-dammed families at 15 DPH (Table S1). In addition, ♀412♂202 specifically expressed five immune-related genes (Table S2). Transposable elements of eukaryotic genomes have also been reported to be activated by stress and trigger an innate immune response [45,46]. Thus, the regulation of immune function may be more active in the ♀412-dammed families. The possible cause of the negative effects of immune function on survivability may be explained by trade-offs between immune and competing physiological functions [47]. Immune function is an energetically costly physiological activity [47]; thus the energy allocated to gastric functions might be lower in ♀412-dammed families.
Another possibility is that the over-activated immune function might raise survival, compensating the gastric dysfunction, especially in the ♀412♂202 family, which maintained its size throughout the 40 days in this study. An acidic gastric environment sterilizes ingested microorganisms, while they survive through the stomachless fish intestine [48]. Gastric dysfunction larvae might, therefore, suffer opportunistic bacterial infection in their digestive tract, because the intestinal epithelium is a potential port of pathogen entry [49]. The over-activated immune function might alleviate bacterial infection in the gastric dysfunction larvae. Since little is known about the microbiome of tuna larvae, further studies are required to investigate the relationship between the gut microbiome and survival of PBT larvae.
In summary, this study indicated that the functional development of the stomach in PBT larvae during the initial feeding of fish larvae is key to the subsequent survival of PBT larvae in the hatchery tank. In addition, expression levels of gastric-function-related genes (ATP4A, ATP4B, PG1, PG2, PG3, MUC5AC, FCGBP, and SLC26A9) in this developmental stage might be useful to distinguish between low- and high-mortality families. Recently, RNA-Seq using high-throughput sequencing technology has been developed as an approach for transcriptome analysis [50]. Unlike microarray technology, RNA-Seq can be applied to non-model organisms with no genomic resources. Thus, although our study focused on the larviculture of PBT and used a microarray platform [15], the strategy highlighted herein can be extended to species of interest in aquaculture to find clues to understand the molecular mechanisms linked to variable family survival. This information may be useful for the optimization of breeding technologies in the aquaculture industry.

Supplementary Materials

The following are available online at https://www.mdpi.com/article/10.3390/cimb43030145/s1.

Author Contributions

Conceptualization, M.Y. and T.S.; methodology, M.Y., K.S. and T.S.; software, M.Y.; validation, M.Y. and T.S.; formal analysis, M.Y. and T.S.; investigation, M.Y. and T.S.; resources, K.K. and Y.T.; data curation, M.Y., K.K. and T.S.; writing—original draft preparation, M.Y.; writing—review and editing, K.S., Y.T., K.K. and T.S.; visualization, M.Y. All authors have read and agreed to the published version of the manuscript.

Funding

This study was supported by a grant for Technological Development for Selection and Secure Stock of Bloodstock for Culture of Bluefin Tuna from the Agriculture, Forestry, and Fisheries Research Council (AFFRC) of Japan.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

All raw and processed microarray data were deposited in the DDBJ Genomic Expression Archive under accession number E-GEAD-447.

Acknowledgments

We thank Takahiro Ohtsuka, Ayako Kondo, and Haruka Ito (Yokohama Field Station, FRA) for assistance with the laboratory work.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Houde, E.D. Patterns and trends in larval-stage growth and mortality of teleost fish. J. Fish Biol. 1997, 51, 52–83. [Google Scholar] [CrossRef]
  2. Robledo, D.; Palaiokostas, C.; Bargelloni, L.; Martínez, P.; Houston, R. Applications of genotyping by sequencing in aquaculture breeding and genetics. Rev. Aquac. 2018, 10, 670–682. [Google Scholar] [CrossRef]
  3. Vandeputte, M.; Haffray, P. Parentage assignment with genomic markers: A major advance for understanding and exploiting genetic variation of quantitative traits in farmed aquatic animals. Front. Genet. 2014, 5, 432. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  4. Sekino, M.; Saitoh, K.; Yamada, T.; Kumagai, A.; Hara, M.; Yamashita, Y. Microsatellite-based pedigree tracing in a Japanese flounder Paralichthys olivaceus hatchery strain: Implications for hatchery management related to stock enhancement program. Aquaculture 2003, 221, 255–263. [Google Scholar] [CrossRef]
  5. Cheung, V.G.; Spielman, R.S. The genetics of variation in gene expression. Nat. Genet. 2002, 32, 522–525. [Google Scholar] [CrossRef] [PubMed]
  6. Oleksiak, M.F. Genomic approaches with natural fish populations. J. Fish Biol. 2010, 76, 1067–1093. [Google Scholar] [CrossRef] [Green Version]
  7. Oleksiak, M.F.; Churchill, G.A.; Crawford, D.L. Variation in gene expression within and among natural populations. Nat. Genet. 2002, 32, 261–266. [Google Scholar] [CrossRef] [PubMed]
  8. Cheung, V.G.; Conlin, L.K.; Weber, T.M.; Arcaro, M.; Jen, K.-Y.; Morley, M.; Spielman, R.S. Natural variation in human gene expression assessed in lymphoblastoid cells. Nat. Genet. 2003, 33, 422–425. [Google Scholar] [CrossRef]
  9. Morley, M.; Molony, C.M.; Weber, T.M.; Devlin, J.L.; Ewens, K.G.; Spielman, R.S.; Cheung, V.G. Genetic analysis of genome-wide variation in human gene expression. Nature 2004, 430, 743–747. [Google Scholar] [CrossRef]
  10. Schadt, E.E.; Monks, S.A.; Drake, T.A.; Lusis, A.J.; Che, N.; Colinayo, V.; Ruff, T.G.; Milligan, S.B.; Lamb, J.R.; Cavet, G.; et al. Genetics of gene expression surveyed in maize, mouse and man. Nat. Cell Biol. 2003, 422, 297–302. [Google Scholar] [CrossRef]
  11. Bicskei, B.; Taggart, J.B.; Glover, K.A.; Bron, J.E. Comparing the transcriptomes of embryos from domesticated and wild Atlantic salmon (Salmo salar L.) stocks and examining factors that influence heritability of gene expression. Genet. Sel. Evol. 2016, 48, 1–16. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  12. Bernatchez, L.; Audet, C.; Bougas, B. The influence of parental effects on transcriptomic landscape during early development in brook charr (Salvelinus fontinalis, Mitchill). Heredity 2013, 110, 484–491. [Google Scholar] [CrossRef]
  13. Sawada, Y.; Okada, T.; Miyashita, S.; Murata, O.; Kumai, H. Completion of the Pacific bluefin tuna Thunnus orientalis (Temminck et Schlegel) life cycle. Aquac. Res. 2005, 36, 413–421. [Google Scholar] [CrossRef]
  14. Tanaka, Y.; Kumon, K.; Ishihi, Y.; Eba, T.; Nishi, A.; Nikaido, H.; Shiozawa, S. Mortality processes of hatchery-reared Pacific bluefin tuna Thunnus orientalis (Temminck et Schlegel) larvae in relation to their piscivory. Aquac. Res. 2017, 49, 11–18. [Google Scholar] [CrossRef] [Green Version]
  15. Yasuike, M.; Fujiwara, A.; Nakamura, Y.; Iwasaki, Y.; Nishiki, I.; Sugaya, T.; Shimizu, A.; Sano, M.; Kobayashi, T.; Ototake, M. A functional genomics tool for the Pacific bluefin tuna: Development of a 44K oligonucleotide microarray from whole-genome sequencing data for global transcriptome analysis. Gene 2016, 576, 603–609. [Google Scholar] [CrossRef] [Green Version]
  16. Saitoh, K.; Suzuki, N.; Ozaki, M.; Ishii, K.; Sado, T.; Morosawa, T.; Tsunagawa, T.; Tsuchiya, M. Natural habitats uncovered?—Genetic structure of known and newly found localities of the endangered bitterling Pseudorhodeus tanago (Cyprinidae). Nat. Conserv. 2017, 17, 19–33. [Google Scholar] [CrossRef] [Green Version]
  17. Uchino, T.; Nakamura, Y.; Sekino, M.; Kai, W.; Fujiwara, A.; Yasuike, M.; Sugaya, T.; Fukuda, H.; Sano, M.; Sakamoto, T. Con-structing genetic linkage maps using the whole genome sequence of Pacific bluefin tuna and a comparison of chromosome structure among teleost species. Adv. Biosci. Biotechnol. 2016, 7, 85–122. [Google Scholar] [CrossRef] [Green Version]
  18. Livak, K.J.; Schmittgen, T.D. Analysis of relative gene expression data using real-time quantitative PCR and the 2−ΔΔCT method. Methods 2001, 25, 402–408. [Google Scholar] [CrossRef] [PubMed]
  19. Maeda, M.; Oshiman, K.; Tamura, S.; Futai, M. Human gastric (H+ + K+)-ATPase gene. Similarity to (Na+ + K+)-ATPase genes in exon/intron organization but difference in control region. J. Biol. Chem. 1990, 265, 9027–9032. [Google Scholar] [CrossRef]
  20. Sachs, G.; Chang, H.H.; Rabon, E.; Schackman, R.; Lewin, M.; Saccomani, G. A nonelectrogenic H+ pump in plasma membranes of hog stomach. J. Biol. Chem. 1976, 251, 7690–7698. [Google Scholar] [CrossRef]
  21. Foltmann, B. Gastric proteinases--structure, function, evolution and mechanism of action. Essays Biochem. 1981, 17, 52–84. [Google Scholar]
  22. Reis, C.A.; David, L.; Nielsen, P.A.; Clausen, H.; Mirgorodskaya, K.; Roepstorff, P.; Sobrinho-Simoes, M. Immunohistochemical study of MUC5AC expression in human gastric carcinomas using a novel monoclonal antibody. Int. J. Cancer 1997, 74, 112–121. [Google Scholar] [CrossRef]
  23. Lang, T.; Klasson, S.; Larsson, E.; Johansson, M.E.; Hansson, G.C.; Samuelsson, T. Searching the Evolutionary origin of epi-thelial mucus protein components-mucins and FCGBP. Mol. Biol. Evol. 2016, 33, 1921–1936. [Google Scholar] [CrossRef]
  24. Yasuike, M.; Iwasaki, Y.; Nishiki, I.; Nakamura, Y.; Matsuura, A.; Yoshida, K.; Noda, T.; Andoh, T.; Fujiwara, A. The yellowtail (Seriola quinqueradiata) genome and transcriptome atlas of the digestive tract. DNA Res. 2018, 25, 547–560. [Google Scholar] [CrossRef] [Green Version]
  25. Manchado, M.; Infante, C.; Asensio, E.; Crespo, A.; Zuasti, E.; Cañavate, J.P. Molecular characterization and gene expression of six trypsinogens in the flatfish Senegalese sole (Solea senegalensis Kaup) during larval development and in tissues. Comp. Biochem. Physiol. Part B Biochem. Mol. Biol. 2008, 149, 334–344. [Google Scholar] [CrossRef]
  26. Xie, J.; Zhu, X.Y.; Liu, L.M.; Meng, Z.Q. Solute carrier transporters: Potential targets for digestive system neoplasms. Cancer Manag. Res. 2018, 10, 153–166. [Google Scholar] [CrossRef] [Green Version]
  27. Xu, J.; Henriksnäs, J.; Barone, S.; Witte, D.; Shull, G.E.; Forte, J.G.; Holm, L.; Soleimani, M. SLC26A9 is expressed in gastric surface epithelial cells, mediates Cl/HCO3 exchange, and is inhibited by NH4+. Am. J. Physiol. Physiol. 2005, 289, C493–C505. [Google Scholar] [CrossRef]
  28. Limoges, M.-A.; Cloutier, M.; Nandi, M.; Ilangumaran, S.; Ramanathan, S. The GIMAP family proteins: An incomplete puzzle. Front. Immunol. 2021, 12, 679739. [Google Scholar] [CrossRef] [PubMed]
  29. Xu, J.; Yang, N.; Xie, T.; Yang, G.; Chang, L.; Yan, D.; Li, T. Summary and comparison of the perforin in teleosts and mammals: A review. Scand. J. Immunol. 2021, 94, e13047. [Google Scholar] [CrossRef] [PubMed]
  30. Efimova, I.; Steinberg, I.; Zvibel, I.; Neumann, A.; Mantelmacher, D.F.; Drucker, D.J.; Fishman, S.; Varol, C. GIPR signaling in immune cells maintains metabolically beneficial type 2 immune responses in the white fat from obese mice. Front. Immunol. 2021, 12, 643144. [Google Scholar] [CrossRef] [PubMed]
  31. Wang, Y.; Bi, X.; Chu, Q.; Xu, T. Discovery of toll-like receptor 13 exists in the teleost fish: Miiuy croaker (Perciformes, Sciaenidae). Dev. Comp. Immunol. 2016, 61, 25–33. [Google Scholar] [CrossRef]
  32. Clark, E.A.; Giltiay, N.V. CD22: A regulator of innate and adaptive B cell responses and autoimmunity. Front. Immunol. 2018, 9, 2235. [Google Scholar] [CrossRef]
  33. Liu, C.; Wang, T.; Zhang, W.; Li, X. Computational identification and analysis of immune-associated nucleotide gene family in Arabidopsis thaliana. J. Plant Physiol. 2008, 165, 777–787. [Google Scholar] [CrossRef] [PubMed]
  34. Liu, P.; Wang, L.; Kwang, J.; Yue, G.H.; Wong, S.-M. Transcriptome analysis of genes responding to NNV infection in Asian seabass epithelial cells. Fish Shellfish Immunol. 2016, 54, 342–352. [Google Scholar] [CrossRef]
  35. Zhang, X.; Liu, Z.; Li, C.; Zhang, Y.; Wang, L.; Wei, J.; Qin, Q. Characterization of orange-spotted grouper (Epinephelus coioides) ASC and caspase-1 involved in extracellular ATP-mediated immune signaling in fish. Fish Shellfish Immunol. 2020, 97, 58–71. [Google Scholar] [CrossRef]
  36. Afzali, B.; Kim, S.; West, E.; Nova-Lamperti, E.; Cheru, N.; Nagashima, H.; Yan, B.; Freiwald, T.; Merle, N.; Chauss, D.; et al. RNF144A shapes the hierarchy of cytokine signaling to provide protective immunity against influenza. Biorxiv 2019, 782680. [Google Scholar] [CrossRef]
  37. Tanji, M.; Yakabe, E.; Kubota, K.; Kageyama, T.; Ichinose, M.; Miki, K.; Ito, H.; Takahashi, K. Structural and phylogenetic comparison of three pepsinogens from Pacific bluefin tuna: Molecular evolution of fish pepsinogens. Comp. Biochem. Physiol. Part B Biochem. Mol. Biol. 2009, 152, 9–19. [Google Scholar] [CrossRef]
  38. Kaji, T.; Tanaka, M.; Takahashi, Y.; Oka, M.; Ishibashi, N. Preliminary observations on development of Pacific bluefin tuna Thunnus thynnus (Scombridae) larvae reared in the laboratory, with special reference to the digestive system. Mar. Freshw. Res. 1996, 47, 261–269. [Google Scholar] [CrossRef]
  39. Miyashita, S.; Kato, K.; Sawada, Y.; Murata, O.; Ishitani, Y.; Shimizu, K.; Yamamoto, S.; Kumai, H. Development of digestive system and digestive enzyme activities of larval and juvenile bluefin tuna, Thunnus thynnus, reared in the laboratory. Aquac. Sci. 1998, 46, 111–120. [Google Scholar]
  40. Murashita, K.; Matsunari, H.; Kumon, K.; Tanaka, Y.; Shiozawa, S.; Furuita, H.; Oku, H.; Yamamoto, T. Characterization and ontogenetic development of digestive enzymes in Pacific bluefin tuna Thunnus orientalis larvae. Fish Physiol. Biochem. 2014, 40, 1741–1755. [Google Scholar] [CrossRef] [PubMed]
  41. Tanaka, Y.; Minami, H.; Ishihi, Y.; Kumon, K.; Higuchi, K.; Eba, T.; Nishi, A.; Nikaido, H.; Shiozawa, S. Differential growth rates related to initiation of piscivory by hatchery-reared larval Pacific bluefin tuna Thunnus orientalis. Fish. Sci. 2014, 80, 1205–1214. [Google Scholar] [CrossRef]
  42. Tanaka, Y.; Minami, H.; Ishihi, Y.; Kumon, K.; Eba, T.; Nishi, A.; Nikaido, H.; Shiozawa, S. Prey utilization by hatchery-reared Pacific bluefin tuna larvae in mass culture tank estimated using stable isotope analysis, with special reference to their growth variation. Aquac. Sci. 2010, 58, 501–508. [Google Scholar]
  43. Tanaka, Y.; Minami, H.; Ishihi, Y.; Kumon, K.; Higuchi, K.; Eba, T.; Nishi, A.; Nikaido, H.; Shiozawa, S. Relationship between prey utilization and growth variation in hatchery-reared Pacific bluefin tuna, Thunnus orientalis (Temminck et Schlegel), larvae estimated using nitrogen stable isotope analysis. Aquac. Res. 2012, 45, 537–545. [Google Scholar] [CrossRef]
  44. Sabate, F.D.L.S.; Sakakura, Y.; Tanaka, Y.; Kumon, K.; Nikaido, H.; Eba, T.; Nishi, A.; Shiozawa, S.; Hagiwara, A.; Masuma, S. Onset and development of cannibalistic and schooling behavior in the early life stages of Pacific bluefin tuna Thunnus orientalis. Aquaculture 2010, 301, 16–21. [Google Scholar] [CrossRef] [Green Version]
  45. Bourque, G.; Burns, K.H.; Gehring, M.; Gorbunova, V.; Seluanov, A.; Hammell, M.; Imbeault, M.; Izsvák, Z.; Levin, H.L.; Macfarlan, T.S.; et al. Ten things you should know about transposable elements. Genome Biol. 2018, 19, 1–12. [Google Scholar] [CrossRef]
  46. Horváth, V.; Merenciano, M.; González, J. Revisiting the relationship between transposable elements and the eukaryotic stress response. Trends Genet. 2017, 33, 832–841. [Google Scholar] [CrossRef]
  47. Urlacher, S.; Ellison, P.T.; Sugiyama, L.S.; Pontzer, H.; Eick, G.; Liebert, M.A.; Cepon-Robins, T.J.; Gildner, T.E.; Snodgrass, J.J. Tradeoffs between immune function and childhood growth among Amazonian forager-horticulturalists. Proc. Natl. Acad. Sci. USA 2018, 115, E3914–E3921. [Google Scholar] [CrossRef] [Green Version]
  48. Miura, T.; Iwata, K.; Zongshe, Z. Biological significance of effects of phytophagous fishes on phytoplankton. Chin. J. Oceanol. Limnol. 1989, 7, 335–338. [Google Scholar] [CrossRef]
  49. Ringø, E.; Olsen, R.E.; Mayhew, T.M.; Myklebust, R. Electron microscopy of the intestinal microflora of fish. Aquaculture 2003, 227, 395–415. [Google Scholar] [CrossRef]
  50. Wang, Z.; Gerstein, M.; Snyder, M. RNA-Seq: A revolutionary tool for transcriptomics. Nat. Rev. Genet. 2009, 10, 57–63. [Google Scholar] [CrossRef]
Figure 2. Parentage assignment for Pacific bluefin tuna larvae at 1, 10, 15, and 40 days post-hatch (DPH) using 11 microsatellite (MS) markers.
Figure 2. Parentage assignment for Pacific bluefin tuna larvae at 1, 10, 15, and 40 days post-hatch (DPH) using 11 microsatellite (MS) markers.
Cimb 43 00145 g002
Figure 3. Hierarchically clustered heat map of 5739 differentially expressed genes (DEGs) among four full-sib families, ♀262♂202, ♀262♂432, ♀412♂387, and ♀412♂202 (four individuals each), of the 15 DPH Pacific bluefin tuna larvae.
Figure 3. Hierarchically clustered heat map of 5739 differentially expressed genes (DEGs) among four full-sib families, ♀262♂202, ♀262♂432, ♀412♂387, and ♀412♂202 (four individuals each), of the 15 DPH Pacific bluefin tuna larvae.
Cimb 43 00145 g003
Figure 4. RT-qPCR analysis of the highly expressed genes, three pepsinogen genes (PG1, PG2, and PG3) and two potassium-transporting ATPase subunit genes (ATP4A and ATP4B), in three families (♀262♂202, ♀262♂432, and♀412♂387) of the 15 DPH Pacific bluefin tuna larvae. The graph represents the relative fold-change values when compared to the ♀412♂202 family, and error bars show the standard error of the mean (SEM). Asterisks indicate a statistically significant difference (* p < 0.05; ** p < 0.01) when compared with the high-mortality family (♀412♂202).
Figure 4. RT-qPCR analysis of the highly expressed genes, three pepsinogen genes (PG1, PG2, and PG3) and two potassium-transporting ATPase subunit genes (ATP4A and ATP4B), in three families (♀262♂202, ♀262♂432, and♀412♂387) of the 15 DPH Pacific bluefin tuna larvae. The graph represents the relative fold-change values when compared to the ♀412♂202 family, and error bars show the standard error of the mean (SEM). Asterisks indicate a statistically significant difference (* p < 0.05; ** p < 0.01) when compared with the high-mortality family (♀412♂202).
Cimb 43 00145 g004
Figure 5. Transcript levels of the PG2 gene in all genotyped 15 DPH progenies (96 individuals) measured by RT-qPCR. Results are expressed as the ratio of the number of copies of the PG2 gene to the number of copies of the β-actin gene. (a) Box plot showing a comparison of the expression level of PG2 between progenies of ♀262 (low mortality, n = 56) and ♀412 (high mortality, n = 21). Dots represent individual expression values. (b) Bean plots showing expression levels of the PG2 gene among full-sib families. Others indicate the full-sib families that include no less than three individuals (<n = 3). The horizontal lines in each bean refer to individual expression values. The thick black horizontal line within each bean is the mean, while the dashed horizontal line is the overall mean across beans. Different letters represent a significant difference (p < 0.05) by one-way analysis of variance with Tukey post hoc analysis.
Figure 5. Transcript levels of the PG2 gene in all genotyped 15 DPH progenies (96 individuals) measured by RT-qPCR. Results are expressed as the ratio of the number of copies of the PG2 gene to the number of copies of the β-actin gene. (a) Box plot showing a comparison of the expression level of PG2 between progenies of ♀262 (low mortality, n = 56) and ♀412 (high mortality, n = 21). Dots represent individual expression values. (b) Bean plots showing expression levels of the PG2 gene among full-sib families. Others indicate the full-sib families that include no less than three individuals (<n = 3). The horizontal lines in each bean refer to individual expression values. The thick black horizontal line within each bean is the mean, while the dashed horizontal line is the overall mean across beans. Different letters represent a significant difference (p < 0.05) by one-way analysis of variance with Tukey post hoc analysis.
Cimb 43 00145 g005
Table 1. Primers and probes used for RT-qPCR assay.
Table 1. Primers and probes used for RT-qPCR assay.
GeneForward PrimerReverse PrimerTaqMan MGB Probe
PG1AACGAGCTGTACTGGCAGATCAAGCCACAACCTGACCATTGACCAGTGGACAGTGTTACC
PG2CCGAGGTCACCTTCACTCTGAGCCAGTTCTGCAACCATAGTAGCTCTGCATCTGCCTACGTC
PG3CCACCTACCTGCCCTCTAGTGAGTGCGGTCGTAGACGGAGTAGTCTCTGTGGATCTTTG
ATP4ATCCTCCAAGAGCCACTGTACCTCACCATGACAACCCTGATACCACAGTGATGAAATGTCG
ATP4BCCATGCCTTGTGTCATCATTAAGATTCTCCTGTCCTTCCAGTATGGTTGAACAGGATCATTGGC
β-actinGAAATCGCCGCACTGGTTGCATCGTCTCCGGCAAATATCCGGAATGTGCAAAG
Table 2. The top 20 highly expressed genes in three families (♀262♂202, ♀262♂432, and ♀412♂387) compared to the family ♀412♂202 at 15 DPH.
Table 2. The top 20 highly expressed genes in three families (♀262♂202, ♀262♂432, and ♀412♂387) compared to the family ♀412♂202 at 15 DPH.
RankingOligo Array Probe ID Putative Gene Product♀262♂202/
♀412♂202
♀262♂432/
♀412♂202
♀412♂387/
♀412♂202
1Ba00000144_g2079Pepsinogen 3 (PG3)26.40 22.88 49.55
2Ba00008844_g23666Hyaluronidase-122.02 26.52 28.03
3Ba00008341_g23272Pepsinogen 3 (PG3)16.78 13.38 27.67
4Ba00001220_g9569Potassium-transporting ATPase alpha (ATP4A)22.39 17.69 16.29
5isotigB109707_cHypothetical protein16.05 16.33 15.60
6Ba00000817_g7404Pepsinogen 1 (PG1)13.95 14.12 14.04
7isotigB46117_nRNA-directed DNA polymerase17.45 15.52 8.50
8Ba00002263_g13715Potassium-transporting ATPase subunit beta (ATP4B)15.41 10.59 15.01
9Ba00000128_g1851Pepsinogen 2 (PG2)13.06 12.98 13.11
10isotigB18854_cHypothetical protein7.91 14.33 6.00
11Ba00000638_g6231Hypothetical protein12.22 6.26 5.80
12Ba00005732_g20720Hypothetical protein6.68 8.80 6.61
13isotigB46118_nProstate stem cell antigen8.78 7.36 4.15
14Ba00000678_g6517Solute carrier family 26, member 9 (SLC26A9)6.47 7.01 5.90
15Ba00000139_g2005Trypsinogen3.81 8.13 6.84
16Ba00003884_g17895.p2.613-1845Solute carrier family 22 member 31 (SLC22A31)5.60 5.27 6.51
17Ba00000257_g3256Fc fragment of IgG binding protein (FCGBP)5.53 6.19 3.57
18Ba00010561_g24778ATP-dependent DNA helicase PIF1 4.77 5.98 4.41
19Ba00004935_g19660Mucin-5AC (MUC5AC)4.06 4.36 5.68
20BaME00000895_g9062Hypothetical protein4.94 4.15 4.68
Values indicate the fold-change of expression when compared with the high-mortality family (♀412♂202).
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Yasuike, M.; Kumon, K.; Tanaka, Y.; Saitoh, K.; Sugaya, T. Linking Pedigree Information to the Gene Expression Phenotype to Understand Differential Family Survival Mechanisms in Highly Fecund Fish: A Case Study in the Larviculture of Pacific Bluefin Tuna. Curr. Issues Mol. Biol. 2021, 43, 2098-2110. https://doi.org/10.3390/cimb43030145

AMA Style

Yasuike M, Kumon K, Tanaka Y, Saitoh K, Sugaya T. Linking Pedigree Information to the Gene Expression Phenotype to Understand Differential Family Survival Mechanisms in Highly Fecund Fish: A Case Study in the Larviculture of Pacific Bluefin Tuna. Current Issues in Molecular Biology. 2021; 43(3):2098-2110. https://doi.org/10.3390/cimb43030145

Chicago/Turabian Style

Yasuike, Motoshige, Kazunori Kumon, Yosuke Tanaka, Kenji Saitoh, and Takuma Sugaya. 2021. "Linking Pedigree Information to the Gene Expression Phenotype to Understand Differential Family Survival Mechanisms in Highly Fecund Fish: A Case Study in the Larviculture of Pacific Bluefin Tuna" Current Issues in Molecular Biology 43, no. 3: 2098-2110. https://doi.org/10.3390/cimb43030145

Article Metrics

Back to TopTop