India continues to be the largest producer of milk in world. Milk production during 2014-15 and 2015-16 was 146.3 million tonnes and 155.5 million tonnes, respectively, showing an annual growth of 6.27% (Anonymous, 2017). Though being the largest milk producing country, the average milk production per cow or buffalo is very low in comparison with advanced dairy countries due to several reasons including low genetic potential of indigenous cows for milk production (Hassan and Khan, 2013). The genetic improvement through grading up and selective breeding including progeny testing project for the indigenous animals is the key to increase the milk production and meet the rising demand for milk and milk products in India. Main objective of any breed improvement projects is the precise technical selection of superior sires and dams at the earliest to bring about faster genetic improvement.

Dairy cattle have traditionally been evaluated on the basis of 305-day lactation yield. Nowadays, as per ICAR (International Committee for Animal Recording) guidelines, dairy cattle are milk recorded at 4-week interval with total minimum 11 test-day records. 305-day lactation yield is calculated using the test interval method as per ICAR guidelines (ICAR, 2016). Random regression model (RRM) is one of the test-day models used widely for genetic analysis of production traits in dairy cattle. Use of test-day records fitting RRM has several benefits including ability to account for environmental factors that affect cows at different stages of lactation at the time of test and account for individual differences in the shape of lactation curves, which includes the persistency of the lactation (Meseret et al., 2015). Field performance recording (FPR) programme is an excellent mean for getting breeding bulls from smallholder dairy farmers under field based progeny testing(PT)project. National Dairy Development Board (NDDB), through Sabarmati Ashram Gaushala (SAG) as the End Implementing Agency (EIA), is carrying out progeny testing project on CBHF cattle under the National Dairy Plan Phase-I, a central sector scheme of the Government of India.  SAG CBHF PT project conducts progeny testing on farmers’ herds in selected villages of Sabarkantha, Panchmahal and Surat districts of Gujarat state. FPR by engaging milk recorder involves huge cost. This may remain a limiting factor in implementing large scale PT project for some agencies. One convenient way of milk recording may be asking farmer to maintain records and using those records for genetic evaluation as it is financially and management point of view feasible. However, the major concern would be the reliability of the records provided by farmers for genetic evaluation. Another strategy could be milk recording on bimonthly recording instead of monthly recording. However, deviations in sire rankings for their genetic merit need to be studied.

So, the present study was conducted with the objectives to investigate the feasibility of different milk recording strategies for genetic evaluation of sires under PT project (i.e. farmer milk recording and bimonthly milk recording), to study probable pattern of error a farmer may make while reporting performance of an animal (test-day milk yield) and to assess the reliability of data provided by farmers for sire evaluation.

Material and Methods

A total of 1,66,426 test-day milk yield (TDMY) records of 18,801CBHFcows recorded during 1993 to 2016 in selected villages of three districts of Gujarat state viz. Sabarkantha, Panchmahal and Surat were retrieved from Information Network for Animal Productivity & Health – Management Information System (INAPH-MIS) developed by NDDB. The historical data of SAG PT project on CBHF cattle collected from INAPH-MIS based on existing milk recording strategy were kept as the base data set. In this milk recording strategy, a designated milk recorder goes to farmer’s doorstep and records the milk yield of animal in morning as well as evening. In case of high milk producing animal, three times milk recording is carried out. Next milk recording is carried out at 4-week interval. In one complete lactation minimum 11 test-day records will be recorded (ICAR, 2016). Animals without sire information were eliminated from the study. Only first lactation records were considered in the study. Therefore, out of 1,66,426 total records of dams and daughters, 63,098 first lactation records of 7,419 daughters sired by 209 sires were used for the study. For better data consistency and technically precise genetic evaluation, the animals with first test-day (minimum days in milk, DIM = 5) to 11^th test-day (maximum days in milk, DIM= 330), 24 to 72 months age at first calving (AFC) and 1 to 40 kg test-day milk yield were retained and others were excluded. In Indian condition, herd means entire village as few numbers of animals kept by a farmer. In present model, the HYMR (Herd * Year of milk recording * Month of milk recording) was kept as random effect, HA (Herd * Age at first calving) and YS (Year of calving * Season of calving) were kept as fixed effects, DIM as fixed regression, animal and permanent environment were kept as random regressions. The minimum observations for animal, HA, YS were 3, 5 and 100, respectively. Finally, 55,703 test-day records of 6112 daughters sired by 192 sires were used for further analysis.

Three data sets were generated from the base data set based on different milk recording strategy: 1) data on existing standard milk recording strategy (i.e. morning – evening milk recording at monthly interval- ST₁); 2) data on bimonthly milk recording strategy (i.e., morning – evening milk recording at bimonthly interval- ST₂); 3)data on milk recording by farmer (i.e., morning – evening milk recording at monthly interval by farmer- ST₃). Data set ST₁ is processed base data set (base data set with above restrictions). For ST₂, out of total 11 test-day records, only 1^st, 3^rd, 5^th, 7^th, 9^th and 11^th records were retrieved from base data set. However, to generate data set for ST₃, information was collected by personally interviewing 300 owners/ milkers of animals from 82 villages under SAG CBHF PT projects in three districts of Gujarat, viz. Sabarkantha, Panchmahal and Surat with the purpose to estimate the magnitude of error a farmer may make while reporting performance of animal. Several factors were examined through a structured questionnaire viz. age and gender of milker /farmer, total number of animals with farmer, total number of calves with farmer, total number of heifers with farmer, percentage of lactating animals out of total animals with farmer, total sale of milk by farmer, possibility of milk recording (MR) at regular interval by farmer, possibility of MR at fixed interval by farmer, possibility of MR twice a day at monthly interval by farmer, possibility of maintenance of milk records by farmer and farmer allowing calf to suckle on the day of MR or not. Out of these factors, total number of animals with farmer, possibility of milk recording at monthly interval by farmer and possibility of milk recording twice a day at monthly interval by farmer were found to have significant effect  (P < 0.05) on error pattern by least squares analysis wherein the above mentioned factors were considered as factors affecting error pattern, which was considered as the trait (error pattern indicates the difference between milk yield and actual milk yield recorded by PT project milk recorder ( the latter being made available from INAPH-MIS)).

The simple arithmetic mean and standard deviation (S.D.) of significant factors affecting error pattern are presented in Table 1. Based on mean and S.D. of three factors, random numbers were generated and those numbers were simulated to ST₁ for generating data set ST₃.

Table 1: Arithmetic mean and standard deviation of factors having significant effect on error pattern

where, A₁ and A₂ are factors for farmer agreed and disagreed for milk recording at regular interval respectively; C₁ and C₂ are factors for farmer agreed and disagreed for milk recording twice a day on monthly basis, respectively; T₁ to T₅ are the groups of farmers having total number of animals 2 to 4, 5 to 6,7 to 9, 10 to 20 and >20, respectively.

Table 2 shows number of records and trait values under three strategies.

Table 2: Information on number of TDMY records, number of daughters, number of sires, minimum TDMY (in kg), maximum TDMY (in kg) and average of TDMY (in kg) in data sets under three different milk recording strategies

Model

Following Random regression model (RRM) was used for analysis of test-day records of dairy cattle:

where,?_?ℎ????is the test-day milk yield of cow k made on day t; HYMR = herd x year of milk recording x month of milk recording as random effect with subclass i; HA = herd x age at first calving as fixed effect with subclass h; YS = year of calving x season of calving as fixed effect with subclass j; ?_?= fixed regression coefficients; ?_?? and ??_??= the l^th random regression for animal additive genetic and permanent environmental effects, respectively, for animal k; ∅_??? = the l^th Legendre polynomial for the test-day record of cow k made on t^th day in milk; nf= the order of polynomials fitted as fixed regressions(Legendre polynomial of order three is used for fixed effects), nr = the order of polynomials for u and pe effects (Legendre polynomial of order three is used for random effects), ?_?ℎ???? = the random residual effect.

The assumptions were:

Var(u) = A ⊗ G; var(pe) = I ⊗ P; var(e) = I ?_?² = R,

where, A is the numerator relationship matrix; ?_?² is error variance; I is identity matrix; ⊗is the Kronecker product; G and P are of the order of polynomial fitted for u and pe effects, respectively.

Genetic analysis was done using WOMBAT software (Meyer, 2011) with AIREML (Average Information Restricted Maximum Likelihood) procedure to estimate variance and covariance and BLUP (Best Linear Unbiased Prediction) was used to estimate breeding values from variance and covariance components estimated using AIREML. Comparisons among different strategies were made using Log likelihood function (logL) (Mrode and Thompson, 2005), Akaike’s information criterion (AIC) (Akaike,1973), Bayesian information criterion (BIC) (Schwarz, 1978), and Percentage Squared Bias (PSB) (Ali and Schaeffer, 1987). These criteria represented with the following equations:

where, b = generalized least square solution

AIC = -2 log L+ 2k and BIC = – 2log L + k log (λ)

where, log L= log-likelihood function, k = number of parameters estimated in model,

λ = n – r(X) where, n = number of test-day milk yield records, r(X) = rank of systematic effects incidence matrix

where, n = number of records, y_hijkt = the observed record, ŷ_hijkt = the record predicted by the particular model

The Spearman’s rank correlation between breeding values of sires derived by various strategies to judge the effectiveness of different strategies was estimated as per Steel and Torrie (1960)-

where, r_s = Rank correlation coefficient, n = Number of sires under evaluation, d_i = Difference of rank between paired items under two models.

The significance of rank correlation was tested by t-test (Kendella and Stuart, 1973) as given below:

It was compared with t-table value with  degree of freedom.

Results and Discussion

Test-day is the day of month on which production performance of animal is recorded. Test-day models are the statistical procedures that consider all genetic and environmental effects directly on a test-day basis.  It improves the accuracy of genetic evaluation, provides better modeling and extending of part lactation is no more needed. Breeding values for 305-day milk yield of total 440 sires (including sire of sire) under three milk recording strategies ST₁, ST₂ and ST₃ were estimated using RRM. The maximum range of breeding value of sires was observed for ST₃ (-559.36 kg to 1478.32 kg) followed by ST₁ (-587.14 kg to 1264.90 kg) and ST₂ (-571.29 kg to 1213.90 kg). Under all three strategies, the average breeding values of sires are positive i.e. higher than average breeding value of the population. In the study, records for population are available on dams and daughters. Thus, on an average, sires are superior to female population with farmer. This is expected in a well-structured PT project as sires are the sons of superior dams.

Since higher values for test-day milk yield is desirable, the sires with positive breeding value were considered as superior and those with negative breeding value were considered to be inferior. The average breeding values of sires under ST₃ were almost double the estimates under ST₁ and ST₂ (Table 3). This result may indicate that farmers have tendency to record higher milk yield of their animals i.e. to inflate the TDMY on higher side.

Table 3:  Breeding value (BV) (in kg) of sires for 305-day milk yield recorded under different milk recording strategies

Reliability (R%) of breeding value of sires is an important criterion indicating the accuracy of prediction of breeding value. It is influenced by the number of daughters per bull as well as variance components of trait under consideration. The maximum range of R% of breeding value was observed for ST₃ (10.00% to 93.20%) followed by ST₂ (11.26% to 91.20%) and ST₁ (12.22% to 91.89%). The average R% of breeding values was found to be highest for ST₁ followed by ST₃ and ST₂ (Table 4). Higher R% for ST₃ as compared to ST₂ does not necessarily mean superiority of ST₃ strategy over ST₂, a further investigation on PSB values for ST₃ and ST₂ revealed clearer picture, as will be discussed in the following. Most of the sires (around 70 %) were observed to have less than 50% reliability under all the strategies. Very less number of sires had greater than 90% reliability (Table5).

It should be noted that though the h²value of the trait was the lowest for ST₃, the number of sires with >70% reliability were observed to be highest for ST₃ which suggested that ST₃ provided poor partitioning of variance, with an indication that part of error variance was being accumulated in to permanent environmental variance, based on the fact that genetic variance exploited was already low. Recent publication on breeding value of sires for 305-DMY for SAG CBHF PT project have shown that the total number of sires evaluated under CBHF PT project were 112 with reliability greater than or equal to 70% with their breeding values ranging from -554 kg to 1204 kg, reliability ranging from 70% to 91%, and number of daughters in analysis per sire ranging from 21 to 143 (NDDB, 2016).

Table 4: Reliability (R%) of breeding values of sires under different milk recording strategies

Table 5: Number of sires with different classes of reliability (R %) of breeding values under different milk recording strategies

The milk recording strategies under the study were further compared based on average estimates of heritability (h²) of test-day milk yield, Log likelihood (logL), Akaike information criterion (AIC), Bayesian information criterion (BIC) and Percentage squared bias (PSB) values to further get an indication of optimum strategy (Table 6).

Table 6: Average heritability of test-day milk yield and comparison criteria for different milk recording strategies

The average estimates of h² for test-day milk yield were almost same for ST₁ and ST₂ and almost double to that for ST₃. This is an indication that ST3 could not exploit as much genetic variation as ST₁ and ST₂. It should also be noted that ST₂ showed slight tendency to decrease heritability as against ST₁, though the difference was negligible. LogL denotes the goodness of fit for models of same data size. Lower the value of logL, better is the model. In the present study, data size for ST₁ and ST₃ was similar while for ST₂, it was lower; as ST₂ represents bimonthly recording. Better logL value revealed by ST₂ compared to ST₁ and ST₃, could be due to reduced data size under ST₂ and hence ST₂ cannot be inferred superior to ST₁ or ST₃ based on logL. ST₁ and ST₃ had almost equal data size but logL value for ST₃ was almost double than that of ST₁.

AIC is a measure of the relative quality of statistical models for a given set of data. BIC is closely related to the AIC and it is based on likelihood function. The highest BIC specifies the best fitting model. Though the same model was used for all strategies, it should be noted that AIC and BIC were compared to check chiefly the goodness of fit which indirectly is an indication of quality of data. For AIC and BIC, the results were similar to that observed for logL.

PSB is an indication of bias in estimation of trait using the model and data specified for each strategy under the study. Lower the value of PSB, better would be the overall estimation of effects and hence lower would be the prediction bias. PSB for ST₁ and ST₂ were same with minor variation due to change in data structure due to removal of alternate month recordings in ST₂. A further investigation on changes in ranks of sires based on breeding values in ST₂ as against ST₁ would clear the picture as to whether there are any variations in estimation of genetic worth of sires in ST₂ relative to ST₁. For ST₃, PSB was much higher (in order of more than ten times) than that for both ST₁ and ST₂. Higher value of PSB for ST₃ compared to ST₁ invariably indicates lower prediction accuracy by ST₃ than ST₁ and ST₂. Thus, ST₃ explained the least variation in dependent variable and provided maximum bias for estimating breeding values for test-day milk yield. Thus, based on the values of comparison criteria such as logL, AIC, BIC and PSB, it was observed that ST₁ and ST₂ provided better accuracy compared to ST₃. Spearman rank correlation is needed to estimate the correlation between ranks of the individuals. Rank correlation coefficients of rankings of sires based on 305-day milk yield breeding values under different strategies are given in Table 7.

Table 7: Rank correlations of sires’ breeding values for 305-day milk yield estimated under different milk recording strategies

These estimates of rank correlations were highly significant (P < 0.001). Rank correlation of ranking of sires by ST₃ with ST₁ and ST₂ are poor than that between ST₁ and ST₂. It should also be noted that though PSB values for ST₁ and ST₂ were almost same, there are changes in rankings of sires in the two strategies. Thus, alternate strategies (ST₂ and ST₃) had poor correlation with the existing strategy (ST₁) which consisted maximum possible information on the population under study with respect to ranking of sires based on their 305-day milk yield breeding values.

NDDB currently publishes breeding values of progeny tested sires with reliability of their breeding values greater than or equal to 70% for cattle. On similar lines, in the present study, rank correlations were also estimated for sires based on their breeding values with greater than or equal to 70 per cent reliability (Table 8). The rank correlations of ranking of sires by ST₃ with ST₁ and ST₂ were poor (0.789 and 0.799, respectively) and lower than that between ST₁ and ST₂ (0.958). These estimates of rank correlation were highly significant (P < 0.001). These results indicated that around 21% and31% changes were observed in ranking of sires based on 305-day milk yield breeding values by ST₂ and ST₃, respectively, as compared to ST₁ for all the sires irrespective of reliability. Around 4% and 21% changes were observed in ranking of sires by ST₂ and ST₃, respectively, as compared to that by ST₁ for sires with R% greater than or equal to 70%. These results indicated higher percentage of changes in sire rankings when ST₃ was followed. This may have a great consequence leading to undesirable reduced rates of genetic gains under large-scale field-based animal breeding programmes.

Table 8: Rank correlations of sires’ breeding values for 305-day milk yield with ≥70% reliability estimated under different milk recording strategies

Conclusion

Among the three milk recording strategies considered for sire evaluation, ST₁ strategy (existing standard test-interval milk recording strategy) is recommended to be the optimum strategy. ST₃ strategy suffers from lower exploitation of genetic variance, biased results and poor conformity of measurement of sire genetic merit relative to ST₁ and ST₂. ST₂ may be opted as an alternative only in case of constraints in deploying ST₁ strategy. Current study suggests changes to the tune of around 21% in ranking of sires by ST₂ relative to ST₁ and subsequent losses. ST₁ is thus the recommended strategy, followed by ST₂ which should be considered with caution, and ST₃ strategy is not recommended for large-scale field-based animal breeding programmes in India.

Acknowledgement

This research was carried out as collaborative work between the National Dairy Development Board (NDDB) and College of Veterinary Science and Animal Husbandry, AAU, Anand. The authors acknowledge the financial and technical support as well as the provision of SAG CBHF PT project related relevant data from NDDB for conducting the research. The authors also acknowledge the field support received from SAG CBHF PT project officials of Surat milk union, Panchmahal milk union, Sabarkantha milk union and Sabarmati Ashram Gaushala for carrying out the research.

References

1. Anonymous (2017). of Animal Husbandry, Dairying & Fisheries, Annual report 2016-17, Ministry of Agriculture & Farmers welfare, Government of India.
2. Kaike, H. (1973). Information theory and an extension of the maximum likelihood principle. In: Petrov, B.N., Csaki, F. (Eds.), Second international symposium on information theory. AkademiaiKaido´, Budapest, Reprinted in Kotz, S., Johnson, N.L., (Eds.), Breakthroughs in statistics, Vol. 1, Springer-Ver- lag, New York.
3. Ali, T.E. and Schaeffer, L.R. (1987). Accounting for covariances among test-day milk yield in dairy cows. J. Anim. Sci.,67: 637–645.
4. Gokhale, S. B. (2013). Rural experiences in implementing cattle and buffalo improvement projects. Paper submitted to International Workshop on “Genetic Improvement of Cattle and Buffalo under NDP-I” held at Anand on 18-19 September, 2013.
5. Hassan and Khan (2013). Performance of crossbred dairy cattle at military dairy farm in Pakistan. The Journal of Animal and Plant Sciences, 23(3):705-714.
6. ICAR, 2016. International Agreement of Recording Practices. International Committee for Animal Recording ed., Interlaken, Switzerland.
7. Kendall, M. G. and Stuart, A. (1973). The Advanced Theory of Statistics, Volume 2: Inference and Relationship. Griffin (Sections 31.19, 31.21).
8. Meyer, K. (2011). WOMBAT-A tool for mixed model analyses in quantitative genetics by restricted maximum likelihoo User notes. Animal Genetics and Breeding Unit, University of New England, Armidale, NSW, Australia.
9. Meseret, S., Tamir, B. and Negussie, E. (2015). Estimation of genetic parameters for first-lactation test-day milk yield in Holstein Friesian cows fitting random regression models. J. Genet. Mol. Biol., 7(4): 25-31.
10. NDDB (2016). National Dairy Development Board, Anand, Report. http://www.nddb.org/services/animalbreeding/genetiimprovement/sireproofs.
11. Schwarz, G. (1978). Estimating the dimension of a model. Stat.,6: 461–464.
12. Steel, R. G. D. and Torrie, J. H. (1960). Principles and procedures of statistics with special reference to the biological sciences. McGrew Hill Book Company Inc. New York, 550.