Instruments
The acceleration data were collected using a wireless sensor system (Ominsense Series 500 Cluster Geolocation System [31]; http://www.omnisense.co.uk/) that includes an embedded tri-axial accelerometer (Xtrinsic MMA8451Q 3-Axis, 14-bit/8-bit Digital Accelerometer with a sensitivity between −8 and +8 g). Accelerometer data was collected at 50 Hz which allowed for effective battery life of approximately 2 days. The wireless sensors contain a 2.4 GHz, IEE 802.15.4a transmitter module to remotely send messages to the CLS-504 location server. The Series 500 sensors can be used to form a wireless mesh sensor-node network that is able to compute relative spatial locations of the sensor nodes using the arrival time of periodic messages sent from each node to its neighbours. In principle, acceleration data could be processed on the sensor in real-time and outputs sent across the network as part of a more general monitoring system. However, in this study, only data from the tri-axial accelerometer were recorded using a 4 GB micro SD flash memory card for posterior data analysis. The sensors were fixed in the same orientation on the right hand side of a neck collar worn by the cows (Fig. 1a). Counterweights (0.5 kg) were used on the neck collars to ensure a stable position of the sensor on the body of the animal. The sensor weighs approximately 0.25 kg in total (including batteries), half the weight of the counterweight. The coordinate frame of the sensors corresponds to X forwards, Y right and Z down as shown in Fig. 1b). At the end of the study, the SD card was removed from the sensor and the accelerometer data was converted from its hexadecimal format to g units (g = 9.81 ms−2).
Study site, animals and observation of behavioural activities
The data collection was carried out on a commercial farm of Holstein dairy cattle located in Essex, UK. The cows where loose housed in a cubicle shed. The herd was milked three times a day at approximately 5 a.m (morning), 1 p.m (afternoon) and 9 p.m (evening). The duration of milking time for each individual cow varied between 1 and 1 ½ h. The herd mean 305-day milk yield was 11,000 litres per cow. Cows were fed a commercial total mixed ration. A total of six cows that had not shown signs of severe lameness, or other disease that might affect their behavioural repertoire, were selected for this study. Cows were selected and collared during morning or afternoon milking and were wearing the collar for a maximum of 2 days (since battery of the sensors could not be guaranteed after this point). Cows were monitored between milking periods; during milking, no visual observations were recorded.
Cow behavioural activities were recorded by observers (ZB and HH) performing a visual focal tracking on each individual cow that was wearing a sensor collar according to the following criteria for each behavioural activity:
-
1.
Feeding (state): cows located at the feeding zone, ingesting food;
-
2.
Lying (state): cows located in a cubicle in a lying down position;
-
3.
Standing (state): cows standing on their four legs;
-
4.
Lying down (transition): cows changing from a standing state to a lying state;
-
5.
Standing up (transition): cows rising from a lying state to a standing state
Drinking, brushing and walking activities were observed less frequently and for short durations and therefore not considered for classification in the algorithm. It should be stressed that these rarer activities and events may still be biologically important in the context of detecting health and welfare status. Hence, although we do not try to classify them here, future studies should also consider methods for detecting these rarer behaviours.
From the data set of visual observations, only the activities of interest for this study were selected to validate the classifier algorithm. The new data set used for validation contains the following observational data:
-
Cow 1: a total of 8 h and 2 min extracted from observation on 28 August 2014, between 08:00–18:00; 29 August 2014, between 08:00–12:40
-
Cow 2: a total of 3 h and 6 min extracted from observations on 28 August 2014, between 08:00–18:00
-
Cow 3: a total of 5 h and 40 min extracted from observations on 3 September 2014, between 15:00–17:30; 4 September 2014, between 08:00–18:00
-
Cow 4: a total of 5 h and 5 min extracted from observations on 15 September 2014, between 15:30–18:30; 16 September 2014, between 07:40–17:40
-
Cow 5: a total of 6 h and 5 min extracted from observations on 3 September 2014, between 15:10–17:10; 4 September 2014, between 08:00–18:10
-
Cow 6: a total of 5 h and 22 min 15 September 2014, between 15:30–18:10; 16 September 2014, between 07:30–17:50
In total, direct visual observations of the cows were completed for 33 h and 20 min, of which 15 h and 30 min were lying, 4 h and 10 min were standing and 13 h and 40 min were feeding. All behavioural observations were entered into a spreadsheet with the start and stop time of every activity and identification of the corresponding cow. Observer and sensor watches were synchronised at the start of the observation period so that observation data could be accurately aligned with the tri-axial accelerometer data retrieved from the sensors in a single database.
Algorithms for behavioural state classification
Raw acceleration data
Figure 3a illustrates example time series of the raw tri-axial accelerometer output for observed periods of lying, standing and feeding behaviour for a single cow. It is clear that there is very little qualitative difference in the acceleration output for the lying and standing behaviours, since for both these behaviours the cow exhibits very little overall movement. When the cow is feeding there is a clear regular shift in the acceleration in the y and z axes that corresponds to the cow moving its head up and down. Figure 3 is only a representative example but similar qualitative patterns in the acceleration output were observed for the other cows in the study. These qualitative observations offer a useful intuitive starting point for determining the most appropriate feature characteristics to include in the classification algorithm.
Feature characteristics
Machine-learning algorithms use feature characteristics (also called summary statistics) calculated from the input data (e.g. the raw accelerometer data) to classify different states (e.g. feeding, lying or standing). The algorithms in this study have been developed using two intuitive and easy to interpret characteristic features based on the biomechanics of the movement behaviour of the cows. These two feature characteristics consist of two different components of the raw acceleration data: a static component caused by the gravity field (SCAY) and a dynamic component caused by the movement of the animal (vectorial dynamic body acceleration, VeDBA [34, 35]). Other studies have used a far larger number of feature characteristics (e.g. 30 or even higher) [5, 6, 8]. In our study, the use of only two features was motivated by the need to reduce computational time and complexity and also to allow more intuitive biological interpretation of the results.
Figure 3b illustrates a typical example time series of running mean in the y-axis and VeDBA output for observed periods of lying, standing and feeding behaviour for a single cow. Low VeDBA output values for lying and standing are caused by the low movement exhibited by cows during these behaviours. In contrast, high VeDBA values obtained for feeding are caused by the upward and downward head movement cows perform during this behaviour. In this figure, it is also possible to observe a small difference in the SCAY outputs between lying and standing. Since the running mean in the y-axis represents the static component caused by the gravity field, output values obtained for this parameter correspond to the orientation of the sensors during the behaviour as seen in Fig. 1c, d. Figure 1c shows an example of the orientation of the sensor when the cow was observed standing, while Fig. 1d shows the orientation of the sensor when the cow was observed lying. The component in the y-axis of the gravity field is given by \( \overrightarrow{{\mathit{\mathsf{g}}}_{\mathit{\mathsf{y}}}}=\mathit{\mathsf{g}}\ast \mathsf{cos}\ \left(\mathsf{180}-\beta \right) \) . Using this expression, a preselected threshold of −0.055 g for the static component in the y-axis corresponds to an angle of β = 86.84° (where an angle of β = 90° can be interpreted as the cow having its neck aligned horizontally). Therefore, the decision-tree classifies standing and lying behaviour if the neck (and therefore sensor) is above or below this threshold. Figure 1c, d are only representative examples, but similar patterns in the static component were found for other cows in this study.
The VeDBA and SCAY feature characteristics are calculated as a mean over a given moving window size centred at the time point of interest (see Additional file 1). This requires a moving window size to be specified before any algorithm is run. A range of moving window sizes was tested for each algorithm and we report results for sizes of 1, 5 and 10 min (Table 2). Results for other moving window sizes are explored for the decision-tree algorithm in Additional file 1.
Machine-learning algorithms
There are a range of different machine-learning algorithms that could be used to classify different animal behaviours. These algorithms can be described as either supervised or unsupervised approaches. A supervised learning algorithm is formed by two processes: training and testing. A supervised learning algorithm uses a known data set to construct a model (training process) that is then used for making predictions on a new data set (testing process). Unsupervised machine-learning algorithms explore the data to find hidden patterns or to cluster the data input in classes with similar statistical properties. In this study, the three following unsupervised algorithms for the classification of the dairy cow behaviours were used: decision-tree, k-means and a HMM. The decision-tree was selected based on its simple structure and low computational cost, making it feasible to be implemented directly in a remote sensor. The selection of the k-means algorithm was based also on the simplicity of its structure and the possibility to compare the decision-tree to methods with similar levels of simplicity (although the k-means may have high computational costs due to a recursive component in the algorithm). The HMM was chosen in order to compare the decision-tree performance with a more sophisticated statistical model that is often used to classify animal behavioural states [23, 24]. Finally, a supervised SVM algorithm was also chosen in order to compare the decision-tree performance to a more complex algorithm that has been used for the classification of accelerometer data to distinguish between different behaviours in dairy cows [8]. The decision-tree and k-means algorithms were custom written by the authors in Matlab [36]. The HMM was applied using the Matlab toolbox for HMMs developed in [37]. The SVM was applied using the machine-learning toolbox provided in [38].
Decision-tree
A full description of the decision-tree algorithm used in this study is available in Additional file 1. We summarise the key features of the algorithm here. The decision-tree algorithm uses two rules with associated thresholds to classify tri-axial acceleration data as either feeding (high activity) or lying or standing (both low activity). The first rule in the decision-tree uses the mean of the VeDBA values and a predefined threshold A to discriminate between cases with high and low energy expenditure activities. Those cases resulting in a high energy expenditure activity are labelled as feeding, and those with low energy expenditure activities are used in the second step of the decision-tree (Fig. 2). The second decision rule of the tree compares the running mean of the acceleration in the y-axis (SCAY) to a predefined threshold B value in order to partition the data into two clusters (mean of static component in the y-axis above or below the threshold value). Cases resulting in values below the threshold are labelled as lying, and those with values above are labelled as standing (Fig. 2). A range of different predefined threshold values were considered (see Additional file 1), and values of A = 0.0413 g and B = 0.055 g were found to give the best performance with this data set. Similarly, to explore the effect of the choice of window size, the performance of the algorithm was investigated using windows ranging from 1–600 s (window sizes above 600 s resulted in too few data points for a fair comparison of performance) and full details are given in Additional file 1.
k-means
Observations for the k-means algorithm are given by the 2-dimensional feature characteristics. The first dimension is represented by the mean of the VeDBA values over the window size, whereas the second dimension is represented by the mean of the acceleration in the y-axis (SCAY). The k-means algorithm discriminates between the observations in one step using both feature characteristics at the same time. This represents a key difference between the decision-tree and the k-means, since the former uses one feature characteristic at each decision rule. A full description of the k-means algorithm is given in Additional file 1.
Hidden Markov model
A sequence of behaviours in dairy cows can be modelled as a first-order HMM with a finite number of hidden states (behaviours) where each activity can be observed through a set of characteristic features (observations). The observations for the HMM correspond to the same characteristic features used for the decision-tree, i.e. mean of VeDBA over the window size and running mean of the acceleration in the y-axis over the window size (SCAY). The hidden Markov model was applied using the Matlab toolbox for hidden Markov models developed in [37]. This toolbox randomly generates an initial transition probability matrix A and an initial probability π. The emission probability distribution B is initialised using a static Gaussian Mixture model. Since the results can depend on the initialisation parameters, we run a total of 100 random initialisations to select the highest scoring model. Further details of the implementation, use and application of the Baum-Welch, the Viterbi and the forward-backward algorithms for HMMs can be found in [39], and further details are given in Additional file 1.
Support vector machines
SVMs are a supervised learning algorithm requiring training and testing processes. In this study, training was performed using k-1 folds and tested in the fold left out. We used a 3-fold cross validation for the implementation of the SVM algorithm. Further details of the SVM algorithm are provided in Additional file 1 and can also be found in [38, 40–42].
Performance of the classification algorithm
Comparison of algorithm classification performance
The performance of the decision-tree classification algorithm was compared across a range of values for the algorithm parameters (window size, thresholds A and B); for details see Additional file 1. The performance of the algorithm was directly compared to alternative classification algorithms such as k-means, HMM and SVM using the same input data set (Table 2). The performance of an automated behavioural classification algorithm can often vary across individuals or breeds of the same species [33]. Hence, we also considered the performance of the decision-tree algorithm at the level of the individual cow. In order to do this, we computed the performance metrics for each individual cow at a window size of 1 min (Table 3). The 1 min window was selected in this context to avoid having only a small number of samples for each individual cow (which can occur at larger window sizes).
Sensitivity and precision
When comparing algorithm classification performance, we considered two performance metrics: the sensitivity of classification and the precision of classification. In standard statistical process control, the sensitivity (Sen) and precision (Pre) are defined as:
$$ \mathsf{S}\mathrm{e}\mathrm{n}=\frac{\mathsf{TP}}{\mathsf{TP}+\mathsf{F}\mathsf{N}}, $$
$$ \mathsf{Pr}\mathrm{e}=\frac{\mathsf{TP}}{\mathsf{TP}+\mathsf{F}\mathsf{P}}. $$
Here, TP (true positive) is the number of instances where the behavioural state of interest that was correctly classified by the algorithm after validation by the visual observer. FN (false negative) is the number of instances where the behavioural state of interest was visually observed in reality but was incorrectly classified as some other behaviour by the algorithm. FP (false positive) is the number of times the behavioural state of interest was (incorrectly) classified by the algorithm but not observed in reality. TN (true negative) is the number of instances where the behavioural state of interest was (correctly) classified as not being observed.
An algorithm for detection of transitions between lying and standing
A further two-step algorithm was developed to detect the transitions between lying and standing (Table 4). The first step of the algorithm (non-specific) uses a threshold over the range of the acceleration in the y-axis to determine if a transition occurs or not. Range in the y-axis represents a good candidate for the threshold due to the biomechanics of the rapid change in this axis when cows exhibit a transition between lying and standing or vice-versa. As described by Martinskainen et al. [8], a cow that lies down bends one front leg, lowers its forequarters then its hindquarters until it settles into a lying position. When a cow stands up, it lunges forward, lifts its hindquarters, then rises to stand up on its four legs. According to this definition and the orientation of the sensors in Fig. 1a, b, a transition movement implies a significant change in the orientation of the sensor in the y-axis (Fig. 4a).
The second step of the transition detection algorithm is performed by applying the decision-tree classification algorithm described previously to infer the anterior and posterior behaviour on either side of the transition and hence discriminate between standing up and lying down. Further details of the transition detection algorithm are given in Additional file 1.
Availability of supporting data
The data collected as part of this study is available in Additional file 2.