Data were collated from five buffalo experimental herds (El-Nattafe El-Gadid, El-Nattafe El-Kadim, Nubariya, Sids and Gimeza) belonging to the Animal Production Research Institute (APRI), Agriculture research center (ARC), Ministry of Agriculture and Land Reclamation, Egypt. A total of 16,370 records were available on 8,327 buffalo calves in which 8,271 and 8,099 records for BW and WW, respectively, were included in the analysis. These calves were born between 1980 and 2018 with pedigree having 188 sires and 2,211 dams (Table 1).

The systematic environmental effects on growth performance traits were examined as fixed effects using least squares methods implemented in general linear model procedure of SAS [9]. These fixed effects included the effects of season of calving (4 seasons), year of calving (39 years), sex (male and female), herd (5 herds), gestation length (GL; 4 levels), weight of dam (WD; 5 levels) and parity (14 parities). The effect of GL had 4 classes; between 295 and 325 days with increment of 15 days. Whereas the WD had 5 classes between 300 and 500 kg with increment of 50 kg. The linear model was fitted as follows:

Where, Y_ijklmnop, the phenotypic record of BW or WW; μ, the effect of the intercept; A_i, the fixed effect of ith season of calving; B_j, the fixed effect of jth year of calving; C_k, the fixed effect of kth sex; D_l, the fixed effect of lth herd; G_m, the fixed effect of mth gestation length; W_n, the fixed effect of nth WD; P_o, the fixed effect of oth parity and e_ijklmnop, random residual assumed to be independent normally distributed with mean zero and variance σ²_e. The significant fixed effects were used to form contemporary groups, which were included in genetic parameters analysis.

Variance components, heritability and breeding values were estimated using the bivariate mixed animal model using Wombat software [10]. The model can be described in the matrix notation follows:

And the assumptions of the variances were as:

Where y, a vector of BW or WW observations; b, a vector of fixed effects with an incidence matrix X; a, a vector of random animal effects with incidence matrix Z; m, a vector of random dam effects with incidence matrix S; c, a vector of random permanent environmental effects of dam with incidence matrix W; and e, a vector of random residual effects. A is the numerator relationship matrix between animals, I is an identity matrix, σ²_a is direct additive genetic variance, σ²_m is maternal additive genetic variance, σ²_am is the direct × maternal genetic covariance, σ²_c is maternal permanent environmental variance, and σ²_e is residual variance. Direct and maternal genetic trends were estimated by regressing yearly mean estimates of breeding value of animal and dam on year of calving.

The principal component (PC) could be used as index to evaluate animals for multiple traits [6]. First, in order to avoid the effect of traits’ different scaling and magnitudes, the estimated direct and maternal breeding values obtained for BW and WW were standardized for mean zero and unit variance. The standardization was performed according to Boligon et al [6] as follows: z_i = (x_i− χ̄)/s_i, where z_i is the standardized value of x_i variable, χ̄ is the mean of ith trait, and s_i is the corresponding standard deviation. Then, the standardized estimated breeding values (SEBVs) were used in the PCA performed using FactoMineR package [11] in R software. Generally, the PCA reduces the information contained in estimated direct and maternal breeding values for BWs and WWs in fewer orthogonal latent variables, PCs, with minimal loss of information or variability [6]. In the current study, the PCs explained the highest percentage of variance (i.e with values greater than one eigenvalues) were used. Using the SEBVs in this analysis, each PC can generate a new value, called a principal component score, which is the result of the sum of the SEBVs of each trait weighed by their respective standardized score coefficient (SSC) [12]. Thus, the PC could be used to evaluate animals for multiple traits. The SSC of each SEBV in each principal component were obtained using the following formula:

Where SSC_ij is the SSC for SEBV of ith trait in the jth PC. The principal component scores were calculated as follows:

Where PCS_jl is the principal component score for the lth animal in the jth principal component, SSC_ij is the SSC for EBVs of ith trait in the jth PC and SEBV_il is the standardised estimated breeding value of the ith trait for the lth animal.

The current study assessed the non-genetic factors affected the growth performance (i.e. BW and WW) in Egyptian buffalo (Table 2). The sex of calf significantly affected (p< 0.0001) BW and WW (p<0.05) where male calves were heavier at birth and weaning than female calves. These differences due to sex effects were in agreement with reported results by Pandya et al [4] in Surti buffalo. These differences are probably due to physiological mechanisms associated with sexual endocrine systems which cause the differences between masculinity and feminist characteristics and the natural buildup of male and female bodies [13]. Season of calving also had a significant (p<0.05) effect on both BW and WW where the heavier BW and WW buffaloes calves born in winter compared to those born in other seasons confirming the significant effect of season on buffalo calves’ growth performance reported in India [4] and Egypt [14]. Such positive effect of winter calving on BW and WW may be due to the appropriate climatic conditions especially temperature, and feeding regimes. Furthermore,, the effect of year of calving on BW had significant effect as in the results stated in other studies [4,14]. The effect of year of calving on BW reflects year to year variability with respect to feeding and management practices and climate changes. Moreover, in the current study, the effect of herd was significant (p< 0.0001) on BW but it did show significant effect on WW. The strong effect of herd may be attributed to differences in management practices among herds during the dam gestation period. Also, Gestation length had a significant effect on BW (p<0.05) but it did not have effect on WW. The calves having gestation length between 260 to 310 days had the heaviest BW, followed by those having gestation length longer than 325 days whereas the calves having gestation length below 296 days were the lightest BW. These results were in agreement with what reported in other studies on beef cattle by Jamrozik and Miller [15], where they found a positive genetic correlation between BW and GL. The current study found a significant (p<0.0001) effect of parity on BW whereas there was no effect for parity on WW. Calf BW was lower for primiparous buffaloes than those of multiparous. Furthermore, BW of calves increased with parity number. This positive relationship has reported in other studies [3,13]. The parity number and age at parturition are reproductive traits which are associated with the physiological status of the female buffaloes and positively effect on calves’ BWs. The effect of parity on BW may be due to the maturity status of the dams buffaloes as in the late parities, they have high body capacity compatible with better development of fetus [16]. In the current study we also found that, the WD had the same trend for BW, confirming the results reported in previous studies [17].

The direct heritability estimates were 0.06 and 0.41 for BW and WW, respectively whereas direct maternal heritability values were 0.03 and 0.14, respectively (Table 3). These heritability values are in the range reported for direct heritability of BW (0.05 to 0.188) and maternal heritability (0.03 to 0.349) on other buffalo populations [4,18,19]. The low direct and maternal heritability values for BW due to maternal environmental variance being more than 7 folds of direct heritability, can be an indication that genetic selection for BW could be effective, when improving dam feeding and management practices during last phase of pregnancy might play an essential role to increase BW of calves [20]. In Brazilian buffaloes calves, Malhado et al [19] reported low heritability estimates 0.09 and 0.03 for direct and maternal heritabilities, respectively. Also, in Colombian buffalo calves, Bolívar et al [21] found low direct and maternal heritability estimates of 0.05 and 0.01 respectively. On the other hand, Pandya et al [4] and Gupta et al [18] obtained heritability values of 0.188 and 0.349, respectively, for BW of Surti Murrah buffalo. Although, the heritability estimates are population specific, the high estimates reported in Surti buffalo were probably attributed to the use of animal model without accounting for maternal effects which may overestimate the direct heritability estimates. In this study, we found a negative correlation between direct and maternal (r_am) effects for BW (−0.16) which was reported by other studies [22]. This negative genetic correlation could be due to another negative correlation between dam and calf resulted from an adverse effect of high nutrition during early growth of calf [22].

In the current study, high (0.41) and moderate (0.14) direct and maternal heritability estimates were obtained for WW, respectively indicating that the genetic selection to improve WW is possible with good response of maternal ability [23]. Such lower estimate of maternal heritability than direct heritability indicated that the direct genetic effect has large influence on WWs in buffaloes [21]. The current heritability estimates for WW are in agreement with estimated reported in Brazilian buffalo. Nonetheless, the current estimates were higher than direct heritability estimate reported by Ashmawy and El-Bramony [24] (0.19) in Egyptian buffalo, and Pandya et al [4] (0.17) in Surti buffalo. Also, in Colombian buffalo, found direct and maternal heritability values of 0.16 and 0.04, respectively [8]. The variation in WW due to maternal environment was moderate but smaller than its value for BW, indicating that the permanent environmental influence of dam on WW has a carryover effect from birth to weaning [25]. The current results showed that the proportion of phenotypic variance due maternal permanent environmental variance (C²) was higher than h²_m for BW and WW, indicating that maternal environmental effect was dominant until weaning [23]. Estimates of C² and h²_m for different beef breeds were close to those obtained in the current study [26]. The direct-maternal genetic correlation (r_am) was negative for WW and higher than its value for BW. This negative r_am was also found in Australian beef cattle [22] and Spanish beef cattle [27]. The negative genetic correlation may be due to a negative correlation between dam and calf resulted from an adverse effect of the high nutrition during early growth of calf [22].

A low genetic correlation between BW and WW (0.05) found in the current study indicated that selection to improve WW in buffalo may not affect BW in a breeding objectives that include both traits. Such relationship is very important to avoid dystocia which may occur for calves with high BW. High estimates of genetic correlation between BW and WW of Egyptian buffalo calves were reported in other studies [18] using animal model which did not account for maternal effects. The estimates genetic correlation between growth traits in Beef cattle are usually moderate to high [3,25,27], but others were low and negative such as those by El-Saied et al [28]. The weak correlation between BW and WW may be a decrease in correlation response because of a decrease in selection intensity and population size [29]. Thus, the low genetic correlation in current study may be attributed to low selection intensity in Egyptian buffalo populations compared to beef cattle. In this study, the maternal genetic r_m and maternal permanent environmental r_c correlations were 0.262 and 0.003 respectively which were lower than estimated values 0.33 and 0.64 reported in Japanese beef cattle [25].

Remarkable fluctuations in breeding values of calves and dams, and negative direct and maternal genetic trend in some years were found (Figure 1). Averages of direct and maternal genetic trend were −0.003, 0.002, −0.020, and 0.001 kg/yr for BW and WW, respectively indicating that there was no specific breeding programs have been applied to improve growth traits in the studied buffalo herds. Therefore, there was no genetic improvement in either trait. Other studies showed improvement in BW and WW where the direct and maternal genetic trends for BW and WW of Brazilian buffalo were 0.006, −0.01, 0.23, and 0.05 kg/yr, respectively [19]. Nonetheless, a negative direct genetic trend was −0.03 kg/yr was reported for water buffalo in Brazil [30].

The total variance explained by PCs was 71.17% in which 45.91% and 25.25% were explained by PC1 and PC2, respectively (Table 4). These results indicated that the first two PCs out of 4 showed eigenvalues greater than 1. Therefore, the first two components were sufficient to explain the most of variation in breeding value estimates for early growth traits in the current buffalo calves herds. The variation in breeding values for growth traits have been studied by Boligon et al [31] where they found the first three PCs explained 93.62% of total additive genetic variation of growth traits in Nellore cattle. In Colombian buffalo, three PCs greater explained 65.87% of the variance of breeding value for growth traits [8]. The two-dimensional graph of PC1 vs PC2 (Figure 2) and the correlation between BV, MV of both traits and PC (Table 5) showed that the direct breeding values of BW were correlated with PC2 and those for direct and maternal breeding values of BW and WW were correlated with PC1.

PC1 explained most of direct and maternal breeding value variation and was correlated with the direct and maternal breeding values of WW and maternal breeding values of BW. Thus, it could be considered as genetic index to WW and also as maternal index for BW and WW. Only 25.25% of direct and maternal breeding value variation was explained by PC2 which was correlated with direct breeding values of BW and considered as genetic index related BW. Therefore, selection to improve WW could be applied using PC1 separate from selection to improve BW in Egyptian buffalo. These results supported the reported genetic correlation between BW and WW in Table 3. The two-dimensional graph of PC1 vs PC2 illustrated in Figure 2, showed a negative correlation between maternal and genetic breeding value of BW and WW supported the direct-maternal genetic correlations reported in Table 3. A study, on Canchim cattle, reported similar trends where the first two principal components explained 73.37% of total additive genetic variance, and the PC1 was considered as genetic index for reproduction traits as well as PC2 was a genetic index correlated with body weight [12].

The PC scores for each animal in the first two principal components were calculated as:

The weights of these indices were SSC for each SEBV. The larger the absolute value of the SSC, the greater the relative importance. The 2 PC scores allowing capture the main information in predicted breeding value for BW and WW. Animals with higher PCS1 could be used to improve WW, and animals with higher PCS2 could be used to improve BWs. Thus, selection with previous index could be used in genetic selection programs of Egyptian buffalo. As previously mentioned regarding the use of PCA is very practical technique to reduce the data dimensionality especially with limited number of observations in the dependent variables under the lacking of national recording system in Egypt. Also, under Egyptian breeding programs, it has been difficult to develop the economic indices and if developed they were rarely used because of the difficulty in economic weights for each trait as the Egyptian market is not stable. Thus, we can consider the PCSs as new traits which can be included in new breeding objectives and used as selection criteria under Egypt condition. Meanwhile, in other studies the PCSs can be used to reduces the genetic variance components accounting for selection bias and correlated response in correlated traits [32].