To explore the different expression patterns among the three types skins, we performed differential gene expression analysis. When comparing PB skin to BP skin and WP skin, the results showed that more DEGs were up-regulated in PB skin, indicating that the formation of black skin is complex and that more genes are needed to participate in the process. Combined with the 10 abundance mRNAs result, there genes rich in melanin synthesis were abundantly expressed in PB skin, while a few of pigment genes were abundantly expressed in WP. It was further suggested that melanin genes were involved in the body color variation of red tilapia during overwintering. Eleven known DEGs were shared by BP_PB, BP_WP and PB_WP comparison groups, of which st3gal1, plxna4, fabp11a and aqp3 plays a vital role in regulating cell proliferation, migration and invasion (Wu et al., 2018; Wang et al., 2020). It has been reported that silencing of st3gal1 gene suppresses melanoma invasion and significantly reduces the survive ability of aggressive melanoma cells of human in metastatic environment (Pietrobono et al., 2020). Plexins family can functionally activate tyrosine kinase receptors in mammalian, such as MET, RON, HER2, and KDR (Swiercz et al., 2008). Overexpression of aqp3 gene can promoted the proliferation and migration of human hepatocytes (Chen et al., 2018). Aqp3 can reduce the differentiation and inhibit the apoptosis of stem cells in human through reducing the expressions of related genes in Wnt/GSK-3 β/β-catenin pathway (Liu et al., 2020). Fabp11a were probably involved in cellular uptake and transport of fatty acids, targeting of fatty acids to transport systems and several signalling pathways in Oryzias latipes (Parmar et al., 2012). In this study, all of these genes expression suggested that skin color variation of red tilapia during overwintering might be related to the proliferation, migration and differentiation of melanocytes.
GO enrichment analysis of DEGs revealed that variations in pigmentation were related to cellular components and biological processes. Most of the DEGs clusters were consistent with previous works on fish, such as zebrafish (Higdon et al., 2013), Midas cichlids (Amphilophus) (Henning et al., 2013) and common carp (Li et al., 2015). KEGG pathway analysis showed that many DEGs were significantly enriched in oxidative phosphorylation, ribosome, ribosome biogenesis in eukaryotes and cardiac muscle contraction in the BP_PB group and PB_WP group. Several studies have shown that high expression of ribosomal protein related genes was associated with the black color in mice (Skarnes et al., 2011). Four of the five highly expressed genes in pigment cells of zebrafish were ribosomal protein (Higdon et al., 2013). There were many DEGs participate in oxidative phosphorylation, cardiac muscle contraction signal pathways in Pristella maxillaris (Bian et al., 2013) and Lutjanus erythropterus (Zhang et al., 2015). In addition, the KEGG results of specific mRNAs showed that ribosome biogenesis in eukaryotes, oxidative phosphorylation were dominant pathways in BP skin. Similar results were found in comparative analysis of BP_PB skin and PB_WP skin. It is suggested that these pathways may play an important role in color variation in red tilapia.
We also found that the DEGs of three skin colors of tilapia were mainly enriched in MAPK signaling pathway, Wnt signaling pathway, tyrosine metabolism, and melanogenesis. Tyrosinase metabolism and melanogenesis pathways have been reported in mammals. And both the Wnt and MAPK signaling pathways are involved in melanophore development in vertebrates (Fujimura et al., 2009, Zhang et al., 2017). For specific mRNA, the tyrosinase metabolism and melanogenesis pathways were dominant in PB skin. It showed that PB skin required more melanin synthesis than BP and WP skins. Meanwhile, we found that some DEGs between the BP and PB skin expressed abundantly in the process of apoptosis and autophagy pathways. The identification of genes in these pigmentation-related term and pathways are informative and are worthy of further study.
Studies have suggested that the mRNA expression levels of genes including tyr, tyrp1, mc1r, mitf, pemel were higher in PB skin. TYR carries out tyrosine hydroxylation to L-DOPA, which is the first step in the biosynthetic pathway of melanin. Under the action of dopachrome tautomerase (DCT) and TYRP1, the dopaquinone (DOPA) chrome was rapidly oxidized and polymerized to form melanin (Braasch et al. 2010, Simon et al., 2009). Therefore, TYR, TYRP1 and DCT are critical enzymes for the formation of melanin. Mutations or dysfunction of tyr or tyrp1 genes lead to melanocyte death or extensive hypopigmentation in zebrafish (Krauss et al, 2015). In our study, compared to WP skin samples, the expression levels of tyr and tyrp1 were the highest in PB skin. This was also consistent with the pigmentation of red crucian carp (Zhang et al, 2017). Mitf is a member transcription factor involved in the development of melanocytes, retinal cells, osteoclasts and mast cells (Minvielle et al., 2010). It has been reported that mitf could directly regulate the expression of multiple genes necessary for the development of melanophores, including tyr, tyrp1, and dct (Cheli et al., 2010). Compared with WP colors, PB skin color was caused by the increase of melanin content, suggesting that mitf may play a potential role in regulating the differentiation and development of melanocytes. Mc1r gene is a key gene in melanogenesis in animals. Alpha-melanocyte stimulating hormone (α-MSH) binds to mc1r, resulting in the decrease of cAMP level. Consequently, melanin biosynthesis process was triggered (Voisey et al., 2001). Previous studies have shown that mc1r mutations were associated with skin color variation in many fish species, such as cavefish, guppy, zebrafish and koi carp (Gross et al., 2009, Tezuka et al., 2011, Richardson et al., 2008, Dong et al., 2020). Similarly, we observed a significant difference in mc1r expression between the red tilapia of three skin colors used in this study. Pmel gene acts as a scaffold in the melanosome by creating a proteolytic fibrillary matrix where melanin is deposited (Solano et al., 2000). Pmel mutations promoted pigment dilution in many animals (Gutierrez et al., 2007). Here pmel was significantly up-regulated in PB skin when compared to WP and BP skin samples. Similarly, we observed the top 10 abundance mRNAs in PB skin, including tyrp1b, tyr, pmelb, pmela, tyrp1a, etc. All genes involved in melanin production, transport and structural proteins for melanin have been verified in red tilapia. In addition, we noticed that the proteins from pigment-related pathways were distinctly integrated together in PPI networks in PB_WP group. It was speculated that those proteins (or genes) could be co-regulated in skin color variation in red tilapia during overwintering. Among them, the most important top 10 genes based on the key nodes in the PPI network, including tyr, mc1r, mitfb, tyrp1b, dct, pmela, pmelb, and mitfa, etc, were consistent with the results of mRNA expression levels.
Regarding the black to pink stage, the mRNA expression levels of asip, tat, hpdb and foxd3 were all up-regulated. Asip gene product blocks melanogenesis by antagonizing the binding of α-MSH to mc1r. Asip mutations were associated with skin color variation in Psetta maxima, zebrafish and medaka (Ceinos et al., 2015, Guillot et al., 2012, Jose´ et al., 2005). In our study, we observed higher expression of asip gene in BP skin transcripts and lower expression of mc1r, which further establishes the role of asip as an antagonistic of the mc1r gene. TAT and HPDB catalyze the substrate tyrosine to form homogeneous acid (HGA). Homogentisate1, 2-dioxygenase catalysis HGA to produce melatonin. Higher level of tat and hpda gene would directly reduce tyrosine level, thereby inhibited the synthesis of melanin (Zhang et al., 2008). In addition, tat gene was the most abundant in WP skin, suggested that it affected the skin variation in red tilapia. Foxd3 is a good candidate for the negative regulator of melanophore development. It can affect the lineage between neural or glial and pigment cells by repressing mitf at the early phase of neural crest migration (Thomas et al., 2009). In addition, overexpression of foxd3 in melb-a mouse melanoblasts blocked the expression of mitf (Lister et al., 2001). Foxd3 was significantly up-regulated in BP skin samples compared to the PB skin, indicated that foxd3 might play a significant role in the black-to-pink color transformation in red tilapia.
In addition, autophagy and apoptotic pathways were able to control the transition from black to pink in red tilapia. In detail, the mRNA level of apoptosis gene, such as baxα, was significantly increased in the body color transformation from the BP skin to the PB and WP skins. Meanwhile, the mRNA levels of beclin1 and autophagy-related genes 7 (atg7), as the autophagy genes, were all upregulated in the BP skin compared with PB skin. Baxα, as one of the homologous proteins of BCL-2, could determine survival or death by an apoptotic stimulus. Overexpression of baxa may accelerate cell death (Oltvai et al., 1993). Beclin1 plays a key role in regulating autophagy and cell death by interacts with either BCL-2 or PI3k class III (Takacs-Vellai et al., 2005). Atg7 activates the ubiquitin-like protein ATGL2, which binds to atg5 and extends the autophagic vesicle membrane. Whole body knock-out of atg7 in mice led to death within 24h after birth (Komatsu et al., 2005). Baxα and atg7 gene were among the top 10 abundance mRNAs in BP skin, further confirmed that appearance of autophagy may lead to melanocyte reduction.