Developmental Stage, Muscle and Genetic Type Modify Muscle Transcriptome in Pigs: Effects on Gene Expression and Regulatory Factors Involved in Growth and Metabolism.

Iberian pig production includes purebred (IB) and Duroc-crossbred (IBxDU) pigs, which show important differences in growth, fattening and tissue composition. This experiment was conducted to investigate the effects of genetic type and muscle (Longissimus dorsi (LD) vs Biceps femoris (BF)) on gene expression and transcriptional regulation at two developmental stages. Nine IB and 10 IBxDU piglets were slaughtered at birth, and seven IB and 10 IBxDU at four months of age (growing period). Carcass traits and LD intramuscular fat (IMF) content were measured. Muscle transcriptome was analyzed on LD samples with RNA-Seq technology. Carcasses were smaller in IB than in IBxDU neonates (p < 0.001), while growing IB pigs showed greater IMF content (p < 0.05). Gene expression was affected (p < 0.01 and Fold change > 1.5) by the developmental stage (5,812 genes), muscle type (135 genes), and genetic type (261 genes at birth and 113 at growth). Newborns transcriptome reflected a highly proliferative developmental stage, while older pigs showed upregulation of catabolic and muscle functioning processes. Regarding the genetic type effect, IBxDU newborns showed enrichment of gene pathways involved in muscle growth, in agreement with the higher prenatal growth observed in these pigs. However, IB growing pigs showed enrichment of pathways involved in protein deposition and cellular growth, supporting the compensatory gain experienced by IB pigs during this period. Moreover, newborn and growing IB pigs showed more active glucose and lipid metabolism than IBxDU pigs. Moreover, LD muscle seems to have more active muscular and cell growth, while BF points towards lipid metabolism and fat deposition. Several regulators controlling transcriptome changes in both genotypes were identified across muscles and ages (SIM1, PVALB, MEFs, TCF7L2 or FOXO1), being strong candidate genes to drive expression and thus, phenotypic differences between IB and IBxDU pigs. Many of the identified regulators were known to be involved in muscle and adipose tissues development, but others not previously associated with pig muscle growth were also identified, as PVALB, KLF1 or IRF2. The present study discloses potential molecular mechanisms underlying phenotypic differences observed between IB and IBxDU pigs and highlights candidate genes implicated in these molecular mechanisms.

Introduction

Modern pig production is mostly based on extensively selected breeds, which show optimized productivity and efficiency [1]. However, in the Mediterranean area these commercial breeds coexist with local breeds, used for the production of unique high-quality traditional pork products. These breeds, usually known as fatty-pigs, are smaller in size, have not undergone intense genetic selection and are less productive than modern breeds [2]. Moreover, due to the traditional rearing system, under extensive free-range conditions, they are exposed to harsh environments and seasonal variations in food availability (associated with the development of a thrifty genotype) [3].

The Iberian pig is the most representative Mediterranean traditional breed, and it has an important commercial value based on high quality dry-cured products in terms of consumers’ health and acceptance [4]. Iberian pig shows special growth, fattening and meat characteristics, mainly as a consequence of its particular high fat deposition and desaturation potential, high food intake associated with leptin resistance [5, 6] and the elevated age at which those pigs are slaughtered [7, 8]. Iberian pigs require a long productive cycle (12–18 months) to develop desirable characteristics expected from Iberian pork products. It has been reported that parameters such as intramuscular fat (IMF) content, fatty acids (FA) profile or redness influence meat quality and improve as pigs mature [8-12].

Due to the low productivity of this breed and to the long time required to obtain Iberian products, the Duroc breed has been introduced as a terminal sire to improve reproductive and growth performances and primal cuts yield. However, the introduction of Duroc genetics is associated with a decrease in meat quality, mainly determined by a decrease in IMF and monounsaturated fatty acids (MUFA) contents [13].

Intramuscular fat content and fatty acid composition are main factors affecting meat quality and strongly depend on genetic type, diet, anatomical location and age [10, 14, 15]. Both the number and size of adipocytes within muscle fibers determine intramuscular fat content. During prenatal development and immediately after birth, preadipocyte differentiation is a very active process that slows down with animal growth [16]. Later in growth, hypertrophy is the most important issue affecting IMF content, although hyperplasia is maintained in the adult animal to a lesser extent [17]. Hypertrophy is driven by triglycerides accumulation in mature adipocytes. The ratio of synthesis (lipogenesis) and degradation (lipolysis) determines the amount of triglycerides stored in these cells. Due to the different timing of those two processes (hyperplasia and hypertrophy), sequential studies across growing stages are required for understanding the set of molecular mechanisms driving fat accretion in Iberian pigs Moreover, fat accumulation also depends on the anatomical location. Several studies reported differences in lipid composition (probably associated with muscle fiber composition) [15, 18], oxidative properties and IMF content [19, 20] across muscles.

Due to the importance of lipid synthesis and accumulation in meat quality [21], new interest has arisen towards the understanding of genetic mechanisms underlying such processes. With this goal, some studies based on microarray technology investigated transcriptome differences between genotypes, such as Iberian vs. Large White or Duroc in endocrine tissues [22]. On the other hand, several studies have addressed the effect of the developmental stage (prenatally and postnatally) on muscle transcriptome [23-25] and the effect of muscle type on transcriptome and proteome, showing important functional differences [26, 27]. For example, 15–30% of proteome has been reported to differ between Longissimus dorsi (LD) and Biceps femoris (BF) muscles [27, 28]. Both LD and BF muscles are of high economic relevance in the Iberian pig industry. So far, LD has been examined in more detail and more frequently but, due to the aforementioned differences, the usefulness of the joint analysis of different muscles is evident, as proposed by Sobol et al. [26].

Previous studies focused on transcriptomic differences between pure Iberian and Duroc-crossbred Iberian pigs using microarray technology in the loin [29] and RNA-Seq technology in BF muscle [30]. However, this is the first RNA-Seq technology-based study assessing genetic differences between genotypes, developmental stages and muscle types aimed at improving the knowledge of the genetic and metabolic basis of meat quality and productive traits in Iberian pigs.

A previous study analyzing BF transcriptome in IB pigs showed phenotypic and transcriptomic differences associated with growth and meat quality [30]. The present study arises from that previous findings aiming to: 1) Address the effects of genetic type, muscle, developmental stage and their interactions on phenotypic parameters; 2) Evaluate changes in muscle gene expression conditional on growth stage, genetic type and muscle, that might be responsible for the observed phenotypic differences and identify pathways and networks in which those genes are involved; 3) Identify transcription factors affecting gene expression in order to establish potential new candidate genes affecting productive parameters and meat quality in purebred and Duroc-crossbred Iberian pigs, which could become targets for future studies as genetic markers for selection.

Materials and Methods

Ethics statement

Animal manipulations were done in compliance with the regulations of the Spanish Policy for Animal Protection RD1201/05, which meets the European Union Directive 86/609 about the protection of animals used in research. The experiment was specifically assessed and approved (report CEEA 2010/003) by the INIA (Instituto Nacional de Investigación y Tecnología Agraria y Alimentaria) Committee of Ethics in Animal Research, which is the named Institutional Animal Care and Use Committee (IACUC) for INIA.

Animals and sample collection

Twenty-two pure Iberian sows raised in the same commercial farm were utilized at their third gestation cycle. All females were managed in the same conditions. Eleven sows were mated to Iberian boars and eleven to Duroc boars. Ten litters were employed for the sampling at birth (5 for purebreds and 5 for crossbreds) and twelve for sampling at growing stage (6 for purebreds and 6 for crossbreds). At birth, nine IB and 10 IBxDU male piglets were randomly selected from the ten litters (two from each litter with the exception of one litter, which provided only a single Iberian male for the study). The piglets were sampled immediately after euthanasia performed by stunning and exsanguination in compliance with RD1201/05 standard procedures. The pigs coming from the remaining litters were subject to standard productive management, being both groups fed ad-libitum the same barley-wheat commercial growing diet up to four months of age (growing period), when seven IB and 10 IBxDU were randomly selected (two from each one of 5 litters and seven litters providing only one male each) were slaughtered as previously described. Blood samples were collected from newborns and growing pigs in sterile heparin blood vacuum tubes (Vacutainer Systems Europe, Meylan, France). Immediately after recovery, the blood was centrifuged at 1500 g for 15 min and the plasma was separated and stored into polypropylene vials at −20°C until assayed for determination of glucose and lipids metabolism-indicating parameters. After blood collection, pigs were slaughtered. Several body measures were obtained with a tape measure: total body length (from the rostral edge of the snout to the tail insertion), ham length (from the anterior edge of the Symphysis pubica to the articulatio tarsi), total length of anterior and posterior limbs (from the distal edge of the hooves to the proximal edge of the scapula or Symphysis pubica, respectively); and thoracic, abdominal and ham circumferences. Carcasses were weighted and samples from LD muscle in newborn and growing pigs were vacuum-packed in low-oxygen permeable film and kept frozen at –20°C until fatty acid composition analysis. Prior to fatty acid analysis, muscle samples were freeze dried for two days in a lyophilizer (Lyoquest, Telstar, Tarrasa, Spain) and ground in a Mixer Mill MM400 (Retsch technology, Haan, Germany) until muscle was completely powdered. For transcriptomic analysis, LD samples were immediately frozen in liquid nitrogen and maintained at –80°C until RNA extraction.

Plasma parameters for glucose (glucose and fructosamine) and lipid metabolism (triglycerides, total cholesterol, high-density lipoprotein cholesterol, HDL-c, and low-density lipoprotein cholesterol, LDL-c) plasmatic levels were measured with a clinical chemistry analyzer (Saturno 300 plus, Crony Instruments s. r. l., Rome, Italy).

Plasma metabolites and carcass data from newborn pigs were analyzed in a previous study [30].

Tissue composition analysis

Longissimus dorsi muscle IMF content was quantified using the method proposed by Segura et al. [31] based on gravimetrical determination of lipid content. Fatty acid methyl esters (FAMEs) were identified by gas chromatography as described by López-Bote et al. [32] using a Hewlett Packard HP-6890 (Avondale, PA, USA) gas chromatograph equipped with a flame ionization detector and a capillary column (HP-Innowax, 30 m × 0.32 mm i.d. and 0.25 μm polyethylene glycol-film thickness). Results were expressed as grams per 100 grams of detected FAMEs.

Statistical analyses of tissue composition

Phenotypic data were analyzed using the general linear model (GLM) procedure using SAS version 9.2 (SAS Inst. Inc., Cary, NC; 2009). The mean and genetic type/muscle/developmental stage were considered as systematic effects, and residual effects as random. Carcass weight was used as covariate when it was significant and removed from the model when it was not. The animal was the experimental unit for all analysis. The results were considered to be significant at p-value < 0.05.

Transcriptomic analysis

RNA extraction

A total of 24 animals were randomly chosen among the available ones, assuring the representation of all available litters, and used to perform transcriptomic analysis (6 animals of each genetic type at each studied age). Total RNA was extracted from 50-100mg samples of LD muscle using the RiboPure TM of High Quality total RNA kit (Ambion, Austin, TX, USA) following the manufacturer’s recommendations. RNA was quantified using a NanoDrop-100 spectrophotometer (NanoDrop Technologies, Wilmington, DE, USA). The quality of the RNA was evaluated using the RNA Integrity Number (RIN) value from the Agilent 2100 Bioanalyzer device (Agilent technologies, Santa Clara, CA, USA). The RIN values ranged from 7.8 to 9.8.

Library construction and RNA sequencing

Sequencing libraries were made using the mRNA-Seq sample preparation kit (Illumina Inc., Cat. # RS-100-0801) according to manufacturer’s protocol. Each library was sequenced using TruSeq SBS Kit v3-HS, in paired end mode with the read length 2x76 bp on a on a HiSeq2000 sequence analyzer (Illumina, Inc). Images from the instrument were processed using the manufacturer’s software to generate FASTQ sequence files.

Mapping and assembly

Sequence reads were analyzed using CLC Bio Genomic workbench software 7.0 (CLC Bio, Aarhus, Denmark). Quality control analysis was performed using the NGS quality control tool, which assesses sequence quality indicators based on the FastQC-project (FastQC-project ). Quality was measured taking into account sequence-read lengths and base-coverage, nucleotide contributions and base ambiguities, quality scores as emitted by the base caller and over-represented sequences [33]. All the samples analyzed passed all the QC parameters having the same length (76 bp), 100% coverage in all bases, 25% of A, T, G and C nucleotide contributions, 50% GC on base content and less than 0.1% over-represented sequences. A hierarchical clustering of the samples was also performed to detect samples that might have been wrongly allocated to a group, samples of unintended or unclean tissue composition or samples for which the processing has gone wrong. Sequence paired-end reads (76bp) were assembled against the annotated Sscrofa10.2 reference genome (http://www.ncbi.nlm.nih.gov/genome/?term=sus+scrofa) using the genome, annotated genes and mRNA tracks. Data was normalized by calculating the ‘fragments per kilo base per million mapped reads’ (FPKM) for each gene [34].

Differential expression analysis

The statistical analysis was performed using the total exon reads as expression values by the Empirical analysis of differential gene expression tool (CLC Bio, Aarhus, Denmark). This tool is based on the EdgeR Bioconductor package [35] and uses count data (i.e. total exon reads) and the same normalization method (trimmed mean of M-values, TMM [35]) for the statistical analysis. Genes were filtered according to two criteria: a minimum mean group expression greater than 0.5 FPKM in at least one group and a Fold-Change (FC) of the expression differences between genotypes, muscles and developmental stages equal or higher to 1.5. Finally, those genes with a p-value ≤ 0.01, corresponding to a false discovery rate (FDR) value ≤ 0.14 in newborns and ≤ 0.20 in growing pigs, were considered as differentially expressed (DE). Gene expression data of BF muscle was previously obtained and analyzed [30] following the same criteria. Three different effects were independently assessed in the differential expression study: developmental stage effect on LD muscle transcriptome, genetic type effect on LD muscle transcriptome (at both ages independently); and muscle type effect (LD vs BF) at birth.

Results validation by RT qPCR

RNA obtained from the 24 animals used in the RNA-Seq assay was utilized to perform the technical validation of the differential expression of some genes that were either affected by the developmental stage, the muscle or the genetic type. Moreover, RNA obtained from all the available animals (16 IB and 20 IBxDU) was used to quantify expression differences by qPCR [30]. The expression of 5 genes: MTUS2, ALDH2, ADAMTS8, HDAC9 and SIM1 (Single-Minded Family BHLH Transcription Factor 1)) was quantified employing methods previously described [30]. The GAPDH, TBP, B2M and ACTB genes were selected as the most stable endogenous genes between the different studied conditions to normalize the data.

The technical validation was performed by studying the Pearson correlation between the expression values obtained from RNAseq data (FPKM) and the normalized gene expression data obtained by RT qPCR, as previously described [30]. To validate the global RNA-Seq results, the concordance correlation coefficient (CCC) [36] was calculated between the FC values estimated from RNA-Seq and qPCR expression measures for the 5 genes analyzed by the two technologies (RNA-Seq and qPCR).

Systems biology study

The biological interpretation of the DE genes was performed using two complementary approaches in order to identify: 1) enriched pathways and networks involving the DE genes, and 2) potential regulators causing the observed changes in gene expression.

Ingenuity Pathway Analysis, (IPA) (Ingenuity Systems, Qiagen, California) software was utilized to identify and characterize biological functions, gene networks and canonical pathways affected by the DE genes. The IPA Canonical Pathways Analysis identified the pathways from the Ingenuity Pathways Analysis library of canonical pathways that were most significant in our dataset. The significance of the association between the dataset and the canonical pathway was measured with Fischer’s exact test, to calculate a P‒value determining the probability that the association between the genes in the dataset and the canonical pathway is explained by chance alone.

Regulatory transcription factors (TRF), which could potentially affect the DE genes in the dataset were also studied by following complementary approaches. First, RIF1 and RIF2 metrics [37, 38] were calculated for the whole set of DE genes obtained conditional on developmental stage (5,812 genes), muscle type (135 genes) and genetic type at birth (261 genes) and at the growing period (113 genes). Candidate TRFs in pigs were obtained from Animal TFDB ("http://www.bioguo.org/AnimalTFDB/BrowseAllTF.php?spe=Susscrofa"). A total of 1,038 TRF were retrieved. Among them, 734, 739, 655 and 738 showed expression values greater than 0.5 FPKM in at least one experimental group when analyzing the developmental stage, genetic type (at birth and at growth) and the muscle type effects, respectively. Those TRFs were thus used in the RIFs metrics approach.

The RIF1 and RIF2 values were computed for the ith TRF as follows: and

Where nde is the number of DE genes, aj and dj the estimated average expression and differential expression of the j DE gene, r1ij and r2ij the co-expression correlation between the i TRF and the j DE gene in each one of the genetic types and being e1j and e2j the expression of the j gene in each genetic type [39]. Both RIF measures for each analyzed TRF were transformed to standardized z-scores by subtracting the mean and dividing by its standard deviation. We identified relevant TRF as those with extreme RIF z-scores according to the corresponding confidence intervals (CI) calculated by bootstrap. In each iteration of bootstrapping, a set of n = number of DE genes were randomly selected from the total expressed genes, and the RIF1 and RIF2 z-scores of the expressed TRF list were calculated for each studied effect. The procedure was repeated 10,000 times for each scenario to obtain the corresponding 95 and 99% CI intervals of both z-scores.

Complementarily, IPA software was used to identify and characterize potential regulators using two different tools, the upstream regulators and the regulators tools. Both of them identify known regulators that may affect the expression of the data set of DE genes. IPA-identified regulators include genes, but also other molecules such as drugs. Thus, out of the identified regulators, only genes that were also included in the RIFs metrics candidate TRF list were considered (genes included in the animal TFDB and with expression values higher than 0.5 FPKM in at least one experimental group).

Using the information obtained from the TRF study, an additional search for enriched pathways and networks was carried out with IPA software considering both, DE genes and TRF conditional on developmental stage, muscle and genetic type.

Results

The transcriptome of LD muscle was studied at two different developmental stages, birth (mean live weight of 1.5 (SD = 0.4) kg) and growth (mean live weight of 64.9 (SD = 10.1) kg). Moreover, the BF muscle transcriptome was previously assessed at birth in the same newborn pigs [30]. The experimental design allowed the study of several main effects regarding muscle transcriptome: developmental stage, genetic type and tissue effects. The biggest transcriptome change was observed between the two studied ages, involving more than 5,800 DE genes (Table 1), followed by the genetic type effect, responsible for differential expression of 261 genes at birth and 113 at the growing stage and the tissue effect, with 135 DE genes between LD and BF muscles (Table 1).

Table 1

Differentially expressed (DE) genes as affected by the three studied main effects.

	Total DE genes	Upregulated genes
		Birth	four months
Developmental stage effect (LD)	5812	3290	2522
		IBxDU1	IB2
Genotype effect at birth (LD)	261	131	130
Genotype effect at four months (LD)	113	25	88
		BF3	LD4
Muscle effect (birth)	135	52	83

1 IBxDU = Iberian x Duroc crossbred pigs

2 IB = Purebred Iberian pigs

3 BF = Biceps femoris muscle

4 LD = Longissimus dorsi muscle

Phenotypic results

Developmental stage and genetic type significantly affected carcass characteristics, plasma biochemical parameters and meat quality traits, such as IMF and fatty acid profile (Table 2). Pure Iberian and crossbred piglets had an average of 1.2 and 1.8 kg birth-weight, respectively (p<0.001; Table 2). Genetic type affected all the carcass phenotypic parameters at this early age as reported in a previous study employing the same newborn animals [30]: IBxDU neonates were bigger and heavier (p < 0.001) than IB newborns. On the other hand, growing IB and IBxDU pigs slaughtered at 59.4 and 68.6 kg live weight, respectively, showed no significant difference in body weight or size, although ham weight and perimeter were higher in IBxDU than in IB pigs (p = 0.039 and p = 0.034, respectively). A significant interaction between developmental stage and genotype was observed for ham weight and circumference measures.

Table 2

Effect of genotype, developmental stage and their interaction on phenotype of pure and crossbred Iberian pigs.

Genetic type (GT)									Stage	GT*Stage
Birth					Growing
Carcass traits	IBxDU1	IB2	SEM3	p-value	IBxDU4	IB5	SEM	p-value	p-value	p-value
Live weight	1.77	1.21	0.07	< .001	68.60	59.40	2.41	0.084	< .001	0.054
Carcass weight	1.41	0.96	0.05	< .001	56.22	48.60	2.17	0.109	< .001	0.074
Carcass length	40.20	35.50	0.54	< .001	121.11	116.17	1.54	0.140	< .001	0.936
Thorax circumference	25.15	22.06	0.39	0.001	89.17	86.00	1.11	0.186	< .001	0.974
Abdomen circumference	18.90	17.28	0.38	0.049	77.28	78.17	1.21	0.724	< .001	0.288
Ham weight	0.16	0.11	0.01	< .001	7.56	6.23	0.29	0.040	< .001	0.017
Ham length	7.45	6.33	0.14	< .001	25.89	24.33	0.36	0.052	< .001	0.541
Ham circumference	12.55	10.89	0.23	0.002	67.83	59.42	1.74	0.034	< .001	0.040
Lipid and glucose metabolism-related plasma indicators
Cholesterol	62.2	102.4	5.60	0.003	107.6	95.1	3.86	0.136	0.144	0.031
LDL6	42.2	45.8	4.40	0.449	67.4	61.5	2.58	0.280	0.002	0.431
HDL7	22.4	41.2	4.25	0.018	30.2	24.0	2.93	0.313	0.452	0.049
TG8	30.0	76.7	5.11	< .001	46.4	53.7	4.03	0.439	0.354	0.026
Fructosamine	169.7	133.7	10.37	0.101	221.3	228.0	8.05	0.691	< .0001	0.135
Glucose	132.4	123.4	10.80	0.684	107.8	92.5	2.91	0.023	0.036	0.804
Longissimus dorsi muscle main fatty acids composition (g/100 g total fatty acids)
IMF9	2.16	2.32	0.50	0.317	2.87	4.05	0.20	0.026	< .0001	0.655
C14:0	2.77	2.56	0.12	0.406	1.14	1.32	0.08	0.295	< .0001	0.230
C15:1	1.51	1.26	0.09	0.171	0.60	0.33	0.04	0.026	< .0001	0.743
C16:0	26.29	26.12	0.18	0.655	24.56	25.44	0.39	0.114	0.001	0.072
C16:1 n-9	2.13	2.03	0.05	0.547	0.20	0.27	0.02	0.099	< .0001	0.232
C16:1 n-7	5.93	5.40	0.22	0.234	3.20	3.34	0.13	0.938	< .0001	0.357
C17:0	1.78	1.50	0.07	0.064	0.47	0.32	0.03	0.039	< .0001	0.197
C17:1	0.94	0.88	0.05	0.564	0.46	0.31	0.03	0.042	< .0001	0.862
C18:0	10.62	9.43	0.32	0.083	12.91	12.69	0.23	0.787	< .0001	0.130
C18:1 n-9	24.06	26.61	0.80	0.129	41.39	42.42	0.49	0.373	< .0001	0.085
C18:1 n-7	6.34	5.82	0.21	0.235	2.67	2.21	0.13	0.051	< .0001	0.623
C18:2 n-6	7.00	8.76	0.56	0.132	7.32	5.36	0.32	0.020	0.040	0.022
C18:3 n-3	0.28	0.25	0.02	0.351	0.18	0.18	0.01	0.839	0.001	0.370
C18:4 n-3	0.13	0.17	0.02	0.351	0.08	0.08	0.00	0.347	0.019	0.222
C20:0	0.36	0.27	0.04	0.320	0.20	0.19	0.01	0.963	0.031	0.532
C20:1 n-9	0.63	0.60	0.05	0.786	0.77	0.78	0.02	0.635	0.007	0.957
C20:2 n-6	0.32	0.33	0.06	0.913	0.61	0.45	0.04	0.040	0.024	0.187
C20:3 n-6	0.63	0.52	0.02	0.017	0.25	0.23	0.02	0.525	< .0001	0.064
C20:4 n-6	5.65	4.97	0.27	0.221	2.02	1.18	0.11	0.007	< .0001	0.578
C22:4 n-6	1.07	0.91	0.05	0.175	0.43	0.54	0.04	0.192	< .0001	0.048
C22:5 n-3	0.36	0.33	0.01	0.254	0.34	0.85	0.04	< .0001	< .0001	< .0001
C22:6 n-3	0.31	0.26	0.04	0.544	0.40	0.99	0.11	0.020	0.003	0.012
∑SFA10	42.57	40.46	0.54	0.068	39.32	40.13	0.55	0.179	0.009	0.012
∑MUFA11	41.41	42.82	0.50	0.174	49.50	49.01	0.64	0.419	< .0001	0.108
∑PUFA12	16.02	16.72	0.44	0.441	11.73	9.71	0.50	0.102	< .0001	0.284
UI13	91.19	91.79	0.44	0.441	79.37	81.39	1.58	0.752	< .0001	0.968
∑n-314	1.35	1.22	0.91	0.746	1.00	2.10	0.16	0.004	0.198	0.001
∑n-615	14.67	15.50	0.07	0.324	10.73	7.65	0.47	0.013	< .0001	0.051
∑n-6/∑n-3	11.61	12.88	0.42	0.344	14.31	3.89	1.22	0.003	0.043	< .0001

1 IBxDU = Iberian x Duroc crossbred pigs (n = 10)

2 IB = Purebred Iberian pigs (n = 9)

3 SEM = Standard error of the mean

4 IBxDU = Iberian x Duroc crossbred pigs (n = 10)

5 IB = Purebred Iberian pigs (n = 7)

6 LDL = Low density lipoproteins

7 HDL = High density lipoproteins

8 TG = Triglycerides

9 IMF = Intramuscular fat

10 ΣSFA = Sum of saturated fatty acids

11 ΣMUFA = Sum of monounsaturated fatty acids

12 ΣPUFA = Sum of polyunsaturated fatty acids

13 UI = Unsaturation index = 1 × (% monoenoics) +2 × (% dienoics) +3 × (% trienoics) +4 × (% tetraenoics) +5 × (% pentaenoics) +6 × (% hexaenoics)

14 Σn3 = Sum of n-3 fatty acids

15 Σn6 = Sum of n-6 fatty acids

Purebred IB newborns had significantly greater plasma levels of total and HDL cholesterol, and triglycerides (TG) than IBxDU neonates (Table 2). However, growing IB and IBxDU pigs showed similar cholesterol and TG levels, resulting in a developmental stage*genetic type interaction effect for these parameters. On the other hand, growing pigs differed in plasma glucose: IBxDU showed higher levels than IB pigs at this stage (p = 0.023).

Regarding loin IMF, total content increased over time (p < 0.0001). Genetic type did not affect IMF content in LD muscle (Table 2) of newborn piglets, while IB showed a 41% increase in IMF content with respect to IBxDU pigs at the growing stage (p = 0.026). Developmental stage also affected IMF composition: MUFA content increased and PUFA content decreased over time (p < 0.0001) in both genetic types. However, the SFA content slightly decreased along growth in IBxDU but not in IB pigs, with a significant interaction between genotype and developmental stage (p = 0.012). Also, C20:3 n-6 concentration (p = 0.017) was higher in newborn IBxDU pigs and SFA content showed a trend (p = 0.068) in the same direction. In the growing period, C15:1, C17:0, C17:1, C18:2-n6, C20:2 n-6 and C20:4 n-6 concentrations were significantly higher in IBxDU than in IB pigs (p = 0.026, 0.039, 0.042, 0.020, 0.040 and 0.007, respectively) while C22:5 n-3 (DPA) and C22:6 n-3 (DHA) proportions were higher in IB pigs (p < 0.0001 and p = 0.020, respectively). As a consequence, total n-3 fatty acids content was higher in IB (p = 0.004) and total n-6 and the ratio n-6/n-3 was higher in the loin of growing IBxDU pigs (p = 0.013 and p = 0.003, respectively).

Data on BF phenotype and transcriptome was available for the same newborn animals from a previous study [30]. Thus, the effect of the studied muscle (BF vs LD) was assessed in IB and IBxDU neonates (Table 3). Biceps femoris and LD muscles showed similar IMF content at this early age (p = 0.158). Differences in total PUFA and n-3 PUFA content were observed, with BF showing higher levels than LD muscle. These differences resulted in a lower n-6/n-3 ratio in BF when compared to LD muscle. When analyzing the global effect of genotype on IMF content and composition of both muscles at birth, higher SFA and lower n-6 PUFA content was observed in IBxDU pigs. No significant interaction was found between the genotype and muscle effects.

Table 3

Longissimus dorsi (LD) and Biceps femoris (BF) muscle characteristics of pure and Duroc-crossbred Iberian piglets at birth.

		Mean			p-value
	Muscle	IBxDU1	IB2	SEM3	Genetic type (GT)	Muscle	GT*Muscle
∑SFA4	LD	42.57	40.79	0.33	0.002	0.182	0.538
	BF	41.66	39.91
∑MUFA5	LD	41.41	42.73	0.32	0.114	0.478	0.653
	BF	41.01	41.77
∑PUFA6	LD	16.02	16.48	0.31	0.057	0.039	0.845
	BF	17.34	18.32
∑n-37	LD	1.35	1.24	0.04	0.145	< .0001	0.942
	BF	1.77	1.69
∑n-68	LD	14.67	15.24	0.30	0.034	0.123	0.833
	BF	15.56	16.63
∑n-6/∑n-3	LD	11.61	12.54	0.36	0.056	0.002	0.892
	BF	8.78	9.96
IMF9	LD	2.16	2.26	0.08	0.060	0.158	0.644
	BF	1.80	2.16

1 IBxDU = Iberian x Duroc crossbred pigs (n = 10)

2 IB = Purebred Iberian pigs (n = 9)

3 SEM = Standard error of the mean

4 ΣSFA = Sum of saturated fatty acids

5 ΣMUFA = Sum of monounsaturated fatty acids

6 ΣPUFA = Sum of polyunsaturated fatty acids

7 Σn3 = Sum of n-3 fatty acids

8 Σn6 = Sum of n-6 fatty acids

9 IMF = Intramuscular fat

Transcriptome analysis

Mapping results

An average of approximately 84 million sequence reads was obtained for each individual sample and were assembled and mapped to the annotated Sscrofa10.2 genome assembly (22,861 genes). In all samples, 67–77% of the reads were categorized as mapped reads to the porcine reference sequence. The FPKM values were used to establish the total number of genes expressed in the muscle transcriptome (> 0.5 FPKM). Approximately 50% of total porcine annotated genes in the Sscrofa10.2 genome assembly were expressed in the studied samples (an average of 11,506 genes out of 22,861 annotated genes).

Effect of developmental stage on Longissimus dorsi transcriptome

A total of 3,290 genes were upregulated in LD muscle at birth when compared to four months-old pigs (p < 0.01, FDR < 0.015), with FC ranging from 1.5 to 219 (IBSP gene showed the largest expression difference). In four months-old pigs, 2,522 genes were upregulated when compared to newborn pigs. Fold changes (FC) ranged from 1.5 to 273, with several immunoglobulin genes showing the biggest expression changes (S1 Table).

Gene expression differences were functionally interpreted using IPA software to detect enriched pathways (S2 Table) and biological functions (Table 4). Moreover, the regulators study gathered 626 genes potentially affecting expression of the DE genes in the dataset (i.e. MSTN, FOXO3, FOXO1, MEF2C or MEF2D). Among them, 270 were identified by IPA software, 142 were identified following the RIFs approach and 344 were found DE between developmental stages (S3 Table). Within the selected TRF, eight (CBX4, HSF1, JARID2, MYC, TCF4, TEAD4, TFEB and YBX3) were identified following the three approaches, being considered as especially relevant regulators.

Table 4

Enriched biological functions in the set of DE genes between LD muscle from newborn and growing Iberian pigs

Upregulated at birth			Upregulated at growing stage
Enriched function	Corrected p-value	z-score1	Enriched function	Corrected p-value	z-score
Invasion of cells	2.66E-10	3.93	Organismal death	4.90E-23	6.95
Size of body	2.23E-06	3.58	Growth Failure	1.08E-04	4.94
Transport of molecule	1.78E-08	3.33	Bleeding	1.31E-03	4.93
Invasion of tumor cell lines	1.67E-06	3.18	Perinatal death	6.00E-06	3.92
Cell movement	1.63E-17	3.09	contractility of skeletal muscle	5.25E-04	3.86
Cell movement of tumor cell lines	3.06E-07	3.05	Hypoplasia of organ	2.73E-03	3.26
Apoptosis of fibroblast cell lines	5.99E-08	3.03	Hypoplasia	1.62E-04	3.15
Adhesion of connective tissue cells	1.40E-03	2.94	Dysgenesis	1.21E-04	3.04
Cancer	9.37E-54	2.94	Congenital anomaly of musculoskeletal system	7.19E-05	2.92
Migration of cells	2.73E-03	2.90	Multiple congenital anomalies	1.13E-07	2.91
Migration of tumor cell lines	9.70E-15	2.77	Polymerization of protein	2.01E-05	2.73
Limb development	1.46E-09	2.77	Fibrosis	1.06E-04	2.39
Invasion of tumor	5.45E-07	2.70	Hypoplasia of thorax	1.41E-03	2.25
Proliferation of fibroblasts	2.37E-04	2.59	Autophagy	1.84E-03	2.18
Adhesion of tumor cell lines	4.38E-09	2.58	Anemia	1.49E-03	1.84
Disassembly of filaments	5.25E-07	2.56	Mass of muscle	1.99E-14	1.77
Synthesis of DNA	4.31E-05	2.56	Adhesion of extracellular matrix	1.75E-05	1.73
Invasion of breast cancer cell lines	9.33E-05	2.56	Autosomal recessive disease	1.61E-03	1.72
Synthesis of carbohydrate	2.27E-03	2.52	Mass of skeletal muscle	2.14E-07	1.69
Tumorogenesis of tissue	2.43E-04	2.46	Blood vessel defect	1.28E-03	1.56
Neoplasia of cells	1.28E-03	2.45	Systemic autoinmune syndrome	7.10E-05	1.53
Invasion of tumor cells	1.95E-03	2.40	Dwarfism	2.68E-08	1.46
Organization of cytoplasm	1.68E-47	2.36	Abnormal bone density	2.78E-03	1.45
Metabolism of carbohydrate	3.24E-04	2.36	Dysplasia of skeleton	6.50E-13	1.38
Endocytosis	1.21E-05	2.33	Aneurysm	3.42E-18	1.38
Proliferation of neuronal cells	2.23E-03	2.30	Vascular tumor	1.09E-20	1.33
Organization of cytoskeleton	3.35E-16	2.28	Hemangioma	4.13E-04	1.33
Development of central nervoussystem	4.95E-06	2.26	Hypertension	2.16E-03	1.33
Development of body trunk	3.87E-05	2.22	Hypoglycemia	3.79E-05	1.24
Outgrowth of cells	5.94E-07	2.19	Adhesion of cell-associated matrix	1.86E-03	1.19
Neoplasia of epitelial tissue	2.20E-15	2.18	Ulcer	1.11E-04	1.13
Formation of celular protrusions	1.88E-04	2.18	Chronic inflammatory disorder	1.15E-03	1.06
Microtubule dynamics	1.50E-08	2.16	Replication of virus	9.27E-05	1.01
Behavior	6.31E-04	2.12	Quantity of muscle cells	7.41E-04	0.99
Cell death of fibroblasts	1.05E-47	2.11	Mass of hind limb muscle	7.62E-09	0.91
Proliferation of cells	2.64E-08	2.10	Function of muscle	1.35E-06	0.89
Invasion of malignant tumor	1.02E-13	2.10	Breast or ovarían cancer	3.28E-13	0.86
Development of cardiovascular system	3.94E-04	2.06	Quantity of muscle	5.41E-05	0.86
Cell surface receptor linked signal transduction	3.70E-04	2.01	Abdominal aortic aneurysm	5.06E-04	0.83
Angiogenesis	6.11E-37	2.00	Stabilization of microtubules	1.79E-03	0.82

1 Z-score: Predicted activation status of biological function. The higher the value, the more activated the functions is predicted to be.

Effect of genetic type on Longissimus dorsi transcriptome

The effect of genetic type on LD muscle transcriptome was analyzed at birth and at growing independently, due to the high influence of developmental stage on gene expression. Genetic type effect on gene expression seemed to be stronger at birth, since 261 genes were DE (p < 0.01, FDR < 0.141) at that young age while less than a half, 113 DE genes were identified between IB and IBxDU growing pigs (p < 0.01, FDR < 0.199). Out of the 261 DE genes observed at birth, 130 were upregulated in IB (FC from 1.8 to 25.6) and 131 in IBxDU pigs (FC from 1.8 to 58.5). At growth, 88 genes were upregulated in IB (FC from 2.1 to 390) and 25 in IBxDU pigs (FC from 2.2 to 88.4) (S4 Table).

Biological interpretation of the DE genes with IPA software retrieved enriched pathways and biological functions in IB and IBxDU pigs at both studied developmental stages (Tables 5 and 6, respectively). Moreover, the regulators analysis identified 122 TRF at birth and 62 TRF at growth, potentially regulating the gene expression changes observed between genetic types. As described in the regulators study performed for the developmental stage effect, we considered TRF those that were either DE, identified by IPA software as regulators or identified in the RIFs study (S5 Table). Sixteen of those TRF (such as ATF4, CEBPA, MYOD1 (Myogenic Differentiation 1), NFE2L2 or REL) were found at both ages. Fifteen regulators (CREB3L1, CREBBP, CSRP3, FOS (FBJ Murine Osteosarcoma Viral Oncogene Homolog), FOXO1, HIRA, HSF1, KLF1 (Kruppel-Like Factor 1), MEF2C, MEF2D, MYOG, NFE2, SF1, SOX4 and TEAD3) were identified following different approaches at birth and 3 (EN1, IRF2 and TCF7L2) at growth. Moreover, a combined analysis performed using IPA software by merging DE genes and TRF datasets revealed the glucocorticoid receptor signaling and adipogenesis as the most enriched pathways at birth (p = 6.14e-11 and p = 1.41e-9, respectively), while the aryl hydrocarbon receptor pathway was the most enriched in growing pigs (p = 2.13e-8). The PPAR signaling and the unfolded protein response pathways were enriched at both stages (S6 Table).

Table 5

Enriched pathways in the set of DE genes conditional on genetic type at birth and at growing.

Birth			Growing
IB1	p-value	Corrected p-value	IB	p-value	Corrected p-value
Role of IL-17A in Psoriasis	2.83E-05	6.60E+00	Serine Biosynthesis	5.57E-09	8.08E-08
PI3K Signaling in B Lymphocytes	6.56E-05	6.60E+00	Superpathway of Serine and Glycine Biosynthesis I	1.29E-08	1.42E-07
ERK5 Signaling	7.65E-05	6.60E+00	Glutathione-mediated Detoxification	1.05E-04	5.14E+00
Melatonin Degradation III	1.87E-04	6.60E+00	LPS/IL-1 Mediated Inhibition of RXR Function	1.09E-04	5.14E+00
Glutamine Biosynthesis I	1.87E-04	6.60E+00	Alanine Degradation III	3.01E-04	5.14E+00
NRF2-mediated Oxidative Stress Response	6.00E-04	4.03E+01	Alanine Biosynthesis II	3.01E-04	5.14E+00
IL-8 Signaling	6.92E-04	4.03E+01	Cell Cycle: G1/S Checkpoint Regulation	2.02E-03	5.14E+00
April Mediated Signaling	1.47E-03	4.03E+01	Cyclins and Cell Cycle Regulation	4.96E-03	5.14E+00
Trehalose Degradation II (Trehalase)	1.81E-03	4.09E+01	Glucose and Glucose-1-phosphate Degradation	6.29E-03	5.14E+00
B Cell Activating Factor Signaling	1.81E-03	4.09E+01	p53 Signaling	1.62E-02	5.14E+00
LXR/RXR Activation	2.80E-03	4.09E+01	Leukotriene Biosynthesis	3.71E-02	5.14E+00
Atherosclerosis Signaling	3.13E-03	6.28E+01	LXR/RXR Activation	4.97E-02	5.14E+00
TNFR1 Signaling	4.68E-03	6.28E+01
Wnt/Ca+ pathway	9.09E-03	6.57E+01
Glucocorticoid Receptor Signaling	1.17E-02	6.57E+01
GDP-glucose Biosynthesis	1.52E-02	2.44E+02
Glucose and Glucose-1-phosphate Degradation	2.27E-02	3.24E+02
Growth Hormone Signaling	2.61E-02	4.33E+02
Melatonin Signaling	2.80E-02	4.41E+02
UDP-N-acetyl-D-galactosamine Biosynthesis II	3.22E-02	4.41E+02
Role of Macrophages, Fibroblasts and Endothelial Cells in Rheumatoid Arthritis	3.99E-02	4.41E+02
Production of Nitric Oxide and Reactive Oxygen Species in Macrophages	4.97E-02	4.41E+02
IBxDU2	p-value	Corrected p-value	IBxDU	p-value	Corrected p-value
AldosteroneSignaling in EpithelialCells	1.11E-10	1.02E-09	Acute Phase Response Signaling	3.90E-05	1.73E-02
ProteinUbiquitinationPathway	1.71E-09	2.09E-08	Retinoate Biosynthesis II	9.68E-05	1.73E-02
Unfoldedprotein response	5.94E-08	1.64E-06	Bupropion Degradation	1.03E-03	4.62E-02
eNOSSignaling	6.94E-06	1.47E-03	Acetone Degradation I (to Methylglyoxal)	1.47E-03	4.62E-02
EndoplasmicReticulum Stress Pathway	1.79E-04	1.91E-01	Retinoate Biosynthesis I	2.02E-03	4.62E-02
PyrimidineRibonucleotidesInterconversion	4.74E-04	1.91E-01	Retinol Biosynthesis	2.80E-03	4.62E-02
Glucocorticoid Receptor Signaling	5.45E-04	1.91E-01	Estrogen Biosynthesis	3.51E-03	4.62E-02
PyrimidineRibonucleotides De Novo Biosynthesis	6.29E-04	1.91E-01	Nicotine Degradation III	5.26E-03	4.62E-02
NeuregulinSignaling	6.29E-04	1.91E-01	Melatonin Degradation I	8.03E-03	4.62E-02
Role of p14/p19ARF in Tumor Suppression	7.27E-04	1.91E-01	Nicotine Degradation II	9.68E-03	4.62E-02
SpermineBiosynthesis	8.41E-04	1.91E-01	IL-17 Signaling	1.10E-02	4.62E-02
AlanineDegradation III	8.41E-04	1.91E-01	GABA Receptor Signaling	1.10E-02	4.62E-02
AlanineBiosynthesis II	8.41E-04	1.91E-01	Superpathway of Melatonin Degradation	1.33E-02	4.62E-02
Hepatic Fibrosis / HepaticStellateCellActivation	2.02E-03	5.58E-01	IL-6 Signaling	4.29E-02	4.62E-02
Spermidine Biosynthesis I	2.02E-03	5.58E-01
PI3K/AKT Signaling	4.96E-03	1.94E+00
Serine Biosynthesis	7.10E-03	2.93E+00
Thioredoxin Pathway	1.17E-02	4.26E+00
Aryl Hydrocarbon Receptor Signaling	1.17E-02	4.26E+00
Epithelial Adherens Junction Signaling	1.62E-02	5.32E+00
Superpathway of Serine and Glycine Biosynthesis I	1.85E-02	5.32E+00
PCP pathway	2.27E-02	6.37E+00
Hypoxia Signaling in the Cardiovascular System	2.80E-02	6.37E+00
Mitotic Roles of Polo-Like Kinase	3.00E-02	6.37E+00
Tight Junction Signaling	4.62E-02	1.03E+01

1 IB = Purebred Iberian pigs

2 IBxDU = Iberian x Duroc crossbred pigs

Table 6

Enriched biological functions in the set of DE genes conditional on genetic type, at birth and at growth.

Birth				Growth
IB1
Function	P-value	Corrected p-value	z-score2	Function	p-Value	Corrected p-value	z-score
Cell death of tumor cell lines	2.49E-06	5.45E-04	-2.82	Size of body	7.37E-03	4.89E-02	-2.356
Apoptosis	2.05E-09	4.03E-06	-2.60	Quantity of connective tissue	3.22E-03	3.47E-02	-1.987
Apoptosis of tumor cell lines	1.16E-07	1.15E-04	-2.57	Lesion formation	2.50E-04	2.92E-02	-1.353
Degradation of protein	2.63E-03	4.00E-02	-2.57	Concentration of lipid	2.45E-04	2.92E-02	-1.137
Necrosis	3.15E-07	2.48E-04	-2.54	Concentration of fatty acid	6.70E-03	4.48E-02	-1.091
Cell death	1.46E-10	5.75E-07	-2.49	Activation of macrophages	1.26E-02	5.59E-02	-1.091
Inflammatory response	2.21E-03	3.55E-02	-2.47	Size of lesion	2.94E-03	3.36E-02	-0.747
Accumulation of myeloid cells	7.89E-03	5.22E-02	-2.42	Size of bone	1.20E-02	5.44E-02	-0.682
Colony formation of tumor cell lines	1.49E-03	2.79E-02	-2.39	Proliferation of cells	1.24E-02	3.36E-02	-0.578
Proliferation of tumor cells	7.12E-03	5.22E-02	-2.35	Concentration of triacylglycerol	1.45E-03	4.93E-02	-0.412
Synthesis of reactive oxygen species	2.91E-04	1.00E-02	-2.35	Quantity of reactive oxygen species	8.63E-03	4.89E-02	-0.283
Colony formation of cells	3.24E-04	1.07E-02	-2.30
Accumulation of neutrophils	5.54E-04	1.45E-02	-2.20
Killing of bacteria	3.49E-04	1.12E-02	-2.19
Generation of reactive oxygen species	9.50E-04	2.04E-02	-2.13
Atrophy of muscle	2.97E-03	4.19E-02	-2.11
Differentiation of connectivetissue	1.24E-03	2.46E-02	-2.05
Metastasis of melanoma celllines	1.47E-04	6.57E-03	-2.00
IBxDU3
Quantity of erythroid progenitor cells	2.14E-03	3.55E-02	2.19	Organismal death	5.13E-03	4.13E-02	2.615
Cytolysis	2.99E-05	2.73E-03	2.18	Necrosis of epitelial tissue	3.53E-03	3.72E-02	1.450
Engulfment of cells	3.67E-03	4.74E-02	2.05	Metabolism of reactive oxygen species	1.09E-02	5.44E-02	1.261
Morphology of cells	7.37E-05	4.53E-03	2.05	Neuronal cell death	8.40E-03	4.93E-02	1.202
Mean corpuscular hemoglobin concentration	4.14E-04	1.22E-02	2.00	Transport of molecule	1.21E-02	5.44E-02	1.154
Cellviability of tumor celllines	4.72E-03	5.22E-02	1.99	Oxidation of lipid	2.74E-04	3.02E-02	1.127
Extension of neurites	3.92E-03	4.94E-02	1.98	Celldeath of epithelialcells	4.45E-03	4.13E-02	1.118
Cell viability	8.69E-05	4.96E-03	1.97	Colonyformation of tumor celllines	1.39E-02	5.96E-02	1.103
Quantity of reticulocytes	1.99E-03	3.36E-02	1.96	Hydrolysis of triacylglycerol	7.45E-06	7.39E-03	1.091
Cytosis	5.19E-03	5.22E-02	1.94	Morbidity or mortality	1.15E-02	5.44E-02	1.067
Cell survival	3.18E-05	2.78E-03	1.86	Inflammation of body cavity	6.07E-03	4.13E-02	1.000
Anemia	1.58E-05	1.81E-03	1.80	Cell death	5.47E-04	3.36E-02	0.980
Phagocytosis of cells	1.49E-03	2.79E-02	1.80	Cancer	7.95E-03	4.93E-02	0.825
Hyperplasia of epidermis	8.90E-05	5.00E-03	1.71	Hypertrophy	1.08E-02	5.44E-02	0.808
Size of cells	8.43E-03	5.22E-02	1.60	Cell death of central nervous system cells	7.74E-03	4.93E-02	0.791
Binding of cells	2.60e-03	3.97E-02	1.49	Mineralization of cells	1.06E-04	1.81E-02	0.762
Binding of granulocytes	2.43e-04	8.61E-03	1.47	Apoptosis	3.24E-03	3.47E-02	0.731
Mass of muscle	1.96e-03	4.53E-02	1.38	Inflammation of organ	2.73E-03	3.36E-02	0.712

1 IB = Purebred Iberian pigs

2 Z-score: Predicted activation status of biological function. The higher the value, the more activated the functions is predicted to be.

3 IBxDU = Iberian x Duroc crossbred pigs

Effect of muscle type, Biceps femoris or Longissimus dorsi, on gene expression

Differences in gene expression were observed between LD and BF muscles, 83 genes showing higher expression levels in LD than in BF muscle (FC from 1.7 to 27.0), while 52 genes were upregulated in BF muscle (FC from 1.8 to 183.2) (p < 0.01, FDR < 0.123). Genes such as HOXA11, PVALB (Parvalbumin) or CXCL13 (upregulated in BF muscle) and IBSP, ZIC1 and MMP13 (upregulated in LD muscle) showed the largest expression differences (S7 Table).

IPA software was used to detect enriched pathways (S8 Table) and biological functions (S9 Table) associated with the DE genes between both muscles. Moreover, the regulators study identified 370 regulatory genes. Among them, 233 were identified by IPA software, 143 were identified following the RIFs approach and 10 were found DE between developmental stages (S10 Table). SIM1 was identified in the three different analysis carried out, which highlights its importance and potential role in regulating muscle gene expression in IB and IBxDU pigs.

In order to validate the results obtained from the RNA-Seq analysis, the relative expression of five DE genes affected by the genetic type, the growing stage or the studied muscle was assessed by qPCR in all the available samples. All of them were validated by qPCR (S11 Table). A concordance correlation coefficient, used to assess technical validation in high throughput transcriptomic studies [36] was calculated (CCC = 0.94) and denoted a high general concordance between RNA-Seq and qPCR expression values. In general good individual correlation values were obtained (S11 Table).

Discussion

Effect of developmental stage and terminal sire line on Iberian pig phenotype

Developmental stage affected carcass characteristics, as expected: all measured parameters increased along time (Table 2). On the other hand, genetic type affected all the carcass phenotypic parameters in newborns: IBxDU were bigger and heavier (p < 0.001) than IB piglets (Table 2). However, growing IB and IBxDU pigs only differed in ham weight and circumference, but no significant difference was found in other body size measures. Differences in ham measures between genotypes observed in juvenile pigs agree with results obtained in pure and Duroc-crossbred Iberian pigs at final slaughter weight [40, 41]. Previous studies also reported that adult Duroc-crossbred Iberian pigs are longer and heavier than their purebred Iberian counterparts [40, 41], but these differences in body weight and size between genotypes are not evident at weaning [29]. Iberian pigs have a lower prenatal growth potential than lean pigs [42]. Thus, the finding of IB newborns being smaller in body size and weight than IBxDU is expected, due to the effect of the leaner and more growth-efficient Duroc sire. Nevertheless, the lack of body weight differences between genotypes at weaning and growing stages is unexpected.

Iberian pigs show some special characteristics, associated with the phenomenon known as thrifty genotype, which might led to differences in voluntary feed intake and energy expenditure between IBxDU and IB pigs [29].Such differences may lead to a compensatory gain in the latter during the suckling period, which would not be reflected on adult pigs due to a much lower growth potential of pure Iberian animals in later stages. However, this hypothesis cannot be tested under the design of the present experiment and further studies specially designed to analyze the evolution of growth and energy balance in different genotypes along early and juvenile growth would be required.

Regarding biochemical parameters, plasma glucose levels were higher in neonates, probably because growing animals were slaughtered after a fasting period, which may deplete plasma glucose levels. At growing, IBxDU pigs showed higher plasma glucose levels than IB pigs, in agreement with the reported decreased plasma glucose levels in obese when compared to lean fasted pigs, probably because of a high rate of glucose utilization for fat synthesis in obese pigs [43]. Fructosamine levels were higher in older pigsin agreement with previous studies [44].

Certain lipid metabolism-related parameters (cholesterol, HDL and triglycerides) were significantly different between IB and IBxDU newborns, but seemed to balance in growing pigs, with a significant interaction observed between developmental stage and genotype for these parameters: newborn IB pigs showed greater plasma cholesterol, HDL and triglycerides than IBxDU piglets [30], but these differences disappeared over time. The evolution over life span in lipid related plasma indicators has not been previously compared between different pig breeds. However, in a study assessing plasma biochemical parameters in wild boars, those under 6 months of age showed slightly higher cholesterol and doubled triglyceride levels than wild boars over 6 months of age [45], while an increase as pigs mature in the concentrations of cholesterol has been reported in lean pigs [46]. Thus, a different regulation of lipid metabolism in swine breeds over life might be responsible for the observed differences in cholesterol, HDL and triglycerides, in agreement with previous findings. Also, the stability of lipid metabolism-related parameters may be associated with a state of ‘healthy’ obesity (without metabolic disorders) that purebred Iberian pigs can maintain under intensive production, whereas crossbred Iberian pigs might accumulate more cholesterol due to the concentrated-based diet provided during the starter and growing production stages.

Intramuscular fat content was also affected by pigs’ growth:total IMF content increased over time (p < 0.0001) in both genetic types (Table 2). In newborns, no differences in loin IMF content were observed between IBxDU and IB pigs, in contrast with results reported for muscle BF of newborn IB piglets, which showed almost 30% more IMF than BF of IBxDU pigs [30]. This supports that regulation of IMF deposition depends, among other factors, on muscle type [47, 48]. However, growing IB pigs showed greater IMF content than growing IBxDU pigs in LD muscle, in concordance with the characteristic difference between genetic types in adult animals [13].

Regarding IMF composition, MUFA content increased over time (p < 0.0001) in both genetic types. The increase in MUFA content was due to an increase in the most abundant fatty acid, oleic acid (C18:1 n-9). Regarding the genetic type effect on IMF composition, a slight difference was observed between genetic types at birth, while a stronger effect was observed at growth. The most remarkable effect was an increase in n-3 fatty acids and a decrease in the ratio n-6/n-3 in IB when compared to IBxDU growing pigs. A decrease in the ratio n-6/n-3 and increase of DPA and DHA fatty acids may promote consumers’ health [49-51]. These results indicate that crossing with Duroc sires decreased meat quality in terms of consumers’ health and IMF concentration in LD muscle of growing pigs, in agreement with differences observed in weaned and adult pigs [13, 29].

Effect of muscle type on intramuscular fat content and composition

We observed higher PUFA content in BF when compared to LD muscle, which was in part due to increased n-3 PUFA content in BF. This also led to a decrease in the ratio n-6/n-3. It is known that the pattern of fatty acids deposition may differ across muscles [15, 18, 52]. In agreement with the results of the present study, the biggest difference previously observed between LD and BF muscle of lean pigs was reported for PUFA content [26], which might be due to the differences in oxidative properties observed between muscles [19, 20]. As previously discussed, a lower ratio n-6/n-3 has been associated with healthier meat. Moreover, when analyzing jointly the data of both muscles, the inclusion of Duroc on Iberian genetics significantly increased SFA content in muscle of IBxDU newborns (Table 3), in accordance with previous results [13]. High SFA consumption increases cardiovascular disease risk in human [53]. Thus, meat from BF muscle, the main muscle within the ham and from pure Iberian pigs, could be considered healthier than meat from LD muscle and from Duroc-crossbred Iberian pigs.

Three main effects (developmental stage, muscle type and genetic type at birth and at the growing stages) were evaluated in the present study.

These effects were assessed in two muscles and time points. Longissimus dorsi and BF muscle were assessed because they represent the most important pork cuts, the loin and the ham. On the other hand, the two time points were selected due to the interest of comparing the initial (birth) and an intermediate (four months of age) stages during the growing period, which in traditional Iberian pig production is considered up to eight months of age [4].

Several genes showed changes in expression according to more than one of the studied effects (Fig 1). Developmental stage affected a high number of genes which were concomitantly affected by muscle or genetic type, probably due to its large impact on muscle transcriptome. However, few common genes were observed between muscle and genetic type effects and no gene was observed to be affected across the three studied conditions. Fourteen common genes (nine of them being known genes) were affected by genetic type in both newborn and growing pigs, and will be further discussed because of their potential phenotypic impact.

Fig 1

Venn diagram showing common differentially expressed genes across studied conditions.

Developmental stage, the biggest studied effect, shows several common differentially expressed genes with other effects, while expression of any gene was affected by the three studied effects.

The developmental stage was the main factor affecting gene expression in pure and Duroc-crossbred Iberian pigs, as 5,800 genes changed their expression levels between both developmental stages. Gene expression has been previously reported to change across age in pigs, especially during early stages of prenatal and postnatal development [54]. Some of the genes showing the highest upregulation in newborn piglets were involved in the development of different tissues, such as bone, cartilage, adipose or muscle tissues (ACTC1, ARHGAP36, IBSP, TNN, ATP6V0D2, COMP, FGF21, DLK1 and several myosin and collagen proteins). Moreover, two genes associated with meat quality were highly upregulated at birth: RETN is associated with adipogenesis and lipogenesis and TNN is involved in the development of the extracellular matrix in pigs [55, 56]. Genes highly upregulated in growing pigs were associated with the immune response (since several immunoglobulin genes such as IGLC, IGLV7, 8, 9 and 10, were identified), but also with protein metabolism (PVALB, UBD). A functional interpretation of the DE genes upregulated in each growing period was carried out. Several metabolic pathways were enriched at both developmental stages, most of them involved in muscle growth (Wnt/β-catenin Signaling, Calcium Signaling, Signaling by Rho Family GTPases or Actin Cytoskeleton Signaling), in agreement with the studied tissue and growth period. Pathways enriched in newborn piglets are involved in cholesterol, triglycerides and other compounds biosynthesis, characteristic of a highly proliferative developmental stage [57, 58]. This is concordant with other enriched functions related to cellular growth and anabolic processes (as invasion of cells, transport of molecule, cell movement, adhesion of connective tissue cells, proliferation of fibro blasts and synthesis of carbohydrate). Contrarily, normal non-proliferating or differentiated cells primarily utilize nutrients to fuel basic cellular processes and predominantly mediate catabolic metabolism to efficiently generate ATP [57, 58], as observed in growing pigs. At this age, enriched pathways were mainly related to catabolic processes (Protein Ubiquitination Pathway, Glycolysis I, Gluconeogenesis I, Glucocorticoid Receptor Signaling and Phospholipase C Signaling) andbiological functions represented a decrease in developmental and growth processes. Moreover, functions such as contractility of skeletal muscle, mass of muscle, quantity of muscle cells and function of muscle, appeared enriched in growing pigs; this suggests a more developed and functional muscular system in older than in newborn piglets [54].

Changes in gene expression and thus, phenotypic consequences, between developmental stages were predicted to be regulated by a number of TRFs. Identified regulators affecting muscle development include the myogenesis regulators MYOD1, MEF2C and MEF2D, essential ones at different stages of muscle growth [54, 59, 60]. The forkhead family members FOXO1 and FOXO3 play a role in muscle but also in adipose tissue development [61-64], while PPARG, PPARGC1B, SIM1, ATF4 and CEBPD regulate adipocyte differentiation, energy homeostasis and lipid metabolism [65-69]. Common genes identified by the different followed procedures were of special interest. Specifically, HSF1, JARID2, TCF4 and YBX3 play known roles in muscle development and protein metabolism [70-74], while TFEB and MYC are involved in adipocyte proliferation, lipid metabolism and fatty acid transport [75, 76]. These TRF potentially regulate muscle and fat accretion in young pigs and are thus, of special interest in the understanding of molecular mechanisms underlying such processes in pigs and other species.

Genetic type significantly affected gene expression over pig´s growth. Nine known genes (GPT2, PSAT1, ART5, ADAMTS8, KCNH2, RASSF9, TP63, ENV and ASB5) were found DE between genotypes at both developmental stages, suggesting an important role of these genes. Interestingly, three of them (GPT2, PSAT1 and ADAMTS8) were differentially regulated at both ages, being upregulated in IBxDU at birth and in IB at four months. GPT2 is involved in gluconeogenesis, fatty acids oxidation and amino acid metabolism [77, 78], while PSAT1 catalyzes serine biosynthesis. ADAMTS8 has an important role in inhibition of angiogenesis [79]. Expression changes observed in these genes, related to muscle growth and metabolism, may suggest a differential muscle growth regulation depending on developmental stage and genetic type, and are in agreement with the observed phenotypic results.

Differentially expressed genes and enriched pathways and functions were of special interest when related to processes that may drive phenotypic differences observed between IB and IBxDU pigs. We found two main processes that seemed to be affected by the genetic type:

Muscle growth:

Genes showing strong upregulation in newborn IBxDU pigs were associated with the extracellular matrix structure (MATN1 and 3 or COL9A1 and 2), connective tissue and muscle development (GDF5, MYH10) and also with protein metabolism and degradation (PVALB, HSPs). MYH10 was also upregulated in IbxDU pigs in two previous studies comparing muscle transcriptome of newborn [30] and weaning [29] IB and IbxDU pigs, indicating a relevant role for this gene in muscle development differences between genotypes. These results are concordant with the greater prenatal development and with the enriched pathways (Table 5) in crossbred newborns.

In growing pigs, pathways involved in the metabolism of non-essential amino acids such as serine, glycine and alanine were enriched in the IB genotype. Those amino acids are necessary for synthesis of proteins and other biomolecules needed for cell proliferation [80], suggesting an active protein synthesis in growing IB pigs. In accordance with this, several genes involved in body size and cell proliferation were upregulated in IB pigs (Fig 2), probably associated with the increased body growth that the smaller IB neonates could have developed, leading to the observed similar body weight in growing IB and IbxDU pigs.

Fig 2

Enriched biological functions related to body growth in growing Iberian (IB) pigs.

The network generated by IPA software shows enriched biological functions in IB pigs (blue color) and genes predicted to be involved in enrichment of these functions. The activation of cell proliferation and body size in growing IB pigs might be related to a compensatory growth occurring in those pigs when compared to crossbred (IBxDU) pigs during the early growing period.

Muscle growth is not only determined by cell proliferation, but also by protein synthesis and degradation and angiogenesis. Protein degradation seems to be an active process in both IB and IbxDU newborn piglets. In the present study, IbxDU pigs showed upregulation of genes (ELANE, MMP9, FBXO32, PVALB, HSPS1, HSPA4L and DNAJA1) and pathways (Table 4) related to protein degradation, although enrichment of muscle degradation or atrophy functions was observed in IB (Fig 3) but not in IbxDU pigs (Table 5), similarly to the results observed in BF muscle [30]. This suggests greater muscle degradation in newborn IB than IbxDU pigs [81]. Accordingly, genes showing high upregulation in four months-old IB pigs were mainly related to protein turnover and degradation (CTSF, ADAMTS8 or CELA2).

Fig 3

Enriched biological functions related to growth and development in newborn Iberian (IB) piglets.

The activation of muscle atrophy and the upregulation of genes involved in protein metabolism and degradation in IB pigs suggests a more active protein turnover in these pigs when compared to crossbred (IBxDU) pigs. Functions associated with animal growth were also enriched in IB newborn pigs.

Energy homeostasis, inflammation and immune system:

Energy homeostasis is tightly regulated in animals due to its importance for normal growth and survival. Pathways involved in cholesterol (LXR/RXR Activation, Atherosclerosis Signaling) and glucose metabolism (Glucocorticoid Receptor Signaling, GDP-glucose Biosynthesis and Glucose and Glucose-1-phosphate Degradation) were enriched in both newborn and growing IB pigs, suggesting an increased energy metabolism in pure Iberian piglets.

A pathway involved in the control of energy homeostasis, cell metabolism and muscle development, the Wnt/Ca+ signaling pathway [82] was enriched (p = 0.009) in newborn IB pigs. This pathway included genes upregulated in IB piglets, which were involved in adipogenesis and in the development of obesity, as NFATC1 and PLCD1 [83, 84]. In agreement, biological functions related to glucose metabolism and lipid accumulation, as concentration of lipid, concentration of fatty acid and concentration of triacylglycerol (Fig 4), were enriched in growing IB pigs. Moreover, the LPS/IL-1 Mediated Inhibition of RXR pathway, enriched in those pigs was reported to positively correlate with fat area [85].

Fig 4

Enriched biological functions related to lipid metabolism in growing Iberian (IB) pigs.

The network generated by IPA software shows metabolism-related enriched biological functions in IB pigs. Concentration of lipid seems to be a more active process in IB pigs, in agreement with phenotypic differences observed in loin intramuscular fat content.

On the other hand, the glucocorticoid receptor signaling pathway, was enriched in both IBxDU and IB newborns, which might be related to the farrowing stress and with the wide range of actions associated to this pathway, from catabolic processes to glucose and energy homeostasis or adipocyte differentiation [86, 87]. However, the upregulated genes in IBxDU pigs associated to this pathway are members of the HSPs family, and thus, its activation might be a consequence of cellular stress, or activated protein catabolism. On the other hand, DE genes involved in this pathway in IB pigs play important roles in protein catabolism (FOXO3A, [64]), but also in lipid metabolism and adipogenesis (NFATC1, FOS, FOXO3A, [83, 88, 89]) and in osteogenesis and glucose uptake (BGLAP, [90]). Thus, although the glucocorticoid receptor signaling pathway was enriched in both IB and IBxDU piglets, a different set of DE genes was involved in each genetic type and thus, different metabolic consequences might be expected.

The juvenile IBxDU pigs showed enrichment of pathways mainly involved in biosynthesis and degradation of retinoids, such as retinoate and retinol or melatonin, involved in the synchronization of the circadian clock, associated with the control of energy homeostasis, among a wide range of functions [91]. The enrichment of these functions is concordant with the upregulation of genes involved in degradation of compounds such as steroids and fatty acids (CYP1A1) [92]. The GABA (Gamma-Aminobutyric Acid Type A) receptor signaling pathway was also found enriched in growing IBxDU pigs. Beyond its function as an inhibitory neurotransmitter [93], the expression of GABA receptor in human muscle was associated with increased resting energy expenditure [94]. Thus, IBxDU pigs might have a greater basal energy expenditure that would decrease their potential for fat accumulation, in agreement with enriched biological functions and with the lower IMF content observed in IBxDU pigs.

Regulators analysis:

We performed a regulatory factors study to investigate the driving molecular mechanisms responsible for the differences in gene expression observed between genetic types. Three different approaches, based on bibliographic (IPA software) and expression and co-expression (RIFs analysis) information were combined as a powerful strategy to identify overlapping TRFs [30].

A total of 122 TRF were identified in newborn and 62 in older pigs. Among them, 16 TRFs were identified at both developmental stages. These common TRF would be expected to have a deeper impact in the final phenotype and thus, should be considered as strong candidate genes driving phenotypic differences in pure and crossbred Iberian pigs. Some of them, such as MYOD1, a well-known myogenic regulator, and some TRF recently associated with its regulation (BHLHE40 and HDAC2; [95, 96]) are necessary for muscle development. Similarly, CTNNB1 (Catenin Beta 1) is a component of the Wnt signaling pathway, related to cell differentiation and metabolism, including myogenesis and muscle regeneration in adult animals [82]. Also, TRF involved in the immune response (IRF9, NFKBIA, REL) were found to affect gene expression in IB and IBxDU newborn and growing pigs. NFKBIA has a critical role in the upregulation of pro-inflammatory factors and is considered a link between immunological stress and obesity [97, 98]. The identification of TRF such as NFKBIA, ATF4 or CEBPA, with direct roles in adipogenesis, and fat accumulation [66, 99] in an obese pig breed such as the Iberian pig [29] is of high interest to understand regulatory mechanisms responsible for fatness regulation.

Regarding the developmental stage-specific TRF identified, 15 regulator genes were identified using different approaches at birth. It is noteworthy that some of them (MEF2C, MEF2D, MYOG, SOX4) are important TRF involved in muscle cell differentiation [59, 100]. In addition, TRF related to protein degradation (CREB3L1, HSF1 and CREBBP) and to blood cells differentiation (KLF1, SF1 and NFE2) were also identified and may play a role in muscle development, due to the need of protein degradation and blood supply in the growing muscle. Moreover, KLF1 was strongly upregulated in IBxDU pigs (25.61), which suggest an important role of this TRF. Other important TRFs identified in the performed analysis were FOS and FOXO1, involved in muscle differentiation and metabolism [63], but also in the regulation of adipogenic genes (as PPARG) expression and IMF accumulation [54, 62]. These TRFs were also previously identified in BF muscle of newborn piglets [30].

In growing pigs, only 3 TRFs (EN1 (Engrailed Homeobox 1), IRF2 (Interferon Regulatory Factor 2) and TCF7L2 (Transcription Factor 7 Like 2)) were identified simultaneously using different approaches. EN1 and IRF2 play important roles in regulating development and cell cycle. In muscle cells, IRF2 activates transcription of the VCAM-1 gene, suggested to play a role in the differentiation of skeletal muscle [101]. EN1 and TCF7L2 interact with the Wnt signaling pathway [102], which negatively regulates adipogenesis [103]. Variants in TCF7L2 gene have been reported to affect insulin secretion and body mass index and to promote type II diabetes [104, 105]. Thus, it seems that lipid metabolism might be more tightly regulated in growing IB and IBxDU pigs.

To better understand the role of the identified regulators on gene expression and on phenotypic differences, information from the DE and regulators analyses was used for biological interpretation, focusing on enriched metabolic pathways (S5 Table). The glucocorticoid receptor signaling and the adipogenesis (Fig 5) were the most enriched pathways at birth. Changes in regulation and function of adipogenesis pathway in newborn piglets may determine the differences in fatness observed in LD muscle of growing IB and IBxDU pigs, in accordance with results observed in BF muscle of newborn piglets [30]. The glucocorticoid receptor signaling may be also implicated in such differences, due to its roles in energy homeostasis and adipocyte differentiation [82, 83].

Fig 5

Adipogenesis pathway.

Several transcription factor identified in the regulators study and one gene upregulated in Iberian pigs are involved in the adipogenesis pathway, most of them regulate PPARg expression or activation.

The aryl hydrocarbon receptor was the most enriched pathway in growing pigs. This pathway is involved in xenobiotics metabolism, although a role in inhibition of lipid biosynthesis and adipocyte differentiation has also been reported [106].

A total of fifty-one pathways were enriched at both stages. Among them, the PPAR signaling, the aryl hydrocarbon receptor, the adipogenesis, the Wnt and the unfolded protein response pathways are particularly interesting because these are specifically involved in adipocyte differentiation and protein degradation. Moreover, finding such pathways at two different growing stages, suggests a long-term activation that may also drive differences in pure and crossbred Iberian pigs in the adulthood.

Effect of muscular tissue, Biceps femoris or Longissimus dorsi, on gene expression at birth

Differences in the transcriptomic profile were observed between BF and LD of newborn IB and IBxDU pigs. However the effect was smaller than that reported when comparing LD and Semimembranosus muscles [27], probably due to the sampling at an early age, when muscle fiber type is still not well determined [107, 108].

In the present study, 83 genes were upregulated in LD muscle, some of them (IBSP, ZIC1 and MMP13) showing large expression differences. ZIC1 and MMP13 are implicated in the early control of myogenesis and in myostatin signaling [109-112]. On the other hand, genes highly upregulated in BF muscle are involved in myogenesis control (HOXA11), by regulating MYOD expression [113], muscle contraction (PVALB), or the immune response and adipocyte differentiation, [114]. PVALB gene expression was also deeply affected by developmental stage (larger expression in growing pigs) and by genetic type in newborn pigs, being upregulated in both BF and LD muscles of IBxDU pigs [30]. The finding of this DE gene across tissues, developmental stages and genetic types suggests an active role of PVALB gene on phenotypic changes.

The biological interpretation retrieved several enriched pathways in LD muscle related to lipid metabolism (Atherosclerosis Signaling and VDR/RXR Activation) muscle development and function (Calcium Signaling, Inhibition of Matrix Metalloproteases and Actin Cytoskeleton Signaling) and to the immune response (Agranulocyte Adhesion and Diapedesis, Leukocyte Extravasation Signaling, Granulocyte Adhesion and Diapedesis and IL-8 Signaling). In agreement with the enriched pathways, several biological functions related to inflammation and immune response were predicted by IPA software to be activated in LD muscle. On the other hand, functions such as proliferation of cells or size of body, related to cell growth and development were enriched in LD muscle of newborn pigs. In agreement, biological functions (neovascularization, development of cardiovascular system and vasculogenesis) and pathways (HIF1α Signalling) related to the supply of oxygen and energy to growing cells, were also enriched in LD muscle suggesting a more active cellular growth.

On the other hand, BF muscle showed enriched pathways involved in adipocyte differentiation and lipid metabolism (LXR/RXR Activation, VDR/RXR Activation, Atherosclerosis Signaling, FXR/RXR Activation and biosynthesis of retinoids, bile acids and thyroid hormone; [115-118]), in agreement with the enrichment of efflux of cholesterol in BF muscle. This suggests a more active lipid metabolism in BF muscle. In agreement, quantity of connective tissue and quantity of carbohydrate functions were also enriched in BF when compared to LD muscle (Fig 6), probably associated to a different energy homeostasis in both muscles. A more active lipid metabolism in BF muscle might allow a deeper impact of the genotype on lipid deposition, leading to the observed significant differences in IMF content between pig genotypes in BF, which are not evident in LD muscle of neonates. Moreover, these differences in lipid metabolism may be associated with the higher IMF content reported in BF than in LD muscle in adult pigs [47].

Fig 6

Enriched biological functions in Biceps femoris (BF) muscle.

The network generated by IPA software shows enriched biological functions in BF muscle, suggesting that lipid metabolism is more active in BF when compared to Longissimus dorsi muscle.

Beyond metabolic pathways and biological functions, we performed a regulatory genes search in the same way as for the other main effects. Among the 370 identified regulators, SIM1 was considered the most significant regulator, as it was highlighted in the three followed approaches. Although no information exists regarding SIM1 gene effects on pigs, the differential expression across tissues, developmental stages and genetic type in newborn piglets, suggest this TRF as a strong candidate gene. Previous information on its effects on adiposity, energy expenditure and food intake in other species such as mouse or humans [119-122], supports its importance as a regulator of fatness traits.

Conclusions

The present study reports the effects of developmental stage, muscle and genetic type on animal phenotype, muscle transcriptome, metabolic pathways and transcriptional regulation, associated with traits of interest. Developmental stage represented the most drastic influence on phenotype and transcriptome. Newborn pigs showed enrichment of anabolic functions, while predominant functions in the growing stage were related to catabolism and muscle functioning, indicating a decrease in developmental and growth processes and a more advanced muscle development in juvenile pigs. Moreover, phenotypic differences regarding body size and plasma biochemical parameters were observed between genetic types at birth but dramatically decreased at growing. This suggests strong differences in early growth patterns and metabolism between them, in spite of the closely related analyzed genotypes. In agreement, IMF differences between genotypes also depended on developmental stage and muscle. Gene expression results support the phenotypic findings, as DE genes and pathways suggest a different timing in growth, proliferative and anabolic processes. Those processes were upregulated in IBxDU newborn pigs (associated with a higher capacity for prenatal growth) and in growing IB pigs (in agreement with a potential compensatory gain during the postnatal period). Differences in metabolism were also observed, and results suggest a more active lipid and glucose metabolism in both newborn and growing IB pigs, in agreement with their greater potential for fat accumulation. An effect of muscle type on muscular metabolic characteristics was observed, as BF muscle showed increased lipid metabolism, while LD was characterized by growth and proliferative processes. The regulatory factors analysis identified several remarkable TRF (as MYOD1, NFKBIA, FOXO1, MEFs, TCF7L2, SIM1 and PVALB), selected due to their identification following different methodological approaches, their identification across different developmental stages and muscles and their potential roles on regulating molecular processes underlying differences in metabolism and productive traits between IB and IBxDU pigs.

Differentially expressed genes conditional on developmental stage (birth vs growth).

(XLSX)

Click here for additional data file.

Enriched pathways in the set of DE genes conditional on developmental stage (birth vs growth).

(XLSX)

Click here for additional data file.

Transcription factors affecting Longissimus dorsi gene expression differences between newborn and growing Iberian pigs.

(XLSX)

Click here for additional data file.

Differentially expressed genes conditional on genetic type (Iberian (IB) vs Duroc X Iberian (IBxDU)) at birth and at growth.

(XLSX)

Click here for additional data file.

Transcription factors affecting Longissimus dorsi gene expression differences between pure and Duroc-crossbred Iberian pigs at birth and at growing stages.

(XLSX)

Click here for additional data file.

Enriched pathways in the set of differentially expressed genes and transcription factors conditional on genetic type at birth and at growth.

(XLSX)

Click here for additional data file.

Differentially expressed genes conditional on muscle (Biceps femoris (BF) vs Longissimus dorsi (LD)).

(XLSX)

Click here for additional data file.

Enriched pathways in the set of differentially expressed genes conditional on muscle type (Longissimus dorsi (LD) vs Biceps femoris (BF)) at birth and growth.

(XLSX)

Click here for additional data file.

Enriched biological functions in the set of differentially expressed genes conditional on muscle type (Longissimus dorsi (LD) vs Biceps femoris (BF)) at birth and growth.

(XLSX)

Click here for additional data file.

Transcription factors affecting differences in gene expression between Longissimus dorsi and Biceps femoris muscles from Iberian pigs at birth.

(XLSX)

Click here for additional data file.

RNA-Seq and qPCR validation results and correlation coefficient (r) between the two used methodologies

(XLSX)

Click here for additional data file.