Skip to main content

Advertisement

Higher plane of nutrition pre-weaning enhances Holstein calf mammary gland development through alterations in the parenchyma and fat pad transcriptome

Article metrics

Abstract

Background

To reduce costs of rearing replacement heifers, researchers have focused on decreasing age at breeding and first calving. To increase returns upon initiation of lactation the focus has been on increasing mammary development prior to onset of first lactation. Enhanced plane of nutrition pre-weaning may benefit the entire replacement heifer operation by promoting mammary gland development and greater future production.

Methods

Twelve Holstein heifer calves (< 1 week old) were reared on 1 of 2 dietary treatments (n = 6/group) for 8 weeks: a control group fed a restricted milk replacer at 0.45 kg/d (R, 20% crude protein, 20% fat), or an accelerated group fed an enhanced milk replacer at 1.13 kg/d (EH, 28% crude protein, 25% fat). At weaning (8 weeks), calves were euthanized and sub-samples of mammary parenchyma (PAR) and mammary fat pad (MFP) were harvested upon removal from the body. Total RNA from both tissues was extracted and sequenced using the Illumina HiSeq2500 platform. The Dynamic Impact Approach (DIA) and Ingenuity Pathway Analysis (IPA) were used for pathway analysis and functions, gene networks, and cross-talk analyses of the two tissues.

Results

When comparing EH vs R 1561 genes (895 upregulated, 666 downregulated) and 970 genes (506 upregulated, 464 downregulated) were differentially expressed in PAR and MFP, respectively. DIA and IPA results highlight a greater proliferation and differentiation activity in both PAR and MFP, supported by an increased metabolic activity. When calves were fed EH, the PAR displayed transcriptional signs of greater overall organ development, with higher ductal growth and branching, together with a supportive blood vessel and nerve network. These activities were mediated by intracellular cascades, such as AKT, SHH, MAPK, and Wnt, probably activated by hormones, growth factors, and endogenous molecules. The analysis also revealed strong communication between MFP and PAR.

Conclusion

The transcriptomics and bioinformatics approach highlighted key mechanisms that mediate the mammary gland response to a higher plane of nutrition in the pre-weaning period.

Background

Replacement heifer rearing accounts for approximately 20% of annual farm costs [1,2,3]. To reduce costs, research has focused on: 1) lowering age at breeding and first calving, and 2) enhancing mammary development before first lactation to increase returns during the productive life of the cow. The first strategy appears the easiest and fastest to implement, being that puberty is tightly associated with body weight (BW), and Holstein heifer pubertal BW is relatively constant (250–280 kg or 40–50% of mature weight [4]). In fact, Tozer and Heinrichs [5] estimated a reduction of 4.3% of heifer rearing costs by decreasing age at first calving by 1 month, mainly due to a shortening of the non-productive portion of a heifer’s life by encouraging earlier herd entry and productivity.

To reach puberty weight earlier, producers must increase feed allowance in the pre-pubertal period. However, when this was done immediately after weaning, mammary parenchymal (PAR) tissue mass and DNA content in Holstein calves decreased by 23 and 32%, respectively [6]. The detrimental effects of increased pre-pubertal nutrient intake (starting post-weaning) on mammary development have long been recognized [7,8,9,10], and recently [6] confirmed in various studies [11,12,13,14,15,16]. Furthermore, the negative effect on mammary development may translate into poorer first lactation performance, as a recent study confirmed that heifers fed at a greater rate produced 14% less milk compared with controls [17].

It is well-accepted that events occurring pre-weaning can have lasting effects throughout a dairy cow’s life [18]. Recent data have indicated that the negative correlation between gain and mammary development may be the opposite in this period to what is observed in the pre-pubertal, post-weaning period [19,20,21,22]. Thus, an enhanced pre-weaning plane of nutrition (targeting greater gains) increases PAR and mammary fat pad (MFP) weight, DNA content of the mammary PAR, and total mammary PAR DNA [23]. Furthermore, increasing pre-weaning average daily gain by 1 kg/d was associated with an increase of 1000 kg or more in milk yield [24].

The mechanisms by which an increase in plane of nutrition pre-weaning enhances mammary development and consequent first lactation performance remain unclear. Therefore, the objective of this study was to obtain a general and holistic view through high-throughput mRNA sequencing of the mammary MFP and PAR transcriptome in Holstein heifer calves fed two distinct planes of nutrition in the pre-weaning period. By identifying genes altered by the dietary treatment, we could uncover molecular pathways and determine possible mediators of the mammary gland response to pre-weaning plane of nutrition. We hypothesized that an enhanced pre-weaning plane of nutrition positively alters mammary gland transcription by upregulating pathways involved in tissue growth and development. Furthermore, we hypothesized that a higher plane of nutrition could promote the inter-tissue communication between PAR and MFP. The effects of an improved plane of nutrition on general body growth and the development of multiple organs, mammary gland included, are reported elsewhere [23, 25].

Methods

Animal handling and experiment design

This experiment was conducted under the review and approval of the Virginia Polytechnic Institute and State University (VT) Institutional Animal Care and Use Committee (#14–045-DASC). The experimental design and animal handling were previously described [25]. Briefly, twelve Holstein heifer calves (6.0 ± 2 d old, 39.0 ± 4.4 kg at the time of arrival, and ≥ 5.5 mg/dL of total serum protein) were purchased from a single commercial producer (located ~ 90 miles from campus), brought to the VT Dairy Farm between May and June of 2014), and randomly assigned to two treatments. The control group was fed 0.45 kg/calf of milk replacer (MR) per day containing 20% crude protein (CP) and 20% fat (a restricted MR, R). The accelerated group was fed 1.13 kg of MR per day at 28% CP and 25% fat (an enhanced MR, EH). Starter feed was introduced at the end of week 4 and kept similar between treatments. Milk replacer was reduced in both treatments to 50% at week 8, to induce weaning. All calves were housed individually with ad libitum access to water. At weaning, calves were euthanized and their whole mammary glands were removed, weighed and dissected. A summary of phenotypic responses obtained and already published [25] is presented in Table 1.

Table 1 Summary of Holstein heifer calf (n = 6/treatment) performance when fed during the pre-weaning period (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/calf/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/head/day, 20.9% crude protein, 19.8% fat, DM basis)

Sample collection and slaughter procedures

A detailed description of sample collection can be found elsewhere [23, 25]. Calves at weaning (8 weeks) were euthanized at VT’s Veterinary Facility (approximately 1 mile from their housing) using a commercial phenobarbital solution administered intravenously (Fatal-Plus, 10 mg/kg of BW, Vortech Pharmaceuticals, Dearborn, MI), and subsequently exsanguinated. Pieces of PAR and MFP (~ 13.0 mg) were sampled from the mammary gland upon removal from the body, frozen by immersion in liquid nitrogen, and stored at − 80 °C.

RNA extraction, library construction, and sequencing

Tissue was weighted (PAR, ~ 0.10 g; MFP, ~ 0.20 g), immediately placed in QIAzol Lysis Reagent (cat#79306, Qiagen) (1.20 mL) and homogenized using a Mini-Beadbeater-24 (cat#112011, Biospec Products Inc.) with two 30 s cycles, and 1 min incubation on ice in between the cycles. Samples were then centrifuged for 10 min at 12,000×g and 4 °C, and the supernatant was transferred to a separate tube and mixed with Chloroform (cat#C298, Fisher Chemical) (0.24 mL). After centrifugation for 15 min at 12,000×g and 4 °C, the aqueous phase was transferred to a new tube, mixed with 100% Ethanol (cat#2701, Decon Labs; 0.90 mL), and total RNA was cleaned using miRNeasy mini kit columns (cat# 217004, Qiagen) following manufacturer’s protocols. During purification, genomic DNA was removed using the RNase-Free DNase Set (cat#79254, Qiagen). Quantity and purity were determined using a NanoDrop ND-1000 (NanoDrop Technologies Inc.), while integrity was assessed via a Fragment Analyzer™ (Advance Analytical). All samples had an RQN score greater than 8.0. RNA samples were stored at − 80 °C until analysis.

RNA-Seq cDNA libraries were constructed using total RNA isolated from both MFP and PAR sampled at weaning (week 8) slaughter. The Illumina TruSeq Stranded mRNA Sample Prep kit was used for single-end read library construction following the manufacturer’s instructions with mRNA enrichment. Complete details of this procedure are available in Additional file 1. Libraries were pooled together and multiplexed across 4 flow cell lanes in the Illumina HiSeq4000 (Illumina Inc., San Diego, CA) platform to obtain an average of 20–30 million reads per sample.

Transcriptome sequencing data processing and statistical analysis

Single-end reads were first filtered using Trimmomatic 0.33 [26] with a minimum quality score of 28 (i.e., base call accuracy of 99.84%) leading and trailing with a minimum length of 30 bp long and subsequently checked using FastQC 0.11.4 (Babraham Institute, Cambridge, UK). Reads were then mapped to the Bos taurus UMD 3.1.1 reference genome (1/29/16 NCBI release) using default settings of STAR 2.5.1b [27] with the quantMode option for gene counts. Further data analysis was conducted using R. 3.2.4 (R Core Team, 2016). Reads uniquely assigned to a gene were used for subsequent analysis. After accounting for high expression genes and library size differences using trimmed mean of M-values normalization in edgeR [28], genes were filtered if 4 samples did not have > 1 count per million mapped reads. Normalization of reads was conducted using the voom variance stabilization function in limma [29]. Differential expression analysis was conducted in limma using a single factor model which included the main effect of diet (2 levels). Raw P-values were adjusted using the false discovery rate (FDR) method [30]. Differences in transcript profiles (differentially expressed genes, DEG) were considered significant at an FDR-adjusted P < 0.05. The focus of this manuscript is on the overall differences in mammary gland transcriptome in response to the main effect of diet; that is, EH vs R.

Bioinformatics analyses

Kyoto encyclopedia of genes and genomes (KEGG) pathways analysis

The dynamic impact approach (DIA) was used for KEGG pathway analysis of DEG. The detailed methodology of DIA is described elsewhere [31]. Briefly, the whole DEG data set with Entrez gene ID, FDR (< 0.05), fold-change (FC), and P-value (< 0.05) was uploaded to DIA. For the analyses, pathways with a minimum of 30% annotated genes in the transcriptome dataset versus the whole genome, and at least 4 genes were selected. Furthermore, pathways related to KEGG category “Human diseases” and Organismal system subcategories “Digestive system”, “Excretory system”, and “Sensory system” were not considered as pertinent to the analyzed tissues.

Function analysis and transcription regulator discovery

Ingenuity pathway analysis (IPA) software (https://www.qiagenbioinformatics.com/products/ingenuity-pathway-analysis/) was used to identify major affected functions and analyze the upstream regulators and their connections with other downstream genes that were differentially expressed. To identify only highly significant functions and upstream regulators, a list of DEG with FC > |1.5| and P-values (< 0.05) was uploaded to IPA, and a core analysis was run using default settings. Functions with p-value < 0.05 and z-score > |2.0| were considered significant, and upstream regulators were considered significant with overlap-P value < 0.05 and z-score > |2.0|.

Tissue interaction and crosstalk

The cross-talk between MFP and PAR was performed using the network capability of IPA as previously described [32]. To highlight how the dietary treatment affects such interactions, DEG (in the EHvsR comparison) considered to code for secreted proteins encompassed the cytokine and growth factor categories, while DEG coding for proteins considered to be receptors that might be able to “sense” the secreted proteins were those in G-protein coupled receptor, ligand-dependent nuclear receptor, transcription regulator, and transmembrane receptor categories. Networks between tissues were built using the IPA Knowledge base. The knowledge base was restricted to tissue and cell line categories “Mammary gland”, “Breast cancer cell lines”, and “Other cell line” when analyzing the effect of MFP secreted molecules on PAR, or to the categories “Adipose”, “Adipocytes”, and “Other cell line” when analyzing the effect of PAR secreted molecules on MFP. Furthermore, when analyzing the involvement of infiltrated immune cells in the tissue crosstalk, the categories were limited to “Immune cells”, “Immune cell lines”, and “Other cell lines”, regardless of the direction of the crosstalk.

Real-time qPCR and statistical analysis

The results from RNA-seq were validated via qPCR for a selected panel of genes representing the top most upregulated (n = 4) or downregulated (n = 4) genes within each tissue. Complete details of primer design, qPCR analysis, and performance are available in Additional file 1. Briefly, the qPCR performed was SYBR Green-based, and results were calculated via the 2-ΔΔCt method. MTG1, PPP1R11, and RPS15A were used as internal controls. Results were subjected to ANOVA and analyzed using repeated measures ANOVA with PROC MIXED within SAS (v9.4). To validate sequencing results, data from each tissue were independently analyzed, and the statistical model included diet as fixed effect. The Kenward-Roger (KR) statement was used for computing the denominator degrees of freedom. Data were considered significant at a P ≤ 0.05 using the PDIFF statement in SAS. The comparison with the RNA sequencing results is reported in Table 2.

Table 2 mRNA expression of target genes selected for validation of sequencing results. Genes were among the top up and downregulated genes in mammary parenchyma or fat pad of Holstein heifer calves (n = 6/treatment) fed pre-weaning (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/head/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/head/day, 20.9% crude protein, 19.8% fat, DM basis)

Results

RNA sequencing and gene expression analyses

A summary of sequencing read alignments and mapping is reported in Additional file 2. Overall, on average, samples had approximately 27 million reads of which 25 million (~ 92%) were uniquely mapped. The statistical analysis identified 15,214 and 14,223 uniquely annotated (EntrezID) genes in PAR and MFP, respectively. Of these, considering an FDR < 0.05, 1561 genes (upregulated 895, downregulated 666) and 970 genes (upregulated 506, downregulated 464) were detected as differentially expressed in PAR and MFP, respectively, when comparing EH to R calves (Fig. 1).

Fig. 1
figure1

Differentially expressed genes (FDR < 0.05) detected in mammary parenchyma (PAR) and fat pad (MFP) of Holstein heifer calves (n = 6/treatment) fed pre-weaning (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/calf/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/head/day, 20.9% crude protein, 19.8% fat, DM basis)

KEGG pathway analysis summary

The DIA analysis yields the impact and flux of all the manually-curated pathways included in the KEGG database. The term “impact” refers to the biological importance of a given pathway as a function of the change in expression of genes composing the pathway (proportion of DEG and their magnitude) in response to a treatment, condition, or change in physiological state [31]. Consequently, the direction of the impact, or flux, characterizes the average change in expression as up-regulation/activation, down-regulation/inhibition, or no change.

Feeding an accelerated (EH) compared to a restricted (R) milk replacer pre-weaning had broad effects on the transcriptome (Fig. 2). All categories, both metabolic (‘Metabolism’) and non-metabolic (‘Genetic information processing’, ‘Environmental information processing’, ‘Cellular processes’, and ‘Organismal systems’), were impacted and up-regulated (positive flux) in PAR of EH vs R calves, with greatest changes in ‘Metabolism’, ‘Environmental information processing’, and ‘Organismal systems’. Similar changes were observed in the MFP response to pre-weaning plane of nutrition; however, ‘Genetic information processing’ was the most impacted category, and contrary to PAR, it was down-regulated (negative flux) in EH vs R calves. Furthermore, despite the ‘Metabolism’ category having a general activation (up-regulation) trend in both tissues, tissue-specific differences could be detected in the fluxes of its subcategories: ‘Energy metabolism’ (downregulated PAR, upregulated MFP), ‘Metabolism of other amino acids’ (upregulated PAR, downregulated MFP) and ‘Metabolism of terpenoids and polyketides’ (upregulated PAR, downregulated MFP). Concerning ‘Metabolism of cofactors and vitamins’, despite different fluxes between the tissues (downregulated PAR, upregulated MFP), these were very close to zero, thus, suggesting the lack of a clear up- or down-regulation in both tissues.

Fig. 2
figure2

Summary of KEGG metabolic subcategories resulting from the DIA analysis in mammary parenchyma (PAR) or fat pad (MFP) of Holstein heifer calves (n = 6/treatment) fed pre-weaning (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/head/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/calf/day, 20.9% crude protein, 19.8% fat, DM basis). For each tissue, the columns represent the effect (impact) and flux responses. The blue bars represent the effect value (0 to 30), and the flux columns represent negative (−) and positive (+) flux (−30 to + 30) based on the direction of the effect. The positive flux (green bars) indicates an upregulation in EH-fed heifer calves, while the negative flux (red bars) indicates an upregulation in R-fed heifer calves

KEGG most impacted pathways and IPA functions

Mammary parenchyma

The top PAR metabolic and non-metabolic pathways are reported in Fig. 3, while the complete dataset of results is available in Additional file 3. The top-20 impacted metabolic pathways are included in the sub-categories ‘Amino acid metabolism’, ‘Energy metabolism’, ‘Glycan biosynthesis and metabolism’, ‘Lipid metabolism’, ‘Metabolism of Cofactors and Vitamins’, and ‘Xenobiotics biodegradation and metabolism’. Among these, all but 5 pathways were up-regulated in EH vs R: ‘Metabolism of xenobiotics by cytochrome P450’, ‘Nicotinate and nicotinamide metabolism’, ‘Nitrogen metabolism’, ‘Retinol metabolism’, and ‘Vitamin B6 metabolism’. Regarding the top-20 impacted non-metabolic pathways, they include the sub-category ‘Circulatory system’, ‘Development’, ‘Endocrine system’, ‘Environmental adaptation’, ‘Immune system’, ‘Nervous system’, ‘Signaling molecules and interaction’, ‘Signal transduction’. All but 3 pathways were up-regulated in EH vs R: ‘Circadian rhythm – mammal’ was down-regulated, while ‘Calcium signaling pathway’ and ‘Neuroactive ligand-receptor interaction’ had null flux, meaning they are biologically important during the pre-weaning physiological period despite plane of nutrition affecting the direction of flux.

Fig. 3
figure3

Dynamic Impact Approach (DIA) results (Impact and Direction of the Impact) for the most impacted metabolic and non-metabolic KEGG pathways (Top 20), grouped by sub-categories of pathways, in mammary parenchyma of Holstein heifer calves (n = 6/treatment) fed pre-weaning (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/head/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/calf/day, 20.9% crude protein, 19.8% fat, DM basis). Blue bars represent the Impact (0–50) of dietary treatment, while red or green bars show the Direction of the Impact (−50 − + 50) (red = upregulation, green = downregulation, in EH vs R)

The IPA analysis identified 16 functions as significantly affected by the pre-weaning plane of nutrition, all predicted to be increased in EH rather than R calves. They encompass activities related to cell proliferation, differentiation, and growth, as well as the involvement of the immune system in PAR pre-pubertal development (Fig. 4).

Fig. 4
figure4

Ingenuity pathway analysis predicted canonical functions in mammary parenchyma of Holstein heifer calves (n = 6/treatment) fed pre-weaning (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/calf/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/head/day, 20.9% crude protein, 19.8% fat, DM basis).

Mammary fat pad

The top MFP metabolic and non-metabolic pathways are reported in Fig. 5, while the complete dataset of results is available in Additional file 3. The top-20 metabolic pathways are represented in the categories ‘Amino acid metabolism’, ‘Carbohydrate metabolism’, ‘Glycan biosynthesis and metabolism’, ‘Lipid metabolism’, ‘Metabolism of cofactors and vitamins’, ‘Metabolism of terpenoids and polyketides’, ‘Nucleotide metabolism’, and ‘Xenobiotics biodegradation and metabolism’. All but 4 pathways were up-regulated in EH vs R: ‘Folate biosynthesis’, Metabolism of xenobiotics by cytochrome p450’, and ‘Terpenoid backbone biosynthesis’ were down-regulated, while ‘Nicotinate and nicotinamide metabolism’ had a null flux. Among the top-20 non metabolic pathways, the sub-categories were ‘Cell communication’, ‘Endocrine system’, ‘Immune system’, ‘Membrane transport’, ‘Replication and repair’, ‘Signal transduction’, ‘Signal molecules and interaction’, ‘Translation’, and ‘Transport and catabolism’. All but 7 of these pathways were up-regulated in EH vs R: ‘Aminoacyl-tRNA biosynthesis’, ‘DNA replication’, ‘Melanogenesis’, ‘Mismatch repair’, ‘Ribosome’, and ‘RNA transport’ were down-regulated, while ‘Cell adhesion molecules (CAMs)’ had a null flux.

Fig. 5
figure5

Dynamic Impact Approach (DIA) results (Impact and Direction of the Impact) for the most impacted metabolic and non-metabolic KEGG pathways (Top 20), grouped by sub-categories of pathways, in mammary fat pad of Holstein heifer calves (n = 6/treatment) fed pre-weaning (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/head/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/calf/day, 20.9% crude protein, 19.8% fat, DM basis). Lines represent the Impact (0–30) of dietary treatment, while red or green bars show the Direction of the Impact (− 30 − + 30) (red = upregulation, green = downregulation, in EH vs R)

The IPA analysis identified 7 significant functions, 6 predicted to be increased and 1 predicted to be decreased. All the functions were related to immune cell trafficking, invasion, and accumulation (Fig. 6).

Fig. 6
figure6

Ingenuity pathway analysis predicted canonical functions in mammary fat pad of Holstein heifer calves (n = 6/treatment) fed pre-weaning (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/calf/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/head/day, 20.9% crude protein, 19.8% fat, DM basis)

Up-stream regulators

The upstream regulator analysis generated 1122 and 685 possible regulators in PAR and MFP, respectively. Of these, only 69 (PAR; 58 predicted activated and 11 predicted inhibited) and 32 (MFP; 28 predicted activated and 4 predicted inhibited) were considered significantly involved in the response to plane of nutrition in EH vs R calves. In MFP, 10 (8 predicted activated, 2 predicted inhibited) of the significant upstream regulator were considered biased by IPA due to the skewness of the down-regulated genes toward up- or down-regulation. The full list of upstream regulators is reported in Tables 3 and 4 for PAR and MFP, respectively. Overall, the predicted regulators in both tissues were part of the molecule categories ‘chemical – endogenous mammalian’, ‘complex’, ‘cytokine’, ‘group’, ‘growth factor’, ‘ligand-dependent nuclear receptor’, ‘transcription regulator’, ‘transmembrane receptor’, and ‘others’. Unique to PAR were the categories ‘kinase’, ‘mature microrna’, and ‘microrna’, while ‘peptidase’, ‘phosphatase’, and ‘transporter’ were unique to MFP.

Table 3 IPA-predicted up-stream regulators in mammary parenchyma (PAR) of Holstein heifer calves (n = 6/treatment) fed pre-weaning (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/calf/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/head/day, 20.9% crude protein, 19.8% fat, DM basis)
Table 4 IPA predicted up-stream regulator in mammary fat pad (MFP) of Holstein heifer calves (n = 6/treatment) performances fed pre-weaning (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/calf/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/head/day, 20.9% crude protein, 19.8% fat, DM basis)

Tissue cross-talk

Our analysis identified 14 and 8 up-regulated differentially expressed genes in EH vs R coding for secreted proteins, and 105 and 64 genes coding for possible receptors in PAR and MFP, respectively (Additional file 4). However, after automatic network construction with the IPA knowledge base, only some of these were deemed relevant to the cross-talk between these specific tissues. PAR differentially expressed (EH vs R) genes encoding 3 potentially secreted proteins (PGF, IL1B, IL1RN) which can influence MFP through 12 possible receptors differentially expressed by the latter (Fig. 7). In contrast, MFP differentially expressed genes encoding 2 potentially secreted proteins (CSF1, TGFB1) that could have influenced PAR through 26 possible differentially expressed receptors (Fig. 8). When assessing the role of immune cells in the cross-talk between PAR and MFP, 4 PAR (CXCL5, IL1B, CXCL10, and CXCL9) (Fig. 9) and 3 MFP (CSF1, NGF, and TGFB1) (Fig. 10) differentially expressed protein coding genes were involved in the tissues cross-talk.

Fig. 7
figure7

Role of mammary parenchyma (blue) in the tissue cross-talk with mammary fat pad (yellow) in Holstein heifer calves (n = 6/treatment) fed pre-weaning (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/head/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/head/day, 20.9% crude protein, 19.8% fat, DM basis). The affected functions in the fat pad are highlighted at the bottom

Fig. 8
figure8

Involvement of the immune system in the tissue cross-talk between mammary parenchyma (blue) and mammary fat pad (yellow) in Holstein heifer calves (n = 6/treatment) fed pre-weaning (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/head/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/head/day, 20.9% crude protein, 19.8% fat, DM basis). The affected functions in the fat pad are highlighted at the bottom

Fig. 9
figure9

Role of mammary fat pad (yellow) in the tissue cross-talk with mammary parenchyma (blue) in Holstein heifer calves (n = 6/treatment) fed pre-weaning (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/head/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/calf/day, 20.9% crude protein, 19.8% fat, DM basis). The affected functions in the parenchyma are highlighted at the bottom

Fig. 10
figure10

Involvement of the immune system in the tissue cross-talk between mammary fat pad (yellow) and mammary parenchyma (blue) in Holstein heifer calves (n = 6/treatment) fed pre-weaning (week 1 to 8 of life) an enhanced (EH) (1.08 kg of powder/calf/day, 28.9% crude protein, 26.2% fat, DM basis) or restricted (R; represents the industry standard or control) milk replacer (0.44 kg of powder/calf/day, 20.9% crude protein, 19.8% fat, DM basis). The affected functions in the parenchyma are highlighted at the bottom

Discussion

Over the past decades, calf nutrition research has focused on developing strategies that facilitate early weaning, delivering a smooth transition from liquid to solid feed. The main focus has been on increasing calf starter intake by reducing milk supply (~ 10% BW), hence facilitating weaning at 7–8 wk., and possibly improving early rumen development and function [33]. However, recent work revealed poor weight gains, higher risk of disease, abnormal behavior, and a reduction in calf welfare when restricted-feeding strategies are used [33]. In contrast, providing greater quantities of milk or MR improves growth and feed efficiency [33]. Furthermore, conventional MR is composed of 20–22% of crude protein, but calves respond better to increased milk allowance when MR contains higher protein and lower fat (up to 30%, with 15–20% fat) [33]. This feeding strategy has been recently described as “accelerated early nutrition”, “accelerated growth”, “enhanced nutrition”, “intensified nutrition”, or “biologically appropriate growth” [34]. The current manuscript compares the effect of a restricted-style management strategy, with an enhanced one on calf parenchyma and mammary fat pad transcriptome.

The mammary parenchyma

The activity of mammary parenchyma

We previously showed that, compared with R, EH calves had greater mammary PAR mass (7.3 fold greater) without changes in PAR composition (e.g. protein or fat concentration) [23]. Furthermore, incorporation of bromo-2′-deoxyuridine indicated greater (~ 60% higher) mammary epithelial cell proliferation in EH than R calves, particularly at the terminal end of the developing ducts (terminally ductal units) [35]. There results were both mirrored at a transcriptome level. The IPA-predicted functions highlighted greater proliferation and differentiation activity of PAR cells (e.g. quantity, proliferation, and differentiation of cells; Fig. 4). The DIA analysis further confirmed the IPA prediction, as greater upregulation of pathways responsible for cell growth, DNA replication and repair, and expression and translation of genetic information was observed. Furthermore, up-regulation of pathways such as ‘cell adhesion molecules’, ‘ECM-receptor interaction’, ‘focal adhesion’, ‘gap junction’, and ‘adherens junction’ illustrate the formation of a cohesive tissue with cell communication and interaction mechanisms in place. Cell-to-cell signaling by contact is a well-known stimulus of cellular proliferation via a PI3K-dependent intracellular signal [36], a cascade that both DIA and IPA analysis predict to be greatly activated in EH calves.

The dissection of the harvested samples focused on the isolation of the mammary PAR. However, the complex nature of the mammary gland tissue will yield a sample mostly composed of PAR cells, together with a heterogeneous ensemble of other tissue types (e.g., nervous, connective). This heterogeneity was detected by the bioinformatics analysis, supporting the notion that the mammary gland is growing as an organ, rather than just due to the proliferation of PAR cells. In fact, both IPA (Fig. 4) and DIA (e.g. ‘axon guidance’ pathway, and all nervous system subcategory pathways) analysis highlighted the growth and development of the nervous system. Furthermore, the up-regulation of pathways such as ‘renin-angiotensin system’ and ‘vascular smooth muscle contraction’ indicate an increased development of the circulatory system, all supported by the activation of the ‘VEGF (vascular endothelial growth factor) signaling pathway’, master regulator of vascular development and of blood and lymphatic vessel function [37].

Nutrient supply and metabolism

A constant and abundant supply of nutrients to allow for the rapid growth of mammary PAR was ensured by the up-regulation in EH calves of the ‘ABC transporters’ (ATP binding cassettes transporters), a large family of macro- and micro-nutrients active transporters [38]. At the same time, the up-regulation of ‘Endocytosis’ and ‘Lysozome’ pathway allow for the use of extracellular matrix protein (an important alternative source of nutrients) when accessible to cells [39]. Furthermore, the up-regulation of the ‘mTOR signaling pathway’ and almost all amino acid metabolism pathways would have helped sustain the protein synthesis machinery of the proliferating tissue, which took advantage of the increase intake of protein provided by the EH MR. Further discussion on amino acid metabolism in the PAR can be found in Additional file 5.

Besides amino acids and proteins as building blocks, proliferating cells require a substantial amount of energy to support neo-synthesis of proteins, and their metabolic activity as a whole [40]. The activation of carbohydrate and lipid metabolism pathways suggests a greater flux of carbon through the citric acid cycle to generate ATP. The increased amount of fatty acid intake via MR likely supplied substrate for beta-oxidation to generate energy, including phosphatidylcholine that was degraded though the alpha-linoleic acid metabolism pathway. Overall, despite the potential increase in the flux through the citric acid cycle, thanks to the supply of acetyl-CoA, the activation of ‘synthesis and degradation of ketone bodies’ suggests a saturation of the cycle’s oxidative capacity. The mammary gland of ruminants is, in fact, capable of ketogenesis even under normal/non-ketogenic conditions [41].

Signaling

The tissue growth observed both phenotypically [23] and at the transcription level was induced by many different, but interrelated, signaling pathways, i.e. MAPK, Hedgehog, Wnt, and Phosphatidylinositol signaling pathway. In all tissues including mammary gland, MAPK plays an important role in complex cellular programs such as proliferation, differentiation, development, transformation, and apoptosis [42]. At least three MAPK families have been characterized: extracellular signal-regulated kinase (ERK), Jun kinase (JNK/SAPK) and p38 MAPK [42]. When feeding a higher plane of nutrition, the pattern of DEG in the pathways and the IPA up-stream regulator analysis revealed up-regulation and activation of all three families. Furthermore, the IPA up-stream regulator analysis predicted the activation of other players in these cascades (i.e., pKa, MEK, Akt; Table 3).

The classical Wnt and Hedgehog signaling pathways, both up-regulated in EH calves, have long been known to direct growth and patterning during embryonic development [43], but their activity in the development of the mammary gland was also observed postnatally [44, 45]. The activation of ‘Phosphatidylinositol signaling pathway’ in EH calves, together with the predicted activation of the PI3K complex and Akt group by IPA, underscore the influence of the PI3K/Akt cascade in the regulation of cell proliferation and differentiation. Akt integrates a plethora of extracellular signals to generate diverse outcomes, including proliferation, motility, growth, glucose homeostasis, survival, and cell death [46].

It is likely that many transcription factors (TF), or transcription regulators, at the end of these signaling cascades were activated. Despite not having measured the activation of specific TF directly in our experiment, bioinformatics analysis allowed us to predict activation of eight TF, all possible down-stream targets of one or more of the three cascades. Each has a role in cell proliferation and differentiation, and mammary gland development: EGR1 [47], FOXM1 [48], MYC [49], NFKB1 [50], NOTCH1 [51], SMAD2 and SMAD3 [52], and SP1 [53].

The DIA and the IPA up-stream regulator analysis revealed a spectrum of possible hormones, growth factors, and endogenous chemicals that could have acted as systemic or local mediators of responses linked to the EH MR (e.g., insulin and GH-IGF-1 axis, GnRH, androgen, ANG, EGF, FGF2, NRG1, hyaluronic acid, and tretinoin). A complete discussion of their role and involvement in mediating the nutritional effect of an EH MR can be found in Additional file 5.

Mammary fat pad

Mammary fat pad is essential for development of the secretory epithelium, providing signals that mediate ductal morphogenesis and, probably, alveolar differentiation [54]. As for all other adipose depots of the animal, the adipocytes residing in the MFP are very sensitive to nutritional changes [55]. In the current study, the mass of the MFP increased 5.9-fold when feeding an EH MR [23]. Opposite to PAR, which usually grows via hyperplasia (increasing cell number though proliferation) instead of hypertrophy (increasing cell dimension), the mass of MFP is controlled by a balance of hyperplasia and hypertrophy, due to the presence of proliferating pre-adipocytes and mature adipocytes [56, 57].

The primary stimulus of adipose tissue growth is dietary energy [58]. In the current study the difference in dietary energy supply between the two groups was substantial, as indicated by the 3-fold greater fat intake of EH vs R calves [25]. Judging by the phenotypic differences in the MFP, this stimulated both tissue hyperplasia (e.g., greater tissue DNA), and hypertrophy (e.g. greater lipid accumulation) in EH heifer calves [23]. Despite some discordance, overall, the transcriptome and bioinformatics analysis support both scenarios. The markedly greater ingestion of lipids would have led to an increase in long-chain fatty acids in the blood. The IPA analysis predicted the involvement of fatty acids and LDL (low density lipoproteins) in the observed patterns of DEG. Fatty acids stimulate preadipocyte proliferation through PPAR-δ (up-regulated in EH vs R) [59], while LPL stimulates adipocyte differentiation and maturation [60]. The activation of preadipocyte proliferation in the EH calf MFP was also supported by the activation of the ‘Insulin signaling pathway’ [61] and nucleotide metabolism (purine and pyrimidine). The IPA predicted activation of AGT [62, 63], PPGF BB [64], TNF [63], and F2 and VEGFA [65, 66], all signaling molecules that have been shown to stimulate this process.

Prior to actively accumulating lipid, the pre-adipocyte must differentiate and mature into adipocytes. Several molecules, besides LDL, that induce and enhance this process were predicted to be activated in EH calves: Pkc(s) [67], Creb (CREB1 in particular in our results) [68, 69], GH1 [70], FGF2 [71], PDGF BB [64], and PPARG [72]. Furthermore, the TF SMAD7, which blocks premature differentiation [73], was predicted to be inhibited in response to the EH MR, and Wnt and the ‘Wnt signaling pathway’, a cascade known to inhibit commitment to the adipogenic lineage and adipogenesis [74] and to be downregulated by PPARγ [75], were predicted to be down-regulated.

Once the adipocytes have matured, their metabolism shifts to that of typical adipose tissue, accumulating lipids in their droplets and releasing them when needed [56]. When fed an EH MR, the lipid content of the MFP increased ~ 6 fold [23]. This accumulation was mainly due to storage of excess dietary lipid from long-chain fatty acid uptake, or because of intracellular metabolism of the fatty acids such as desaturation as indicated by the activation of ‘Biosynthesis of unsaturated fatty acids’ pathway. Previous data [76] showed no effect of pre-weaning plane of nutrition on blood insulin concentration, but the accumulation of lipids in the current study was probably driven by the activation of the insulin cascade, as suggested by the up-regulation of the ‘Insulin signaling pathway’. Further discussion on the MPF lipid and energy metabolism can be found in Additional file 5.

The mammary gland as an organ: Tissues interaction and cross-talk

The transcriptome analysis of both PAR and MFP provided great insights on the effect of pre-weaning nutrition in pre-pubertal mammary development. However, despite gaining a holistic perspective on each of its component, the same perspective is lost when considering the mammary gland as a whole organ. In fact, the two tissues are not separate entities, rather collaborative players in the development of the gland. The PAR expansion is of primary interest for the producer in the context of maximizing future milk production. However, the MFP plays a fundamental role in PAR development [77]. In fact, the MFP, together with its various constituents, facilitates many functions of the gland: (i) it houses a vascular and lymphatic system, (ii) it provides a three-dimensional matrix for growth and expansion of the PAR, and (iii) it functions as a local site for hormone action, the provision of lipids, and growth factor synthesis [77].

Regarding growth factors potentially released, the MFP from EH heifers tended to secrete more IGF-1 (1.4x; FDR = 0.11), IGFALS (1.7x; FDR = 0.04), FGF-1 (1.8x, FDR = 0.09), and FGF-2 (1.3x; FDR = 0.09). Of these, only FGF-2, a confirmed mitogen in bovine mammary gland [77], was predicted to be activated by the IPA analysis in both PAR and MFP. Despite IGF-1 not being predicted to be involved in the observed changes in PAR due to EH MR, the up-regulation of its binding proteins IGFALS and IGFBP2, in MFP and PAR, respectively, suggests a possible prolonged effect that might have contributed to the enhanced mammary gland development in EH heifer calves.

Immune cell involvement

Recent studies have illustrated the importance of immune cells and their mediators during the various stages of mammary gland development [78].

Information about the involvement of immune cells in mammary development using non-ruminant models has been extensively reviewed elsewhere [78]. In these models, immune cells such as eosinophils, mast cells, macrophages and T-cell occupy unique sites during the various stages of mammary gland development where they contribute to a diversity of effector functions. Bovine specific data [79] suggest these cells cluster closely with the terminal ductal units, clusters of ductules arising from large ducts typical of bovine, in the developing bovine mammary gland. The DIA ‘Immune system’ category, and the IPA canonical pathways (Fig. 6), upstream regulator analysis (Tables 3 and 4), and immune-specific networks (Figs. 7, 8, 9, 10) displayed a substantial involvement of immune-related functions in both tissues. As it is hard to precicely separate the PAR from its surrounding stroma, it is not surprising that expression of immune cell markers for all aformentioned cell-types were detected in both MFP and PAR: CD4 and CD32 (general immune cells marker), CD68 (macrophage), CD301 (M2-macrophage), EMR1 (eosinophil), and FCER1G (mast cell). Furthermore, CD11d (M1-macrophage), FCER1A (mast cell), and CXCR3 (Th1 T cell) were expressed in PAR only. Their expression not only confirms the previously determined involvement of immune cells in dairy cow mammary gland development [79], but results indicated an effect of plane of nutrition on immune cell quantity within each tissue.

According to a model based on non-ruminant species [78], mast cells first localize to the stromal regions surrounding the terminal end buds (TEB), bulb-shaped structures that direct the growth of the ducts throughout the rest of the fat pad. Next, macrophages are recruited and migrate and localize to the neck of TEBs, and lastly, eosinophils are recruited to the head of TEBs. Bovine specific data confirm the recuitment of eosinophils at latter stage of early development, while no clear priority between macrophages and mast cells seemed to emerge [79]. As previously shown by differences in marker expression (FCER1A and FCER1G), a greater number of mast cells can be hypothesized to reside both in the MFP and around the ducts in the PAR of EH calves. The absence of mast cells was shown to reduce TEBs and number of ducts, and ductal length in mouse, suggesting a fundamental role of this immune cell type in mammary gland development [80]. The mechanism of action is still unclear, but mediators contained in their granules, such as histamine, and serine proteases and their activating enzyme DPPI (dipeptidyl peptidase I), seem to be involved [80]. Bioinformatics analysis indicated that histamine synthesis might have been enhanced in PAR of EH calves, supporting the role of mast cells in the mediation of the greater mammary gland mass in these calves. Furtermore, the mast cell specific triptase beta 2 (a serine protease) was upregulated (1.8x, FDR = 0.05) in PAR of EH heifer calves.

Together with mast cell, macrophages seem to be the primary immune cell type recruited to the bovine mammary gland [79]. CD68 expression suggested a higher amount in both the MFP and PAR of EH calves, with a tendency for a greater population of anti-inflammatory M2 macrophages in the MFP. Despite being infiltrated in both tissues, the MFP seemed to be primarily responsible for their recruitment. In fact, no changes in the chemokines CCL2 and CX3CL1 were detected in PAR, while they were upregulated (CCL2, 1.9x, FDR = 0.12; CXC3CL1, 1.3x, FDR = 0.05) in MFP of EH heifer calves.

In mice, once recruited by the MFP, macrophages move closer to the neck of the TEBs following the PAR signal through CSF1 (Colony stimulating factor-1), a key player in macrophage proliferation and survival [81]. Here they stimulate ductal elongation through reorganization of surrounding collagen [78]. Indeed, feeding the EH MR increased CSF1 expression in the MFP (1.4x, FDR = 0.01), but not in the PAR, and this molecule was predicted to be part of the tissue cross-talk network as an MFP signal to the PAR. Despite parenchymal tissue being extraneous to the recruitment of monocytes to the mammary gland, EH feeding increases PAR expression of CXCL10 (C-X-C Motif Chemokine Ligand 10), a chemokine able to recruit monocytes [82]. As suggested by the generated tissue crosstalk network, this molecule could have affected the MFP, enhancing monocyte recruitment in the surrounding tissue, or could have worked as a local attractant in the PAR to enhance the migration of macrophages from the MFP to the PAR.

Lastly, as previously reported [79], the recruitment of eosinophils to the mammary gland can also be hypothesized in the current animals. In fact, in the MFP, but not in PAR, we were able to detect the expression of EMR1, an eosinophil marker. Despite the lack of difference due to plane of nutrition in the expression of this marker, the MFP of EH calves had an upregulation of eotaxin CCL24 (1.2x, FDR = 0.05), an eosinophil chemoattractant. Furthermore, the IPA predicted the activation of IL5 in the MFP of EH calves. IL5 is a potent eosinophil chemoattractant, which can also increase their activity. In mice lacking either eotaxins or IL5, TEB numbers and ductal branching were impaired [83, 84]. Together with macrophages, eosinophils are a source of TGFβ, a growth factor that not only was predicted to be activated by the IPA in both MFP and PAR, but also was predicted to participate in the crosstalk between MFP and PAR in EH calves. TGFβ might be the primary mediator of eosinophil action on mammary development. This growth factor has an important role in coordinating patterning and growth of the mammary ductal tree, reducing excessive elongation of the single ducts, favoring instead the formation of a complex 3D structure [85]. Besides TGFβ, eosinophils can also release epidermal growth factor (EGF), which supports epithelial cell maintenance during processes of tissue development, repair, and remodeling [78]. EGF was among the predicted activated growth factors involved in the PAR response to an EH MR.

Signs of T-cell infiltration were observed in both PAR and MFP, but differences between nutritional groups were observed only in the PAR, as EH calves had greater expression of CXCL9, a T-cell chemoattractant, and CD3 (subunits ε, 1.3x, FDR = 0.05; and γ, 1.4x, FDR = 0.04), a marker of the T-cell lineage. T helper cell 1 (Th1) cytokines, including IL-2 and interferon-γ (IFNγ), are responsible for directing cell-mediated immunity, including activation of macrophages, while T helper cell 2 (Th2) cytokines such as IL-4, IL-5, IL-6, IL-10, and IL-13 are involved in processes such as matrix deposition and tissue remodeling, all functions needed for mammary gland development [78]. Among these, IL-2, IL-4, IL-5, and INFγ all were predicted to be activated by the IPA, suggesting a role of T-cells in the immune-cell mediated development of the mammary gland induced by an EH MR. Their putative action in the current study was mainly as recruiter and enhancer of the activity of other immune cell types.

Conclusions

The current bioinformatics analyses of the transcriptome highlight the physiological pathways and signaling cascades that mediate the heifer response to an enhanced MR plane of nutrition in the pre-weaning period. Furthermore, they suggest an increased degree of cross-talk between PAR tissue and the surrounding MFP, with the latter increasing the local release of mammogenic signals and orchestrating the infiltration of immune cells with key roles in the development of the mammary gland.

In light of the current and previous results, enhancing nutrient intake in the pre-weaning period of Holstein heifer calves appears to be a feasible strategy to stimulate mammary development prior to onset of first lactation. Future studies will need to assess lactation performance of calves fed an accelerate milk replacer prior to weaning.

Abbreviations

ADF:

Acid detergent fiber

AGT:

Angiotensinogen

Akt:

protein kinase B

ANGPT2:

Angiopoietin 2

AR:

Androgen receptor

BW:

Body weight;

CD11d:

Cluster of differentiation 11d

CD301:

Cluster of differentiation 301

CD32:

Cluster of differentiation 32

CD4:

Cluster of differentiation 4

CD68:

Cluster of differentiation 68

CP:

Crude protein

Creb:

cAMP response element binding protein

CSF1:

Colony stimulating factor 1

CXCL10:

C-X-C motif chemokine ligand 10

CXCL5:

C-X-C motif chemokine ligand 5

CXCL9:

C-X-C motif chemokine ligand 9

CXCR3:

C-X-C motif chemokine receptor 3

DEG:

Differentially expressed genes

DIA:

Dynamic impact approach

DM:

Dry matter

DPPI:

Dipeptidyl-peptidase I

EGF:

Epidermal growth factor

EGR1:

Early growth response 1

EH:

Enhanced milk replacer

EMR1:

EGF-like module-containing mucin-like hormone receptor-like 1

ERK:

Extracellular signal-regulated kinase

F2:

coagulant factor II (prothrombin)

FC:

Fold change

FCER1A:

Fc fragment of IgE receptor Ia

FCER1G:

Fc fragment of IgE receptor Ig

FDR:

False discovery rate

FGF-1:

Fibroblast growth factor 1

FGF2:

Fibroblast growth factor 2

FOXM1:

Forkhead box M1

FSH:

Follicle-stimulating hormone

GH:

Growth hormone

GH1:

Growth hormone 1

GnRH:

Gonadotrophin releasing hormone

IGF-1:

Inusulin-like growth factor

IGFALS:

Insulin like growth factor binding protein acid labile subunit

IL-10:

Interleukin 10

IL-13:

Interleukin 13

IL1B:

Interleukin 1β

IL-2:

Interleukin 2

IL-4:

Interleukin 4

IL-5:

Interleukin 5

IL-6:

Interleukin 6

INFγ:

Interferon γ

IPA:

Ingenuity pathway analysis

JAK/STAT:

Janus kinase/signal transducer and activator of transcription

JNK/SAPK:

stress-activated protein kinase/c-Jun NH(2)-terminal kinase

LDL:

Low density lioprotein

MAPK:

Mitogen-activated protein kinase

MFP:

Mammary fat pad

MR:

Milk replacer

MTG1:

Mitochondrial ribosome associated GTPase 1

MYC:

MYC proto-oncogene

NDF:

Neutral deternget fiber

NFKB1:

Nuclear factor κ B Subunit 1

NGF:

Nerve growth factor

NOTCH1:

Notch 1

NRG1:

Neuregulin 1

PAR:

mammary parenchyma

PDGF BB:

Platelet-derived growth factor BB

PI3K:

Phosphoinositide 3-kinase

Pkc:

Protein kinase C

PPAR:

Peroxisome proliferator-activated receptor

PPARG:

Peroxisome proliferator-activated receptor γ.

PPP1R11:

Protein phosphatase 1 regulatory inhibitor subunit 11

R:

Restricted milk replacer

RPS15A:

Ribosomal protein S15a

SMAD2:

SMAD family member 2

SMAD3:

SMAD family member 3

SP1:

transcription factor Sp1

TDU:

Terminal duct unit

TEB:

Terminal end bud

TF:

Trasncription factor

TGFB1:

Transforming growth factor beta 1

TGFβ:

Transforming growth factor β

Th1:

Type 1 T helper

Th2:

Type 2 T helper

TNF:

Tumor necrosis factor

VEGFA:

Vascular endothelial growth factor A

References

  1. 1.

    Radcliff RP, VandeHaar MJ, Skidmore AL, Chapin LT, Radke BR, Lloyd JW, Stanisiewski EP, Tucker HA. Effects of diet and bovine somatotropin on heifer growth and mammary development. J Dairy Sci. 1997;80(9):1996–2003.

  2. 2.

    Heinrichs AJ. Raising dairy replacements to meet the needs of the 21st century. J Dairy Sci. 1993;76(10):3179–87.

  3. 3.

    Davis Rincker LE, Weber Nielsen MS, Chapin LT, Liesman JS, Daniels KM, Akers RM, Vandehaar MJ. Effects of feeding prepubertal heifers a high-energy diet for three, six, or twelve weeks on mammary growth and composition. J Dairy Sci. 2008;91(5):1926–35.

  4. 4.

    Sejrsen K. Relationships between nutrition, puberty and mammary development in cattle. Proc Nutr Soc. 1994;53(1):103–11.

  5. 5.

    Tozer PR, Heinrichs AJ. What affects the costs of raising replacement dairy heifers: a multiple-component analysis. J Dairy Sci. 2001;84(8):1836–44.

  6. 6.

    Sejrsen K, Huber JT, Tucker HA, Akers RM. Influence of nutrition of mammary development in pre- and postpubertal heifers. J Dairy Sci. 1982;65(5):793–800.

  7. 7.

    Swanson EW. Effect of rapid growth with fattening of dairy heifers on their lactation ability. J Dairy Sci. 1960;43(3):377–87.

  8. 8.

    Little W, Kay RM. The effect of rapid rearing and early calving on subsequent performance of dairy heifers. Anim Prod. 1979;29:131.

  9. 9.

    Gardner RW, Schuh JD, Vargus LG. Accelerated growth and early breeding of Holstein heifers. J Dairy Sci. 1977;60:1941–8.

  10. 10.

    Sejrsen K. Mammary development and milk yield in relation to growth rate in dairy and dual-purpose heifers. Acta Agric Scand. 1978;28:41–6.

  11. 11.

    Bettenay RA. Effect of growth rate and mating age of dairy heifers on subsequen production over four years. Aust J Exp Agric. 1985;25:263–9.

  12. 12.

    Harrison RD, Reynolds IP, Little W. A quantitative analysis of mammary glands of dairy heifers reared at different rates of live weight gain. J Dairy Res. 1983;50(4):405–12.

  13. 13.

    Hohenboken WD, Foldager J, Jensen J, Madens P, Andersen BB. Breed and nutritional effects and interactions on energy intake, production, and efficiency of nutrient utilization in young bulls, heifers, and lactating cows. Acta Agri Scand. 1995;45:92–8.

  14. 14.

    Peri I, Gertler A, Bruckental I, Barash H. The effect of manipulation in energy allowance during the rearing period of heifers on hormone concentrations and milk production in first lactation cows. J Dairy Sci. 1993;76(3):742–51.

  15. 15.

    Petitclerc D, Chapin LT, Tucker HA. Carcass composition and mammary development responses to photoperiod and plane of nutrition in Holstein heifers. J Anim Sci. 1984;58(4):913–9.

  16. 16.

    Sejrsen K, Foldager J. Mammary growth and milk production capacity of replacement heifers in relation to diet energy concentration and plasma hormone levels. Acta Agric Scand. 1992;42:99–105.

  17. 17.

    Radcliff RP, Vandehaar MJ, Chapin LT, Pilbeam TE, Beede DK, Stanisiewski EP, Tucker HA. Effects of diet and injection of bovine somatotropin on prepubertal growth and first-lactation milk yields of Holstein cows. J Dairy Sci. 2000;83(1):23–9.

  18. 18.

    Heinrichs AJ, Heinrichs BS. A prospective study of calf factors affecting first-lactation and lifetime milk production and age of cows when removed from the herd. J Dairy Sci. 2011;94(1):336–41.

  19. 19.

    Brown EG, Vandehaar MJ, Daniels KM, Liesman JS, Chapin LT, Forrest JW, Akers RM, Pearson RE, Nielsen MS. Effect of increasing energy and protein intake on mammary development in heifer calves. J Dairy Sci. 2005;88(2):595–603.

  20. 20.

    Daniels KM, Capuco AV, McGilliard ML, James RE, Akers RM. Effects of milk replacer formulation on measures of mammary growth and composition in Holstein heifers. J Dairy Sci. 2009;92(12):5937–50.

  21. 21.

    Meyer MJ, Capuco AV, Ross DA, Lintault LM, Van Amburgh ME. Developmental and nutritional regulation of the prepubertal bovine mammary gland: II. Epithelial cell proliferation, parenchymal accretion rate, and allometric growth. J Dairy Sci. 2006;89(11):4298–304.

  22. 22.

    Meyer MJ, Capuco AV, Ross DA, Lintault LM, Van Amburgh ME. Developmental and nutritional regulation of the prepubertal heifer mammary gland: I. parenchyma and fat pad mass and composition. J Dairy Sci. 2006;89(11):4289–97.

  23. 23.

    Geiger AJ, Parsons CL, Akers RM. Feeding a higher plane of nutrition and providing exogenous estrogen increases mammary gland development in Holstein heifer calves. J Dairy Sci. 2016;99(9):7642–53.

  24. 24.

    Soberon F, Raffrenato E, Everett RW, Van Amburgh ME. Preweaning milk replacer intake and effects on long-term productivity of dairy calves. J Dairy Sci. 2012;95(2):783–93.

  25. 25.

    Geiger AJ, Parsons CL, James RE, Akers RM. Growth, intake, and health of Holstein heifer calves fed an enhanced preweaning diet with or without postweaning exogenous estrogen. J Dairy Sci. 2016;99(5):3995–4004.

  26. 26.

    Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30(15):2114–20.

  27. 27.

    Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, Batut P, Chaisson M, Gingeras TR. STAR: ultrafast universal RNA-seq aligner. Bioinformatics. 2013;29(1):15–21.

  28. 28.

    Zhou Y, Lin N, Zhang B. An iteration normalization and test method for differential expression analysis of RNA-seq data. BioData Min. 2014;7:15.

  29. 29.

    Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, Smyth GK. Limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 2015;43(7):e47.

  30. 30.

    Storey JD, Tibshirani R. Statistical significance for genomewide studies. Proc Natl Acad Sci U S A. 2003;100(16):9440–5.

  31. 31.

    Bionaz M, Periasamy K, Rodriguez-Zas SL, Hurley WL, Loor JJ. A novel dynamic impact approach (DIA) for functional analysis of time-course omics studies: validation using the bovine mammary transcriptome. PLoS One. 2012;7(3):e32455.

  32. 32.

    Moyes KM, Sorensen P, Bionaz M. The impact of Intramammary Escherichia coli challenge on liver and mammary transcriptome and cross-talk in dairy cows during early lactation using RNAseq. PLoS One. 2016;11(6):e0157480.

  33. 33.

    Khan MA, Weary DM, von Keyserlingk MA. Invited review: effects of milk ration on solid feed intake, weaning, and performance in dairy heifers. J Dairy Sci. 2011;94(3):1071–81.

  34. 34.

    Drackley JK. Calf nutrition from birth to breeding. Vet Clin North Am Food Anim Pract. 2008;24(1):55–86.

  35. 35.

    Geiger AJ, Parsons CLM, Akers RM. Feeding an enhanced diet to Holstein heifers during the preweaning period alters steroid receptor expression and increases cellular proliferation. J Dairy Sci. 2017;100(10):8534–43.

  36. 36.

    Nelson CM, Chen CS. Cell-cell signaling by direct contact increases cell proliferation via a PI3K-dependent signal. FEBS Lett. 2002;514(2–3):238–42.

  37. 37.

    Koch S, Claesson-Welsh L. Signal transduction by vascular endothelial growth factor receptors. Cold Spring Harb Perspect Med. 2012;2(7):a006502.

  38. 38.

    Wilkens S. Structure and mechanism of ABC transporters. F1000Prime Rep. 2015;7:14.

  39. 39.

    Palm W, Park Y, Wright K, Pavlova NN, Tuveson DA, Thompson CB. The utilization of extracellular proteins as nutrients is suppressed by mTORC1. Cell. 2015;162(2):259–70.

  40. 40.

    DeBerardinis RJ, Lum JJ, Hatzivassiliou G, Thompson CB. The biology of cancer: metabolic reprogramming fuels cell growth and proliferation. Cell Metab. 2008;7(1):11–20.

  41. 41.

    Baird GD: Ketogenesis in slices of liver and lactating mammary gland from normal and ketotic cows. Biochim Biophys Acta 1965, 111(1):339–341.

  42. 42.

    Zhang W, Liu HT. MAPK signal pathways in the regulation of cell proliferation in mammalian cells. Cell Res. 2002;12(1):9–18.

  43. 43.

    Taipale J, Beachy PA. The hedgehog and Wnt signalling pathways in cancer. Nature. 2001;411(6835):349–54.

  44. 44.

    Hatsell S, Frost AR. Hedgehog signaling in mammary gland development and breast cancer. J Mammary Gland Biol Neoplasia. 2007;12(2–3):163–73.

  45. 45.

    Yu QC, Verheyen EM, Zeng YA. Mammary development and breast Cancer: a Wnt perspective. Cancers (Basel). 2016;8(7).

  46. 46.

    Wickenden JA, Watson CJ. Key signalling nodes in mammary gland development and cancer. Signalling downstream of PI3 kinase in mammary epithelium: a play in 3 Akts. Breast Cancer Res. 2010;12(2):202.

  47. 47.

    Tarcic G, Avraham R, Pines G, Amit I, Shay T, Lu Y, Zwang Y, Katz M, Ben-Chetrit N, Jacob-Hirsch J, et al. EGR1 and the ERK-ERF axis drive mammary cell migration in response to EGF. FASEB J. 2012;26(4):1582–92.

  48. 48.

    Wierstra I, Alves J. FOXM1, a typical proliferation-associated transcription factor. Biol Chem. 2007;388(12):1257–74.

  49. 49.

    Hynes NE, Stoelzle T. Key signalling nodes in mammary gland development and cancer: Myc. Breast Cancer Res. 2009;11(5):210.

  50. 50.

    Brantley DM, Chen CL, Muraoka RS, Bushdid PB, Bradberry JL, Kittrell F, Medina D, Matrisian LM, Kerr LD, Yull FE. Nuclear factor-kappaB (NF-kappaB) regulates proliferation and branching in mouse mammary epithelium. Mol Biol Cell. 2001;12(5):1445–55.

  51. 51.

    Bouras T, Pal B, Vaillant F, Harburg G, Asselin-Labat ML, Oakes SR, Lindeman GJ, Visvader JE. Notch signaling regulates mammary stem cell function and luminal cell-fate commitment. Cell Stem Cell. 2008;3(4):429–41.

  52. 52.

    Gupta V, Thiagalingam A, Maheswaran S. Smad-signaling in mammary gland development and tumorigenesis. Curr Signal Transd T. 2007;2(2):111–20.

  53. 53.

    Lu S, Archer MC. Sp1 coordinately regulates de novo lipogenesis and proliferation in cancer cells. Int J Cancer. 2010;126(2):416–25.

  54. 54.

    Hovey RC, Aimo L. Diverse and active roles for adipocytes during mammary gland growth and function. J Mammary Gland Biol Neoplasia. 2010;15(3):279–90.

  55. 55.

    Lohakare JD, Sudekum KH, Pattanaik AK. Nutrition-induced changes of growth from birth to first calving and its impact on mammary development and first-lactation Milk yield in dairy heifers: a review. Asian Australas J Anim Sci. 2012;25(9):1338–50.

  56. 56.

    Roh SG, Hishikawa D, Hong YH, Sasaki S. Control of adipogenesis in ruminants. Anim Sci J. 2006;77(5):472–7.

  57. 57.

    Neville MC, Medina D, Monks J, Hovey RC. The mammary fat pad. J Mammary Gland Biol Neoplasia. 1998;3(2):109–16.

  58. 58.

    Rutkowski JM, Stern JH, Scherer PE. The cell biology of fat expansion. J Cell Biol. 2015;208(5):501–12.

  59. 59.

    Jehl-Pietri C, Bastie C, Gillot I, Luquet S, Grimaldi PA. Peroxisome-proliferator-activated receptor delta mediates the effects of long-chain fatty acids on post-confluent cell proliferation. Biochem J. 2000;350(Pt 1):93–8.

  60. 60.

    Stanton LA, van de Venter M, Litthauer D, Oelofsen W. Effect of lipoproteins on the differentiation of 3T3-L1 and human preadipocytes in cell culture. Comp Biochem Physiol Part B Biochem Mol Biol. 1997;116(1):65–73.

  61. 61.

    Geloen A, Collet AJ, Guay G, Bukowiecki LJ. Insulin stimulates in vivo cell proliferation in white adipose tissue. Am J Phys. 1989;256(1 Pt 1):C190–6.

  62. 62.

    Crandall DL, Armellino DC, Busler DE, McHendry-Rinde B, Kral JG. Angiotensin II receptors in human preadipocytes: role in cell cycle regulation. Endocrinology. 1999;140(1):154–8.

  63. 63.

    Hausman DB, DiGirolamo M, Bartness TJ, Hausman GJ, Martin RJ. The biology of white adipocyte proliferation. Obes Rev. 2001;2(4):239–54.

  64. 64.

    Virakul S, Dalm VA, Paridaens D, van den Bosch WA, Mulder MT, Hirankarn N, van Hagen PM, Dik WA. Platelet-derived growth factor-BB enhances Adipogenesis in orbital fibroblasts. Invest Ophthalmol Vis Sci. 2015;56(9):5457–64.

  65. 65.

    Strande JL, Phillips SA. Thrombin increases inflammatory cytokine and angiogenic growth factor secretion in human adipose cells in vitro. J Inflamm (Lond). 2009;6:4.

  66. 66.

    Bagchi M, Kim LA, Boucher J, Walshe TE, Kahn CR, D'Amore PA. Vascular endothelial growth factor is important for brown adipose tissue development and maintenance. FASEB J. 2013;27(8):3257–71.

  67. 67.

    Fleming I, MacKenzie SJ, Vernon RG, Anderson NG, Houslay MD, Kilgour E. Protein kinase C isoforms play differential roles in the regulation of adipocyte differentiation. Biochem J. 1998;333(Pt 3):719–27.

  68. 68.

    Zhang JW, Klemm DJ, Vinson C, Lane MD. Role of CREB in transcriptional regulation of CCAAT/enhancer-binding protein beta gene during adipogenesis. J Biol Chem. 2004;279(6):4471–8.

  69. 69.

    Reusch JE, Colton LA, Klemm DJ. CREB activation induces adipogenesis in 3T3-L1 cells. Mol Cell Biol. 2000;20(3):1008–20.

  70. 70.

    Chaves VE, Junior FM, Bertolini GL. The metabolic effects of growth hormone in adipose tissue. Endocrine. 2013;44(2):293–302.

  71. 71.

    Kakudo N, Shimotsuma A, Kusumoto K. Fibroblast growth factor-2 stimulates adipogenic differentiation of human adipose-derived stem cells. Biochem Biophys Res Commun. 2007;359(2):239–44.

  72. 72.

    Bionaz M, Chen S, Khan MJ, Loor JJ. Functional role of PPARs in ruminants: potential targets for fine-tuning metabolism during growth and lactation. PPAR Res. 2013;2013:684159.

  73. 73.

    Choy L, Derynck R. Transforming growth factor-beta inhibits adipocyte differentiation by Smad3 interacting with CCAAT/enhancer-binding protein (C/EBP) and repressing C/EBP transactivation function. J Biol Chem. 2003;278(11):9609–19.

  74. 74.

    Christodoulides C, Lagathu C, Sethi JK, Vidal-Puig A. Adipogenesis and WNT signalling. Trends Endocrinol Metab. 2009;20(1):16–24.

  75. 75.

    Moldes M, Zuo Y, Morrison RF, Silva D, Park BH, Liu J, Farmer SR. Peroxisome-proliferator-activated receptor gamma suppresses Wnt/beta-catenin signalling during adipogenesis. Biochem J. 2003;376(Pt 3):607–13.

  76. 76.

    Daniels KM, Hill SR, Knowlton KF, James RE, McGilliard ML, Akers RM. Effects of milk replacer composition on selected blood metabolites and hormones in preweaned Holstein heifers. J Dairy Sci. 2008;91(7):2628–40.

  77. 77.

    Hovey RC, McFadden TB, Akers RM. Regulation of mammary gland growth and morphogenesis by the mammary fat pad: a species comparison. J Mammary Gland Biol Neoplasia. 1999;4(1):53–68.

  78. 78.

    Reed JR, Schwertfeger KL. Immune cell location and function during post-natal mammary gland development. J Mammary Gland Biol Neoplasia. 2010;15(3):329–39.

  79. 79.

    Beaudry KL, Parsons CL, Ellis SE, Akers RM. Localization and quantitation of macrophages, mast cells, and eosinophils in the developing bovine mammary gland. J Dairy Sci. 2016;99(1):796–804.

  80. 80.

    Lilla JN, Werb Z. Mast cells contribute to the stromal microenvironment in mammary gland branching morphogenesis. Dev Biol. 2010;337(1):124–33.

  81. 81.

    Hamilton JA. Colony-stimulating factors in inflammation and autoimmunity. Nat Rev Immunol. 2008;8(7):533–44.

  82. 82.

    Petrovic-Djergovic D, Popovic M, Chittiprol S, Cortado H, Ransom RF, Partida-Sanchez S. CXCL10 induces the recruitment of monocyte-derived macrophages into kidney, which aggravate puromycin aminonucleoside nephrosis. Clin Exp Immunol. 2015;180(2):305–15.

  83. 83.

    Colbert DC, McGarry MP, O'Neill K, Lee NA, Lee JJ. Decreased size and survival of weanling mice in litters of IL-5−/ -mice are a consequence of the IL-5 deficiency in nursing dams. Contemp Top Lab Anim Sci. 2005;44(3):53–5.

  84. 84.

    Gouon-Evans V, Rothenberg ME, Pollard JW. Postnatal mammary gland development requires macrophages and eosinophils. Development. 2000;127(11):2269–82.

  85. 85.

    Daniel CW, Robinson S, Silberstein GB. The role of TGF-beta in patterning and growth of the mammary ductal tree. J Mammary Gland Biol Neoplasia. 1996;1(4):331–41.

Download references

Acknowledgements

The authors would like to acknowledge Land O’ Lakes, Inc. (St. Paul, MN) for providing milk replacer for the trial. The authors would like thank the Virginia Tech farm crew and veterinary staff for assistance throughout the trial.

Funding

This research was supported by USDA Agricultural and Food Research Initiative (AFRI) competitive grants program, no. 2016–67015-24565 “Impact of Pre-weaning Nutrition on Endocrine Induction of Mammary Development in Dairy Heifers”, and 2016–67011-24703, predoctoral fellowship awarded to A. J. Geiger.

Availability of data and materials

The sequencing data have been submitted to Gene Expression Omnibus database (GEO Accession: GSE102435).

Author information

RMA and AG designed and coordinated the animal trial. JCM and SZ performed statistical analysis of the sequencing reads, and MVR, with help from SZ, performed the bioinformatics analysis. REB and MVR performed qPCR. MVR, wrote the main draft of the manuscript, with inputs from JJL and RMA. All authors read and approved the final manuscript.

Correspondence to J. J. Loor.

Ethics declarations

Ethics approval

This experiment was conducted under the review and approval of the Virginia Polytechnic Institute and State University Institutional Animal Care and Use Committee (#14–045-DASC).

Consent for publication

Non applicable.

Competing interests

The authors declare that they have no competing interests.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Additional files

Additional file 1:

Complete details of the protocols used for library preparation, primer design and sequences, and qPCR analysis and performance. (DOCX 28 kb)

Additional file 2:

Summary of read counts per sample and details of the alignment procedure performed with the STAR 2.5.1b package (DOCX 14 kb)

Additional file 3:

Complete result datasets generated by the DIA analysis. The KEGG pathways are sorted by category and subcategory. Blue bars represent the Impact along the transition period, while red and green bars show the Direction of the Impact (red = upregulation, green = downregulation). (XLSX 60 kb)

Additional file 4:

Complete list of genes considered in the development of the tissue cross-talk networks between mammary fat pad and parenchyma. (XLSX 37 kb)

Additional file 5:

Extended discussion on Parenchymal metabolism and molecular signaling, and Fat Pad lipid and energy metabolism. (DOCX 75 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Vailati-Riboni, M., Bucktrout, R.E., Zhan, S. et al. Higher plane of nutrition pre-weaning enhances Holstein calf mammary gland development through alterations in the parenchyma and fat pad transcriptome. BMC Genomics 19, 900 (2018) doi:10.1186/s12864-018-5303-8

Download citation

Keywords

  • Mammary gland
  • Milk replacer
  • Transcriptome
  • Calf