Gene expression studies of developing bovine longissimusmuscle from two different beef cattle breeds

Temporal changes in gene expression

Tables 1 and 2 list expression values for differentially expressed (DE) genes for which reliable annotations could be obtained using the Basic Local Alignment Search Tool (BLASTN [15]). These genes were differentially expressed across developmental time, regardless of breed. Twenty genes showed decreasing expression over time and 35 genes showed increasing gene expression over time. All genes in Tables 1 and 2 showed expression changes of more than 2-fold over the time course examined. Gene Ontology annotation of the DE genes using EASE [16] showed that genes involved with the biological processes, morphogenesis and cellular growth were overrepresented in the genes that show decreasing expression over time (Figure 1A). When annotations with respect to cellular component were examined, 69% of the genes with decreasing expression could be classified as belonging to "extracellular matrix" (data not shown). Collagen genes alone represent 68% of the DE genes in Table 1. On the other hand, in the list of genes with increasing expression over time, genes with biological process annotations belonging to muscle cell contraction and metabolism are overrepresented (Figure 1B). The main observation from the examination of the developmental time course is therefore that transcripts associated with extracellular matrix synthesis are increasingly less represented as a proportion of total tissue RNA as muscle development proceeds, while at the same time transcripts of muscle fiber structural, contractile and metabolism genes are gaining in prevalence when maturation of muscle fibres takes place in preparation for postnatal life. In most instances, the most marked increase or decrease in gene expression can be observed between d 195 of fetal development and birth (Tables 1 and 2).

Gene2	GenBank Accession3	Number of array elements4	Gene expression ratio5 135 d/60 d	Gene expression ratio5 195 d/60 d	Gene expression ratio5 birth/60 d
cadherin 11, type 2, OB-cadherin (osteoblast) (CDH11)	DW521798	1	0.26**	0.23**	0.10**
collagen, type I, alpha 1 (COL1A1)	CF613531	14	0.87	0.74	0.16**
collagen, type I, alpha 2 (COL1A2)	CF613950	9	0.75	0.64**	0.15**
collagen, type III, alpha 1 (COL3A1)	CF613798	29	0.72	0.78	0.15**
collagen, type XII, alpha 1 (COL12A1)	CF614200	1	0.48**	0.23**	0.12**
collagen, type XV, alpha 1(COL15A1)	DW521810	1	0.58**	0.41**	0.18**
fibrillin 1 (FBN1)	CF765332	4	0.57**	0.60**	0.23**
fibronectin 1 (FN1)	CF614211	7	0.39**	0.36**	0.08**
follistatin-like 1 (FSTL1)	CF613847	1	0.61**	0.35**	0.13**
guanine nucleotide binding protein (G protein), alpha inhibiting activity polypeptide 1 (GNAI1)	DW521797	1	0.58**	0.50**	0.20**
glypican 3 (GPC3)	DW521793	1	0.91	0.84	0.11**
Insulin-like growth factor binding protein 5 (IGFBP5)	CF614176	1	0.71	0.72	0.21**
lumican (LUM)	CF613654	2	0.79	0.64**	0.16**
osteoglycin (OGN)	CF614268	1	0.59**	0.67*	0.19**
plastin 3 (T isoform) (PLS3)	DW521782	1	0.49**	0.43**	0.20**
osteonectin (SPARC)	CF613734	8	0.69*	0.68*	0.21**
osteopontin (SPP1)	DW521777	3	0.53**	0.19**	0.15**
stathmin 1/oncoprotein 18 (STMN1)	DW521774	1	0.54**	0.38**	0.10**
tubulin, alpha 1 (TUBA1)	CF613655	2	0.54**	0.46**	0.17**
vimentin (VIM)	CF613639	7	0.51**	0.49**	0.16**

¹* indicates a statistically significant decrease in gene expression at P < 0.05 level; ** indicates a statistically significant decrease in gene expression at P < 0.01 level.

²Sequence annotation was carried out using approved gene symbols (human genome nomenclature) via the IBISS4 database [17]. BLAST [15] scores > 110 and e-values < 2 E-23 were accepted for annotation purposes.

³in cases where more than 1 microarray element was used to calculate the expression data, the accession number of 1 representative element is shown.

⁴Number of microarray elements representing the same gene that were used to calculate the expression data.

⁵To obtain the ratios shown, averaged absolute microarray signal intensity values for d 135, 195 and birth were divided by the signal intensities measured at d 60.

Gene2	GenBank Accession3	Number of array elements4	Gene expression ratio5 135 d/60 d	Gene expression ratio5 195 d/60 d	Gene expression ratio5 birth/60 d
actin, alpha 1, skeletal muscle (ACTA1)	CF614694	6	0.80	2.49**	4.09**
actinin, alpha 3 (ACTN3)	CF614924	3	0.62	3.79**	5.63**
adenylate kinase 1 (AK1)	CF614671	1	1.38**	4.25**	4.60**
aldolase A, fructose-biphosphate (ALDOA)	CF615006	1	1.28*	3.27**	5.82**
ankyrin repeat domain 1 (cardiac muscle) (ANKRD1)	CF614403	2	0.83	2.62**	6.62**
mitochondrially encoded ATP synthase 6 (MT-ATP6)	CF613499	1	1.10	2.69**	5.48**
chaperone, ABC1 activity of bc1 complex like (S. pombe) (CABC1)	CF615094	1	0.51	1.18	6.43**
calsequestrin 1 (fast-twitch, skeletal muscle) (CASQ1)	CO729180	1	1.08	3.93**	5.00**
creatine kinase, muscle (CKM)	CF614901	15	2.15**	3.04**	7.49**
creatine kinase, mitochondrial 2 (sarcomeric) (CKMT2)	CF614479	1	0.60	1.76**	10.08**
cardiomyopathy associated 5 (CMYA5)	CF615342	4	0.77	1.63**	4.79**
cold shock domain protein A (CSDA)	CF615315	3	0.64	1.74**	4.71**
cytochrome C (CYC1)	CF614629	1	2.14**	2.53**	6.90**
enolase 3 (beta, muscle) (ENO3)	CF614534	16	1.93**	3.58**	5.25**
fructose-1,6-bisphosphatase 2 (FBP2)	CF614569	2	0.53	2.39**	6.73**
heat shock 70 kDa protein 1B (HSPA1B)	CF614728	1	1.03	2.17**	4.02**
heat shock 90 kDa protein 1, alpha (HSPCA)	DW521445	1	1.30	3.63**	6.10**
lactate dehydrogenase A (LDHA)	CF614400	5	0.70	2.40**	4.31**
similar to dysferlin interacting protein 1, transcript variant 1 (LOC616223), mRNA	CO729195	2	1.57**	3.62**	4.97**
myoglobin (MB)	CF615014	18	1.79**	2.51**	6.68**
Myosin binding protein C, fast type (MYBPC2)	CF615254	14	0.89	2.73**	9.83**
myosin, heavy polypeptide 1, skeletal muscle, adult (MYH1)	CF615300	3	0.59	2.70**	5.85**
myosin, heavy polypeptide 2, skeletal muscle, adult (MYH2)	CF614936	13	1.40**	2.70**	5.37**
myozenin 1 (MYOZ1)	CF614412	18	2.22**	2.89**	7.66**
ncRNA orthologous with NEAT1	CF614858	12	0.91	1.97**	6.75**
PDZ and LIM domain 3 (PDLIM3)	CF615153	10	0.92	3.88**	5.05**
6-phosphofructo-2-kinase/fructose-2,6-biphosphatase 3 (PFKFB3)	CF615087	1	0.36	0.51	6.27**
phosphofructokinase, muscle (PFKM)	DW521807	1	1.01	4.12**	6.71**
phosphoglycerate mutase 2 (muscle) (PGAM2)	CF614969	6	1.12	4.03**	4.65**
phosphoglucomutase 1 (PGM1)	CF615240	7	1.06	3.50**	5.47**
Phosphorylase, glycogen; muscle (PYGM)	CF614749	6	1.06	3.17**	7.27**
solute carrier family 25 (mitochondrial carrier; adenine nucleotide translocator), member 4 (SLC25A4)	CF614535	3	0.80	1.93**	5.14**
titin-cap (telethonin) (TCAP)	CF614598	5	2.81**	4.15**	8.78**
tropomodulin 4 (muscle) (TMOD4)	CF614407	2	0.79	3.62**	6.53**
Titin immunoglobulin domain protein (myotitilin) (TTID)	CF614674	5	0.58	1.63**	4.61**

¹* indicates a statistically significant increase in gene expression at P < 0.05 level; ** indicates a statistically significant increase in gene expression at P < 0.01 level.

²Sequence annotation was carried out using approved gene symbols (human genome nomenclature) via the IBISS4 database [17]. BLAST [15] scores > 110 and e-values < 2 E-23 were accepted for annotation purposes.

³In cases where more than 1 microarray element was used to calculate the expression data, the accession number of 1 representative element is shown.

⁴Number of microarray elements representing the same gene that were used to calculate the expression data.

⁵To obtain the ratios shown, averaged absolute microarray signal intensity values for d 135, 195 and birth were divided by the signal intensities measured at d 60.

Two genes, FSTL1 and IGFBP5 (Table 1), which are implicated in the regulation of cell growth and differentiation, are significantly downregulated in developing bovine muscle. These genes were chosen for quantitative reverse transcription-PCR (qRT-PCR) validation (see below).

Twelve microarray probes that reported increasing expression levels during development all matched to a single consensus sequence (btcn26022) in the IBISS4 database [17]. Using the NETBLAST function in IBISS4, a match with a non-coding (nc) RNA transcript orthologous with the recently described NEAT1 transcript was identified [18] (Table 2).

Gene expression differences correlated with sire breed

A number of genes showed smaller, yet still significant differences in gene expression between breeds in at least one time point [see Additional file 2]. Table 3 provides a summary for some of these classes of genes, for example the collagen genes, which were more highly represented in the transcriptome of late-stage Wagyu fetuses. All differentially expressed ribosomal protein genes and cytochrome c oxidase genes were detected as significantly up-regulated in the muscle of Piedmontese crossbred calves at birth (Table 3 and [Additional file 2]). Ribosomal protein and cytochrome c oxidase genes were not detected as differentially expressed over developmental time (Tables 1 and 2).

qRT-PCR validation of microarray gene expression patterns

Using the MATR3 qRT-PCR assay, we were unable to confirm the >2 fold higher gene expression detected by microarray in d60 muscle RNA samples from Wagyu-sired fetuses (Table 5). In order to clarify what may have caused this discrepancy, qRT-PCR assays for MATR3 were performed on the amplified RNA samples that had been used for the microarray analysis. This analysis showed that the RNA samples amplified from the d 60 longissimus muscle of Piedmontese-sired fetuses contained higher amounts of MATR3 RNA than the non-amplified RNA sample (Figure 2). Specifically, this elevated expression was due to one individual RNA sample (P663, data not shown). qRT-PCR testing of MATR3 expression therefore confirmed that the microarray analysis, which was carried out with amplified RNA, had accurately identified the elevated MATR3 expression levels in the d 60 sample. The results did, however uncover a discrepancy between amplified and total RNA for at least one particular individual and one gene.

Temporal changes in gene expression

This microarray study of fetal bovine muscle development detected a large number of differentially expressed genes across developmental time. The expression pattern of muscle fiber structural and metabolism genes observed in this study is correlated to the timing of contractile differentiation of muscle fibers during fetal development. The overall trend of decreasing representation of extracellular matrix gene expression in fetal bovine muscle over time is in accordance with histological observations that show the proportion of extracellular matrix material on a muscle section decreasing from 70% to 27% over the time course studied here [7]. The increasing representation of muscle fiber-specific gene expression over time may therefore be a reflection of the proliferation and differentiation of muscle fibers within a connective tissue "scaffold" established early during development [20].

An earlier transcriptome study of bovine muscles during ontogenesis [10] did not identify connective tissue gene expression as a significant contributor to differential gene expression across developmental time. This may be a consequence of the skeletal muscle focused macroarray that was used in their study. The bovine fat/muscle microarray used in this study was constructed from muscle and subcutaneous fat cDNA libraries [21]. The differentially expressed microarray elements that showed decreasing gene expression over time, mainly had their origin in the subcutaneous fat cDNA library.

The observation that the most marked gene expression changes occur between d 195 and birth is in accordance with other studies that have shown the dramatic adaptations of skeletal muscle in preparation for birth [10]. As the postnatal muscle samples studied here were collected within 24 hours after birth, our study may, in addition, detect gene expression differences that reflect birth trauma.

We were unable to detect transcriptional changes that reflect the replacement of fetal with adult isoforms of contractile protein transcripts during the final trimester of pregnancy, as we are using a microarray based on adult skeletal muscle gene expression [21].

The vast majority of temporal changes in gene expression in fetal longissimus muscle comprise structural and metabolic components of extracellular matrix and muscle fibers. Two exceptions are FSTL1and IGFBP5, which have been implicated in muscle cell growth and differentiation. The expression patterns for these two genes were confirmed by qRT-PCR and their reduced expression over the developmental time course is in accord with postulated roles in the inhibition of muscle cell differentiation [22, 23]. Both genes have been reported as early targets of MyoD expression in a mouse model of MyoD-induced differentiation [24]. In this model of myogenesis, Fst1 was down-regulated, whereas Igfbp5 was strongly up-regulated. This discrepancy with our results may indicate species-specific differences in IGFBP5 action in bovine muscle development, be due to the fact that our results are derived from an in vivo, not an in vitro model of myogenesis, or hint at the proposed dual role of this gene in the regulation of muscle cell development [25]. Our results lend further support to the notion that FSTL1, while its precise function is still unknown, has a crucial role in muscle development.

The NEAT1orthologue, a ncRNA species, is found at increasing levels during development of cattle muscle. NEAT1 was one of only three large abundant ncRNAs found enriched in the nucleus and highly conserved across mammalian species [18]. This, together with its association with SC35 splicing domains, suggests a fundamental function in mRNA metabolism for the NEAT1 transcript. The increasing representation of this ncRNA species in developing cattle muscle may be a reflection of the increasing number of nuclei in developing myofibres. Myonuclei accumulate in muscle cross-sections at an increasing rate up to about 180 days post conception [26]. The 6.8 fold higher representation of the NEAT1 orthologue transcript in the muscle of newborn cattle, compared with the muscle at d60, may therefore indicate a correspondingly higher proportion of nuclei per cell volume sampled. The observations of Ansay [27] on the DNA content of developing cattle muscle lend support to this view.

Individual gene expression differences correlated with sire breed

While the microarray analysis detected a large number of genes that were significantly differentially expressed in at least one of the time points studied, only 2 of those genes showed a more than two-fold difference between the two sire breeds. We therefore examined the expression of myostatin (GDF8) itself by qRT-PCR. GDF8 RNA was expressed at significantly higher levels in the muscle of Piedmontese-sired fetuses at d 60 and d 135, in accordance with the findings of Forbes et al. [28] that the non-functional GDF8 molecule of the Piedmontese breed failed to repress the GDF8 promoter and therefore led to higher GDF8 transcript levels in the fetal muscle. Our observations therefore confirm that the GDF8 mutation did affect gene regulation in the fetal longissimus muscle as predicted.

One of the genes showing more than two-fold elevated expression in the muscle of newborn Wagyu-sired calves, FABP5, may indicate a higher level of adipogenic differentiation in the longissimus muscle of these animals, even at this early stage of their development. FABP4 showed a similar pattern of breed-specific differences, and the differential expression of both genes was confirmed by qRT-PCR (Tables 3 and 5). Fatty acid binding proteins are expressed at very high levels by adipocytes and are themselves, by interaction with peroxisome proliferator-activated receptors (PPAR) involved in mediating the effects of fatty acids on gene regulation [29]. Genetic variation in the bovine FABP4 gene has been found to be associated with intramuscular fat content and subcutaneous fat depth [30]. The possibility that FABP5 and FABP4 may be involved in defining the early onset of intramuscular fat deposition in the Wagyu breed remains to be explored. On the other hand, this finding may be an indication of relatively lower expression of FABP4 and FABP5 in Piedmontese-sired fetuses and therefore a consequence of the biochemical changes initiated by the myostatin mutation. In support of this latter possibility, a recent study of bovine fetuses carrying a myostatin mutation found a gene marker of adipocyte differentiation down-regulated in semitendinosus muscle at d 260 of gestation [11].

Global gene expression differences correlated with sire breed

The microarray used in these experiments was based on adult muscle gene expression and is therefore less likely to detect expression differences in developmental regulatory genes involved with determining the double-muscled phenotype. However, the comprehensive view of developing muscle biochemistry that these experiments describes, points to subtle differences in timing of gene expression, as well as differences in the translational capacity of the Piedmontese muscle and an altered biochemical profile as the main consequences of the GDF8 mutation. Our study provides evidence for gene expression timing differences by the large number of ribosomal genes, which tend to be expressed at a relatively lower level in the longissimus muscle of Piedmontese-sired fetuses at the earlier time points and are then relatively more highly expressed at birth. This confirms the findings by Steelman et al. [31], who found genes involved in translation overrepresented among the up-regulated genes in 5-week-old myostatin null mice. A study of d 30- d 32 bovine embryos had identified a number of ribosomal protein genes as differentially expressed in embryos carrying the GDF8 mutation [32]. Cassar-Malek et al. [11], in a gene expression study in bovine semitendinosus muscle from d 260 fetuses carrying myostatin mutations also described higher levels of expression of ribosomal protein genes at this time point. Support for the view that timing differences characterize the development of GDF8 mutants comes from the work of Deveaux et al[33], which showed evidence for a delay in contractile differentiation of muscle from Belgian Blue fetuses, and Cagnazzo et al. [9], who showed that the fetuses from a more highly muscled breed of pigs showed a later peak of myogenesis-related gene expression during fetal development.

Our study confirms the findings by Cassar-Malek et al. [11] of decreased collagen gene expression in the muscle of late stage bovine fetuses carrying a myostatin mutation. This finding lends support to our interpretation of the temporal changes of collagen gene expression (above). As the GDF8 mutation results in larger numbers of muscle fibers that have larger cross-sectional areas on average, connective tissue fibroblasts (with collagen as the main part of their expression signature) should contribute relatively less to the overall gene expression profile.

GDF8 mutants are known to have a more glycolytic adult muscle phenotype than wild type animals [33] and this is clearly reflected in the gene expression comparisons between breeds at this time point. For example, the Piedmontese-sired fetuses show higher expression of genes that reflect fast, glycolytic fiber structural differentiation such as ACTN3, MYH1 and MYH2, and MYBPC2, coupled with higher expression of genes for glycolytic metabolic enzymes such as ALDOA and LDHA. These findings confirm the observations of Cassar-Malek et al[11] on bovine muscle and also of Steelman et al. [31], who documented a shift away from slow muscle fiber isoforms in myostatin null mice.

At birth, a suite of gene expression changes suggesting an increased ATP supply to the muscle of the Piedmontese-sired fetuses is observed. In fact, there is a coordinated up-regulation of genes encoding enzymes that determine mitochondrial activity at multiple steps in cellular respiration. These include key oxidative enzymes, in particular isoforms of cytochrome c oxidase (COX5B, 6C, 7B, 7C), the terminal enzyme in the citric acid cycle. Cagnazzo et al[9] made a similar observation in their comparison of transcription profiles of developing muscles from two contrasting breeds of pigs. They found that a highly muscled breed of pigs (Pietrain) showed more active transcription of ATP metabolism-related genes during late fetal development.

This surge in ATP supply may reflect increased energetic demand in the immediate postnatal period in animals that experience very rapid postnatal growth, such as Piedmontese-sired calves or Pietrain piglets. Two major ATP consuming cellular processes are the maintenance of transmembrane ion gradients and protein turnover [34]. The demand for sarcolemmal ion gradients is unlikely to increase at this time. In fact, the expression of ATPase genes involved in this process (ATP1A2, ATP2A1, ATP2A2) is relatively lower in the muscle of Piedmontese-sired calves.

Therefore, taken together with the observation that a suite of ribosomal proteins is upregulated in Piedmontese-sired calves, we suggest that the extra ATP might be supporting more active ribosomal synthetic machinery. Consequently, the gene expression changes captured by this microarray appear to describe the energetic modifications necessary for the rapid muscular growth of newborn GDF8 mutants.

Animals

All animal experimentation complied with the Animal Ethics requirements of NSW Agriculture/GARAS. Hereford cows were artificially inseminated or mated to one of 5 different Wagyu sires or one of 6 different Piedmontese sires. All Piedmontese sires were homozygous for the GDF8 missense mutation in exon 3 and none of the Wagyu sires carried the mutation. 3 Fetuses from each sire breed were recovered by caesarean section at around 60, 130 and 195 of gestation. 3 newborn calves per sire breed were euthanased by lethal injection within 24 hours of birth. Details of the sampled fetuses and calves were recorded [see Additional file 4]. Longissimus muscle was dissected immediately after death of the fetuses or calves and tissue samples snap frozen in liquid nitrogen.

Microarray methods

The process for RNA extraction and purification was performed as described in Lehnert et al. [21]. Anti-sense RNA (aRNA) amplification was performed using the MessageAmp aRNA Kit (Ambion, Austin, TX). Indirect labeling procedures and hybridizations were performed as described in Lehnert et al. [21] with the following modifications. The cDNA was generated using a modified amino C6 dT random primer [35] (Geneworks, Adelaide, SA) and the reverse transcriptase Superscript III (Invitrogen, Carlsbad, CA). The cDNA (pre-labeling) and the post-labeled cDNA were purified using the QIAquick PCR purification columns (Qiagen, Valencia, CA). The hybridization mixture used a reduced detergent concentration of 0.2%SDS.

We utilized a bovine cDNA microarray constructed from two cattle cDNA libraries; longissimus muscle and subcutaneous fat tissue derived from a 24 mo-old grass-fed Angus steer [21]. The array contained 9,600 bovine cDNA probes, printed in duplicate.

Experimental design and data acquisition criteria

The weight of the individual fetuses [see Additional file 4], was analyzed by fitting an ANOVA model that contained the effects of breed, time, sex, and breed by time interaction. The ANOVA model accounted for 97.8% of the total variation in fetus weight and significant effects included time (P < 0.001), breed by time (P < 0.05) and breed (P < 0.10). The design of the microarray experiment therefore focused on the developmental aspect of the study, as well as exploring sire breed.

Factors considered in designing the experiment included the availability of RNA for each sample of interest as well as the cost of the arrays themselves. Figure 3 illustrates the configuration layout for the experimental design that was developed with a total of three biological replicates, two breeds and four time points that were compared using 25 microarray hybridizations. The design was chosen to allow a focus on the developmental aspects of the study, but to also permit a breed comparison to be carried out. The entire set of expression data was deposited on the GEO database [14] and can be accessed using accession numbers GPL4196 and GPL4197.

We used the GenePix 4000A optical scanner (Molecular Devices, Sunnyvale, CA) and the GenePixPro 5.1 image analysis software (Molecular Devices) to quantify the gene expression level intensities, employing editing criteria detailed in [36]. Editing criteria were applied separately for the red (Cy5) and for the green (Cy3) intensity channels so that a different number of observations for each channel were obtained. These resulted in a total of 887,706 gene expression intensity readings (443,292 red and 444,414 green) that were background corrected and base-2 log transformed. The arithmetic mean and standard deviation (in brackets) for the red and green intensities were 11.93 (2.02) and 11.97 (1.95), respectively.

Mixed-model equations and identification of differentially expressed genes

The following linear mixed-effect model was fitted to the data:

Y _ijktmn = μ + C _ijk + G _m + AG _ijm + DG _km + TG _tm + ε _ijktmn

where Y _ijktmn represents the n-th background-adjusted, normalized base-2 log-intensity from the m-th gene (probe) at the t-th treatment (fetal age and breed sample) from the i-th array, j-th printing block and k-th dye channel; C represents a comparison group fixed effect defined as those intensity measurements from the same array slide, printing block and dye channel; G represent the random gene (probe) effects with 8,845 levels; AG, DG, and TG are the random interaction effects of array × gene, dye × gene, and treatment × gene, respectively. Finally, ε is the random error term.

Variance components for random effects were estimated using restricted maximum likelihood. Differentially expressed (DE) genes were identified after processing the appropriate linear combination of the solutions (best linear unbiased predictions) of TG via model-based clustering, addressing the developmental as well as the breed comparison contrasts of interest [37, 38].

A total of six contrasts were considered in the identification of DE genes. These included the two within-breed over-time contrasts and the four within-time across-breed contrasts. For each contrast, a two-component normal mixture model was fitted and posterior probabilities of belonging to the non-null component used to identify DE genes for an estimated experiment-wise false discovery rate of < 1% were computed as described by McLachlan et al. [39].

Synthesis of complementary DNA (cDNA) for quantitative reverse transcription PCR (qRT-PCR) validation

TRIzol-extracted (Invitrogen) total RNA from longissimus muscle samples was further purified on RNeasy mini columns and treated on-column with DNase I (Qiagen) according to the manufacturer's instructions. Total RNA concentrations and relative purity were determined spectrophotometrically by measurement of UV absorbance at 260 and 280 nm (GeneQuant; GE Healthcare, Little Chalfont, UK). Reverse transcription of 1.5 μg of total RNA was performed with the Omniscript cDNA synthesis kit (Qiagen) using a mixture of oligonucleotides [100 μM OligodTVN (5'-TTTTTTTTTTTTTTVN-3' where V = A,C,G and N = A,C,G,T) and 1 μM 18Sprimer (5'-CACACGCTGAGCCAGTCAGT-3')]. The cDNA was diluted 1:10 in 10 mM Tris-HCl (pH8.0) for target gene measurements, whereas the 18S rRNA reference gene measurements were performed on a 1:200 dilution of cDNA.

Total RNA from d 60 fetal muscle samples and amplified RNA from the d 60 samples that were used to probe the DNA microarray were separately reverse transcribed. 1.5 μg of total or amplified RNA was reverse transcribed with the Superscript III first strand cDNA synthesis kit (Invitrogen). The cDNA synthesis used random hexamers and was performed at 37°C for 1 hr in the presence of RNase inhibitor (RNaseOUT, Invitrogen). The cDNA was diluted 1:10 in 10 mM Tris-HCl (pH8.0) and was used to measure the mRNA levels of MATR3 and RPLP0.

qRT-PCR analyses

Statistical assessment of qRT-PCR gene expression data

Statistical analyses were performed using the Procedure GLM of SAS (SAS v9.1, SAS Institute Inc., Cary NC, USA) with an ANOVA model that included the expression of the genes (both target and reference genes), as measured by their threshold cycle (C_t) of the PCR reaction, as dependent (response) variable and biological replicate (with two levels), breed (with two levels), age (with four levels), breed by age interaction (with eight levels) and residual as independent variables. Biological replicates were fitted as nested within breed by age interaction. Least square means of Ct for each gene and at each level of breed, age, and breed by age interaction were compared for differential expression. To ensure overall protection level, only probabilities associated with pre-planned comparisons were used.