A transcriptomic approach to study the effect of long-term starvation and diet composition on the expression of mitochondrial oxidative phosphorylation genes in gilthead sea bream (Sparus aurata)

Animals, feeding trial, sampling and growth performance

RNA extraction, cDNA library construction, normalisation and 454 sequencing

Total RNA was extracted from 30 mg of liver or white skeletal muscle using the RNeasy tissue and RNeasy fibrous tissue, respectively, mini kits (Qiagen, Hilden, Germany) according to the manufacturer’s instructions. RNA concentration and purity was determined spectrophotometrically at 260/280 nm using Nanodrop ND-1000 (Thermo Fischer Scientific, Waltham, MA, USA). RNA integrity was determined with an Agilent 2100 bioanalyzer (Agilent Technologies, Santa Clara, CA, USA). Only samples with RNA Integrity Number (RIN) > 9.2 were used for subsequent studies. Pools consisting of 1 μg of total RNA isolated from 6 individuals per condition (fasting and feeding with diets HLL, MHL, MLH, LHH and LLH; 36 fish and RNA samples in total) were used to construct liver and white skeletal muscle cDNA libraries. The dsDNA synthesis was performed using a MINT-Universal cDNA synthesis kit (Evrogen, Moscow, Russia). To increase the presence of rare transcripts, the cDNAs libraries were normalised using TRIMMER cDNA normalization kit (Evrogen, Moscow, Russia) following manufacturer’s instructions. Sequencing of cDNAs libraries was performed using a GS FLX 454 platform (Roche, Basel, Switzerland) at the CCiTUB of the Universitat de Barcelona (Barcelona, Spain).

Transcriptome assembly and annotation

Pre-processing of raw reads to remove low quality bases, primers and adapters and quality control were performed using Cutadapt [35], Prinseq [36] and FastQC (Babraham Bioinformatics, Cambridge, UK). The reads were de novo assembled using the GSAssembler software (Roche, Basel, Switzerland). Three assemblies were performed: one for sequence data from liver, another for white skeletal muscle, and an incremental assembly using both liver and white skeletal muscle data (hybrid assembly). Gene annotation of unique sequences for the three assemblies was performed with local BLASTx and BLASTn searches against non-redundant protein and nucleotide sequence databases of the NCBI’s QBLAST using GPRO software [37] with a significant threshold of E-value <1e-4. Gene names and putative functions were assigned to each sequence based on the highest alignment score among BLAST best 10 matches. GPRO software was also used to perform functional analysis by mapping annotated unique sequences to the Gene Ontology (GO) database. The number of sequences associated to GO terms was calculated under the categories of biological process, molecular function and cellular component.

Microarray design, hybridization and data analysis

An Agilent custom high-density oligonucleotide microarray (8 × 60 k; ID 079501; Agilent Technologies, Santa Clara, CA, USA) was designed to contain 2 different 60-mer probes for each of the 25,392 assembled unique sequences present in the S. aurata hybrid transcriptome. Labelling, hybridisation and scanning were performed using the Two-Color Microarray-Based Gene Expression Analysis v. 6.5 protocol (Agilent Technologies, Santa Clara, CA, USA) according to the manufacturer’s instructions. The microarray analysis was performed using total RNA isolated from liver and skeletal muscle from four fish per condition (starved fish and fish fed diets HLL, MHL and LLH). Briefly, for each sample 200 ng of total RNA was labelled with Cy3 or Cy5 using Low Input Quick Amp Labeling Kit. The Two-Color and RNA Spike-In Kit, Two-Color was used to monitor microarray workflow for linearity, sensitivity and accuracy. Labelled cRNA was purified with the RNeasy Mini Kit (Qiagen, Hilden, Germany) and quantified using a NanoDrop spectrophotometer (Thermo Scientific, Waltham, MA, USA). After assessing successful dye incorporation and sample integrity, 2.5 μg of each labelled sample was hybridised to the custom-made oligo-microarray at 65 °C for 17 h following Gene Expression Hybridization Kit instructions. A double loop hybridisation with dye swap experimental design was adopted with n = 4 fish per condition (starvation and feeding with diets HLL, MHL and LLH) and a total of 16 hybridisations (8 per tissue) [38]. Scanning was performed using an Agilent Microarray Scanner G2565BA. Outlier spots and spot intensity for Cy3 and Cy5 channels were extracted using Agilent Feature Extraction software version 10.7. Loess and Aquantile normalisation for within-and inter-array normalisation, respectively, was applied using the R-Bioconductor package [39, 40]. Data analysis was only considered for transcriptomic annotations of OXPHOS genes with E-value <1e-10 and HSP/hit >30. A linear model analysis using Limma [41] was conducted to select differentially expressed genes between conditions.

Quantitative real-time RT-PCR

Statistics

Data obtained by performing qPCR were analysed by one-way ANOVA using the SPSS software Version 22 (IBM, Armonk, NY, USA) and are represented as mean ± standard deviation. When statistical significance was found for ANOVA, the Student-Newman-Keuls post hoc test was used to determine differences among treatments.

454 pyrosequencing and assembly

To obtain a deep-coverage database useful for developing nutritional studies in S. aurata, total RNA pools were made from the liver and skeletal muscle of 23-day starved fish or fish fed with five different diets for the same time period. The experimental diet HLL (High protein, Low lipids, Low carbohydrates) has a macronutrient composition similar to the diet of wild S. aurata. The other diets presented partial substitution of protein by lipids and/or carbohydrates. Liver and white skeletal muscle cDNA libraries were constructed, normalised to increase the presence of rare transcripts, and subsequently sequenced using the 454 FLX technology. Two sequencing runs of liver and skeletal muscle libraries yielded 812,770 and 691,433 reads, respectively. After clean-up and removal of primers and adapter sequences, sequencing of liver and white skeletal muscle libraries generated 699,991 and 555,031 high-quality reads with an average length of 374 and 441 bp and a total of 262 and 245 Mbp sequence data, respectively. Three assemblies were performed: liver, white skeletal muscle and an incremental assembly by integrating data from liver and white skeletal muscle (hybrid transcriptome). A summary of the 454 assemblies is listed in Table 3. Unique reads represented 96.8 and 96.4% of total reads of liver and skeletal muscle sequencing results, respectively (Fig. 1a). Assembly of the hybrid transcriptome was performed by aligning 1,119,088 reads (89.2% of high quality reads), which generated 23,956 contigs with average and N50 lengths of 901 bp and 1322 bp, respectively, and 102,569 singletons. Contigs were grouped into 19,530 isotigs with average and N50 lengths of 1330 and 1536 bp, respectively, and isotigs were further grouped into 15,954 isogroups. For the liver and skeletal muscle, aligned reads (88.6% and 91.2% of high quality reads, respectively) generated 13,313 and 12,333 contigs with an average length of 934 bp (N50 = 1195 bp) and 1056 bp (N50 = 1364 bp), and 43,460 and 30,302 singletons, respectively. The number of isotigs constructed was 11,545 with average and N50 lengths of 1208 bp and 1315 bp, respectively, for the liver, and 10,599 isotigs with average and N50 lengths of 1390 bp and 1491 bp, respectively, for the skeletal muscle. Isogroups in the liver and skeletal muscle transcriptomes were 10,033 and 9189, respectively (Table 3). Length distribution of assembled unique sequences is shown in Fig. 1b. Liver and white skeletal muscle transcriptomes presented 76.1 and 86.4% of unique sequences above 500 bp.

	Liver	Muscle	Hybrid
Data generation
Number of reads	812,770	691,433	1,504,203
Total Megabases	298.9	275.9	574.8
Average length of reads (bp)	368	399	382
After clean-up processing
Number of high quality reads	699,991	555,031	1,255,022
Total high quality Megabases	261.6	244.7	506.3
Average length of high quality reads (bp)	374	441	404
Assembly
Number of contigs	13,313	12,333	23,956
Average length contigs	934	1056	901
N50 contig length	1195	1364	1322
Number of Isotigs	11,545	10,599	19,530
Average length of isotigs	1208	1390	1330
N50 isotig length	1315	1491	1536
Number of Isogroups	10,033	9189	15,954

Annotation and gene ontology

In the present study, 64.3%, 73.9% and 83.1% unique sequences of liver, skeletal muscle and hybrid transcriptomes, respectively, were annotated using BLASTx and BLASTn against the NCBI non-redundant protein and nucleotide collection database. To categorise gene products, 6162, 5581 and 4793 unique sequences of liver, skeletal muscle and hybrid transcriptomes, respectively, were annotated with specific GO terms. Figure 2 summarises GO terms for the assembled transcriptomes at fourth level for the categories of biological process and cellular components, and at third level for molecular function GO terms. GO terms distribution was similar for the liver, skeletal muscle and hybrid transcriptomes. Terms describing biological process were most abundant for nucleobase-containing compound metabolic process (GO:0006139; 17–19%), gene expression (GO:0010467; 13–14%), transport (GO:0006810; 11–12%), and signal transduction (GO:0007165; 10–11%) (Fig. 2a). GO terms describing molecular function were highest for nucleotide binding (GO:0000166; 12–14%), hydrolase activity (GO:0016787; 13%), protein binding (GO:0005515; 10–12%), transferase activity (GO:0016740; 12%), and nucleic acid binding (GO:0003676; 11–13%) (Fig. 2b). Most abundant GO terms describing the results for cellular components were intracellular (GO:0005622; 50–52%), membrane (GO:0016020; 37–39%), plasma membrane (GO:0005886; 3–4%), extracellular space (GO:0005615; 2%), and endomembrane system (GO:0012505; 2%) (Fig. 2c).

GO analysis was also conducted on annotated sequences that were exclusively present either in the liver transcriptome (2457 sequences out of 19,494) or in the skeletal muscle transcriptome (2174 sequences out of 17,935). GO terms in the categories of biological process (level 4) and molecular function (level 3) are shown in Fig. 3. The more remarkable biological process GO terms showing different percentage of annotated sequences in liver versus skeletal muscle were transport (GO:0006810; 23.6 vs. 20.7%), carbohydrate metabolic process (GO:0005975; 7.1 vs. 5.0%), lipid metabolic process (GO:0006629; 6.5 vs. 3.8%), defense response (GO:0006952; 4.1 vs. 0.7%), cellular amino acid metabolic process (GO:0006520; 3.6 vs. 2.1%), cell differentiation (GO:0030154; 4.5 vs. 7.6%), embryo development (GO:0009790; 1.7 vs. 4.8%), actin cytoskeleton organisation (GO:0030036; 1.5 vs. 3.1%), pattern specification process (GO:0007389; 1.1 vs. 2.1%), embryonic morphogenesis (GO:0048598; 0.9 vs. 3.1%), and mitotic cell cycle (GO:0000278; 0.6 vs. 1.9%). With regard to molecular function, GO terms with different percentage of annotated sequences in liver and skeletal muscle included transferase activity (GO:0016740; 15.8 vs. 12.3%), oxidoreductase activity (GO:0016491; 7.5 vs. 5.7%), ligase activity (GO:0016874; 4.3 vs. 2.3%), lipid binding (GO:0008289; 2.8 vs. 0.7%), carbohydrate binding (GO:0030246; 2.1 vs. 0.3%), nucleotide binding (GO:0000166; 13.0 vs. 20.0%), and nucleic acid binding (GO:0003676; 9.9 vs. 13.4%).

A BLAST top-hit species distribution of gene annotations in the hybrid S. aurata transcriptome showed the highest homology with NCBI non-redundant database sequences from Stegastes partitus (23.5%), followed by Larimichthys crocea (16.3%), Neolamprologus brichardi (7.5%), Sparus aurata (3.6%), Takifugu rubripes (3.5%), Cynoglossus semilaevis (3.5%) and Kryptolebias marmoratus (3.4%) (Fig. 4).

Microarray and gene expression analysis

With the aim of identifying new biomarkers related to OXPHOS of interest for nutritional studies in S. aurata, we used the obtained transcriptomic data for designing an oligonucleotide microarray to further analyse changes in gene expression due to long-term starvation or feeding diets differing in macronutrient composition. To this end, three groups of fish were fed 23 days on diets HLL (High protein, Low lipids, Low carbohydrates; with a macronutrient composition similar to the diet of wild S. aurata), MHL (Medium protein, High lipids, Low carbohydrates; with a composition similar to commercial diets for S. aurata culture) and LLH (Low protein, Low lipids, High carbohydrates; with partial substitution of protein by carbohydrate compared to the HLL diet). A fourth group of fish was deprived of food for the same period. Growth performance of starved fish and fish fed diets HLL, MHL and LLH, expressed as mean SGR (%/day) ± SD, was −0.89^a ± 0.54, 2.54^c ± 0.21, 2.07^bc ± 0.33, and 1.43^b ± 0.51, respectively (different superscript letters indicate significant differences, p < 0.05; n = 3 groups of 8 fish each). The expression of a total of 129 genes involved in mitochondrial respiratory chain, oxidative phosphorylation and ATP translocation across the mitochondrial inner membrane were analysed in the liver and skeletal muscle of S. aurata (Additional file 1). Despite recent efforts to increase transcriptomic data concerning genes involved in the OXPHOS pathway in S. aurata [29], 38 out of 129 genes analysed in the present study are newly described genes for S. aurata that have not been previously reported in public databases. Long-term starvation deeply affected the expression of an important part of analysed genes. Considering statistical significance with an adjusted P value <0.05, 24–38 genes (depending on the diet supplied) out of 129 were significantly upregulated in the liver of starved fish, while food restriction downregulated 10–14 genes. Differentially expressed genes in the liver included 8 upregulated and 5 downregulated subunits out of 42 genes associated to NADH:ubiquinone oxidoreductase (complex I of the electron transfer chain), 1 upregulated and 2 downregulated subunits out of 7 from succinate dehydrogenase (complex II), 4 upregulated subunits out of 14 from ubiquinol-cytochrome c reductase (complex III), 15 upregulated and 4 downregulated subunits out of 34 from cytochrome c oxidase (complex IV), 5 upregulated and 2 downregulated subunits out of 19 from F1F0-ATP synthase (OXPHOS complex V), and 2 upregulated genes out of 3 ADP/ATP translocases (Additional file 2). Food deprivation resulted in 17–40 upregulated and 21–28 downregulated OXPHOS genes in the skeletal muscle (Additional file 3).

Figure 5 shows a heat map hierarchical clustering of differentially expressed genes (fish fed with diets HLL, MHL or LLH versus starved fish) with an adjusted P value <0.05 and with a difference of at least 2-fold in the normalised intensity ratio (Cy5/Cy3 or Cy3/Cy5) for one or more feeding conditions. With a fold change greater than 2, all differentially expressed genes were upregulated by starvation in the liver of S. aurata, and 6 genes corresponded to cytochrome c oxidase subunits (COX6A2, COX5B1, COX4I2, COX7A2, COX8B and COX6B1). Indeed, 4 cytochrome c oxidase components (COX6A2, COX5B1, COX4I2 and COX7A2) showed a fold change greater than 4 regardless of the diet supplied (fold change ranging from 4.9 to 86.7). In addition to cytochrome c oxidase subunits, SLC25A6 and COQ10 were also found among upregulated genes with greater fold change expression in the liver of starved fish (SLC25A6: 44.9 to 60.1 fold change depending on the diet supplied; and COQ10: 14.4 to 19.6 fold change). In contrast, starvation resulted in both upregulated and downregulated genes, 18 and 10 respectively, with a fold change greater than 2 in the skeletal muscle. The magnitude of expressional changes due to starvation was lower in the skeletal muscle than in the liver. No genes in the skeletal muscle showed a fold change greater than 4 when compared starved fish versus fish fed diets HLL, MHL and LLH. Remarkably, diet composition affected the expression of OXPHOS genes in the skeletal muscle to a greater extent than in the liver. As a general trend, fish fed LLH exhibited higher mRNA levels than fish fed HLL for the majority of significantly upregulated genes in the skeletal muscle of starved fish (Additional file 3).

To validate microarray data, the mRNA levels of several genes involved in OXPHOS were determined by RT-qPCR in liver and skeletal muscle samples from S. aurata submitted to long-term starvation or fed diets differing in macronutrient composition. The mRNA abundance of eigth genes involved in mitochondrial electron transfer chain (NDUFB8, NDUFS1, COQ10, UQCR11A, COX5B1 and COX6A2), ATP synthesis (ATP5B) and ADP/ATP translocation (SLC25A6) were analysed. The mRNA levels of NDUFB8 and NDUFS1, components of NADH:ubiquinone oxidoreductase, were not affected by the nutritional condition in the liver and the skeletal muscle. (Fig. 6a and b). The expression of COQ10, a protein required for coenzyme Q activity in the electron transport chain [42], was also determined in tissue samples of treated fish. In the liver, COQ10 mRNA significantly increased after starvation, reaching values 17-fold higher than fish fed diet HLL, and 9.3 and 10.2-fold higher than fish fed diets LLH and MHL, respectively. A different regulation of COQ10 expression was observed in the skeletal muscle by the nutritional condition: a significant 3-fold decrease was found in starved fish compared to those fed a high protein/low carbohydrate diet (HLL) (Fig. 6c).

The mRNA abundance of UQCR11A and COX5B1, components of ubiquinol-cytochrome c reductase and cytochrome c oxidase, respectively, of the respiratory chain, showed a similar behaviour as a result of starvation or feeding different diets. Starvation significantly increased the expression of both proteins in the liver (UQCR11A: 1.8 and 1.9-fold compared to fish fed HLL and MHL, respectively; and COX5B1: 3.2 to 4.5-fold compared to fed fish regardless of the diet) and the skeletal muscle (1.8 to 2.0-fold for UQCR11A and 2.4 to 9.2-fold for COX5B1 compared to fed fish regardless of the diet). Diet composition did not affect UQCR11A expression, while COX5B1 mRNA levels significantly increased 2.9 and 3.9-fold in the skeletal muscle of fish fed diet LLH compared to those fed HLL and MHL, respectively (Fig. 6d and e).

Another component of cytochrome c oxidase, COX6A2, was also greatly affected by food deprivation in the liver. Starvation significantly increased the hepatic mRNA levels of COX6A2 23.1 to 34.6-fold, depending on the diet supplied. In contrast, COX6A2 expression was not affected by the nutritional condition in the skeletal muscle (Fig. 6f). Neither starvation nor diet composition significantly affected the expression of ATP5B, a subunit of the catalytical core of F1F0-ATP synthase, in the skeletal muscle (Fig. 6g). SLC25A6, also known as ADP/ATP translocase 3, is a member of the solute carrier family 25 that facilitates the exchange of ADP and ATP across the mitochondrial inner membrane by an antiport mechanism [43]. In the liver, starvation significantly increased 14.3 to 19.3-fold the expression of SLC25A6 compared to fed fish. The effect of starvation depended on the diet composition in the skeletal muscle: fish fed diet HLL presented the lower SLC25A6 mRNA levels, which were 7.6 and 8.2-fold higher in fish fed diet MHL and starved fish, respectively (Fig. 6h).

Transcriptomic analysis

The goal of the present study was to increase currently available transcriptomic data for S. aurata and identify OXPHOS genes affected by the nutritional status in this species. 454 sequencing of liver and white skeletal muscle cDNA normalised libraries from S. aurata juveniles submitted to starvation or feeding on five selected diets allowed us to generate transcriptomes for liver, skeletal muscle and a hybrid transcriptome integrating data from both tissues. The number of isotigs (19,530), contigs (23,956) and average isotig and contig length (1330 bp and 901 bp, respectively) of the hybrid transcriptome indicates that effective data assemblies were carried out compared to values reported for other fish species of interest in aquaculture [44]. In contrast to previous studies, prior to 454 pyrosequencing we constructed normalised libraries using tissue samples from long-term starved fish and fish fed five diets differing in macronutrient composition to allow maximal representation of transcripts abundance and diversity in regard to dietary condition. Furthermore, the present study reports the first liver transcriptome for S. aurata obtained using NGS technology. In regard of the skeletal muscle, assembly of five cDNA libraries from the fast skeletal muscle of fed and short-term fasted fish submitted to various rearing temperatures allowed Garcia de la Serrana et al. [21] to identify 50,515 isotigs with an average length of 454 bp (11,545 isotigs with an average length of 1208 bp in the present study), while Calduch-Giner et al. [20] obtained 7808 contigs with an average length of 968 bp for the same tissue (12,333 contigs with an average length of 1056 bp in the present study). The length of assembled sequences of the present study allowed annotation of a significant part of unique sequences, improving transcriptome depth.

Functional annotation allowed us to analyse frequency of GO terms in assembled transcriptomes and within annotated sequences exclusively expressed in the liver or the skeletal muscle. Concerning GO terms, liver, skeletal muscle and hybrid transcriptomes exhibited an overall similar functional profiling. However, after filtering common sequences present in the liver and skeletal muscle transcriptomes, marked differences in GO term frequency were found between both transcriptomes for biological process and molecular function categories. Consistent with the plasticity and complexity of the hepatic metabolism, the liver transcriptome presented a higher number of GO terms describing biological processes such as transport, carbohydrate metabolic process, lipid metabolic process, defense response and cellular amino acid metabolic process, and molecular functions related to enzyme activity. During growth skeletal muscle develops an adaptive response that includes tissue repair, cytoskeleton organization and remodeling processes [45]. Consistent with a higher growth rate in S. aurata juveniles than in adults [46], GO terms associated with biological processes involved in development and morphogenesis, such as cell differentiation, embryo development and actin cytoskeleton organisation, were more abundant in the skeletal muscle transcriptome from S. aurata juveniles.

BLAST similarity analysis of unique sequences in the S. aurata hybrid transcriptome against NCBI non-redundant database showed higher homology to genes from ten fish species belonging to the Percomorphaceae subdivision (68.3% of annotations), as it is the case for S. aurata [47]. The low number of annotations matching sequences from S. aurata highlights the limited genomic information available for this species in the NCBI database.

Nutritional regulation of OXPHOS in S. aurata

Transcriptomic data obtained in the present study was further used to identify biomarker genes of interest for nutritional studies by analysing the effect of nutritional status and diet composition on the expression of genes involved in OXPHOS and ATP translocation across the mitochondrial inner membrane in S. aurata. Given that starvation periods can induce compensatory growth in commercial fish farming [48], we studied the effect of long-term starvation on the expression of OXPHOS genes in S. aurata. Microarray analysis and validation of microarray data by qPCR pointed to a marked increase in the mRNA levels of a high proportion of OXPHOS genes as a result of long-term starvation in the liver of S. aurata. In the skeletal muscle, starvation upregulated the majority of OXPHOS genes displaying changes in expression >2-fold change, although at a considerably smaller magnitude than in the liver. In contrast to the liver, food deprivation downregulated also a significant amount of OXPHOS genes in the skeletal muscle of S. aurata. Thus, our findings indicate that the OXPHOS process in the liver is more sensitive to long-term starvation. In support of this hypothesis, mitochondrial-related changes in fasted mice are also greater in the liver than in the skeletal muscle [49]. Indeed, it was reported that fasting does not affect the expression of OXPHOS genes in the skeletal muscle of humans [50].

Bearing in mind that cytochrome c oxidase is considered to be the rate-limiting reaction of the electron transfer chain [23, 24], the fact that starvation markedly increased the expression of multiple components of cytochrome c oxidase is consistent with upregulation of OXPHOS in the liver of food-deprived S. aurata (15 out of 35 genes significantly upregulated with an adjusted P value <0.05; 6 of them with a fold change >2). Our findings are in agreement with previous evidences highlighting increased OXPHOS activity as an important physiological adaptation in the liver during starvation. Similarly as in S. aurata, starvation increases cytochrome c oxidase activity and the transcription of OXPHOS genes in the liver of mice [26, 27], and enhances the activity of cytochrome c oxidase in human fibroblasts [28]. Moreover, a microarray analysis conducted on liver samples of mice indicated that 75% of diet restriction results in a modest but significant increased expression of an important number of OXPHOS genes [51]. In contrast to our findings, a recent report addressing the effect of 10 days of starvation on the expression of 88 genes of the OXPHOS pathway in S. aurata juveniles concluded that food-deprivation mostly downregulated the expression of OXPHOS genes in the liver, while starvation upregulated a number of OXPHOS genes in the white skeletal muscle and cardiac muscle, in most cases with a fold change <2 [29]. The starvation period (10 days versus 23 days in the present study) may explain striking differences among studies. Indeed, long periods of starvation are frequent during the life cycle of many fish species as a result of migration, reproduction and food availability [52]. Hence, metabolic adaptation to food deprivation includes a gradual depletion of hepatic glycogen levels in S. aurata, and starvation periods longer than 18 days are required to decrease the liver glycogen content to barely detectable levels [53, 54]. Therefore, it is not surprising that the expression of OXPHOS genes in S. aurata tissues could exhibit different patterns depending on the starvation period. In agreement with this hypothesis, three weeks of starvation were necessary to promote a significant increase in the protein content of succinate dehydrogenase, cytochrome c oxidase and F1F0-ATP synthase in the skeletal muscle of fine flounder (Paralichthys adspersus) [55]. Thus, stimulation of OXPHOS in the liver of S. aurata by long-term starvation may result from a metabolic adaptation that would involve increased lipolysis, enhanced hepatic fatty acid oxidation and ketogenesis, and reduced glucose uptake and oxidation in peripheral tissues. Consistently, starvation stimulates ATP-dependent processes in the liver, such as gluconeogenesis and ureagenesis, and β-oxidation of fatty acids to provide ketogenic substrates [56]. Therefore, increased OXPHOS activity may be essential to facilitate substrate oxidation and supply ATP in the liver of long-term starved S. aurata.

Among cytochrome c oxidase components significantly upregulated in the liver of starved S. aurata, COX6A2 and COX5B1 were more sensitive to food deprivation. Both components may exert an important role in cytochrome c oxidase activity during starvation. In fact, the expression of COX5B, which enhances cytochrome c oxidase activity [57], is repressed by high concentrations of glucose in Saccharomyces cerevisiae [58]. The hypothesis that COX6A2 expression may be relevant under food deprivation is supported by the fact that COX6A2^−/− mice are very sensitive to food deprivation and loose more weight during starvation than wild type mice [59]. Consistent with an increased expression of cytochrome c oxidase, the mRNA levels of cytochrome c (CYCS), a small one-electron carrier that shuttles electrons from ubiquinol-cytochrome c reductase to cytochrome c oxidase, also increased markedly in the liver of starved fish.

In addition to CYCS and cytochrome c oxidase components, upregulation of two critical proteins in mitochondrial respiration such as COQ10 and SLC25A6 closely fits with an increased OXPHOS function in the liver of starved S. aurata. COQ10 is a member of the steroidogenic acute regulatory protein (StAR)-related lipid transfer (START) domain superfamily that is located in the mitochondrial inner membrane, where binds ubiquinone in its hydrophobic pocket [60]. Deletion of COQ10 gene causes respiratory deficiency and the inability to oxidise NADH and succinate in yeast, while overexpression of COQ10 restores respiratory electron transport in the yeast COQ10 null mutant [42, 61–63]. These observations led to hypothesise that COQ10 is essential for mitochondrial respiration by facilitating ubiquinone biosynthesis and acting as a chaperone for the transport of ubiquinone between electron transfer chain complexes [42, 60, 63]. SLC25A6 is an ADP/ATP translocase that is ubiquitously expressed at levels that are proportional to the respiratory activity of the tissue [64]. As a core component of the mitochondrial permeability transition pore, SLC25A6 can strongly induce apoptosis in cultured human-derived cells [65, 66]. ADP/ATP translocases are integral proteins that supply cellular energy by coupling mitochondrial respiration with ADP/ATP exchange across the mitochondrial inner membrane. Presence of different ADP/ATP translocases in eukaryotic organisms may enable control of expression levels in response to a variety of stimulus and energy requirements, rather than displaying functional differences [43]. Thus, since long-term starvation caused modest effects on other ADP/ATP translocases (SLC25A4 and SLC25A5) in the liver of S. aurata, our findings suggest that increased mRNA levels of SLC25A6 may be critical to support ATP-dependent processes in the liver during starvation. The fact that starvation markedly upregulated COQ10, COX6A2 and SLC25A6 in the liver, while they were barely affected in the skeletal muscle, indicates that the hepatic expression of COQ10, COX6A2 and SLC25A6 can therefore be used as sensitive markers of the nutritional status in S. aurata.

In contrast with sensitivity of the expression of OXPHOS genes in the liver of long-term starved S. aurata, diet composition did not cause significant changes in the expression pattern of OXPHOS genes. Nevertheless, the expression of OXPHOS genes in the skeletal muscle was more affected by macronutrient dietary composition. Compared to fish fed HLL, the supply of a low-protein/high-carbohydrate diet (LLH) promoted a trend to upregulate mRNA levels of most of the genes whose expression significantly increased in the skeletal muscle of starved fish, such as COX5B1. Given that macronutrient composition of LLH is far from the dietary composition of wild S. aurata and resulted in the lowest SGR values among fed fish, our findings are consistent with the fact that partial substitution of dietary protein by carbohydrates reduces weight gain and results in expression levels of appetite-regulating peptides close to the values observed in starved S. aurata [67].