# Weakest link pruning of a dendrogram

Jiacheng Ge 111Department of Statistics, Stanford Univ., and Robert
Tibshirani222Department of Departments of Biomedical Data Science and of
Statistics, Stanford Univ.

###### Abstract

Keywords: unsupervised learning, clustering, dendrogram, pruning

Hierarchical clustering is a popular method for identifying distinct groups in
a dataset. The most commonly used method for pruning a dendrogram is via a
single horizontal cut. In this paper, we propose a new technique “weakest link
optimal pruning”. We prove its superiority over horizontal pruning and provide
some examples illustrating how the two methods can behave quite differently.

##  1 Introduction

Hierarchical clustering is a widely used method for identifying distinct
groups in data in an unsupervised setting. There are two main types of
hierarchical clustering strategies: agglomerative (bottom-up) and divisive
(top-down). Both strategies use a greedy approach and we focus on
agglomerative clustering in this paper. One benefit of hierarchical clustering
is that any legitimate distance metric can be applied. Common distance metrics
include average linkage, single linkage, and complete linkage. The hierarchy
of clusters is usually presented in a dendrogram with desired clusters as its
branches. The height of a node is proportional to its within-cluster
dispersion. “Pruning” the dendrogram involves the removal of some of its lower
branches.

The most popular dendrogram pruning method is the horizontal cut: we merge the
leaves into several clusters either by specifying the desired number of
clusters or the cut height. Every contiguous branch of objects below the cut
is grouped together. This technique is a straightforward, simple method with
many appealing qualities.

Recall that a dendrogram is constructed from the bottom up, by merging the
closest pair of nodes at each stage, and it is displayed with the height of
each node equal to the distance at which the merge occurs. From this, one
might expect that the horizontal cut yields the tightest set of clusters for a
given number of leaves. Due to the greedy nature of the merges, however, this
turns out not to be the case, as we show here.

Dynamic Tree Cut (Langfelder et al.) [6] is a dendrogram pruning method based
on investigating the dendrogram’s branch shapes. It has two variants: “Dynamic
Tree” cut and “Dynamic Hybrid” cut. The Dynamic Tree cut is a top-down
approach that involves adaptively breaking down and recombining clusters until
a stable number of clusters is reached. The Dynamic Hybrid cut, on the other
hand, is a bottom-up process that involves identifying branches that meet
certain criteria and then assigning unassigned objects to the clusters found
in the first step if they are sufficiently close. Dynamic tree cut is flexible
and has been found to work well on bioinformatic tasks. However, it is
computationally intensive and hard to interpret. Also, it is still unclear how
to select the height and shape parameters and how to estimate the number of
clusters.

Dynamic programming, which was created by Richard Bellman [2], is a technique
that is used in both mathematical optimization and computer programming to
solve complex problems by breaking them down into simpler sub-problems and
solving them recursively. We aim to minimize a chosen loss when pruning a
dendrogram into $k$ groups. While it is possible to enumerate all
combinations, a more efficient method is to identify the optimal subtrees by
dynamic programming [1]. To optimally cut a tree into $k$ clusters, it
equivalently determines the optimal number of groups in the left and right
branches of the root ($k_{l}$ and $k_{r}$, respectively), with the condition
that $k_{l}+k_{r}=k$. This process is then repeated for each branch,
recursively until it reaches the leaves of the dendrogram or all remaining
nodes fall into one cluster. An important drawback of this method is that it
does not necessarily produce a nested sequence of subtrees.

In Section 2, we propose weakest link optimal pruning and prove theoretically
its advantage over the horizontal pruning method. In this paper, we mainly
focus on the comparison between the horizontal cut and our optimal cut. We
provide some examples in Section 3, including simulated examples with and
without clusters, and some popular real datasets, illustrating how the two
methods behave differently. In Section 3.4, we embed the two pruning methods
into a classification task using a real dataset and compare their
performances. We combine the optimal pruning sequence with the Gap test in
Section 4 to estimate the ideal number of clusters. Finally, Section 5
contains some discussion.

##  2 Our proposed method: weakest link pruning

###  2.1 Set up

Our approach follows closely the weakest link pruning method used in
classification and regression trees (CART) [3], a supervised learning method.
In our unsupervised clustering problem, the only change we make here is in the
loss function, where we replace misclassification loss with the within-cluster
dispersion loss.

Let $T$ be a dendrogram and $t$ be a node in $T$. If there is only one
observation in a node $t$, we call it a terminal node or leaf. Let $\tilde{T}$
denote all terminal nodes in $T$. Denote by $D(x_{i},x_{i}^{\prime})$ the
distance between points $x_{,}x_{i}^{\prime}$: most commonly we use squared
error loss $||x_{i}-x_{i}^{\prime}||^{2}$. The loss function $R(t)$ is defined
as the sum of the distances between all possible pairs of data points in node
$t$, i.e.

$R(t)=\underset{x_{i},x_{i^{\prime}}\in t}{\sum}D(x_{i},x_{i^{\prime}}).$

The loss function $R$ for $T$ is defined as the total dispersion loss within
all leaves, i.e.

$R(T)=\underset{t\in\tilde{T}}{\sum}R(t)=\underset{t\in\tilde{T}}{\sum}\,\underset{x_{i},x_{i^{\prime}}\in
t}{\sum}D(x_{i},x_{i^{\prime}}).$

Let $T_{max}$ denote an unpruned dendrogram with only one observation in each
terminal node. If $T$ is a subtree of $T_{max}$, we denote $T\subseteq
T_{max}$. Our goal is to find the “optimal” dendrogram $T\subseteq T_{max}$.
However, $R(T)$ alone is not a useful metric for finding the “optimal” subtree
because it always prefers larger trees. It is easy to see that $R(T)=0$ when
each leaf only contains one case. So, we introduce a complexity penalty to
control the balance between dispersion and tree size. For any subtree
$T\subseteq T_{max}$, we define its complexity as $|\tilde{T}|$, the number of
terminal nodes in $T$. Let $\alpha\geqslant 0$ be a real number and called the
complexity parameter. We define the cost-complexity measure $R_{\alpha}(T)$
as:

$R_{\alpha}(T)=R(T)+\alpha|\tilde{T}|.$

Our general approach is to (a) for each $\alpha$, find the subtree that
minimizes $R_{\alpha}(T)$ and then (b) estimate the best value $\hat{\alpha}$
and hence the optimal subtree $T(\hat{\alpha})$ using a method for estimating
the best number of clusters. In Section 4 below we illustrate this using the
Gap Test, but emphasize that other methods for estimating the number of
clusters could be used.

For (a), our goal is to find subtree $T(\alpha)\subseteq T_{max}$ which
minimizes $R_{\alpha}(T)$; $T(\alpha)$ is called the minimizing subtree for
$\alpha$, i.e.

$R_{\alpha}(T(\alpha))=\underset{T\subseteq
T_{max}}{\text{min}}R_{\alpha}(T).$

The smallest such $T(\alpha)$, in terms of complexity, is defined as the
smallest minimizing subtree for the complexity parameter $\alpha$. By its
definition, the smallest minimizing subtree for a given $\alpha$ is unique.

In this paper, we use squared distance as the dispersion loss, i.e.

$D(x_{i},x_{i^{\prime}})=\|x_{i}-x_{i^{\prime}}\|^{2},$

then

$R(t)=\underset{x_{i},x_{i^{\prime}}\in
t}{\sum}\|x_{i}-x_{i^{\prime}}\|^{2}\quad{\rm and}\quad
R(T)=\underset{t\in\tilde{T}}{\sum}R(t)=\underset{t_{j}\in\tilde{T}}{\sum}\,\underset{x_{i},x_{i^{\prime}}\in
t_{j}}{\sum}\|x_{i}-x_{i^{\prime}}\|^{2}.$

###  2.2 Weakest link optimal pruning

To find the smallest minimizing tree, we adapt weakest link pruning: we
successively collapse the internal node that causes the smallest per-node rise
in $R(T)$ until a single-node (root) dendrogram is produced. Given a
dendrogram of $n$ leaves, choosing to cut the internal node with the smallest
per-node increase in $R(T)$ is equivalent to finding the subtree with fewer
than $n$ leaves which has the smallest value of $R$. By performing weakest
link pruning, we generate a sequence of minimal cost-complexity subtrees.
Notably, the subtrees are decreasing in their complexity, but their number of
leaves may not be consecutive. In other words, some tree sizes could be
skipped.

On pages 66-71 of CART book [3], the following facts are established:

  1. 1.

For every value of $\alpha$, there exists a (unique) smallest minimizing
subtree.

  2. 2.

The weakest link sequence of minimal cost-complexity trees contains the
smallest minimizing subtree.

The proofs of these claims follow directly from the CART results [3].

###  2.3 Optimal pruning versus horizontal pruning

First, we claim that weakest link pruning always yields a total dispersion
loss equal to or less than the horizontal cut for a given number of terminal
nodes. We will prove this by contradiction. Let $T_{max}$ be a dendrogram and
$\\{T_{s_{i}}\\}$ be its sequence of minimal cost-complexity subtrees. Assume
that there exists a horizontally pruned dendrogram $T_{h}\subseteq T_{max}$
that is not in $\\{T_{s_{i}}\\}$ and the dendrogram has a smaller loss than
any of the optimal pruned subtrees, i.e.
$R(T_{1})\leqslant\underset{s_{i}}{{\rm min}}R(T_{s_{i}})$. The horizontal cut
is denoted by $h_{1}$. All the nodes equal to or below $h_{1}$ are already
collapsed while all nodes above are not. Denote the pruned dendrogram before
$h_{1}$ happening by $T_{h^{\prime}}$. Hence:

  * •

If all the nodes above the $h_{1}$ have a larger per-node rise in
$R(T_{h^{\prime}})$ than all the nodes equal to or below $h_{1}$ when being
collapsed, then by the definition of weakest link pruning, $T_{1}$ must be in
$\\{T_{s_{i}}\\}$.

  * •

If there exists a node $t_{1}$ above $h_{1}$ and a node $t_{2}$ equal to or
below $h_{1}$ such that per-node rise of $R(T_{h^{\prime}})$ when cutting
$t_{1}$ is smaller than cutting $t_{2}$, then we can reduce $R(T_{h})$ by
collapsing $t_{1}$ and expanding $t_{2}$.

In both cases, there’s a contradiction to our assumption. Therefore, the
horizontal cut solution is either in the weakest link merge sequence or the
resultant subtree has a higher loss, implying that weakest link pruning is
never worse than horizontal cut for a fixed number of leaves.

Second, recall that a sequence of minimal cost-complexity subtrees may skip
some tree sizes. Without loss of generality, let’s assume the tree
complexities in the optimal subtree sequence are
($n,n-1,...,k+c+1,k+c,...,k,...,1$). When we have a subtree of $k+c+1$ nodes,
we have the option to prune it into k+c, $k+c-1,...,k$ leaves. Since we choose
the next merge by comparing the smallest per-node increase in $R$, then the
best subtree with $k+c-1,k+c-2,...,k+1$ leaves must have bigger losses than
the best subtree with $k+c$ leaves or $k$ leaves. By transitivity, the
horizontal pruned trees with $k+c-1,k+c-2,...,k+1$ leaves will have even
larger losses than the best minimal cost-complexity subtree with $k+c$ leaves
or $k$ leaves.

To summarize, given a horizontally pruned dendrogram $T_{h}\subseteq T_{max}$,
we can always find a $T$ in the sequence of minimal cost-complexity subtrees
of $T_{max}$ such that $R(T_{h})\geqslant R(T)$. This implies that horizontal
pruning can never do better than weakest-link optimal pruning.

##  3 Some examples

###  3.1 A simple example

Here we illustrate that the optimal pruning and horizontal cuts can behave
differently. Consider a dataset of five (scalar) data points: 13, 0, 10, 1, 3,
indexed by 1, 2, 3, 4, 5. We created a dendrogram for the dataset via average
linkage, which is shown in Figure 1. The left plot labels the leaves by their
data values and the right one shows the node numbers 1–5.

(a) Leaves labelled by their data values

(b) Leaves labelled by the node numbers 1–5

Figure 1: Cluster dendrograms for a simple example

The left plot of Figure 2 shows the cuts in the horizontal pruning process
using the cutree function in the R language package stats [7]. This process
collapsed the lowest internal node at each step from bottom to top. The right
plot labels the merges in the optimal pruning process. Notice that the cuts
were not increasing in height at each step. (It is worth mentioning that it’s
also common for some internal nodes to be skipped.)

(a) Cuts in horizontal pruning

(b) Cuts in optimal pruning

Figure 2: The cuts in two different pruning processes

Figure 3 details the subtrees produced by two different pruning methods.
Notably, resulting subtrees with 5, 4, 2, and 1 terminal nodes were the same
for both methods. However, optimal pruning merged node 1 and node 3 before
merging node (2 4) and node 5. And its resulting within-cluster sums of
squares (WSS) is 10 when there were 3 leaves, smaller than the value 14 for
the horizontal cuts.

(a) Horizontal pruning process

(b) Optimal pruning process

Figure 3: The sequence of subtrees from two different pruning strategies.

###  3.2 Experiments with simulated datasets

####  3.2.1 Random datasets with no clustering

We randomly generated 200 datasets $\mathit{X_{i}}$ with number of
observations $n_{i}$ between 30 and 100, and number of features $p_{i}$
between 1 and 50. We chose $n_{i}$ and $p_{i}$ uniformly sampled and $x_{ij}$
$\overset{\mathrm{iid}}{\sim}$ $\mathcal{N}(0,\mathit{I_{p_{i}}})$. Then, we
ran the two pruning methods, horizontal and optimal weakest-link pruning, on
the 200 sets of data for comparison. Figure 4 shows the log of WSS when the
dendrogram was pruned down to between 2 and 25 terminal nodes. Each colored
circle represents the mean of the (log) WSS and the error bar represents one
standard error. For optimal pruning, if the subtree with a certain number of
leaves was not in the sequence of optimal pruned subtrees, we increased the
number of leaves to the closest available number. For the optimal_skip
pruning, if the subtree with a certain number of leaves was not in the
sequence of optimal pruned subtrees, we skipped this subtree, thus the loss is
defined as NA.

We see that the optimal cut has consistently smaller WSS, sometimes by a large
margin. For some specific numbers of clusters, the variance of the
optimal_skip WSS is higher because some number of clusters are skipped and the
amounts of losses are reduced.

Figure 5 shows the relative improvement in the loss for weakest-link pruning
relative to horizontal cut. In this plot, the relative differences were only
calculated when the subtree with the clusters was a member of the sequence of
minimizing subtrees. The loss reduction is quite remarkable and the median of
loss reduction (after removing NAs due to skipping numbers of clusters) was
0.51. As demonstrated by both plots, optimal pruning shows better performance
than the horizontal cut when the data is completely random.

Figure 4: Log WSS mean and SE comparison of the two pruning methods on 200
random datasets with no clusters. “Optimal pruning” increases the number of
leaves to the closest available number if the subtree with a certain number of
leaves is not in the sequence of optimal pruned subtrees, while “Optimal_skip”
skips the subtree. Figure 5: Relative reduction in WSS, for optimal pruning
versus horizontal cut on 200 random datasets with no clustering.

####  3.2.2 Random datasets with well-separated clusters

We randomly generated 200 datasets with number of observations $n_{i}$ between
30 and 100, number of features $p_{i}$ between 1 and 50, and the numbers of
clusters $c_{i}$ between 3 and 15, all uniformly sampled. We generated
$x_{ij}$ $\overset{\mathrm{iid}}{\sim}$ $\mathcal{N}(0,\mathit{I_{p_{i}}})$
first, and then shifted the data from each cluster by a vector of length
$p_{i}$ whose components repeat an integer uniformly sampled from 1 to $c_{i}$
without replacement. The sizes of all clusters were equal, except for the last
cluster which is adjusted for indivisibility. The 2D projection of an example
dataset produced by Rtsne R package [5] is illustrated in figure 6. We can see
that the clustered data are well separated with a small amount of overlap.

Figure 6: TSNE projection of a random dataset with well-separated clusters:
100 points, 15 dimensions, 10 clusters

We ran the two pruning methods, horizontal and optimal weakest-link pruning,
on the 200 datasets. Figure 7 shows the log of WSS when the dendrogram was
pruned down to between 2 and 25 terminal nodes. The colored circles are the
mean of log WSS, and the error bar represents one standard error. For optimal
pruning, if the subtree with a certain number of leaves was not in the
sequence of optimal pruned subtrees, we increased the number of leaves to the
closest available number. For the optimal_skip pruning, if the subtree with a
certain number of leaves was not in the sequence of optimal pruned subtrees,
we skipped this subtree and thus the loss was set to NA.

Figure 8 shows the relative reduction in WSS for weakest-link pruning relative
to horizontal cut. In this plot, the relative differences were only calculated
when the subtree with the input number of clusters was a member of the
sequence of minimizing subtrees. The results were similar to the previous
results for the no-cluster setting, but the improvement is less dramatic. The
median of relative reduction in loss (after removing NAs due to skipping
numbers of clusters) was 0.085.

Figure 7: Log WSS mean and SE comparison of the two pruning methods on 200
random datasets with well-separated clusters. Optimal pruning increases the
number of leaves to the closest available number if the subtree with a certain
number of leaves is not in the sequence of optimal pruned subtrees, while
Optimal_skip skips the subtree. Figure 8: Relative reduction in WSS, for
optimal pruning versus horizontal cut on 200 random datasets with well-
separated clusters.

###  3.3 Some real data examples

Here we applied the pruning methods to six popular datasets from the book
“Introduction to Statistical Learning”, Second Edition (ISLR2) [4]: Boston,
Khan(train), Khan (test), NCI60, Portfolio, and USArrests. Details can be
found in table 1.

Dataset | Number of observations | Number of features |  |  |
---|---|---|---|---|---
Boston | 506 | 13 |  |  |
Khan (train) | 63 | 2308 |  |  |
Khan (test) | 20 | 2308 |  |  |
NCI60 | 64 | 6831 |  |  |
Portfolio | 100 | 2 |  |  |
USArrests | 50 | 4 |  |  |
Table 1: Six real datasets from the ISLR2 text.

We ran the two pruning methods on the six real datasets for comparison.
Combining results for all six datasets together, Figure 9 shows the log of WSS
when each dendrogram was pruned to between 2 and 19 terminal nodes. The
colored circles represent the mean of (log) WSS and the error bars represent
one standard error. If the subtree with a certain number of leaves was not in
the sequence of optimal pruned subtrees, we increased the number of leaves to
the closest available number. Figure 10 shows the relative reduction in loss
when we replaced the horizontal cut with weakest link cut. In this plot, the
relative differences were only calculated when the subtree with the input
number of clusters was a member of the sequence of minimizing subtrees. The
reduction in loss for optimal pruning versus horizontal cut was quite
substantial.

Figure 9: Log WSS mean and SE comparison of two pruning methods on six real
datasets. Optimal pruning increases the number of leaves to the closest
available number if the subtree with a certain number of leaves is not in the
sequence of optimal pruned subtrees, while Optimal_skip skips the subtree.
Figure 10: Relative reduction in WSS, for optimal pruning versus horizontal
cut on 6 real datasets

###  3.4 NCI60 classification Problem

The NCI60 dataset comes from cancer cell line microarrays, consisting of 6,830
gene expression measurements on 64 cancer cell lines. There are 14 different
cancer types in this dataset: ”BREAST”, ”CNS”, ”COLON”, ”K562A-repro”,
”K562B-repro”, ”LEUKEMIA”, ”MCF7A-repro”, ”MCF7D-repro”, ”MELANOMA”, ”NSCLC”,
”OVARIAN”, ”PROSTATE”, ”RENAL”, ”UNKNOWN”. First ignoring the cancer labels,
we built a dendrogram using average linkage, and performed optimal and
horizontal pruning methods, pruning down to 14 terminal nodes in each case.

Figure 11 shows the sequence of horizontal cuts. We label every 5 cuts with
transparent thin blue lines and the final cut with the solid thicker blue
line. Figure 12 shows the sequence of optimal cuts. We label cuts with
transparent blue numbers and the last cut with the bigger solid blue “46”.

Figure 11: The horizontal pruning steps resulting in 14 leaves for NCI-60
Figure 12: The optimal pruning steps resulting in 14 leaves for NCI-60

The final dendrogram produced by horizontal pruning is shown in figure 13 and
the optimal pruned one is shown in figure 14. The WSS for the horizontal
pruned tree was $2,544,265.78$ while the final WSS for the optimal pruned tree
was $915,484.12$, around only one-third of the value from the horizontal
counterpart. The median relative reduction in WSS (after removing NAs due to
skipping numbers of clusters) was 0.36. As we can see from Figure 15, optimal
cut dominates horizontal cut for all numbers of clusters and the discrepancy
is more significant when the number of clusters is between 5 and 40.

Figure 13: NCI-60: The horizontal pruned tree with 14 leaves  Figure 14:
NCI-60: The optimal pruned tree with 14 leaves.

(a) WSS comparison

(b) Relative reduction in WSS

Figure 15: Comparison of two pruning methods NCI-60

Finally, for both pruning methods, we focused on the resulting trees with 14
terminal nodes (the number of cancer classes). We used a majority vote of the
cell lines in each cluster to determine a predicted cancer label for that
cluster. We then compared these predictions to the true sample labels. The
confusion matrices are shown in Figure 16. The classification error rate was
0.41 for the horizontal pruned tree, and 0.28 for the optimal pruned tree. The
results suggest that our optimal pruning has yielded a more meaningful
grouping of the observations.

(a) Using horizontal pruning

(b) Using optimal pruning

Figure 16: Confusion matrices for classifying 14 cancer cell lines

##  4 Choosing the best-sized subtree

With the sequence of optimal subtrees in hand, we can estimate the best-sized
subtree among the sequence. This is equivalent to choosing the optimal value
of the cost-complexity parameter $\alpha$. In CART, this is done via cross-
validation, but this approach is not available in the unsupervised setting. So
instead we apply a method for choosing the optimal number of clusters. There
are many such methods, with two popular ones being the silhouette test
(Rousseeuw) [8] and the Gap test (Tibshirani et al.) [9]. For illustration, we
use the Gap test here.

The Gap test compares the changes in within-cluster dispersion to what would
be expected under a suitable reference null distribution using the output of
any clustering algorithm, such as hierarchical clustering. More specifically,
it compares the curve of log within-cluster dissimilarity against the number
of clusters to the curve achieved from data uniformly distributed throughout a
rectangular space. The optimal number of groups can be chosen based on many
criteria, here we chose the place where the gap between the two curves is most
significant.

We randomly generated 20 datasets $\mathit{X_{i}}$ with the number of
observations $n_{i}$ between 20 and 30, the number of features $p_{i}$ between
1 and 30, and the number of well-separated clusters equal to 4. We sampled
both $n_{i}$ and $p_{i}$ uniformly. Then we used the monte function in the R
package fungible [10] to generate the $x_{ij}$. The sizes of all clusters were
equal, except for the last cluster which adjusted for indivisibility. Then we
ran the Gap test with the optimal pruning method on the 20 datasets. If the
subtree with a certain number of leaves was not in the sequence of optimal
pruned subtrees, we increased the number of leaves to the closest available
number. The Gap test chose the correct number of clusters (4) 18 out of 20
times. Figure 17 shows an example that the optimal number of clusters was 4.
Of course, the Gap test (and other tests) will not perform as well when the
clusters are not well separated.

Figure 17: A Gap test example of 44 data points and 25 dimensions. The
estimated number of clusters is 4, the true number.

##  5 Discussion

In this paper, we have proposed the pruning of a dendrogram by a weakest-link
approach. This gives a sequence of the tightest clusters, for each fixed
number of groups. This sequence can be used in an exploratory way or can be
combined with a procedure such as the Gap test to estimate the best set of
clusters.

A public domain R language package pruneClust will soon be available on the
CRAN repository.

Acknowledgments. We thank Jacob Bien for helpful comments, and Tal Galili for
help with the dendextend package R package. RT is funded by the National
Institutes of Health (5R01 EB001988-16) and the National Science Foundation
(19 DMS1208164)

## References

  * [1] [clustering] how and where should you cut a dendrogram?  https://gmarti.gitlab.io/ml/2017/05/12/cut-a-dendrogram.html, author=Marti, Gautier, year=2017, month=May.
  * [2] Richard Bellman.  Dynamic Programming.  Dover Publications, 1957.
  * [3] Leo Breiman, Jerome Friedman, Charles J. Stone, and R. A. Olshen.  Classification and Regression Trees.  Taylor & Francis, Andover, England, UK, 1984.
  * [4] Gareth James, Daniela Witten, Trevor Hastie, and Robert Tibshirani.  An Introduction to Statistical Learning: with Applications in R.  Springer, 2013.
  * [5] Jesse H. Krijthe.  Rtsne: T-Distributed Stochastic Neighbor Embedding using Barnes-Hut Implementation, 2015.  R package version 0.16.
  * [6] Peter Langfelder, Bin Zhang, and Steve Horvath.  Defining clusters from a hierarchical cluster tree: the Dynamic Tree Cut package for R.  Bioinformatics, 24(5):719–720, 11 2007.
  * [7] R Core Team.  R: A Language and Environment for Statistical Computing.  R Foundation for Statistical Computing, Vienna, Austria, 2021.
  * [8] Peter J. Rousseeuw.  Silhouettes: A graphical aid to the interpretation and validation of cluster analysis.  Journal of Computational and Applied Mathematics, 20:53–65, 1987\.
  * [9] Robert Tibshirani, Guenther Walther, and Trevor Hastie.  Estimating the number of clusters in a data set via the gap statistic.  Journal of the Royal Statistical Society: Series B (Statistical Methodology), 63(2):411–423, 2001.
  * [10] Niels G. Waller.  fungible: Psychometric Functions from the Waller Lab., 2022.  version 2.2.1.


# Cattle Identification Using Muzzle Images and Deep Learning Techniques

Geoffrey Kimani
CyLab-Africa / Upanzi Network
Kigali
<EMAIL_ADDRESS>
Oluwadara Adedeji
CyLab-Africa / Upanzi Network
Kigali
<EMAIL_ADDRESS>
Patrick Fashingabo
CyLab-Africa / Upanzi Network
Kigali
<EMAIL_ADDRESS>
Moise Busogi
Carnegie Mellon University Africa
<EMAIL_ADDRESS>
Edith Luhanga
Carnegie Mellon University Africa
<EMAIL_ADDRESS>
Karen Sowon
Carnegie Mellon University
<EMAIL_ADDRESS>
Lenah Chacha
CyLab-Africa / Upanzi Network
<EMAIL_ADDRESS>

###### Abstract

Traditional animal identification methods such as ear-tagging, ear notching,
and branding have been effective but pose risks to the animal and have
scalability issues. Electrical methods offer better tracking and monitoring
but require specialized equipment and are susceptible to attacks. Biometric
identification using time-immutable dermatoglyphic features such as muzzle
prints and iris patterns is a promising solution. This project explores cattle
identification using 4923 muzzle images collected from 268 beef cattle. Two
deep learning classification models are implemented - wide ResNet50 and
VGG16_BN and image compression is done to lower the image quality and adapt
the models to work for the African context. From the experiments run, a
maximum accuracy of 99.5% is achieved while using the wide ResNet50 model with
a compression retaining 25% of the original image. From the study, it is noted
that the time required by the models to train and converge as well as
recognition time are dependent on the machine used to run the model.

Code files 111Code files: https://github.com/peter716/Animal_Biometrics_System
and data files. 222Data files: https://zenodo.org/record/6324361#.ZAsiTNJBzb3

_K_ eywords Animal Biometrics $\cdot$ Deep Learning $\cdot$ Image compression
$\cdot$ Muzzle images

##  1 Introduction

###  1.1 Background

In the realm of animal identification, classical and electrical methods have
long been employed, but they possess certain shortcomings. Classical methods
involve permanent or temporal identification techniques. Permanent methods,
such as ear notching, freeze branding, and hot-iron branding, leave permanent
marks on animals but are often invasive and cause pain. Temporal methods, like
ear-tagging, are less painful but susceptible to fraud and manipulation. These
methods require manual work, animal restriction, and can be unreliable in
preventing fraud or duplication [1, 2, 3, 4] .

Electrical methods, such as ear tags, injectable transponders, and rumen
boluses, rely on electronic devices to identify and track cattle. While these
methods are less invasive and easy to apply, they can cause discomfort to the
animals and require additional equipment, maintenance, and calibration.
Moreover, they can be susceptible to fraud and attacks [1, 2, 4].

To overcome the limitations of classical and electrical methods, biometric
techniques have gained prominence in animal identification. Biometric
modalities utilize unique physiological or behavioral features to assign a
distinct identity to each animal [1]. For bovine animals, these features can
include retinal vascular patterns, iris patterns, or muzzle prints.

Retinal vascular patterns are based on the unique blood vessel patterns in the
retina. They provide accurate identification, are relatively unaffected by
corneal injuries, and can be used in various animals. However, certain
inherited retinal diseases and congenital ocular disorders may impact the
recognition process. Acquiring high-quality retinal images using techniques
like optical coherence tomography (OCT) allows the extraction of features such
as vessel diameter, branching patterns, and curvature. These features are then
used to create a unique template for matching and identification [1, 5, 6].

Iris patterns, found in the colored part of the eye, are another biometric
feature. They are determined by the distribution of pigments in the iris,
resulting in a stable and distinctive identifier. However, iris patterns can
be affected by diseases and aging, which may reduce accuracy. Extracting iris
texture, color, and intensity using techniques like Gabor filters enables the
creation of unique templates for recognition [1, 7, 8, 9].

Muzzle prints, analogous to human fingerprints, offer an accurate and time-
immutable identifier. High-quality images of the muzzle can be used to extract
unique patterns of beads and ridges. While the effect of diseases on muzzle
prints is not extensively studied, conditions altering the textual structure
of the muzzle may affect identification. Templates generated from muzzle
prints are compared to verify an animal’s identity [1].

Various approaches have been proposed to implement biometric identification
systems. Some employ deep learning techniques, such as convolutional neural
networks (CNNs) or deep belief networks (DBNs), to process muzzle images.
These models achieve high accuracies, but they often require large amounts of
data and computational resources. Other approaches utilize template matching
by creating databases of unique templates extracted from biometric features.
These templates are compared with incoming data to determine a match or
similarity score, enabling identification [4, 10, 11, 12, 13].

Biometric methods offer non-invasive, reliable, and accurate means of
identifying animals. Techniques such as retinal vascular patterns, iris
patterns, and muzzle prints provide distinct features for identification. Deep
learning models and template matching approaches have shown promising results
in animal identification systems, offering improved accuracy and non-
invasiveness compared to classical and electrical methods.

###  1.2 Problem statement

Animal identification has been a problem long addressed in various ways with
the aim of accurately recognizing individual animals. The identification
techniques have been broadly categorized into four classes: permanent methods,
temporary methods, electrical methods, and animal biometrics [1]. In the
present study, permanent and temporary methods are amalgamated into a single
category termed classical methods. This category encompasses practices such as
ear-tagging, ear notching, ear tattooing, hot iron branding, and freeze
branding. However, these methods have serious shortcomings, including fraud,
infections, animal distress, equipment loss (ear tags), and limited
scalability [1, 2]. Electrical methods, characterized by rumen boluses, Radio
Frequency Identification (RFID) ear tags, and injectable glass tags [1], offer
advantages over classical methods by addressing tracking and animal monitoring
beyond identification requirements. However, they come with high installation
costs due to specialized equipment and personnel. Moreover, electrical methods
are vulnerable to Denial of Service (DoS) attacks, tag modification, and
system spoofing [1, 2, 4]. Biometric identification, leveraging time-immutable
dermatoglyphic features, offers harmless animal identification. These methods
provide automated authentication, high security, reliability, and zero memory
load on humans [1]. Muzzle prints, retinal vascular patterns, and iris
patterns are among the various techniques used for biometric identification.

In the current study, biometric systems are examined for the purpose of cattle
identification, specifically focusing on the unique identification of cattle
through muzzle images. Previous research [1] has explored three biometric
identifiers for cattle: muzzle images, iris patterns, and retinal vascular
patterns. The dataset employed in this study is identical to that used by
[13], comprising 4923 muzzle images from 268 beef cattle 333Beef Cattle
Muzzle/Noseprint database for individual identification
https://zenodo.org/record/6324361#.ZAsiTNJBzb3. The study further narrows its
focus to applications in resource-constrained environments where high-quality
digital cameras may not be readily accessible.

###  1.3 Use cases discussion

Animal identification is very important in establishing ownership, preventing
fraud, and identifying theft. In this section, various use cases are presented
that could benefit from the implementation of a biometric identification
system. These use cases are particularly challenging to address using
traditional classical and electrical methods of identification. There are
various applications where such a system could be useful including disease
tracking, insurance claims, and collateral for financial loans among others.
With the biometric systems being able to accurately identify individual
cattle, banks, and insurance companies could find this of interest to them
since they would be able to establish ownership and proof of claims for
insurance or collateral for bank loans. This would suppress insurance fraud
and allow farmers to use their cattle as security for bank loans. Government
and research institutions dealing with animal production and animal health
might find the system useful as well. A fully connected system with a front-
end dashboard for user interaction and a database can be built on top of the
proposed deep-learning model to keep track of vaccination programs as well as
track animal disease trajectories. This information might also be useful to
public health officials who need to certify meat that is sold for consumption.
Security agencies might also find this system useful in cases of identifying
missing or swapped cattle and establishing ownership in cases of dispute. It
should be noted that the biometric system alone may not adequately address all
challenges inherent to a specific use case. As such, it may require
integration with other systems to achieve full utility. For example, in
disease tracking, the system would ideally rely on an existing system that
keeps track of ownership, region, vaccinations administered, animal movements
and interactions(electrical methods are useful in this case), etc. whether it
is an independent system or extension of the biometric system.

##  2 Methods

###  2.1 Dataset

####  2.1.1 Data collection

The dataset employed in the current study was sourced from previous research
[13], comprising 4923 muzzle images from 268 beef cattle. The dataset is
available at 3 under Creative Commons license 4.0 and is provided as a
directory of cattle images for each of the 268 cattle. Each cattle had at
least 4 images representing it while some cattle had up to 70 images. These
images were obtained from cattle in their natural environment while grazing
under natural lighting conditions and a high-quality mirrorless digital camera
with a 70–300 mm F4-5.6 focal lens. The muzzle areas in the images were first
horizontally aligned and then cropped out to focus the muzzle area resulting
in images uniformly resized to 300 by 300 pixels. It is important to note that
though a cursory inspection of the muzzle prints may show that the color and
texture of the prints look similar in different animals the patterns of beads,
grooves, and ridges are different upon a comprehensive analysis.

####  2.1.2 Data splitting

For data preparation, the dataset was partitioned into 70% for training and
30% for evaluation. To guarantee the representation of each individual cattle
in both the training and testing sets, an iteration was performed over the
image folders, selecting a train and test set for each individual animal. This
approach ensures that even the least represented animals are included in both
the training and testing sets.

####  2.1.3 Data augmentation

In this section, various data augmentation techniques were applied, including
image flipping, color jittering, blurring, and rotations. Images were resized
to dimensions of 300 by 300 pixels. Random horizontal flips were introduced to
augment the dataset by generating additional training examples with varied
orientations, thereby enhancing the model’s generalization performance. Color
jittering was employed to adjust the brightness, contrast, saturation, and hue
of the images, further contributing to the diversity of the training examples.
Gaussian blur was applied to smooth out noise and pixelation, aiding in
feature extraction and classification. Data augmentation was conducted using
PyTorch modules.

###  2.2 Modeling

In the present study, two models were selected from the 59 deep-learning
models proposed in the baseline research. The VGG16_BN model was chosen due to
its superior performance in the baseline paper, and wide ResNet50 was selected
as an alternative, representing a newer architecture. A transfer learning
approach was employed, utilizing models pre-trained on ImageNet. These models
are readily available through
PyTorch444https://pytorch.org/vision/stable/models.html.

####  2.2.1 Architecture Discussion

The VGG16 model is a convolutional neural network architecture that was
developed by the Visual Geometry Group at the University of Oxford [14]. The
model has 16 layers, consisting of 13 convolutional layers and 3 fully
connected layers, and is known for its simplicity and effectiveness in image
classification tasks. The VGG16 model architecture is characterized by a
series of convolutional layers, and each convolutional layer applies multiple
3x3 convolutional kernels to the input feature maps to produce a set of output
feature maps. The number of kernels in each convolutional layer increases as
the network goes deeper, which allows the model to learn increasingly complex
and abstract features. The max-pooling layers are used to reduce the spatial
dimension of the feature maps. The last few layers of the network are fully
connected layers that map the extracted features to a final set of output
classes. In PyTorch, the VGG16 model can be imported and used as a pre-trained
model for transfer learning, as described in the documentation 555
https://pytorch.org/vision/stable/models/generated/torchvision.models.vgg16.html.
The pre-trained VGG16 model in PyTorch has been trained on the large-scale
ImageNet dataset, which contains over 1 million images and 1000 classes [15].
This pre-training allows the model to learn a set of general features that can
be transferred to other image classification tasks with relatively little
additional training. This is particularly useful when working with smaller
datasets, as transfer learning can help to improve the performance of the
model and reduce overfitting.

The Wide ResNet architecture is based on the ResNet architecture, which is a
deep convolutional neural network that uses residual connections to improve
the performance of the network. It is similar to the ResNet architecture but
uses wider layers with more filters to improve the accuracy of the network.
The model starts with an input layer followed by a convolution layer with 64
filters and a kernel size of 7x7 which are then followed by batch
normalization and ReLU activation layers finally a max pooling layer. This
section is then followed by at least four residual blocks each containing
multiple convolution layers with 128, 256 or 512 filters. A global average
pooling layer follows to reduce the spatial dimensions of the feature maps and
finally followed by a fully connected layer and a softmax activation layer
[16]. In PyTorch, the pre-trained Wide ResNet model is available as part of
the torchvision library 666
https://pytorch.org/vision/stable/models/generated/torchvision.models.wide_resnet50_2.html#torchvision.models.wide_resnet50_2
which provides a number of pre-trained models for use in computer vision
tasks. Figure 1 shows the architectures side by side

Figure 1: Pretrained Models Architecture
(Source: Deep Transfer Learning for Land Use and Land Cover Classification: A
Comparative Study, Raoof Naushad et. al.)

In the current study, modifications were made to the architecture to include a
custom fully connected layer tailored to the specific number of classes, in
this case, 268. A linear layer was added to the pre-trained model, accepting
$n$ inputs—where $n$ is the number of classes output by the model, 1000 for
IMAGENET1K_V1—and outputting 256 classes. This was followed by a ReLU
activation layer to aid the model in learning the decision boundaries for
different classes, as well as a dropout layer to mitigate overfitting.
Subsequently, another linear layer was introduced, taking the 256 classes and
outputting the 268 classes specific to this study. The architecture concludes
with a softmax activation layer, providing a probability distribution over the
predicted classes. Figure 2 shows our new architecture skeleton.

Figure 2: Modified Architecture

####  2.2.2 Image compression

A primary objective of this study was to enhance model robustness in the
context of low-quality images. This consideration is particularly relevant for
potential applications in Africa, where acquiring high-quality images akin to
those used in the baseline paper may be challenging. Furthermore, in real-
world implementations, lower-grade cameras may be employed as opposed to high-
quality ones, thereby affecting image quality. Consequently, the study sought
to explore techniques that would facilitate experimentation with low-quality
images, given the unavailability of a secondary dataset for this purpose.

The study explored the use of image compression techniques. Compression in the
context of information theory is the process of reducing the size of a data
file by removing redundant or unnecessary information. This is an important
process as it reduces the amount of information being transmitted. There are
two main schemes in compression; lossy and lossless compression. Each scheme
is determined by the data rate against the distortion level as shown in figure
3. The data rate refers to the amount of data needed to represent the
compressed image usually measured in bits per second (bps) while distortion
refers to the difference between the original and resulting data. On the
x-axis is the distortion level and the y-axis is the data rate. If you have a
high data rate and zero distortion, it is a lossless compression scheme,
conversely, if you have a lower data rate in the resulting data than in the
original and a distortion greater than 0 you have a lossy compression scheme.
Our interest was the lossy compression scheme that allows us to lose some
image information consequently losing some quality in the images to be
compressed.

Figure 3: Compression scheme theory
(Source: BITMOVIN https://bitmovin.com/lossy-compression-algorithms/)

Lossy image compression has two steps - quantization and transform coding. The
quantization step involves dividing the image into small, fixed-size blocks of
pixels, and then mapping the color values of each block to a limited range of
possible values. The mapping process discards some of the original color
information, resulting in a loss of detail, but allows for smaller file sizes
[16]. In the coding step, the quantized data is compressed using algorithms
that reduce the amount of information needed to represent the image. A common
technique is to use transform coding, in which a mathematical transform is
applied to the image data to reduce redundancy before encoding it [16]. The
information lost in this step is permanently lost and cannot be recovered. A
good example of a lossy algorithm popular in image compression is the discrete
cosine transform(DCT) algorithm.

The DCT algorithm transforms the image data from the spatial domain (pixel
values) to the frequency domain (frequency components that make up the image).
In the frequency domain, the DCT separates the image data into a set of
coefficients that describe how much different frequency components contribute
to the image. The DCT algorithm is very efficient and is used in JPEG(Joint
Photographic Experts Group) compression, a popular industry-grade image
compression scheme.

In addressing the issue, the study employed image compression as a strategy.
While a custom implementation of the Discrete Cosine Transform (DCT) algorithm
was initially considered, it proved to be time-inefficient. As a result, the
study opted for the use of existing, optimized libraries. Specifically, the
Pillow library was utilized for its image compression capabilities. This
library allows for the specification of image quality during the save
operation. For instance, specifying a quality of 50 results in an image
retaining 50% of its original quality. The library also offers an optimize
parameter for further file size reduction. Pillow is a fork of the Python
Image Library (PIL). In this study, images were compressed to retain three
different quality levels—100%, 50%, and 25%—and the model’s performance was
evaluated at each of these quality levels. Figure 4 shows sample images
compressed at 50% of the original quality. The images on the top row are the
original while the bottom are the compressed ones. A visual cursory
examination of the images may not show huge differences but the crispness,
sharpness, and color differentiation of object edges within the image are
reduced (possibly not apparent to the eye) by compression.

Figure 4: Original and Compressed Images Sample

##  3 Results

A series of experiments were conducted to assess the performance of the
selected models. It was observed that the wide ResNet50 model achieved the
highest accuracy, approximately 99.4%, outperforming the baseline’s highest
accuracy of 98.7% achieved with the VGG16_BN model. Early stopping mechanisms
were implemented to facilitate quicker model convergence. Notably, both
selected models converged in fewer than 40 epochs, which is more efficient
than the 50 epochs reported in the baseline study.

It is noted that the runtime of the models is dependent on the machine spec as
seen in tables 3 and 3. These experiments were run on two different machines
both provided under the Google Colab platform. Results in 3 were obtained
using the free account which comes with an Intel Xeon CPU @2.20 GHz, 13 GB
RAM, Tesla K80 accelerator, and 12 GB GDDR5 VRAM. Results in 3 were obtained
using a paid version of colab featuring more powerful specs Intel(R) Xeon(R)
CPU @ 2.30GHz Tesla P100-PCIE-16GB accelerator and 12 GB GDDR5 VRAM with a
high-RAM option available. The difference in these machines explains the
difference in the total time taken as well as time per epoch between the two
experiments. With more powerful machines the time per epoch and total time
decreases considerably. For example in a separate experiment using an NVIDIA
A100-SXM4-40GB machine, it took 25 seconds per epoch with a total of 30 epochs
run bringing the total time to around 12.5 minutes.

The performance of the models actually improved after the lossy compression
but it was noted that the epochs required to converge are more. Upon examining
Table 3, it is observed that when retaining 100% of the original image
quality, the model converges after 28 epochs, achieving an accuracy rate of
99.2%. When 50% of the image quality is retained, it takes 8 more epochs to
converge but a higher accuracy of 99.4% is attained while retaining 25% of the
original quality the model converges after 30 epochs with an accuracy of
99.5%. The same is true for the VGG16_BN model as a similar trend is observed.
Given that the images are of lower quality, the models take longer to learn
the features and this could explain the additional epochs.

In terms of accuracy measurement, both the initial accuracy achieved after the
first epoch and the final accuracy upon model convergence are reported. These
values are presented in the accuracy column of the respective tables, with the
top value indicating the starting accuracy and the bottom value representing
the final accuracy. The VGG16_BN model starts with a very low accuracy but
also attains a lower accuracy than the wide ResNet50 model. While quality is
reduced, the images are quite large so the models are still able to extract
useful features.

The time taken by each of the models was also documented.Experiments were
conducted on two different platforms: Colab and Colab Pro. While the exact
machine specifications are not available due to Colab’s dynamic assignment, it
is noted that Colab Pro generally offers more powerful machines. Wide ResNet50
required a longer total training time compared to VGG16_BN, which was run on
Colab Pro. Interestingly, the total time taken was nearly similar across
different levels of image compression. However, higher compression levels did
result in longer training times. Regarding time per epoch, compressed images
led to shorter epoch times for both models. It was also observed that the more
powerful machine allocated for VGG16_BN on Colab Pro resulted in significantly
shorter epoch times compared to those for Wide ResNet50.

Tables 3 and 3 summarize the results from the two models comparing the
performance based on image quality retained.

Experiment |  File sizes |  Percent of the Original Quality |  Accuracy |  Epochs |  |
---|---|---|---|---|---|---
(Early stopping) |  Total Time |  Time per epoch (mins) |  |  |  |
Normal run |  622MB (4923) |  100% |  9% $\sim$99.2% |  28 |  1 hr 12 min |  $\sim$3.2
Pill JPEG Compression |  76MB (4923) |  50% |  7% $\sim$99.4% |  36 |  1 hr 29 min |  $\sim$2.4
Pill JPEG Compression |  45MB (4923) |  25% |  7% $\sim$99.5% |  30 |  1hr 10 min |  $\sim$2.3

Table 1: Wide ResNet50 Results

Experiment |  File sizes |  Percent of the Original Quality |  Accuracy |  Epochs |  |
---|---|---|---|---|---|---
(Early stopping) |  Total Time |  Time per epoch (mins) |  |  |  |
Normal run |  622MB (4923) |  100% |  0.7% $\sim$97.7% |  38 |  53 min |  $\sim$1.4
Pill JPEG Compression |  76MB (4923) |  50% |  0.6% $\sim$97.97% |  39 |  53 min |  $\sim$1.38
Pill JPEG Compression |  45MB (4923) |  25% |  0.4% $\sim$98.04% |  43 |  58 min |  $\sim$1.35

Table 2: VGG16_BN Results

##  4 Conclusion

###  4.1 Significance

In the present study, a biometric system for cattle identification was
developed, utilizing VGG16 and Wide ResNet50 models. Notably, the Wide
ResNet50 model achieved a peak accuracy of 99.5% following image compression.
These results are significant for application in a context where high-quality
images cannot be obtained owing to high equipment cost as well as field
operation costs. In the African context, for example, most farmers are small-
holder farmers owning a few tens of cattle as well as herding nomadic
communities owning hundreds of cattle. However, it is hard for the farmers or
herders to acquire high-quality camera in order to use a biometric
identification system, additionally, while agricultural extension officers can
obtain a high-definition camera for a larger community, using such cameras
require expertise and need careful handling especially when using them in the
field. This project shows that even low-quality images can be used to a high
degree of accuracy for the biometric identification system therefore requiring
no specialized equipment.

A significant limitation of this approach is scalability. Adding new cattle
would require retraining the model in order to learn the new cattle classes,
if the number of cattle is many, this could be a significant effort. There are
two alternatives to scale the model; 1) scale-up and 2) scale-out. The scale-
up approach entails having a single model with all the individual cattle
enrolled and as seen, this requires retraining the model in cases of new
cattle classes. An alternative strategy proposed in this study involves a
scale-out approach, utilizing an ensemble of models to accommodate a larger
number of cattle classes. In this configuration, multiple models operate in
unison to perform identification tasks. Ensemble models are a common industry
practice and that approach can be used to address the scalability issues
arising from this methodology.

###  4.2 Future Work

To improve on this work, more models can be considered as implemented in the
baseline paper. In addition, an interface can be implemented, such that muzzle
images can be uploaded and inference can be made to identify the cattle. This
can also be interlaid with a database system that contains more information
about the cattle such as the vaccine and disease history, and more. This
interface will facilitate a use case such as disease tracking. Some issues
related to scalability were also raised during the research. To scale the
work, as suggested, an ensemble of models can be used, such that retraining is
done on just the newly added cattle classes.

## Acknowledgments

This work was done in collaboration with CyLab Africa/Upanzi Network and
Carnegie Mellon University Africa.

## References

  * [1] Ali Ismail Awad.  From classical methods to animal biometrics: A review on cattle identification and tracking.  Computers and Electronics in Agriculture, 123:423–435, April 2016\.
  * [2] Santosh Kumar, Sanjay Kumar Singh, Ravi Shankar Singh, Amit Kumar Singh, and Shrikant Tiwari.  Real-time recognition of cattle using animal biometrics.  Journal of Real-Time Image Processing, 13(3):505–526, September 2017.
  * [3] Michael Neary and Ann Yager.  Methods of Livestock Identification.  https://www.extension.purdue.edu/extmedia/as/as-556-w.pdf, 2002\.
  * [4] Santosh Kumar and Sanjay Singh.  Cattle Recognition: A New Frontier in Visual Animal Biometrics Research.  Proceedings of the National Academy of Sciences, India Section A: Physical Sciences, 90, May 2019.
  * [5] Michalis Georgiou, Kaoru Fujinami, and Michel Michaelides.  Retinal imaging in inherited retinal diseases.  Annals of eye science, 5:25, September 2020.
  * [6] Stewart Muchuchuti and Serestina Viriri.  Retinal Disease Detection Using Deep Learning Techniques: A Comprehensive Review.  Journal of Imaging, 9(4):84, April 2023.  Number: 4 Publisher: Multidisciplinary Digital Publishing Institute.
  * [7] Shengnan Sun, Shicai Yang, and Lindu Zhao.  Noncooperative bovine iris recognition via SIFT.  Neurocomputing, 120:310–317, November 2013.
  * [8] Tariq Mehmood Aslam, Shi Zhuan Tan, and Baljean Dhillon.  Iris recognition in the presence of ocular disease.  Journal of the Royal Society, Interface, 6(34):489–493, May 2009\.
  * [9] Mateusz Trokielewicz, Adam Czajka, and Piotr Maciejewicz.  Iris Recognition Under Biologically Troublesome Conditions \- Effects of Aging, Diseases and Post-mortem Changes.  In Proceedings of the 10th International Joint Conference on Biomedical Engineering Systems and Technologies, pages 253–258, 2017\.  arXiv:1809.00182 [cs].
  * [10] Yongliang Qiao, Daobilige Su, He Kong, Salah Sukkarieh, Sabrina Lomax, and Cameron Clark.  Individual Cattle Identification Using a Deep Learning Based Framework**The authors acknowledge the support of the Meat Livestock Australia Donor Company through the project: Objective, robust, real-time animal welfare measures for the Australian red meat industry.  IFAC-PapersOnLine, 52(30):318–323, January 2019.
  * [11] Santosh Kumar, Amit Pandey, K. Sai Ram Satwik, Sunil Kumar, Sanjay Kumar Singh, Amit Kumar Singh, and Anand Mohan.  Deep learning framework for recognition of cattle using muzzle point image pattern.  Measurement, 116:1–17, February 2018.
  * [12] Rotimi-Williams Bello, Abdullah Talib, and Ahmad Sufril Azlan Mohamed.  Deep Learning-Based Architectures for Recognition of Cow Using Cow Nose Image Pattern.  Gazi University Journal of Science, September 2020.
  * [13] Guoming Li, Galen E. Erickson, and Yijie Xiong.  Individual Beef Cattle Identification Using Muzzle Images and Deep Learning Techniques.  Animals : an Open Access Journal from MDPI, 12(11):1453, June 2022\.
  * [14] Karen Simonyan and Andrew Zisserman.  Very Deep Convolutional Networks for Large-Scale Image Recognition, April 2015.  arXiv:1409.1556 [cs].
  * [15] Olga Russakovsky, Jia Deng, Hao Su, Jonathan Krause, Sanjeev Satheesh, Sean Ma, Zhiheng Huang, Andrej Karpathy, Aditya Khosla, Michael Bernstein, Alexander C. Berg, and Li Fei-Fei.  ImageNet Large Scale Visual Recognition Challenge.  International Journal of Computer Vision, 115(3):211–252, December 2015.
  * [16] Lossy Compression Algorithms: Everything a Developer Needs to Know \- Bitmovin, March 2020.  Section: Developers.