Ocean warming combined with lower omega-3 nutritional availability impairs the cardio-respiratory function of a marine fish

The highly unsaturated fatty acids of the n-3 series (also called long-chain omega-3 HUFA, and simply HUFA through this paper) are critical structural and functional components of the cell membranes of most animal tissues, making them vitally important for numerous key physiological functions in aquatic and terrestrial organisms. HUFA are, however, weakly synthesised de novo by animals, and must be obtained through food (Glencross, 2009; Tocher, 2010). In natural food webs, they are mainly supplied by marine microalgae (Crawford and Broadhurst, 2012; Pahl et al., 2010). This HUFA production, sitting at the base of aquatic food webs, determines HUFA availability for the upper trophic levels, as the fatty acid (FA) content in the tissues of animals is a reflection of their diet (Arts et al., 2001; Crawford and Broadhurst, 2012; Gladyshev et al., 2009). In marine fish, membranous HUFA are mainly represented by ecosapentaenoic acid (EPA, 20:5n-3) and docosahexaenoic acid (DHA, 22:6n-3; Glencross, 2009; Tocher, 2010), mostly supplied by diatoms. Fish constitute one of the main vectors of HUFA from aquatic ecosystems to human populations (Arts et al., 2001).

Global change, through warming, acidification or hypoxia, is leading to the lowered production of HUFA at the base of food webs, due to the alteration of the FA biosynthesis pathways and species assemblages of microalgae (Galloway and Winder, 2015; Hixson and Arts, 2016). Reductions of 8.2% for EPA and 27.8% for DHA are expected as a consequence of a temperature increase of 2.5°C (Hixson and Arts, 2016). This lowered HUFA production may decrease the HUFA availability for higher trophic levels (Gladyshev et al., 2009, 2011; Kang, 2011). For example, a 1°C increase in sea surface temperature in the SW Pacific Ocean is expected to decrease EPA and DHA content in tuna by 3% and 1.5%, respectively, as a result of a decreasing biomass of HUFA-rich diatoms in favour of HUFA-poor dinoflagellates (Pethybridge et al., 2015). Although a few studies have described the ability of marine fish to stimulate, at the molecular level, the FA desaturation pathways in response to a HUFA dietary depletion, this process has never been shown to compensate for HUFA depletion in tissues (Glencross, 2009; Sargent et al., 1997; Tocher, 2010; Vagner et al., 2007, 2009).

Such changes in the HUFA composition of consumers may have cascading detrimental effects throughout the world's ecosystems, affecting the physiology of many species, and consequently their role within ecosystems, and culminating in an overall decline in the global availability of these nutrients for human health (Arts et al., 2001; Gladyshev et al., 2009; Hixson and Arts, 2016). Despite this, the various consequences of changing HUFA availability remain largely understudied.

Previous studies have reported that the degree of lipid unsaturation in cell membranes is correlated with the activity of membrane-bound enzymes and rates of transmembrane proton leakage, with consequences for enzyme activities and whole-animal metabolism (Hulbert and Else, 2000; Kraffe et al., 2007; Martin et al., 2012). For example, Wagner et al. (2004) suggested that a low dietary HUFA/saturated FA ratio caused a decline in Atlantic salmon swimming performance tested at 10–12°C; however, the opposite relationship was found in Atlantic salmon exposed to 9°C (McKenzie et al., 1998). These contrasting results for the same species suggested potential interactions between HUFA food content and acclimation temperature on fish physiology. This interaction was recently demonstrated in a species highly tolerant to temperature variations, the golden grey mullet, Chelon auratus: the detrimental effects of a low HUFA diet were greater in warmer conditions (20°C versus 12°C), resulting in greater energy expenditure to reach a given swimming speed (Vagner et al., 2015). Highlighting these combined effects is of particular importance given the strong seasonal temperature variations that occur in temperate coastal areas, and the predicted scenario of a 3°C increase in ocean water temperature by 2100 (IPCC, 2013).

The underlying mechanisms governing the metabolic response observed by Vagner et al. (2015) have never been explored, but they are key to further understanding and predicting fish responses to the nutritional and thermal challenges induced under climate change. In the present study, we hypothesised that temperature and HUFA food content affect fish metabolism through a modulation of cardiac function, which is key to supplying the required amount of oxygen and nutrients to tissues (Claireaux, 2005; Clark et al., 2005). HUFA depletion in cardiomyocyte membranes combined with a warmer environment could modify the conformation of the ionic pumps and channels embedded in cell membranes (Raynard and Cossins, 1991; Turner et al., 2003), with consequences for the cardiomyocyte contraction capacity. This contraction is normally initiated after depolarisation of the cell, and graded by an increase in the concentration of free intracellular calcium (Ca²⁺), supplied through: (1) an external Ca²⁺ influx across the sarcolemma, via either L-type Ca²⁺ channels or a Na⁺/Ca²⁺ exchanger, and (2) internal Ca²⁺ mobilisation from storage in the sarcoplasmic reticulum (SR). In fish, Ca²⁺ influx derives mainly from L-type Ca²⁺ channel activity (Chatelier et al., 2006; Driedzic and Gesser, 1994; Vornanen, Shiels and Farrell, 2002). Relaxation is then rapidly regulated by a decrease in free cytoplasmic Ca²⁺ largely governed via the Na⁺/Ca²⁺ exchanger (in reverse mode), but also via Ca²⁺-ATPase pumps (SERCA), and sometimes via sarcolemmal (SL) Ca²⁺-ATPase pumps (PMCA).

In this context, this study evaluates the cardiac response of fish as a first step towards understanding the combined effects of environmental temperature and lowering HUFA food content, by integrating responses at molecular (gene expression), tissue (cardiac performance) and individual (oxygen consumption during swimming effort) levels. This work aims to improve our understanding of the integrative response of fish to environmental changes, which is a challenge burdened with ecological and socio-economic consequences, as well as public health concerns.

We used golden grey mullets, Chelon auratus (Risso 1810), from the same experimental groups as those used by Vagner et al. (2015): they were exposed to two highly contrasting HUFA diets (a high and a low HUFA diet), chosen according to the needs defined for several fish species, and two contrasting environmental acclimation temperatures that are representative, respectively, of mean winter (12°C) and mean summer (20°C) temperatures in the Pertuis Charentais area (Rochebonne, France). Their cardio-respiratory performance was then measured at the individual level in terms of metabolic rates and associated swimming effort. The physiological mechanisms underlying this response were studied: (1) on isolated hearts in order to evaluate cardiac performance in terms of contraction and relaxation capacity and the relative involvement of SL channels and/or exchangers and the SR in excitation–contraction coupling, and (2) at the transcriptional level to evaluate differences in gene expression in heart tissue according to the experimental conditions tested.

All fish manipulations were performed according to the Animal Care Committee of France (ACCF). No specific permissions were required from the Departmental Service of Fisheries to collect C. auratus in their natural environment, as this is not a protected or endangered species, and as a relatively small number of fish were collected for this research (n=160 individuals) compared with their abundance in this fishing area. The protocol was approved by the ACCF (approval number: 17-300-2). All fish manipulations were performed under anaesthesia (tricaine methane sulphonate MS-222; 0.1 g l⁻¹, Sigma-Aldrich, St Quentin-Fallavier, France), and all efforts were made to minimise suffering.

Juvenile golden grey mullets (initial mean±s.e.m. mass: 6.1±0.2 g; initial mean±s.e.m. standard length: 7.0±0.1 cm) were caught in the marine coastal area of La Rochelle (France) in 2012 and transported in an aerated tank to the laboratory, where all experiments were conducted. Upon arrival, fish were transferred into four indoor tanks (volume: 400 l; n=40 fish per tank) that were individually supplied with aerated recirculated sand-filtered natural seawater and equipped with an external biological filter (Eheim professional 3 2080, Eheim, Deizisau, Germany). Fish were progressively acclimated to the water tank temperature, which was kept constant by a recirculating water system (TECO TR20, Conselice, Italy), and maintained in a temperature-controlled room (20°C) exposed to a 12 h light:12 h dark photoperiod cycle. Temperature (18.7±0.04°C), salinity (29.3±0.5) and oxygen (87.9±2.6% air saturation) were monitored daily using a conductimeter (WTW model Oxi 340i, Weildeim, Germany). After a few days of acclimation to the experimental structure, each of the four initial tanks was separated in two tanks (n=8 tanks of 400 l each with n=20 fish per tank), and water temperature was progressively modified to reach 12°C (mimicking winter temperature in a temperate coastal area: 12.5±0.1°C; n=4 tanks) or 20°C (mimicking summer temperature in a temperate coastal area: 19.8±0.1°C; n=4 tanks). Fish were initially fed with a commercial diet (Le Gouessant^® aquaculture, Lamballe, France) once a day for 1 month.

They were then progressively fed with two experimental isolipidic and isoproteic diets differing in their HUFA content (the experimental diet to commercial diet ratio was increased every 2 days from 25%:75% to 50%:50%, 75%:25% and finally to 100%:0%). The two experimental diets were formulated and made at the Ifremer-PFOM unit, UMR 6539 LEMAR (Plouzané, France): a high HUFA diet (HH: 4.8% EPA+DHA on a dry matter basis), and a low HUFA diet (LH: 0.2% EPA+DHA on a dry matter basis). The LH diet was obtained by replacing the fish oil present in the HH diet with soybean oil. Diets were formulated taking a standard diet as the reference (1.2% EPA+DHA on a dry matter basis), as previously formulated for the same species (Vagner et al., 2014). The composition and FA content of the two experimental HH and LH diets are summarised in Vagner et al. (2015). For each of these two diets, the fish in two tanks per experimental temperature were fed for at least 4 months (2% of biomass per day) to ensure a good incorporation of dietary lipids into their cell membranes (Vagner et al., 2007, 2014, 2015). The four experimental conditions were created in duplicate and called HH20, LH20, HH12 and LH12.

Fish swimming and metabolic performance were assessed using a swim-tunnel respirometer (Loligo Systems, Tjele, Denmark) in June 2013 for the 20°C-acclimated individuals (i.e. after 4 months of exposure to the experimental diet) and in September 2013 for the 12°C-acclimated individuals (i.e. after 7 months of exposure to the experimental diet). This lag allowed the fish to at least double their body mass from the beginning of the feeding period, ensuring an accurate oxygen consumption signal detection in the swim-tunnel respirometer.

The experimental setup of the swim-tunnel respirometer was the same as described in Vagner et al. (2014, 2015). Briefly, the respirometer (volume: 10 l) was composed of: (i) a swim chamber with a square working section (40 cm length, 10 cm height, 10 cm width) and (ii) a hydraulic system placed upstream to promote a laminar flow in the swim chamber. No correction was made for the solid blocking effects of the fish in the working section, as the calculated fractional error was <5% of the working section area (Webb, 1975). The flow in the respirometer was generated by an electric motor with a propeller. It was calibrated before the start of experiments and the speed ranged between 0 and 150 cm s⁻¹. A flush pump allowed water exchange between the respirometer and the external bath, in which water temperature and oxygenation were controlled. Temperature was kept constant (12 or 20°C) by a recirculating water system from the external bath (TECO TR20).

[O₂] versus time was plotted, and for each Ṁ_O2 measurement, a linear regression was adjusted (Graphical Analysis 3.4, Beaverton, OR, USA). This allowed the determination of Δ[O₂]/Δt and a regression coefficient illustrating Ṁ_O2 measurement accuracy. Only those measurements with a regression coefficient above 0.9 were considered. The bacterial Ṁ_O2 was measured for half an hour before and after each experiment, and the mean of both was subtracted from the Ṁ_O2 measured.

Fish starved for 96 h (HH20 group: mean±s.e.m. mass: 36.7±0.6 g, standard length: 12.72±0.3 cm, n=14; LH20 group: 33.1±2.4 g, 12.3±0.3 cm, n=13; HH12 group: 11.6±0.6 g, 8.9±0.2 cm, n=14; LH12 group: 12.3±0.8 g, 8.9±0.2 cm, n=12) were randomly sampled in tanks using a net, and individually tested in the swim-tunnel respirometer at their acclimation temperature. Fish were subjected to a step-protocol based on progressive swimming speed increments following the protocol fully described in Vagner et al. (2014). At each swimming step, the velocity was maintained for 20 min and Ṁ_O2 was measured. Oxygen saturation never fell below 75% of air saturation during Ṁ_O2 measurements. The water in the respirometer was renewed through a flush pump during the swimming speed increment phase (i.e. between two consecutive swimming steps). This allowed the oxygen saturation to return to >85% of air saturation. The speed increments were repeated until fish exhaustion. Fish were considered to be exhausted when: (1) they were maximally laterally bent on the grid placed at the rear of the swim chamber, forming a C-shape, or (2) they were unable to swim away from this grid for more than 10 s of stimulation through the tunnel window at the rear of the swim chamber. The stimulation involved lighting (with a 150 W spotlight) the rear part of the chamber, which created a high contrast with the shadow area at the front.

Muscle was collected for lipid analysis above the lateral line in the caudal part of the fish (Martinez, 2003). It was immediately frozen in liquid nitrogen and stored at −80°C for further lipid analysis. The heart was stored in liquid nitrogen before being transferred to a −80°C freezer for storage pending transcriptomic analyses.

The following parameters were calculated for each fish: U_crit, standard metabolic rate (SMR), active metabolic rate (MMR), aerobic scope (AS) and net cost of transport (NCOT) per unit distance.

MMR (mg O₂ kg⁻¹ h⁻¹) was considered as the highest Ṁ_O2 recorded during the U_crit test.

AS (mg O₂ kg⁻¹ h⁻¹) was determined for each fish as the difference between MMR and SMR, and represents the aerobic scope of the individual to perform its energy-demanding activities (Brett, 1964).

The lipid analysis methodology and the resulting data used in the present study were first reported in Vagner et al. (2015). Briefly, lipid analysis was performed on fish muscle at the beginning (n=10 fish as control; t₀) and after the period of exposure to the experimental diet (n=13 fish per experimental condition). Total lipid was extracted and fractionated into neutral (lipids of reserves) and polar lipids (lipids in membranes). FA analyses were performed on both of these lipid classes by gas chromatography. The muscle tissue was chosen because it is involved in anaerobic/aerobic metabolism during swimming performance. It was always taken on the right side, above the lateral line and below the spinal column, just behind the operculum. For technical reasons, FA analysis was not performed on cardiac tissue (the amount of tissue was not sufficient to perform both lipid and transcriptomic analyses), but numerous mammalian and fish studies have established that the HUFA composition of cardiac membrane phospholipids is directly influenced by dietary fat intake, and similar levels of HUFA were measured in muscle and hearts (Bell et al., 1993; Bell et al., 1991; Lopez-Jimena et al., 2015; Owen et al., 2004; Thomassen and Røsjø, 1989). The results of individual FA composition are expressed as a percentage of total identified fatty acid methyl esters (FAME).

Chemical analyses of the experimental diets were performed in duplicate for each sample and both the method and the resulting data are also reported in Vagner et al. (2015).

Another set of 24 h-starved fish (HH20 group: mean±s.e.m. mass: 49.3±2.2 g, standard length: 13.9±0.2 cm, n=15; LH20 group: 50.3±2.6 g, 14.1±0.2 cm, n=17; HH12 group: 14.3±0.7 g, 9.4±0.2 cm, n=20; LH12 group: 14.4±1.1 g, 9.4±0.3 cm, n=7) were randomly sampled in September 2013 (i.e. after 7 months of exposure to the experimental diet) in tanks using a net, and individually anaesthetised (MS-222; 0.1 g l⁻¹). The preparation of cardiac tissue for measuring isometric tension and adrenergic sensitivity was as detailed in Shiels and Farrell (1997). Briefly, the heart was rapidly excised and placed in physiological saline (pH 7.8 with Tris base). For each heart (HH20 group: mean±s.e.m. heart mass: 41.8±2.3 mg; LH20 group: 45.0±1.9 mg; HH12 group: 20.8±1.0 mg; LH12 group: 21.4±1.3 mg), two 1 mm-wide ventricular strips were dissected out and tied at either end with a spring and stretched. Four strips (i.e. n=2 per fish) were tested simultaneously. Each strip was lowered into a water-jacketed organ bath containing 20 ml of oxygenated physiological saline thermoregulated at 20 or 12°C (Ministat, Huber, Edison, NJ, USA) corresponding to the temperature at which hearts of the 20°C- and 12°C-acclimated fish were tested, respectively. They were stretched until the tension reached a peak and were then allowed to equilibrate for 30 min under basal stimulation (0.2 Hz). Strips were then subjected to the protocols described below, during which they were triggered using a stimulator and two flattened platinum electrodes. Signals from the isometric force transducer (MLT0201/184, ADInstruments, Dunedin, New Zealand) were amplified (Power Lab, ADInstruments) and displayed on a chart recorder (LabChart 7 pro, ADInstruments).

Adrenaline solution (10 µmol l⁻¹, Sigma-Aldrich) was used to investigate the adrenergic sensitivity of ventricular tissue. As a β-type adrenergic receptor agonist, adrenaline increased the density of the L-type calcium current by elevating the open-probability of the L-type calcium channels (Hove-Madsen, 1992; Vornanen, 1996, 1998) after phosphorylation (Shiels et al., 1998). This increases the force of contraction. Adrenaline also stimulates the Ca²⁺-ATPase pump of the SR, accelerating relaxation. As adrenaline is naturally present at nanomolar concentrations in resting fish (Gamperl, 2004), a low adrenaline concentration (10 nmol l⁻¹) was used as the control for the subsequent tests on (1) the effects of maximal adrenergic solution (10 µmol l⁻¹ adrenaline) and (2) the effects of ryanodine (see below). The adrenaline doses used in this study (10 nmol l⁻¹ and 10 µmol l⁻¹) were, respectively, the minimum and maximum values of the dose–response curve (Keen et al., 1994). Under conditions of increased energetic demand, such as stress, the adrenaline concentration in circulation can indeed reach 10 µmol l⁻¹ (Mc Donald and Milligan, 1992). The control solution of adrenaline (10 nmol l⁻¹) was added to each of the four baths. After a 10 min equilibration period, the four ventricular strips were exposed to the force–frequency trials where the pacing frequency was successively 0.2 (i.e. 12 beats min⁻¹), 0.5, 0.8, 1, 1.6 and 2 Hz (i.e. 120 beats min⁻¹). At the end of the trials, the frequency was returned to the basal level of 0.2 Hz. The adrenergic sensitivity was then evaluated by adding a saline physiological solution containing adrenaline (10 µmol l⁻¹) to two of the four ventricle strips. The two other strips were used as controls, by adding the tonic concentration of adrenaline (10 nmol l⁻¹) in order to account for muscle fatigue when adrenaline (10 µmol l⁻¹) was added (Shiels and Farrell, 2000). After a 10 min equilibration period, the strips were exposed to the different stimulation frequencies listed above (from 0.2 to 2 Hz).

To test the relative contribution of the SR to the isometric contraction of the ventricle, 10 µmol l⁻¹ of ryanodine (Shiels and Farrell, 1997) was added to each of the four baths containing the two ventricular strips pre-treated either with the tonic control adrenaline solution (10 nmol l⁻¹) or with 10 µmol l⁻¹ adrenaline. After a 30 min period of stabilisation, the strips were exposed to successive increasing stimulations from 0.2 to 2 Hz. During the stabilisation period, ryanodine is assumed to be bound irreversibly to the Ca²⁺ channel on the SR (Rousseau et al., 1987), rendering it ineffective and unable to contribute to force development.

The data were analysed using the software LabChart 7 pro (ADInstruments). The contraction rate (mN mm⁻² min⁻¹), half-relaxation rate (mN mm⁻² min⁻¹), time to peak tension (min) and peak tension (mN mm⁻²) were calculated at each pacing frequency. Pumping capacity (mN mm⁻² min⁻¹), the product of peak tension and heart frequency responses, served as an index of the heart power output (Matikainen and Vornanen, 1992).

Total RNA was isolated from heart tissues stored at −80°C (see ‘Fish sampling after U_crit test’, above) using the RNeasy Plus Mini Kit (Qiagen, Hilden, Germany) containing a genomic DNA eliminator solution, following the manufacturer's protocol. To obtain enough RNA for library preparation, at least 10 µg of heart tissue from two to three individuals exposed to the same experimental condition were pooled. A total of five RNA extracts (sample replicates) per experimental condition (LH20, LH12, HH20, HH12) were obtained. RNA integrity was estimated using an Agilent 2100 Bioanalyzer with Agilent RNA 6000 nanochips and reagents (Agilent Technologies, Les Ulis, France). RNA integrity number ranged from 6.6 to 8.7 (mean±s.d. 7.65±0.6), indicating the high and consistent quality of RNA samples. RNA concentration was measured using a NanoDrop ND-1000 spectrophotometer (Wilmington, DE, USA).

Libraries were constructed at the Dyneco lab (Ifremer, Plouzané, France) based on 500 ng of total RNA and using the TruSeq stranded RNA kit (Illumina^®, San Diego, CA, USA) with reagent quantities adapted by Dr M. Le Gac (Dyneco) from the GeT-PlaGe (Toulouse, France) based on Illumina (Illumina^®) protocols. mRNA was isolated, and the poly-A tails purified using poly-T oligo-attached magnetic beads via two rounds of purification. During the second elution of the poly-A mRNA, the mRNA was fragmented and primed for cDNA synthesis. The cleaved mRNA fragments were then copied into first-strand cDNA using reverse transcriptase and random hexamer primers. This step was followed by second-strand cDNA synthesis, which removed the RNA template and synthesised a replacement strand, incorporating dUTP in place of dTTP to generate double-strand cDNA. The incorporation of dUTP quenches the second strand during amplification, because the polymerase does not incorporate past this nucleotide. AMPure XP beads were used to separate the double-strand cDNA from the second strand reaction mix. At the end of the process, the cDNA was blunt-ended. Then, for each sample, a unique mRNA simple adapter index was ligated at the end of the cDNA fragments. The libraries obtained were then amplified by PCR performed on a Peqstar Peqlab thermocycler (Peqlab Biotechnologie GmbH, Erlangen, Germany) with a PCR primer cocktail that anneals to the ends of the adapters, and using the following protocol: 98°C for 30 s; 16 cycles of 98°C for 10 s, 60°C for 30 s, 72°C for 30 s; 72°C for 5 min. The number of cycles was minimised to avoid skewing the representation of the library. After cleaning of PCR products, library quality was assessed on an Agilent 2100 Bioanalyzer using Agilent high-sensitivity DNA analysis chips and reagents (Agilent Technologies). Libraries were quantified using the Kapa Library Quantification Kit (Kapa Biosystems, Wilmington, MA, USA) with ROX Low qPCR master mix, on a QuantStudio 3 thermocycler (Applied Biosystems, Foster city, CA, USA) according to the manufacturer's instructions: 95°C for 5 min, 35 cycles of 95°C for 30 s, 60°C for 45 s.

The 20-indexed cDNA libraries (mean±s.d. insert length 280.9±38 bp, range 240–373 bp) were normalised at 3 nmol l⁻¹ and pooled in equal volume (5 µl) before being sent for sequencing (2×150 bp) on one Illumina HiSeq3000 lane at the GeT-PlaGe France Genomics Sequencing Platform (Toulouse, France).

Read quality was assessed using FastQC v.0.11.5 after each step of the pipeline until assembly. Raw data were trimmed using Trimmomatic v.0.36 (Bolger et al., 2014) and de-duplicated (exact duplicate reads were removed) using the Dedupe function within Geneious v.10.1.3 (Kearse et al., 2012) on the default settings. Reads were filtered using Deconseq v.0.4.3 (Schmieder and Edwards, 2011) to remove potential contaminant sequences from rRNA of eukaryotes, viruses, bacteria and archaea (SILVA database; Quast et al., 2012) and from human transcripts (Deconseq database). The de novo transcriptome assembly was performed with Trinity v.2.4.0 (Grabherr et al., 2011) on the default settings.

Assembly quality was evaluated using four criteria (as recommended in Honaas et al., 2016): (1) the proportion of reads mapping back to the assembly (determined using Bowtie v. 2.1.0; Langmead and Salzberg, 2012), (2) the recovery of genes conserved across eukaryotes, vertebrates and Actinopterygii (assessed using BUSCO v.3; Simão et al., 2015; Waterhouse et al., 2018), (3) the contig length distribution (N50) and (4) the number of genes. The assembly statistics are provided in Table S1.

TransDecoder (Haas et al., 2013) was used to detect putative coding regions within transcripts (open reading frames ≥100 amino acids long). Homology to known proteins was determined using blastp and pfam searches following the pipeline provided on the TransDecoder website (https://github.com/TransDecoder/TransDecoder/wiki). Transcripts were functionally annotated using Trinotate within Trinity. Transcript quantification was performed on the trimmed reads, followed by an analysis of differential gene expression among the four experimental conditions using Salmon (Patro et al., 2017) within Trinity. Finally, we used gene ontology (GO) terms to describe gene functions and their interactions (Gene Ontology Consortium, 2004). We performed a GO term enrichment analysis using GOseq (Young et al., 2010) as implemented in Trinity. In this analysis, we used the GO terms assigned to differentially expressed transcripts in order to test whether specific biological processes, cellular components or molecular functions (Gene Ontology Consortium, 2004) were over-represented among these differentially expressed transcripts. For example, the biological process ‘regulation of ion transmembrane transport’ (Table 1) is found to be ‘enriched’ if this process is statistically over-represented among the annotated differentially expressed transcripts (which can be under- or over-expressed). We used REVIGO (Supek et al., 2011) to reduce the functional redundancy between these GO sets and summarise the functional response to temperature and HUFA concentration in the diet. The analysis pipeline is shown in Fig. S1.

Data are presented as means±s.e.m. All statistical analyses were performed with R v.3.3.1 (2016). The level of significance for statistical analyses was always set at α=0.05.

Slopes of growth regression curves were determined for each fish. The effects of diet and temperature on these slopes, as well as on TGC, U_crit, SMR, MMR, AS, NCOT, haematocrit, HSI, CSI and FI, were tested using a two-way ANOVA (Table 2). The fish muscle/food content FA ratio was calculated for DHA, EPA and two of their precursors, i.e. the stearidonic acid 18:4n-3 and the alpha-linolenic acid 18:3n-3. The effects of diet and temperature on these ratios were also tested using a two-way ANOVA. When two-way ANOVA indicated significant interactions between diet and temperature, the differences between means were compared using a Tukey HSD test. Prior to ANOVA analyses, homoscedasticity was tested using Levene's test. Normality was tested on the residuals using the Shapiro test.

The effects of acclimation temperature, food, stimulation frequency and their interactions on cardiac variables (i.e. contraction rate, peak tension, pumping capacity, relaxation rate) were measured with the idea of comparing overall trends in heart tissue function in each experimental group. They were tested by applying linear mixed-effects models. Individual identity was considered as a random factor in order to take into account the repeated exposure of the heart of each individual to the successive frequencies and chemical stimulations. All models were computed using the lmer fitting routine in R (lme4 package; R Core Team 2013).

The pumping capacity curve was fitted using a polynomial model for each individual to determine the stimulation frequency at which optimum pumping capacity occurs. These frequency values were then compared between groups using non-parametric Wilcoxon tests.

The FA composition in fish muscle from the four experimental groups (HH20, LH20, HH12 and LH12) is fully presented in Vagner et al. (2015). Briefly, and as expected, the final fish FA composition reflected that of the diet, and HH fish presented a higher HUFA content than LH fish. In addition, EPA content was significantly affected by an interaction between diet and temperature, and was significantly less incorporated in LH20 fish.

The fish muscle/food content ratios for 18:3n-3, 18:4n-3, DHA and EPA are presented in Fig. 1. Values of this ratio above 1 suggest the selective incorporation or synthesis of this FA in fish tissue, to cope with potential deficiencies. Values around or below 1 suggest that the FA does not need to be retained, because its supply from food is close to the fish's physiological needs. Here, the 18:3n-3 ratio was below 1 for all experimental groups, suggesting that the supply of this FA from both diets corresponded to the needs of all groups, independently of their acclimation temperature. This FA ratio was significantly affected by temperature (ANOVA: F_1,46=17.11, P=1.48×10⁻⁴), diet (ANOVA: F_1,46=4.54, P=3.80×10⁻²) and an interaction between the two factors (ANOVA: F_1,46=6.43, P=0.01; Fig. 1), and was significantly highest for the HH12 group. This indicates the better incorporation of this FA in this group.

The 18:4n-3 and EPA ratios were also significantly affected by temperature (ANOVA 18:4n-3: F_1,46=49.68, P=7.71×10⁻⁹; EPA: F_1,46=8.36, P=5.8×10⁻³), diet (ANOVA 18:4n-3: F_1,46=125.49, P=9.75×10⁻¹⁵; EPA: F_1,46=160.47, P=2×10⁻¹⁶) and the interaction between the two factors (ANOVA 18:4n-3: F_1,46=23.22, P=1.61×10⁻⁵; EPA: F_1,46=13.00, P=7.64×10⁻⁴). These ratios were higher in LH groups than in HH fish, and were significantly highest in the LH12 group.

The DHA ratio was significantly affected only by diet, and reached a ratio value of about 8 in LH groups, which was about 4 times higher than for the HH groups (ANOVA temperature: F_1,46=2.90, P=0.09; food: F_1,46=25.12, P=8.46×10⁻⁶, food×diet: F_1,46=0.27, P=0.60).

The TGC was significantly affected by temperature (ANOVA: F_1,104=140.62, P=∼2×10⁻¹⁶; Table 2), and was higher at 20°C than at 12°C; however, it was not affected by diet (ANOVA: F_1,104=0.03, P=0.85) or by an interaction between temperature and food (ANOVA: F_1,104=1.19, P=0.28). The fish growth curves according to time are presented in Fig. S2.

The three body condition indices (FI, HSI and CSI) were significantly affected by temperature. The condition index FI was significantly higher for fish reared at 20°C (ANOVA: F_1,53=24.31, P=8.48×10⁻⁶ ; Table 2), while HSI and CSI were significantly higher for the 12°C groups (ANOVA HSI: F_1,50=5.93, P=1.85×10⁻²; CSI: F_1,50=4.70, P=0.03; Table 2). None of these parameters were affected by food (ANOVA HSI: F_1,50=2.68, P=0.11; FI: F_1,53=0.24, P=0.62; CSI: F_1,50=0.08, P=0.77) or by an interaction between food and temperature (ANOVA HSI: F_1,50=3.11, P=0.08; FI: F_1,53=1.08, P=0.30; CSI: F_1,50=2.90, P=0.09).

As expected, Ṁ_O2 increased exponentially with swimming speed (data not shown). The effects of diet and temperature on U_crit, SMR, MMR, AS, NCOT and haematocrit are summarised in Table 2.

U_crit was not affected by diet (ANOVA: F_1,49=0.05, P=0.83), temperature (ANOVA: F_1,49=2.34, P=0.13) or an interaction between the two (ANOVA: F_1,49=1.34, P=0.25).

SMR and NCOT were both affected by temperature: the 20°C groups showed an approximately 33% lower SMR (ANOVA: F_1,48=11.92, P=0.1×10⁻²) and a NCOT almost twice as high as that of the 12°C groups (ANOVA: F_1,44=65.50, P=2.97×10⁻¹⁰). SMR and NCOT were not affected by food (ANOVA SMR: F_1,48=0.28, P=0.60; NCOT: F_1,44=0.09, P=0.76) or by an interaction between food and temperature (ANOVA SMR: F_1,48=3.99, P=0.05; NCOT: F_1,44=0.81, P=0.37).

Both MMR and AS were affected by temperature (ANOVA MMR: F_1,45=47.88, P=1.34×10⁻⁰⁸; AS: F_1,43=68.76, P=1.84×10⁻¹⁰), and were globally higher at 20°C than at 12°C. They were also affected by an interaction between temperature and diet, and were significantly highest in the LH20 group and lowest in the LH12 group (ANOVA MMR: F_1,45=6.31, P=0.02; AS: F_1,43=5.90, P=0.02). Neither MMR nor AS was affected by food (ANOVA MMR: F_1,45=0.03, P=0.90; AS: F_1,43=0.02, P=0.90).

Haematocrit was 20% higher in the 20°C groups than in the 12°C groups (ANOVA: F_1,45=9.01, P=0.4×10⁻²), and 13% higher in HH than in LH groups (ANOVA: F_1,45=5.96, P=0.02), but there was no interaction between temperature and food (ANOVA: F_1,45=0.54, P=0.46).

Under control conditions (10 nmol l⁻¹ adrenaline; Fig. 2D), the pumping capacity increased to its maximum value before slightly decreasing at high frequencies, whatever the group (Lmer: frequency effect: χ²₁=101.8, P=2.2×10⁻¹⁶). The 12°C groups exhibited a significantly higher pumping capacity (Lmer: χ²₁=7.5, P=0.6×10⁻²; Fig. 2D), which reached a maximal value at lower stimulation frequencies than for the 20°C groups (1.08±0.40 Hz for 12°C groups versus 1.43±0.18 Hz for 20°C groups; median±IQR; Fig. 3; ANOVA: F₁=14.7, P=0.3×10⁻³). Neither an effect of diet (Lmer: χ²₁=1.2; P=0.28) nor any interactions between the factors tested was observed on pumping capacity (P>0.05).

Conversely, peak tension and half-relaxation rates were particularly affected by an interaction between temperature, diet and frequency, and were highest in LH12 fish (Fig. 2B,C). Both of these variables tended to be lowest at high frequencies (peak tension, interaction temperature×frequency×food: χ²₁=6.83, P=0.9×10⁻²; half-relaxation rate, temperature×frequency×food: χ²₁=7.1, P=0.8×10⁻²).

The contraction rate was affected by an interaction between food and frequency (lmer: food×frequency: χ²₁=21.3, P=3.9×10⁻⁶), as well as by an interaction between temperature and frequency (temperature×frequency: χ²₁=12.5, P=0.4×10⁻³; Fig. 2A). It tended to be highest in the 12°C groups, especially at 0.5 and 1 Hz, and highest in the LH groups, especially under 1 Hz, for each temperature considered.

After the addition of a high adrenaline concentration (10 µmol l⁻¹; Fig. 4), the four variables varied greatly depending on the pacing frequency (lmer: half-relaxation rate χ²₁=5.9, P=0.02; peak tension χ²₁=36.2, P=1.7×10⁻⁹; pumping capacity χ²₁=32.7, P=1.1×10⁻⁸; contraction rate χ²₁=24.6, P=7.2×10⁻⁷). The variation in contraction rate was significantly affected by the interaction between pacing frequency and food (lmer χ²₁=5.9; P=0.02). No effect of temperature was observed on these variables (P>0.05).

The blocking of ryanodine receptors by adding ryanodine induced a decrease in peak tension in all groups compared with control conditions (i.e. after basal 10 nmol adrenaline level injection; Fig. 5A), whatever the pacing frequency (all curves below the 0% variation line); however, this decrease was more marked in the 12°C groups than in the 20°C groups, and particularly below 1 Hz (lmer: temperature×frequency: χ²₁=14.06; P=1.77×10⁻⁴).

Following 10 µmol l⁻¹ adrenaline injection (Fig. 5B), the addition of ryanodine significantly reduced the peak tension in all groups, except for the LH20 group, which demonstrated a relatively stable peak tension when stimulated from 1 to 2 Hz (points near the 0% variation; lmer: temperature×food: χ²₁=3.95, P=0.04; temperature×frequency: χ²₁=15.96, P=6.46×10⁻⁵).

De novo transcriptome assembly yielded 328,877 transcripts with a contig N50 of 1474 bp. Most (92.1%) quality-filtered reads mapped back to the assembly. Transcriptome assembly was successful as most genes conserved among eukaryotes (95.7%), metazoans (95%) and actinopterygians (74.4%) were recovered in our assembly (BUSCO v3; Simão et al., 2015; Waterhouse et al., 2018).

Transcription levels differed mostly between temperature treatments (Fig. 6). At high HUFA levels [false discovery rate (FDR) threshold of 5%], 3.77% (n=2319) and 1.4% (n=860) of transcripts were respectively up- and down-regulated at 20°C versus 12°C. A similar pattern was observed at low HUFA levels (4.68% and 1.88% transcripts up- and down-regulated, respectively). Expression changes between temperature treatments (measured as log₂ fold-change, logFC) varied between −10 and 8.9 (median: −3.7) at high HUFA levels, and between −9.1 and 10.1 (median: −2.6) at low HUFA levels. No differentially expressed transcripts were detected between HUFA treatments, regardless of the temperature.

However, we observed an interaction between temperature and diet, as 4.63% and 1.85% of transcripts were up- and down-regulated at LH20, compared with HH12. Similarly, 3.42% and 1.34% of transcripts were up- and down-regulated at HH20, compared with LH12. In these cases, expression changes were of the order of −9.4 and 9.5 logFC (median: −1.1). When considering all differentially expressed genes found between temperature treatments, we observed that the median logFC was significantly stronger for up-regulated genes (−4) than for down-regulated ones (1.5; Wilcoxon rank sum test on absolute values: P<0.001). In summary, an increase in temperature resulted in more up-regulated transcripts with stronger fold-change than down-regulated transcripts. A large proportion of the transcripts mobilised in this response to increased temperature were the same, whether diet varied or not, and ranged between 27.5% of transcripts in common between the HH20–LH12 and HH12–LH20 comparisons and 51.1% of transcripts in common between the HH12–LH20 and HH12–HH20 comparisons (Fig. 7; for logFC>2).

Among all comparisons between different temperatures (HH12–HH20, LH12–LH20 and HH12–LH20; ignoring the LH12–HH20 comparison as it is ecologically improbable), identified genes up-regulated at 20°C are involved in muscle contraction (e.g. Myosin-7 and -10 isoforms; logFC between 2.44 and 7.32). At this temperature, we also detected up-regulation of a Ras-related gene (RAB44; logFC between 6.55 and 10.02), coding for a protein belonging to the plasma membrane, involved in intra-cellular membrane transport and enabling calcium ion binding.

At low HUFA, a total of 486 GO terms were significantly enriched (FDR threshold of 5%; Tables S2, S3) among differentially expressed transcripts between the 12 and 20°C temperature treatments (LH12–LH20 comparison). At high HUFA, 99 GO terms were enriched among differentially expressed transcripts between 12 and 20°C (HH12–HH20 comparison).

There were 73 common GO terms between the HH12–HH20 and LH12–LH20 comparisons (Fig. 7). These 73 GO terms represented 74% (i.e. 73/99 GO terms) of the total GO terms for HH12–HH20, but only 15% (i.e. 73/486 GO terms) of the total GO terms for LH12–LH20. In other words, 85% (413/486) of the enriched GO terms detected for the LH12–LH20 comparison were only detected in that specific comparison, while only 36% (26/99) of the enriched GO terms detected for the HH12–HH20 comparison were only detected in that specific comparison. Therefore, the functional response to temperature elevation was more specific at low HUFA than at high HUFA (Fig. 7A). Similarly, there were more significantly enriched GO terms associated with the HH12–LH20 comparison (n=311) than with the HH12–HH20 comparisons (99 GOs). The functional response for HH12–LH20 was more specific (24% of enriched GO terms in common with HH12–HH20) than that for HH12–HH20 (77% of enriched GO terms in common with HH12–LH20; Fig. 7). Depleted GO terms were never detected among these comparisons.

Four classes of biological processes recurred in all comparisons: cytoskeleton organisation, transmembrane transport (including ion transport and metal ion transport), actin filament processes and regulation of heart and muscle contraction. In all comparisons, cellular components included contractile fibre part, sarcolemma and cell junctions. Finally, molecular functions mainly included motor activity, actin or protein binding, potassium channel regulatory activity and ion transmembrane transport activity. The 73 enriched GO terms found in common between the HH12–HH20 and LH12–LH20 comparisons were reduced to key biological processes summarised as the regulation of heart rate (including action potential, positive regulation of ion transmembrane transport including calcium, ATPase activity and muscle contraction), cytoskeletal organisation and membrane assembly, and actin filament-based processes (including striated muscle cell development and cytokinesis). Cellular components included the contractile fibre part, myosin filaments, the cytoskeleton and sarcolemma, as well as cell contact zones. Molecular functions mainly involved protein binding (calmodulin, actin, ankyrin and cytoskeletal protein binding), motor activity and sodium transmembrane transporter activity.

We then looked at the 235 enriched GOs that were specific to the HH12–LH20 comparison, when compared with HH12–HH20; these GO terms are linked to the combined effect of elevated temperature and lowered concentration of HUFA in the diet (Fig. 8, Table 1). These GO terms pertained mainly to processes involving muscle contraction and development, cytoskeletal dynamics and transmembrane dynamics (Fig. 8, Table 1). More specifically, among the 165 enriched GO terms pertaining to biological processes (including heart muscle contraction, the regulation of ion transmembrane transport and cytoskeleton organisation), 14 involved calcium-mediated processes, such as the regulation of voltage-gated calcium channel activity (GO:1901385, including 11 differentially expressed transcripts), the positive regulation of calcium-transporting ATPase activity (GO:1901896, including seven differentially expressed transcripts), calcium ion transmembrane transport via high voltage-gated calcium channel (GO:0061577, including four differentially expressed transcripts), regulation of cardiac muscle contraction by regulation of the release of sequestered calcium ion (GO:0010881, including 16 differentially expressed transcripts) and regulation of the release of sequestered calcium ion into cytosol through sarcoplasmic reticulum (GO:0010880, including 21 differentially expressed transcripts). Among the 33 differentially expressed transcripts involved in these processes (26 transcripts are involved in more than one GO term), isoforms of the Ryanodine receptor 2 (RYR2; logFC ranging between 2.95 and 7.71) and the voltage-dependent L-type calcium channel subunit alpha 1C (CAC1C; logFC 4.11-6.64) were the most strongly up-regulated at 20°C. We also detected isoforms of Ankyrin (ANK2 and ANK3; logFC 2.58-5.71), which modulates the attachment of actin to the cytoskeleton (differentially expressed transcripts are presented in Table S4). REVIGO also detected many GO terms associated with morphogenesis (anatomical structure, neuron projection, hypochord, cardiac muscle tissue morphogenesis, adult heart development). Similarly, for cellular components, GO terms were associated with the cell membrane (T-tubules, polysynaptic membrane, plasma membrane region, cell–cell junctions, cell surface furrow, calcium and sodium channel complexes). A large proportion of GO terms involved in molecular function were linked to binding (e.g. cadherin binding, S100 protein binding, ion channel binding, cell adhesion molecule binding), structural molecule activity conferring elasticity and the structural constituent of the cytoskeleton.

Using an integrative approach, this study offers insights into the underlying mechanisms governing the fish cardio-respiratory phenotype in response to temperature and lowering HUFA food content. It reveals that the metabolic effects observed at the individual level may result from the modification of heart tissue functioning, particularly at high levels of effort.

The selective retention of EPA and DHA in LH fish shows that the HUFA content supplied by the LH diet was not enough to cover the requirements in fish tissue, confirming the HUFA deficiency in the LH diet. Conversely, the HH diet supplied enough HUFA to cover the requirements of fish tissue, as shown by the HUFA ratio between muscle and diet, which was always near one or above.

Temperature explained most of the variation in individual, physiological, systemic and molecular responses observed in fish acclimated at 12 and 20°C (e.g. Fig. 6). Thus, and not surprisingly, the 20°C acclimation temperature increases the metabolism (MMR and AS) and energy allocation to physiological performance, such as growth (TGC and FI index) (Claireaux and Lagardère, 1999; Fry, 1971).

Interestingly, fish acclimated at 12°C showed (i) a higher SMR than those acclimated at 20°C and, concomitantly, (ii) improved cardiac performance under the basal adrenaline level. This was shown by a higher pumping capacity, peak tension, relaxation rate and ability to reach a maximum pumping capacity at lower frequency than the 20°C groups. This effect of temperature on cardiac performance is in accordance with what has already been reported in rainbow trout (Keen et al., 1994). Similarly, as is usually found (Driedzic and Gesser, 1994; Gamperl, 2004; Keen et al., 2017), the better performance at a colder temperature was associated with a larger heart (CSI) in the present study. This would compensate for the higher blood density at low temperature by increasing the pumping capacity in order to supply blood through the entire organism (Keen et al., 2017). This increase in heart mass at low temperature is usually attributed to hyperplasia (increased cell numbers) and/or hypertrophy (increased myocyte size) resulting from a more efficient protein and/or FA synthesis or degradation (Driedzic and Gesser, 1994; Gamperl, 2004; Keen et al., 2017). Accordingly, in our study, biological processes involved in muscle cell development (GO:0055001) and striated muscle development (GO:0055002) were significantly enriched in all comparisons involving temperature (HH12 versus HH20, HH12 versus LH20, LH12 versus HH20, LH12 versus LH20; Table S3).

This cardiac remodelling at low temperature was also associated with a decreasing peak tension in the presence of ryanodine. From a functional standpoint, this indicates that the cardiac contraction at 12°C depends more on the mobilisation efficiency of intracellular Ca²⁺ stored in the SR than on an extracellular Ca²⁺ source. This increase in the SR channel contribution to Ca²⁺ mobilisation, as usually found in cold-acclimated fish (Imbert-Auvray et al., 2013; Keen et al., 1994; Shiels et al., 2002; Vornanen, 1996, 1998) could be explained by a SR proliferation in heart cells (Bowler and Tirri, 1990; Keen et al., 1994, 2017; Shiels and Farrell, 1997; Shiels et al., 2011), a higher ryanodine sensitivity (Vornanen et al., 2005) and/or a partial (or complete) inactivation, due to cold, of the SL channels/exchanger allowing extracellular Ca²⁺ influx into the cells (Vornanen, 1998; Vornanen et al., 2002).

Cardiac remodelling at low temperature is thus demonstrated in C. auratus, but this remodelling may have a cost, as suggested by the higher SMR measured at 12°C versus 20°C. However, the causative relationship between cardiac remodelling and elevated SMR remains to be fully investigated.

Moreover, higher cardiac performance and oxygen consumption were observed in LH12 fish versus HH12 fish under physiological conditions (low adrenaline level). This indicates that the assumptions made above to explain a higher cardiac performance in cold-acclimated fish compared with 20°C-acclimated fish cannot be the only processes involved, and that diet interacts with temperature to affect performance. The higher cardiac performance of LH12 versus HH12 could be the result of the use of saturated and monounsaturated FA as major metabolic fuels, preferred over HUFA, and highly present in the LH diet (Chatelier et al., 2006; Egginton, 1996; Henderson and Sargent, 1985; McKenzie, 2001). Conversely, there is no indication of the preferential use of different metabolic fuels in LH20 versus HH20 groups, as fish exposed to these two 20°C conditions demonstrated similar heart performance.

This study clearly demonstrates that the significant changes observed among the physiological, cardiac and molecular variables measured between 12 and 20°C were exacerbated by a lower HUFA availability (LH12–LH20 treatments). In particular, the mechanisms involved in Ca²⁺ mobilisation capacity during excitation–contraction coupling appear specific to the LH20 group compared with the other groups, and especially compared with the LH12 group.

In more detail, LH20 fish demonstrated a reduction of peak tension due to blocked ryanodine receptors, which was compensated for by the application of a high adrenaline dose (10 µmol l⁻¹) from 1 Hz. Firstly, this could reflect the increased adrenergic sensitivity of the ventricular tissue through an increase of β-adrenergic receptor density (Gamperl, 2004; Keen et al., 1994). Secondly, this suggests that the LH20 ventricle can increase its force-generating ability by increasing extracellular Ca²⁺ influx via sarcolemmal channels. A differential expression of genes involved in Ca²⁺ membrane transport via high voltage-gated calcium channels (such as L-type channels) supported this specific cardiac muscle contractility process in LH20 [GO:0061577 detected among the 235 enriched GO terms that were specific to the HH12–LH20 comparison (Fig. 7B, 8) and therefore indicative of the combined effect of high temperature and low HUFA dietary content]. Complementarily, a higher number of transcripts encoding L-type and ryanodine receptors were found in that group compared with the others (Table S4). Moreover, GO terms involving the SR were only significantly enriched (i.e. process over-represented among differentially expressed transcripts) in comparisons involving LH20 (LH12–LH20 and HH12–LH20 comparisons). In particular, the biological process named regulation of release of sequestered Ca²⁺ into cytosol by sarcoplasmic reticulum (GO:0010880) was detected for both comparisons (Table S2, S3). Thus, combining data from physiological tests (exposure to ryanodine and adrenaline under cardiac stimulation) and transcriptomics (expression levels of Ca²⁺ channels) suggests that, under stressful conditions (10 µmol l⁻¹ adrenaline), the low HUFA diet promotes heart contraction in warm-acclimated fish through a Ca²⁺ mobilisation depending probably on both the adrenergic pathway via SL channels and SR Ca²⁺ stores.

This strategy of LH20 fish was associated with the highest MMR and AS in our study, but also a relatively higher NCOT to swim at U_crit. From a metabolic standpoint, this suggests that fish reared in LH20 conditions utilise a more expensive cardiac strategy. In their natural environment, fish face several energy-demanding activities concomitantly to swimming, such as digestion, growth or escape from predators (Fry, 1971). All these activities are mainly governed by the heart’s capacity to supply oxygen to tissues (Claireaux, 2005; Clark et al., 2005). Even though our results were obtained experimentally, without any constraints on foraging or predator risk related to the natural environment, they suggest that cardiac strategies for Ca²⁺ mobilisation set up by LH20 fish may not be optimal and may have ecological consequences, such as a reduction in further long-term growth.

Conversely, the LH12 hearts showed the highest peak tension and half-relaxation rate, especially at low pacing frequencies, under basal metabolic conditions (Fig. 2). This suggests a more efficient Ca²⁺ mobilisation during excitation–contraction coupling in these hearts. Concomitantly, the decreasing peak tension under ryanodine treatment added to the basal adrenaline dose (10 nmol l⁻¹) and the enrichment of GO:0010880 (regulation of release of sequestered Ca²⁺ ion into cytosol by sarcoplasmic reticulum) in differentially expressed transcripts in the LH12–LH20 and HH12–LH20 comparisons support a high SR involvement in Ca²⁺ mobilisation in LH12 hearts (Tables S2 and S3).

Interestingly, at all pacing frequencies, LH12 cardiac performance was not significantly modulated under stressful conditions (10 µmol l⁻¹ adrenaline) compared with that observed under 10 nmol l⁻¹ adrenaline (i.e. no significant variation in the contraction rate, peak tension relaxation rate and pumping capacity; Fig. 4). This suggests a limited capacity to mobilise extracellular Ca²⁺ to maintain ventricle contraction, and thus a low adrenergic sensitivity in this group. The SR contribution also seems to be reinforced under stressful conditions in LH12, as shown by the more intense reduction of peak tension by ryanodine in 10 µmol l⁻¹ adrenaline pre-treated ventricles. Even if not significant, this SR contribution would tend to compensate for the low adrenergic sensitivity of the LH12 group at all pacing frequencies, shown by stable cardiac performance under adrenaline stimulation (10 µmol l⁻¹).

Taken together, these results suggest that lowering HUFA availability combined with high temperature induces (1) a proliferation of SL and SR Ca²⁺ channels and (2) a higher force-generating ability by increasing extracellular Ca²⁺ influx via sarcolemmal channels when the heart has to sustain excessive effort due to stress and/or exercise. The same diet given at low temperature, however, would induce a lower adrenergic sensitivity with a higher SR contribution, increasing under stressful conditions. In addition, while no enriched GO term corresponding to transmembrane receptors was detected, 11 biological process GO terms pertaining to the regulation of ion transmembrane transport (GO:0034765) were enriched in LH20 versus HH12 and LH12. Similarly, the plasma membrane cellular component (GO:0005886) was significantly involved. Finally, multiple GO terms involved in ion channel regulation (GO:0015459, GO:0016247, GO:0022838, GO:0022836, GO:0008514) were significantly enriched in the comparisons involving LH20 (LH12–LH20 and HH12–LH20; Table S2; Fig. 8).

One hypothesis to explain such results could be that the lesser dietary supply in EPA and DHA in LH groups versus HH groups, reflected in muscle tissue (Vagner et al., 2015), could also have impacted the membrane composition of cardiomyocytes (in which the sarcolemmal channels are embedded) and that of sub-cellular organelle membranes, such as the SR (in which the ryanodine receptors are embedded) (Porter et al., 1996). This lesser HUFA content may have then impacted the functionality of the embedded channels and receptors, as shown in other species (Turner et al., 2003; Ushio et al., 1997). The transcriptomic results support this hypothesis, showing the significant enrichment of GO:0044325, involved in ion channel regulation, in the comparisons involving the LH diet (i.e. HH12–LH20, LH12–LH20, LH12–HH20; Table S3).

The differences in cardiac performance observed between LH20 and LH12 fish could result from the different contribution to the Ca²⁺ pool of the two membrane sub-types (sarcolemmal and organellar), the HUFA composition of which varies with both acclimation temperature and energetic demand (Bell et al., 1986). While the DHA content was similarly low in the two LH groups, lipid analysis previously performed on fish muscle demonstrated that the LH diet induced a lower EPA content in muscle membranes at 20°C than at 12°C (Vagner et al., 2015). This lower EPA content in LH20 fish could be a response to increasing temperature with a decreasing degree of unsaturation known to preserve the membrane fluidity (Person-Le Ruyet et al., 2004; Raynard and Cossins, 1991; Wodtke, 1991).

According to the abundant literature reporting that dietary FA shapes heart FA composition in fish (Bell et al., 1991, 1993; Lopez-Jimena et al., 2015; Thomassen and Røsjø, 1989), we assume that the cardiomyocyte cell and sub-cellular membranes of LH20 fish also demonstrated a low EPA content, as found in their muscle. The different receptor/channel transcription and functionality demonstrated between LH20 and LH12 groups could thus be explained by the lower EPA content in the LH20 groups. Modulation of the calcium pathway by EPA and DHA has been observed in the ventricular myocytes of mammals (O'Neill et al., 2002; Swan, 2003; Turner et al., 2003; Wu et al., 2003; Xiao et al., 1997). For example, in rat ventricular myocytes, an acute treatment with EPA induced inhibition of the cardiac L-type calcium channel (Xiao et al., 1997) and a modulation of SR function in order to participate in this way in the positive antiarrhythmic effects on heart attributed to the SR (O'Neill et al., 2002; Swan, 2003). According to the literature, and to our results, we offer the hypothesis that in C. auratus, dietary EPA could be an efficient driver in the functioning of cardiomyocyte membrane sub-types, and that its deficiency in fish cardiomyocytes would result in a higher energetic cost in order to maintain an efficient cardiomyocyte functionality.

This physiological response will probably have detrimental consequences for the fish itself and for the functioning of the whole ecosystem. Our understanding of these processes is still limited, because the lessening HUFA content in food webs remains an understudied factor in the context of global change. Through this original and integrative approach, we hope that this study will spark interest in considering the effects of HUFA availability on ectotherms, given its significant functional interaction with temperature.

The authors are very grateful to M. Durollet for her help with fish rearing, and to M. Le Gac and J. Quéré from the DYNECO lab of the Ifremer Institute of Plouzané (France) for lab space and their great help with the preparation of DNA libraries. The authors also thank the two reviewers for their constructive criticism on the manuscript. We thank the bioinformatics platform Toulouse Midi-Pyrenees (Bioinfo GenoToul) for providing help and computing resources.