Detection of difficult calvings in the Polish Holstein-Friesian Black-and-White heifers

Abstract

The aim of this study was to check the effectiveness of classification and regression trees (CART), support vector machines (SVM) and generalized linear model (GLZ) in the detection of difficult calvings in the Polish Holstein-Friesian Black-and-White heifers and to identify the most influential predictors of such calvings with special regard to oestrogen receptor alpha (ERα) and cytochrome P450 (CYP19) aromatase genotypes. A total of 529 calving records containing four continuous and four categorical predictors were analysed. Calving category (difficult vs. easy) was a dependent variable. The most significant predictors were as follows: gestation length, body condition score index and calving age (for CART), as well as ERα genotype (for GLZ). Sensitivity (percentage of correctly diagnosed difficult calvings) and specificity (percentage of correctly diagnosed easy calvings) on a test set were 70.0% and 94.7% (for CART) as well as 50.0% and 99.0% (for SVM and GLZ), respectively. The studied models were characterized by good detection performance and high specificity and reliability of predictions.

Keywords:

• dystocia
• data mining
• classification and regression trees
• support vector machines
• cattle
• diagnosis

1. Introduction

The incidence of difficult calvings (dystocia) in cattle is highest from among all farm animal species (Purohit et al. 2011). Moreover, heifers are characterized by up to three times higher dystocia rate than multiparous cows (Meyer et al. 2001). The frequency of dystocia in dairy breeds (Holstein-Friesian, Norwegian Red, Swedish Red and White) from different countries of Europe, the USA, Australia and New Zealand was estimated at 3–23% (Mee 2008). This condition results in many adverse effects, such as increased mortality and culling rates, higher risk of calving difficulties at subsequent calving and susceptibility to diseases (metritis, retained placenta, milk fever), reduced fertility (longer calving-to-conception interval, increased number of artificial inseminations per conception and delayed first service), production and dry matter intake, body weight and body condition loss during lactation and disadvantageous metabolic changes (Mee 2004; Barrier 2012). Dystocia also considerably affects the health and welfare of calves (increased perinatal mortality, reduced immunity to diseases, higher level of physiological stress; Barrier et al. 2012; Barrier et al. 2013). All these factors are associated with considerable financial losses, estimated in the USA at approx. 28 USD per average heifer (Dematawewa & Berger 1997), although the total cost of particularly difficult calving may be much higher. There are many direct risk factors of dystocia, with the most important being mismatch between the size of a foetus and its dam (more than 50% of obstetrical cases), foetal malpresentation (20–40% of obstetrical cases) and maternal causes such as uterine torsion (less than 10% of cases) and hypocalcaemia (less than 17% of cases; Schuenemann et al. 2011).

For heifers, the most significant factors are foetopelvic incompatibility, abnormal foetal disposition and incomplete dilation of the vulva and cervix (Mee 2004). They are in turn affected by many indirect causes (calf birth weight, gestation length, parity, calf sex, breed, nutrition, climate, heifer body weight, body condition, calving age, etc.; Mee 2008), including genetic ones. Genes potentially influencing calving difficulty are, among others, oestrogen receptor alpha (ERα) gene, mediating the action of oestrogens on the maintenance of normal physiological processes associated with female reproduction and playing a key role in the pathology of female sex organs (Szreder et al. 2011), as well as cytochrome P450 aromatase (CYP19) gene, responsible for oestrogen synthesis (Kowalewska-Łuczak 2010).

That being so, it seems justified to make an attempt at earlier detection of difficult calvings in dairy heifers based on available information on potential risk factors. Such a possibility is given by, among others, classification and regression trees (CART; Breiman 1984) and support vector machines (SVM) known from the field of data mining as well as a more traditional statistical method, i.e. generalized linear model (GLZ). The first one consists in asking a series of simple questions. An answer given at each stage determines the next (if any) question, which ultimately leads to a network of questions resembling a tree structure (Nisbet et al. 2009). The second method (SVM) is based on the application of linear models to the creation of non-linear boundaries separating classes through the transformation of the input pattern space to a new space (by means of non-linear mapping), where the constructed model can represent non-linear decision boundaries in original space (Witten et al. 2011). The last aforementioned method (GLZ) allows for linking the values of a dependent variable following a specific distribution to a set of independent variables (both continuous and categorical) via an appropriate link function, whose parameters are estimated (using the maximum likelihood method) in such a way so that the values observed in a training data-set are most probable (Sokal & Rohlf 1995). The most recent applications of CART to classification tasks in animal farming include modelling of calving ease and stillbirths in cattle (Piwczyński et al. 2013), somatic cell count in cow milk (Piwczyński & Sitkowska 2012) and lamb mortality (Piwczyński et al. 2012), detection of cows with conception difficulties (Grzesiak et al. 2011), classification of lactation curves in cattle (Pietersma et al. 2002; Pietersma et al. 2003), whereas SVM were used for, among others, identification of cattle breeds (Rodero et al. 2012), classification of bull semen samples (Aguiar et al. 2012), mastitis (Miekley et al. 2013) and lameness (Alsaaod et al. 2012) detection in dairy cows, classification of cow behavioural patterns (Martiskainen et al. 2009), detection of heifers and steers, in which prohibited growth stimulators had been used (Nebbia et al. 2011) and the classification of sleep and lying periods in calves (Hokkanen et al. 2011). On the other hand, GLZ with a binomial distribution of dependent variable and the logit link function was used, among other things, to analyse risk factors of dystocia (Mee et al. 2011) and stillbirths (Gulliksen et al. 2009) in dairy cows.

Therefore, the aim of this study was to compare the quality of the CART, SVM and GLZ models used for dystocia detection in the Polish Holstein-Friesian Black-and-White heifers, verify detection performance of these models and indicate the most influential predictors of calving difficulty with special regard to the aforementioned genotypes.

2. Material and methods

A total of 529 calving records from the Polish Holstein-Friesian Black-and-White heifers collected between 2004 and 2009 were used in this study. The initial data-set comprised 606 cases and was subsequently reduced to 529 cases after deleting erroneous or incomplete records without information on genotypes, calving age and other predictors, which was impossible to impute. Heifers were managed in a free-stall system in barns with outside runs on a farm located in the West Pomeranian province and fed a total mixed ration (TMR). The mean milk yield in the analysed period (305 d of lactation) was 9952 kg. Only heifers without any serious illnesses before calving were included. Four continuous (Table 1) and four categorical (Table 2) explanatory variables were used in total. These included proportion of Holstein-Friesian genes in heifer genotype (in %) – %HF (X₁), gestation length (in days) – GEST (X₂), body condition score index – BCSI (X₃), which is a difference between the mean BCS before calving (BCS_pre) and BCS at calving (BCS_post), calving age (in months) – AGE (X₄); ERα genotype (GG or AG) determined using the BglI restriction enzyme – BGLI (X₅), ERα genotype (AA or AG) determined using the SnaBI restriction enzyme – SNABI (X₆), cytochrome P450 aromatase genotype (AB or AA) determined using the PvuII restriction enzyme – CYP19 (X₇) and calving season [autum–winter (AW) and spring–summer (SS)] – SEASON (X₈).

Table 1. Means and standard deviations (SD) of continuous predictors in the training (n = 397), and test (n = 132) data-sets and both data-sets combined (n = 529).

Variable	Mean/SD	Train	Test	Total
%HF (%)	Mean	88.44	90.08	88.85
	SD	7.75	6.43	7.47
GEST (days)	Mean	277.81	277.65	277.77
	SD	2.86	3.14	2.93
BCSpre (scores)	Mean	2.77	2.76	2.77
	SD	0.34	0.34	0.34
BCSpost (scores)	Mean	2.74	2.75	2.74
	SD	0.37	0.34	0.36
BCSI (scores)	Mean	0.02	0.04	0.03
	SD	0.16	0.19	0.17
AGE (months)	Mean	26.60	25.95	26.43
	SD	3.85	2.42	3.56

%HF, proportion of Holstein-Friesian genes; GEST, gestation length; BCSI, body condition score index, which is a difference between the mean BCS before calving (BCS_pre) and BCS at calving (BCS_post); AGE, heifer calving age.

Table 2. Distributions of categorical predictors and the dependent variable in the training (n = 397) and test (n = 132) data-sets and both data-sets combined (n = 529).

	Training set		Test set		Total
Variable category	n	%	n	%	n	%
BGLI
GG	363	91.44	122	92.42	485	91.68
AG	34	8.56	10	7.58	44	8.32
SNABI
AA	370	93.20	119	90.15	489	92.44
AG	27	6.80	13	9.85	40	7.56
CYP19
AB	46	11.59	19	14.39	65	12.29
AA	351	88.41	113	85.61	464	87.71
SEASON
SS	269	67.76	84	63.64	353	66.73
AW	128	32.24	48	36.36	176	33.27
DYST – dependent variable
B	350	88.16	112	84.85	462	87.33
A	47	11.84	20	15.15	67	12.67

BGLI, ERα genotype (GG or AG) determined with PCR-RFLP using the BglI restriction enzyme; SNABI, ERα genotype (AA or AG) determined with PCR-RFLP using the SnaBI restriction enzyme; CYP19, cytochrome P450 aromatase genotype (AB or AA) determined with PCR-RFLP using the PvuII restriction enzyme; SEASON, calving season (autumn-winter – AW and spring-summer – SS); DYST, calving difficulty category (A – difficult or B – easy).

BCS was assessed on a 5-point scale (Ferguson et al. 1994; Flamenbaum et al. 1995) and the optimum was set at 3.50 scores with higher values deducted from 3.50, i.e. 3.75 was converted to 3.25, etc. This was necessary to appropriately reflect the optimum BCS in the modelling process. All the analysed polymorphisms were the A/G transitions. The first A/G transition recognized by the BglI nuclease was located at position 151 within the 5′-region of ERα gene (GenBank Acc. No. AY340597; Szreder and Zwierzchowski 2004). The second one recognized by the SnaBI nuclease was located at position –1213 relative to +85 splicing acceptor site in exon 1 of ERα (GenBank Acc. No. AY332655; Szreder et al. 2007). The last polymorphism recognized by the PvuII restriction enzyme was located at position 1044 within the promoter 1.1 of the CYP19 gene (EMBL Acc. No. Z69241; Vanselow et al. 1999).

Heifer calving class (A – difficult, B – easy) was a target (dependent) variable (Zaborski and Grzesiak 2011). Initially, calving difficulty was assessed by one experienced zootechnician on the following rank scale:

1. An easy, spontaneous calving;

2. An easy calving with little help from man;

3. A difficult calving, with substantial help from man and veterinarian's assistance; and

4. A difficult calving requiring veterinarians' intervention with potential complications, such as placental abruption, etc.

The rank scale was subsequently converted into nominal scale with two categories: A – difficult (ranks 3 and 4) and B – easy (ranks 1 and 2).

The whole data-set of information records was randomly divided into two subsets: a training set (397 records), used for training the CART and SVM models and estimating parameters of GLZ, and a test set (132 records), used for evaluating the detection performance of the models on new data, not utilized previously during model building. The division was made according to the 3:1 rule (75% records for training and 25% for testing) proposed by Ivakhnenko and Yurachkovskii (1987). In this way, it was possible to assess the effectiveness of the models in their practical application. Statistica® 10 PL software (StatSoft Inc., Tulsa, OK, USA) was used for model construction and evaluation. In each analysis, misclassification costs for both classes were equal and the a-priori probabilities of class membership were calculated from the proportions of records in the training data-set.

When building the CART models, three different splitting criteria were used (Gini index, chi-square and G-square measures; Hastie et al. 2009; StatSoft 2013). The tree pruning was based on a minimum cost and complexity of 10-fold cross-validation, in which the training data-set was randomly divided into approximately 10 equal subsets, from which nine were used for model training and one for assessing its generalization ability. The whole process was repeated 10 times, so that each subset served as a validation subset once and the error on this subset was averaged (Haykin 2009).

In the SVM model construction, its C-SVM type was used, in which the following error function (F_E) was minimized during the training process (Hill & Lewicki 2006):

under the following conditions:

where: F_E – error function, C – capacity constant determining an error reduction rate, w – vector of coefficients, b – constant, ξ_i – parameter determining the treatment of inseparable records, N – number of training records, y – class label (1 or –1), x_i – independent variable, φ – kernel used for the transformation of data from the input pattern space to the feature space.

In the present study, the C constant was determined based on a 10-fold cross-validation, allowing us for finding its optimal value, which was equal to 10 (StatSoft 2013), whereas the selected kernel type was a linear one (Hsu et al. 2003):

where: x_i, x_j – vectors.

Before SVM training, all predictors were scaled to the interval between 0 and 1, in order to reduce the number of support vectors and the risk of overtraining.

In the case of GLZ, the binomial distribution of the dependent variable and the logit link function were used (Koronacki & Ćwik 2005):

where: p – conditional probability of assuming the value of 1 by the dichotomous-dependent variable for the given predictor values.

After evaluating all effects for the model built based on the training data-set, the method of the best subset was applied to the selection of predictors. As a result, the following variables were retained in the model: GEST, BCSI, AGE and SNABI. The normal distribution of residuals was also verified.

After the stage of model building and training, their quality on the training data-set was assessed using the following probabilities: sensitivity (proportion of correctly diagnosed difficult calvings), specificity (proportion of correctly diagnosed easy calvings) and accuracy (proportion of correctly diagnosed difficult and easy calvings). In addition, statistically significant differences were determined using the test for proportions. Akaike information criterion (AIC) was also calculated in order to compare model quality (CART, SVM and GLZ) and its lower value indicated a better model.

The next stage of the present study was an identification of the most influential predictors (those with the greatest effect on the dependent variable class) for the CART and GLZ methods, based on the so-called importance coefficient (StatSoft 2013) and the Wald test statistic (Koronacki & Ćwik 2005), respectively. In the case of SVM, this identification was impossible using the default settings of Statistica® 10 PL software.

After evaluating the quality of the CART, SVM and GLZ models and identifying significant predictors, their ability to detect heifers with calving difficulties was verified based on the previously created test set, whose records had not been used for model building. In this way, it was possible to estimate approximately model performance during its potential use in breeding practice. The following measures were calculated: sensitivity (Se), specificity (Sp), accuracy (Acc), false alarm rate [P(FP); proportion of heifers calving easily that were indicated as those with difficult calvings by the model], false negative rate [P(FN); proportion of heifers with difficult calvings erroneously indicated as easy calvers] as well as the a-posteriori probabilities of true positive [P(PSTP)] and true negative [P(PSTN)] responses, showing the reliability of predictions made by the models. Statistically significant differences between probabilities were tested using the test for proportions. Moreover, the receiver operating characteristic (ROC) curves were plotted and the area under the curves (AUC) was calculated. The ROC curve may be used for classifiers with a dichotomous-dependent variable. It shows the relationship between sensitivity and specificity for all possible cut-off values. An ideal classifier is characterized by the ROC curve with coordinates (0;0), (0;1) and (1;1) (AUC = 1.0), whereas one whose performance does not differ from a completely random classification has a diagonal ROC curve (AUC = 0.5; Witten et al. 2011).

3. Results

The layout of CART obtained in this study for the classification of calvings in dairy heifers is shown in Figure 1. As can be seen, the predictor allowing for the best division of calving records in the tree root (ID = 1) was gestation length and the value of 279.5 days was selected as a threshold. On this basis, 54 animals (node ID = 3) were classified as dystotic ones (with a confidence of 66.7%) and 343 as easy calvers (with a confidence of 96.8%, ID = 2). In order to further divide the records in the second node (ID = 2), the CART algorithm chose BCSI, assuming 0.125 points as a threshold value. Animals whose BCSI was above this value were assigned to the ‘difficult calving’ category (with a confidence of 83.3%; ID = 5), whereas those with BCSI lower than or equal to the threshold were classified to the ‘easy calving’ category (with a confidence of 98.2%; ID = 4).

Figure 1. Schematic diagram of CART for dystocia prediction in dairy heifers.

All node impurity measures (Gini index, chi-square and G-square statistics) resulted in the same tree structure and identical classification outcomes. The results of this classification (on the training set) are presented in Table 3.

Table 3. Sensitivity (Se), specificity (Sp), accuracy (Acc) and akaike information criterion (AIC) for individual models based on the training set (n = 397).

Model	Se (%)	Sp (%)	Acc (%)	AIC
CART	87.23^a	94.57^a	93.70	19.97
SVM	70.21^b	97.71^b	94.46	24.13
GLZ	72.34^b	98.00^b	94.96	15.11

Note: Values in the same column without a common superscript letter are significantly different (P < 0.05).

In the case of SVM, the final model included 67 support vectors, including 32 for the ‘difficult calving’ class and 35 for the ‘easy calving’ class. The decision constant was 4.48. The estimated parameters for GLZ are shown in Table 4. In the verification of the GLZ model, a significant deviation of residuals from the normal distribution was observed.

Table 4. The GLZ model parameters.

Parameter	Estimate	Standard error	Wald test statistic	P value
Intercept	−193.7755	25.6389	57.1214	0.0000
GEST	0.7304	0.0947	59.4423	0.0000
BCSI	24.0010	5.6539	18.0206	0.0000
AGE	−0.4544	0.1632	7.7567	0.0054
SNABI	−1.0733	0.3744	8.2158	0.0042

GEST, gestation length; BCSI, body condition score index; AGE, calving age; SNABI, ERα genotype determined using the SnaBI restriction enzyme.

Classification results for CART, SVM and GLZ on the training set are given in Table 3. Statistically significant differences in sensitivity and specificity were noticed (P < 0.05) with the CART model characterized by the highest sensitivity and the GLZ model by the highest specificity. However, the AIC value depending also on the number of model parameters favoured GLZ. The most significant predictors for CART were as follows: GEST, BCSI and AGE (Figure 2), whereas those for GLZ included: GEST, BCSI, SNABI and AGE (Table 4). Probabilities for individual models calculated based on the test set are presented in Table 5.

Figure 2. Predictor importance for CART.

Table 5. Probabilities (in %) for individual models calculated based on the test set (n = 132).

Model	Se	Sp	Acc	P(FP)	P(FN)	P(PSTP)	P(PSTN)
CART	70.00^a	94.64^a	90.91	5.36^a	30.00^a	70.00^a	94.64
SVM	50.00^b	99.11^b	91.67	0.89^b	50.00^b	90.91^b	91.74
GLZ	50.00^b	99.11^b	91.67	0.89^b	50.00^b	90.91^b	91.74

Note: Values in the same column without a common superscript letter are significantly different (P < 0.05).

CART, classification and regression trees; SVM, support vector machines; GLZ, generalized linear model; Se, sensitivity; Sp, specificity; Acc, accuracy; P(FP), probability of false positive response; P(FN), probability of false negative response; P(PSTP), a-posteriori probability of true positive response; P(PSTN), a-posteriori probability of true negative response.

Statistically significant differences (P < 0.05) between CART and the two remaining models were found in sensitivity (higher for CART), specificity (lower for CART), probability of false positive results (higher in the case of CART) and false negative results (lower in the case of CART) as well as the a-posteriori probability of true positive response (lower for CART).

The ROC curves plotted based on the test set and their corresponding AUC values for CART, SVM and GLZ are shown in Figure 3. The highest AUC value was characteristic of GLZ. A slightly smaller AUC was found for CART, whereas SVM had the lowest AUC value.

Figure 3. The ROC curves and their corresponding AUC values (in parentheses) for individual models.

4. Discussion

In the case of CART model obtained in the present study, the first division of the whole data-set in the tree root node was based on gestation length. Somewhat different result was reported by Piwczyński et al. (2013), who applied CART and chi-squared automatic interaction detector (CHAID) trees to the analysis of significant predictors of calving difficulty in dairy cows (Polish Holstein-Friesian Black-and-White breed). In that study, only the fourth partition was based on pregnancy length and this variable was used for splitting twice (at the threshold values of 282 and 284 days, respectively, which were similar to those in the present work).

The use of different node impurity measures for the CART model construction did not significantly affect classification results, similarly as in the study by Piwczyński et al. (2013), where the misclassification rate was the same for Gini index and entropy function. In the case of GLZ, a significant deviation from the normal distribution of residuals was observed in the present study, which could negatively influence the value of test statistic (Stanisz 2007).

The values of sensitivity, specificity and accuracy obtained on the training set and showing the model quality can be regarded as high (Table 3). Only in the case of the first measure, there was quite a large difference in favour of CART, which was also confirmed statistically. In the work by Piwczyński et al. (2013), the overall accuracy of dystocia detection in cows was 61.5% for both CART and CHAID, and so it was much lower than that in the present study. On the other hand, Zaborski and Grzesiak (2011), in the classification of difficult calvings in dairy heifers using artificial neural networks, obtained higher sensitivity (89.1–90.9%), especially in comparison with SVM and GLZ in the present study, at lower specificity (88.3–89.6%) and accuracy (88.6–89.5%). The method of classification functions used for the same purpose (Zaborski 2010) gave much worse results in terms of sensitivity (60%) with comparable specificity (98%) and identical accuracy (93%). The application of discriminant function analysis to the prediction of difficult calvings in beef heifers was investigated by Basarab et al. (1993), who obtained sensitivity in the range of 26–40% depending on heifer breed and the moment of prediction, which was much lower than that in the present study. However, specificity was similar (97%) and accuracy was also lower (85–89%) compared with those in the present work. In another study on dystocia prediction in beef heifers and cows using discriminant function analysis, Morrison et al. (1985a,b) recorded similar sensitivity (78–87%) with lower specificity (81–86%) and accuracy (82–85%). In other applications of CART to animal farming, Piwczyński et al. (2013) obtained the accuracy on a training set of 94% when investigating stillbirths in dairy cows.

The next stage of the present study was the identification of variables with the greatest contribution to the determination of calving difficulty class. It should be emphasized that in the present study, only predictors available before calving were included. Therefore, factors with a documented effect on dystocia such as calf sex, its birth weight, the ratio of cow body weight to calf birth weight, etc., were not investigated. At present, despite quite commonly used ultrasonography for pregnancy diagnosis in dairy cattle, the sex of calves was not determined using this method, so this variable could not have been used either. The most significant predictor for CART and GLZ was pregnancy length (Figure 2; Table 4). This result is consistent with those of sensitivity analysis of the artificial neural networks applied to dystocia detection in heifers (Zaborski & Grzesiak 2011), and with the order of independent variables in discriminant function analysis used for the same purpose (Zaborski 2010). Also, in the study by Piwczyński et al. (2013), in which the CART and CHAID trees were utilized, gestation length was an important predictor of calving category (third to fourth position according to the importance criterion applied). In general, both too short and too long gestation result in more frequent calving difficulties. Fiedlerová et al. (2008) observed a significant non-linear relationship between gestation length in the Czech Holstein cows and calving difficulty, assuming the value of 275.8 days as an optimum number ensuring the easiest labour, which was very close to the threshold found in the present study. More difficult calvings in animals with shorter gestation, may result, among other things, from a higher calf perinatal mortality rate associated with an increased risk of labour complications (Johanson & Berger 2003), whereas in heifers with longer pregnancy duration, they may be a consequence of excessive calf birth weight (Mee 2008). According to the latest estimates (Nogalski & Piwczynski 2012), an optimal gestation length in the population of Polish Holstein-Friesian Black-and-White cattle in terms of calving ease and stillbirth risk was 275–277 days, whereas this value (in reference to calving ease) ranged between 274 and 279 days in US Holsteins (Norman et al. 2011). At a genetic level, only a weak genetic correlation (0.18–0.38) was found between gestation length and dystocia rate in the Danish Holstein cattle (Hansen et al. 2004).

The next significant factor for both CART and GLZ was BCSI, i.e. the difference between heifer's BCS assessed before calving and BCS at calving. This result is in line with the sensitivity analysis of artificial neural networks and the sequence of variables determined by discriminant function analysis (Zaborski 2010; Zaborski & Grzesiak 2011), as well as with the results obtained by Basarab et al. (1993) in dystocia detection in beef heifers. An optimal BCS during the periparturient period ensures an easy calving. The currently recommended BCS in heifers is 2.75–3.00 points on a 5-point scale, which positively affects gestation length and reduces the risk of perinatal calf mortality (Lorenz et al. 2011). On the other hand, overconditioned heifers show a greater risk of dystocia and metabolic diseases in early lactation, whereas too thin ones do not reach an appropriate size and body weight at first calving (Schröder and Staufenbiel 2006). Nevertheless, a higher BCS does not always lead to a difficult labour, which may result from the difference in the growth of overfed heifers. This means that heifers utilizing additional energy for increasing the size of their frame and, to a lesser extent, calf birth weight may have relatively easy calvings, whereas those with excessive amounts of fat in the pelvic region, decreased lumen of the birth canal and weaker musculature of the uterus caused by overfeeding have more difficult calvings (Kolkman 2010). However, it should be added that not all authors are agreed that BCS significantly affects the course of parturition. For instance, Berry et al. (2007) did not find such an association.

Calving age of heifers was also a significant factor contributing to the determination of calving difficulty category. In the aforementioned work on dystocia detection using artificial neural networks (Zaborski & Grzesiak 2011), this variable was less important, whereas discriminant function analysis confirmed its high significance (Zaborski 2010). Similar conclusions were drawn by Basarab et al. (1993) in their study on dystocia detection in beef heifers. This variable played an important role in discriminant function models including variables available before breeding and during pregnancy. According to Mee (2008), a higher risk of dystocia generally concerns very young and relatively older heifers, aged 2–3 years, although some authors (King et al. 1993; Ettema & Santos 2004) showed lack of a direct relationship between calving age and difficulty.

The last considered factor, indicated by GLZ, is the ERα-SnaBI genotype. However, its role in the CART model was negligible. As already mentioned, oestrogen receptors, including ERα, as transcription factors, enable the action of these hormones on target cells and tissues (Szreder & Zwierzchowski 2004). Moreover, it has been shown that ERα, together with glucocorticoids, may be involved in the labour-initiating mechanisms (Kimmins et al. 2003). However, it should be emphasized that despite the unquestionable role of oestrogens and their receptors in the reproduction process regulation, the latest molecular research has not revealed any significant effect of the polymorphism within the ERα gene investigated by the authors on calving difficulty, although it significantly affected calf sex (Szreder et al. 2011). From among other genes with a possible or confirmed effect on dystocia, the following should be mentioned: hedgehog interacting protein (HHIP), affecting direct calving difficulty and calf size; microRNA mir-1256 (MIR1256), associated with maternal calving difficulty, sialic acid binding Ig-like lectin 5 (Siglec-5), affecting sire and daughter calving ease, rump width, stature and body depth, genes coding for zinc finger proteins (e.g. zinc finger protein 827), which play an important role in the regulation of bone and skeletal development in mammals; genes coding for the insulin-like growth factor (IGF)-binding proteins (e.g. IGF-binding protein 2 precursor gene – IGFBP-2, IGF-binding protein-3 gene – IGFBP-3), which inhibit IGFs with IGFBP-2 hypothesized to act as a growth-promoting hormone during gestation (Cole et al. 2011; Purfield et al. 2014); the non-SMC condensin I complex, subunit G (NCAPG) gene, modulating foetal growth in cattle (Eberlein et al. 2009); the solute carrier family 44, member 5 (SLC44A5) gene, associated with dystocia rate (Sugimoto et al. 2012); the secreted phosphoprotein 1 (SPP1) gene, also known as osteopontin (OPN), the integrin-binding sialoprotein (IBSP) gene as well as the matrix extracellular phosphoglycoprotein (MEPE) gene included in a cluster of bone-tooth mineral extracellular matrix phosphoglycoproteins associated with bone and cartilage morphogenesis (Olsen et al. 2010); the cluster of differentiation 37 (CD37) gene, affecting direct calving difficulty (Cole et al. 2009); the IGF1 gene, influencing the incidence of twin pregnancies in cattle (Kim et al. 2009), the ribosomal protein S20 (RPS20) gene and the small nucleolar RNA U54 (SNORD54) gene, both of which affect foetal growth (Pausch et al. 2011) as well as the bone morphogenetic protein receptor type 1A (BMPRA1), growth differentiation factor 2 (GDF2), peptidoglycan recognition protein 1 (PGLYRP1) and IGF-like family member 1 (IGFL1) genes which are candidates for calving ease (Cole et al. 2011; Macciotta et al. 2014). In addition, polymorphisms within the growth hormone, growth hormone receptor and leptin genes were found to affect growth traits and gestation length in dairy breeds (Biswas et al. 2003; Giblin et al. 2010; Maltecca et al. 2011).

Other candidate genes for birth weight include phosphatidylinositol glycan anchor biosynthesis, class N (PIGN), KH domain containing, RNA binding, signal transduction associated 3 (KHDRBS3), ribosomal protein L22-like 1 (RPL22L1), LIM homeobox 4 (LHX4), zinc finger protein 75a (ZNF75A), ATP-binding cassette, sub-family A, member 12 (ABCA12), netrin G2 (NTNG1), regulatory factor X, 6 (RFX6), fibronectin leucine rich transmembrane protein 2 (FLRT2), mitogen-activated protein kinase kinase kinase 5 (MAP3K5), chromosome 21 open reading frame, human C14orf180 (NRAC), AMME, dystrophin (DMD), PP2R1A, ataxia telangiectasia mutated (ATM), BTG, CD33 and LIH1 (Cole et al. 2014). Moreover, genes related to both regulation of growth and response to hypoxic stress, i.e. protein tyrosine kinase 2 beta (PTK2B) and hypoxia inducible factor 1, alpha subunit (HIF1A), as well as the ghrelin/obestatin prepropeptide (GHRL) gene, which has been linked to foetal growth in other species, are candidates for gestation length and dystocia (Maltecca et al. 2011). Finally, polymorphisms within the guanine nucleotide-binding protein Gs subunit alpha (GNAS), paternally expressed gene 3 (PEG3), pregnancy-associated plasma protein-A2 (PAPP-A2) and placental growth factor (PGF) genes have recently been shown to affect direct and maternal calving difficulty (Magee et al. 2010; Seidenspinner et al. 2011; Sikora et al. 2011; Wickramasinghe et al. 2011). The remaining predictors (calving season and genotypes) did not play any important role in calving class determination.

After building and training the CART and SVM models and estimating parameters of GLZ, their performance was verified on new data that were not used previously for model construction (the so-called test set). Sensitivity (percentage of correctly diagnosed difficult calvings) was acceptable only in the case of CART (Table 5). SVM and GLZ were characterized by a low value of this measure, taking into account the fact that the dependent variable had only two categories. However, it should be emphasized that the number of difficult calvings in the training data-set was much smaller compared with the second class, which undoubtedly influenced detection results, whereas specificity (percentage of correctly diagnosed easy calvings) and accuracy (percentage of correctly recognized heifers from both classes) were high. In the study on dystocia detection in heifers using artificial neural networks, Zaborski and Grzesiak (2011) obtained much higher sensitivity (75–83%) at lower specificity (82–88%) and accuracy (82–86%), while the application of discriminant function analysis and classification functions to this issue (Zaborski 2010) resulted in lower sensitivity (42%) at equally high specificity (97%) and accuracy (90%). The use of the aforementioned method for dystocia detection in Angus heifers was investigated by Arthur et al. (1999), who obtained lower sensitivity (0–40%) and specificity (73–90%) than those in the present study. On the other hand, discriminant function analysis applied to the diagnosis of difficult calvings in heifers and cows with different proportion of Chianina genes (Morrison et al. 1985b) yielded the sensitivity of 57% (higher than that in the present study, except for CART), specificity of 90% and accuracy of 84% (lower than the values in the present work). However, one should bear in mind that the most important measure of detection performance of the models on new, previously unseen data, is sensitivity, which was not sufficiently high in our study (except for CART), although it was comparable with the above-mentioned results.

Additional detection measures calculated from the test set were the probability of FP and FN. The former is associated with an erroneous diagnosis of normal heifers as dystotic ones by the model, whereas the latter indicates the proportion of heifers with difficult calvings diagnosed as ‘easy calvers’. It is obvious that the second situation results in more serious consequences. Hence, it important to minimize P(FN). In the present work, P(FP) was very low (Table 5) and P(FN) assumed relatively high values; however, this results directly from sensitivity and specificity obtained in the present study, since these measures are closely connected with each other, i.e. Se + P(FN) = 1 and Sp + P(FP) = 1. The above-mentioned authors (Morrison et al. 1985b; Basarab et al. 1993; Arthur et al. 1999) obtained similar or higher P(FP) (1–27%) and P(FN) was in quite a broad range of values (43–100%). The probabilities describing the reliability of diagnosing animals from both categories were P(PSTP) and P(PSTN), which assumed high values in the present study. They were higher or comparable with those of other authors (10–57% and 88–94%, respectively; Morrison et al. 1985b; Basarab et al. 1993).

In order to more comprehensively evaluate the detection performance of the CART, SVM and GLZ models, the ROC curves were plotted and the AUC values were computed (Figure 3). These values in the present study were slightly lower than those in the detection of difficult calvings in heifers using artificial neural networks and classification functions (0.85–0.89; Zaborski 2010; Zaborski & Grzesiak 2011). For comparison, in the application of the CART and CHAID decision trees to the analysis of factors affecting dystocia and stillbirths in dairy cows (Piwczyński et al. 2013), the AUC values ranged from 0.61 to 0.71. However, it should be emphasized that the aforementioned values were calculated based on a training set or cross-validation.

One of the important aspects of developing statistical and data mining models is their possible use for minimizing the incidence of dystocia. As already mentioned, the course of parturition in cattle is affected by many factors. Some of them are difficult to capture and may occur without clear reasons. Others have been well described in the available literature. The use of factors that can be identified and that affect calving difficulty may, to some extent, determine its final outcome. The results of the present study and those of other authors show that the inclusion of these predictors in various statistical models may result in a greater precision of dystocia detection in heifers. This, in turn, may help to alleviate its disadvantageous consequences resulting in substantial economic losses. Additionally, when using data mining methods, it is possible to adjust prognosis depending on changing conditions. Moreover, the detection of heifers with dystocia, even if biased, is much more advantageous (the possibility of taking some countermeasures – greater supervision or veterinarian's assistance) than the situation in which difficult calving occurs without any prior information on it.

Finally, it should be stated that the present study is not free from limitations. One of them is a small sample size used for analysis. Although it was statistically justified, a larger number of training cases would be necessary to definitely confirm the results obtained in the present work. Moreover, calving records of heifers were collected on only one farm. The farm was large as on domestic scale (approx. 700 cows) and was one of many similar farms located in the province. It was situated in the lowland area; however, climatic conditions are similar within the country. Also, the practical application of the models in breeding practice would require their higher effectiveness in diagnosing dystotic animals in order to limit the costs associated with unnecessary help provision or intervention.

4. Conclusion

In conclusion, it should be stated that the best quality was characteristic of the CART model, which had also the best detection performance evaluated on new data that were not used previously for classification model construction. In the case of SVM and GLZ, percentage of correctly diagnosed animals with calving difficulties was relatively low. Therefore, the practical use of these models in animal farming would not be advisable. However, all presented methods were characterized by high specificity and reliability of the predictions being made. The most important variables determining calving category were gestation length, body condition score index, calving age and, to a lesser extent, ERα genotype. The two remaining genotypes (and calving season) did not significantly contribute to prediction accuracy. Future research should focus on the possibility of using different predictors and other data mining methods in order to achieve even better detection performance.

Funding

This work was supported by the Polish Ministry of Science and Higher Education [grant number 517-01-028-3962/17].

Additional information

Funding

Funding: This work was supported by the Polish Ministry of Science and Higher Education [grant number 517-01-028-3962/17].