The ALLELEMATCH Package for R: Supplementary Documentation with Tutorials

P. Galpern, M. Manseau, P. Hettinga, K. Smith, and P. Wilson

2026-07-23

1 Introduction

This R vignette can also be found as a PDF document supplement to the publication describing the allelematch package. Here, we describe the operation of the package in more detail and illustrate its use in a tutorial format.

The allelematch package was developed to identify unique genotype profiles in situations where there are likely to be multiple samples of each individual present. Such conditions arise in wildlife conservation genetics where non-invasive sampling of feather, hair, and fecal material can produce large datasets containing an unknown number of animals.

Finding unique individuals appears on the surface to be a trivial exercise; a matter of sorting genotypes into identical groups. When there is unlikely to be error in profiling and missing data is absent, identity analysis is straightforward and software for this purpose has long been available (Kalinowski et al., 2007; Peakall & Smouse, 2006; Wilberg & Dreher, 2004). The task becomes much more complicated if data sets contain missing information at some loci, or are subject to even a modest amount of genotyping error. In these cases samples from the same individual may differ because of allele dropouts or false alleles, or they may be identical except for some missing information. These issues make it ambiguous which samples belong to which individual, or how many individuals are present.

This is the niche for allelematch: applications of genotype profiling where identity must be established in suboptimal circumstances. In conservation genetics these conditions are common. Collecting high quality samples can often be challenging, and as data sets grow larger in size, quality control may be uneven. allelematch may also be useful in laboratory quality control, i.e., to visually highlight genotyping errors in data sets where the same individual is profiled multiple times for confirmation.

2 Approach

allelematch operates by comparing rows of a multilocus genotype matrix, where rows are samples, and columns for each locus give the names of the alleles that are present. For diploid codominant data, there are typically two columns for each locus. It finds the similarities between the samples using a metric which is a form of the Hamming distance (Hamming, 1950) modified to account for missing data. It then uses hierarchical clustering and a dynamic method for identifying clusters on a dendrogram. This step is done using the Dynamic Tree Cut package for R (Langfelder et al., 2008) which defines groups of similar samples at a specified threshold of dissimilarity.

2.1 Algorithm

Specifically, amUnique(), the workhorse function of the package, performs the following operations:

  1. Determine a similarity score for each pair of samples. Let \(A\) be a matrix \([a_{ij}]_{m \times n}\) with elements that are allele names or missing data, where \(m\) is the number of samples, and \(n\) is the total number of columns containing the allele calls. Find the similarity score,

    \[s_{p,q} = \frac{N_{\text{match}_{p,q}}}{n} + \frac{N_{\text{missing}_{p,q}}}{2n}, \forall p,q \in \{1, 2, \dots, m\} \quad (1)\]

    where \(N_{\text{match}_{p,q}}\) is the number of elements (alleles or missing data) that match in rows \(a_{p,*}\) and \(a_{q,*}\), and \(N_{\text{missing}_{p,q}}\) is the number of elements that are missing in either row \(a_{p,*}\) or row \(a_{q,*}\). A score of \(s_{p,q} = 1\) means that the elements of both rows are identical. Scores less than one have been penalized \(\frac{1}{n}\) for each mismatching element, and penalized a smaller amount, \(\frac{1}{2n}\), where an element is missing in one row but not missing in the other row.

  2. Produce a symmetric dissimilarity matrix from the similarity scores, \(s_{p,q}\).

\[D_{m,m} = 1 - \begin{bmatrix} 1 & s_{1,2} & \dots & s_{1,m} \\ s_{2,1} & 1 & \dots & s_{2,m} \\ \vdots & \vdots & \ddots & \vdots \\ s_{m,1} & s_{m,2} & \dots & 1 \end{bmatrix}\]

  1. Perform agglomerative hierarchical clustering on the dissimilarity matrix, \(D\), using the complete linkage method.

  2. Use a hybrid dynamic tree cutting method (Langfelder et al., 2008) to identify clusters, where a cut height parameter \(\hat{d}\) gives the maximum dissimilarity between any two samples in a cluster. This parameter is related to the similarity score as \(\hat{d} = 1 - \hat{s}\).

  3. Find the consensus sample for each cluster, if a cluster contains more than one sample at this \(\hat{d}\). allelematch provides two supported options1 for declaring consensus samples:

    1. The consensus is the sample that is most similar to other samples in the cluster (calculated using Equation 1).

    2. The consensus is the sample with the least amount of missing data.

  4. Declare as singletons samples which are the only member of their cluster.

  5. Produce a new matrix \(A' = [a'_{ij}]_{m' \times n}\) consisting of the consensus samples if they exist, and the singleton samples.

  6. IF \(A'\) contains only singletons continue to step 9. ELSE let \(A = A'\) and return to step 1.

  7. Declare as unique all rows of the matrix \(A'\).

  8. Compare unique matrix \(A'\) to all original samples. For each row \(a'_{q,*}\) of the unique matrix \(A'\), find all rows \(a_{p,*}\) of the original matrix \(A\) that have a similarity score \(s_{p,q} \ge \hat{s}\). For example,

\[s_{p,q} = \frac{N_{\text{match}_{p,q}}}{n} + \frac{N_{\text{missing}_{p,q}}}{2n}, \forall p \in \{1, 2, \dots, m\}, q \in \{1, 2, \dots, m'\} \quad (2)\]

  1. Declare as matches rows \(a'_{q,*}\) that have a score \(s_{p,q} \ge \hat{s}\) for exactly one row \(a_{p,*}\).

  2. Declare as multiple matches rows \(a'_{q,*}\) that have a score \(s_{p,q} \ge \hat{s}\) for two or more rows \(a_{p,*}\).

  3. Declare as unclassified rows \(a'_{q,*}\) that do not have a score \(s_{p,q} \ge \hat{s}\).

2.2 Notes

3 Tutorial

The tutorial is organized around example data sets. Five examples review the use of allelematch for identifying samples representing unique individuals in datasets where these individuals may be sampled multiple times.

The tutorials provide sample code, and the graphical output from R has been inserted as figures. However, the most important output from allelematch is in the form of HTML documents. These illustrate in colour how samples match and mismatch. Note that within one example multiple HTML files may be generated and these are numbered sequentially (e.g. example1_1.html, example1_2.html, etc.).

It is also possible to reproduce these examples directly in R. All example data sets are provided with the package and the tutorials demonstrate how to access this data.

3.1 Quantifying Missing Data

In datasets where missing indices are prevalent or individual loci suffer from a high drop-out rate, allelematch calculations may approach mathematical stability limits. The first suggested step in the workflow is to pre-screen your dataset for missing allele calls.

Use amPreCheck() to generate a comprehensive summary of missing data across loci and samples, alongside the minimum number of shared loci discovered between any pair of individuals. If a dataset features samples that share zero overlapping loci, amPreCheck() will issue a critical warning to help you identify these constraints before they trigger runtime errors during clustering.

Users are strongly encouraged to enforce quality control filtering before executing full allelematch optimization routines. For example, a reliable threshold for stable clustering is ensuring that at least 75% of the total locus panel is called for every individual sample. While lower thresholds can sometimes be tolerated, severe data sparsity may introduce instability or bias into downstream genetic analyses. We will demonstrate the utility of this function in each of the three examples below.

3.2 Identifying unique individuals

Example 1: High quality data set

The data in this example are simulated2 to represent a high quality data set that might result from a laboratory protocol where samples were run multiple times to confirm their identity. It has no genotyping error, a near-zero missing data load, and approximately 60% of the individuals have been artificially resampled more than once.

In typical usage data would be imported into R in the most convenient manner. This may be the read.csv() function for example. In this and the following examples we load the data supplied with the package using the data() function.

library(allelematch)

data(amExample1)

Next, we create an amDataset object to prepare the data for use with other allelematch functions. The sample number serves as the index, missing data is coded as “-99”, and we explicitly exclude the knownIndividual column. We then pre-screen the data using amPreCheck().

example1 <- amDataset(
  amExample1, 
  indexColumn = "sampleId",
  ignoreColumn = "knownIndividual", 
  missingCode = "-99"
)

amPreCheck(example1)

The report for this dataset indicates high data quality: no single locus exceeds 5% missing data, the average individual missing data load is only 0.5%, and the individual with the highest missing data load is at 10%. Furthermore, every pair of samples shares at least 9 out of 10 loci. Given these metrics, the function provides a safe verdict, indicating the dataset is well-suited for clustering.

The next step is to determine the optimal dissimilarity criterion for identifying unique individuals. This criterion is best understood in terms of the number of allele mismatches allowed (alleleMismatch parameter; or \(\hat{m}\) (see Section 2.2 Notes). This routine evaluates the number of unique individuals identified across a range of parameter values to suggest an optimal threshold.

amUniqueProfile(example1, doPlot = TRUE)

The output profile (Figure 1) shows that an optimum for this data is allowing 2 alleles to mismatch (alleleMismatch = 2). Note how multipleMatch, the number of samples that match more than one unique individual, is zero at the optimum. This tells us that every sample has been declared as a unique individual or as a match to a unique individual.

Figure 1: High quality dataset.
Figure 1: High quality dataset.

This is described as a ZeroSecondMinimum profile, indicating that multipleMatch reaches zero at a parameter setting other than alleleMismatch = 0. Simulations indicate that with this profile morphology the software can independently determine the alleleMismatch parameter required to identify the correct number of unique genotypes (\(\pm 3\%\)).3

The final step is to use the alleleMismatch criterion to identify which samples are unique and which samples are matches.

uniqueExample1 <- amUnique(example1, alleleMismatch = 2)

This should be reviewed carefully before proceeding. The results can be viewed in an HTML format as follows.

## Save to disk
summary(uniqueExample1, html = "example1_1.html")

## View in default browser
summary(uniqueExample1, html = TRUE)

In this example, 12 unique genotypes were identified among 20 samples. The alleleMismatch parameter we set at an earlier step allows up to two alleles to be different for samples to be declared identical. Examining the HTML output, there are no cases where mismatches occurred. The mismatching alleles would be highlighted in red. There is one case with missing data at one locus (unique genotype 12), and pink highlighting is used to indicate missing data. Recall how a missing allele in one genotype but not the other penalizes the score by \(\frac{1}{2n}\) (by Equation 1). For this reason sample 19 (also unique genotype 12) matches itself with a score of 1, but matches sample 20 with a score of 0.95.

The same results can also be saved in a CSV spreadsheet format. This format lacks some of the supplementary information as well as the colour formatting that can be helpful for reviewing the results.

summary(uniqueExample1, csv = "example1_1.csv")

An analysis-ready dataset can also be prepared from the CSV output. We encourage caution with this step. Typically the output of the analysis must be reviewed in detail before an analysis-ready data set can be used uncritically. For example, checking that \(P_{\text{sib}}\) probabilities are all below a threshold value is a minimum requirement as allelematch does not use this as a criterion.

summary(uniqueExample1, csv = "example1_2.csv", uniqueOnly = TRUE)

When the multiple sampling was simulated in this data set, we kept track of which individual each sample came from in the knownIndividual column. For typical non-invasive sampling applications, however, such information would not be known. We will repeat the analysis using this column as the meta-data to demonstrate that allelematch performed correctly.

example1chk <- amDataset(
  amExample1, 
  indexColumn = "sampleId",
  metaDataColumn = "knownIndividual", 
  missingCode = "-99"
)

uniqueExample1chk <- amUnique(example1chk, alleleMismatch = 2)

summary(uniqueExample1chk, html = "example1_2.html")

Notice how the individual identifier (a three-letter code) is consistent between the genotype declared unique and the genotypes that match it. Also note how each unique identifier appears only once in the list of unique genotypes, indicating that the analysis has not overestimated the number of individuals.

Example 2: Good quality data set

The data in this example have also been simulated4, this time to reflect the qualities of a good quality data set, where genotyping error and missing data exist, but these can be confidently handled by allelematch without manual intervention. At each locus, a random 4% of heterozygotes lost their second allele to simulate an allele dropout, and a random 4% of samples at each locus had alleles set to missing.

Load the data and create an amDataset object. The knownIndividual column is again kept as meta-data for instructional purposes. Then, we pre-screen the data using the amPreCheck() function.

data(amExample2)

example2 <- amDataset(
  amExample2, 
  indexColumn = "sampleId",
  metaDataColumn = "knownIndividual", 
  missingCode = "-99"
)

amPreCheck(example2)

The report for this dataset indicates good data quality: no single locus exceeds 4.7% missing data, the average individual missing data load is 4.6%, but the individual with the highest missing data load is at 20%. Furthermore, every pair of samples shares at least 6 out of 10 loci. Given these metrics, the function provides a safe verdict, indicating the dataset is well-suited for clustering.

Find the optimal alleleMismatch parameter for this data set.

amUniqueProfile(example2, doPlot = TRUE)

The output profile (Figure 2) demonstrates a ZeroSecondMinimum profile suggesting that the identified alleleMismatch parameter can be applied with confidence.

Figure 2: Good quality dataset.
Figure 2: Good quality dataset.

Conduct the unique analysis with the recommended alleleMismatch setting, and review the HTML output.

uniqueExample2 <- amUnique(example2, alleleMismatch = 3)

summary(uniqueExample2, html = "example2_1.html")

allelematch identified 100 unique genotypes from these 148 samples, and did so with no errors (as illustrated by the knownIndividual identifier). There are a number of examples of allelematch correctly matching samples despite allele mismatches (the mismatching alleles are highlighted in red).

\(P_{\text{sib}}\) is not given in these mismatching cases by default because partial matching samples are not “identical”. However, if we treat mismatching alleles as if they were missing, \(P_{\text{sib}}\) can be calculated for the non-missing loci as follows:

uniqueExample2 <- amUnique(example2, alleleMismatch = 3, doPsib = "all")

summary(uniqueExample2, html = "example2_2.html")

It is important to stress again that allelematch does not identify unique genotypes using the \(P_{\text{sib}}\) criterion, but rather presents this value to allow the user to assess the probability that:

Matching samples represent siblings of unique genotypes rather than duplicate samples of the same individual.

Unique genotypes represent siblings of other unique genotypes.

Example 3: Marginal quality data set

The data in this example have been simulated5 to represent a data set of marginal quality where the use of allelematch combined with careful manual review of the results is required to achieve a confident assessment of the unique genotypes. At each locus, a random 4% of heterozygotes lost their second allele to simulate an allele dropout, and a random 10% of samples at each locus had alleles set to missing.

Load the data and create an amDataset object. Once again, we retain the knownIndividual for instructional purposes. Then, we pre-screen the data using the amPreCheck() function.

data(amExample3)

example3 <- amDataset(
  amExample3, 
  indexColumn = "sampleId",
  metaDataColumn = "knownIndividual", 
  missingCode = "-99"
)

amPreCheck(example3)

The report for this dataset indicates marginal data quality: no single locus exceeds 9.7% missing data, the average individual missing data load is 9.7%, but the individual with the highest missing data load is at 50%. Here, at least one pair of sample shares only 2 out of 10 loci. Still, given these metrics, the function provides a safe verdict, indicating the dataset is well-suited for clustering.

Find the optimal alleleMismatch parameter for this data set.

amUniqueProfile(example3, doPlot = TRUE)

The output profile (Figure 3) demonstrates a NonZeroSecondMinimum profile, which is the first sign that we are dealing with a marginal data set. When this sort of profile is produced, amUniqueProfile() can be somewhat error-prone in finding the optimal value for the alleleMismatch parameter, and our first concern is making sure that it chose the best value.

Figure 3: Marginal quality dataset.
Figure 3: Marginal quality dataset.

We are looking for the second minimum in the multipleMatch variable (the first minimum being at alleleMismatch = 0), and the function appears to have identified this correctly at alleleMismatch = 6.

Conduct the unique analysis with the recommended alleleMismatch setting, and review the HTML output.

uniqueExample3 <- amUnique(example3, alleleMismatch = 6)

summary(uniqueExample3, html = "example3_1.html")

Two new issues appear in this example that did not exist for the previous two data sets. The header of the HTML output informs us that there are 2 unclassified samples and 4 multipleMatch samples. We will consider these in turn.

Unclassified samples are those which just exceed the criterion of similarity to be declared a match, but as an artifact of the dynamic tree cutting algorithm are not sufficiently different to be declared unique (i.e., to form their own clusters). Often this is because of missing data rather than mismatching data. If the analysis protocol does not permit their exclusion, the most robust approach is to reprofile these samples at their missing loci. However, this may not always be possible. Here we demonstrate an approach to simplify the manual classification of these samples.

The task is to determine whether these unclassified samples should be declared a unique genotype or as a match to an existing unique genotype. We use a function which conducts a pairwise analysis, comparing all the rows in the dataset of the first argument with all rows in the dataset of the second argument and returning those that exceed a similarity score (i.e., by Equation 2). In this case, the first dataset consists of the two unclassified samples from the previous step, and the second consists of all the samples declared as unique in the previous step. These datasets can be accessed from the uniqueExample3 object produced by amUnique().

unclassifiedExample3 <- amPairwise(
  uniqueExample3$unclassified,
  uniqueExample3$unique, 
  alleleMismatch = 7
)

Note how we used the next highest criterion, alleleMismatch = 7, because we want to see if these samples are only slightly more different than existing unique genotypes. We can then examine this in HTML format.

summary(unclassifiedExample3, html = "example3_2.html")

Using this output, we must now make a judgement call about what types of evidence are sufficient for declaring differences among individuals. Here, the rows with yellow highlighting are the unclassified samples and the rows without highlighting are unique genotypes.

For the first unclassified sample, it is mostly missing data and one mismatch that is driving the borderline status. If the five missing loci were reprofiled, it is likely that these two rows would match at a much lower criterion, and therefore sample 208 should be declared a match of sample 204. Thanks to the knownIndividual information retained from when the data were simulated, we can confirm this supposition. Indeed, the two samples did come from the same individual, ACN.

The second unclassified sample achieves its borderline status with four missing loci, but it also differs from its closest unique genotype because of three mismatching alleles. In this case, it is more likely that sample 251 is a unique genotype. Again, this conclusion is supported by the knownIndividual column with samples representing individuals ADA and ABN.

MultipleMatch samples are those that match more than one unique genotype. These occur when samples do not all sort into clearly defined groups representing unique genotypes. As the profiling shows, the numbers of these uncertain samples varies with the alleleMismatch criterion.

If we set the criterion too low, unique genotypes will not be sufficiently different and samples will appear to match multiple unique individuals. If we set a high criterion, it might not matter as long as unique genotypes are well-differentiated (e.g., Example 1; alleleMismatch = 3 to alleleMismatch = 7).

However, in cases where there is low allelic diversity and heterozygosity in the dataset, there could be insufficient information to distinguish unique genotypes when the criterion is set too high, and once again unique genotypes may be improperly determined, resulting in samples that will match multiple unique genotypes.

We could examine the earlier output from the unique analysis to resolve these samples (the questionable genotypes are flagged CHECK UNIQUE and the samples that cause this are flagged MULTIPLE MATCH). However, it is simpler to again use a pairwise approach to bring together just the relevant samples. Here we compare the four multipleMatch samples against the unique genotypes. The criterion is not changed this time from the original unique analysis.

multipleMatchExample3 <- amPairwise(
  uniqueExample3$multipleMatch,
  uniqueExample3$unique, 
  alleleMismatch = 6
)

summary(multipleMatchExample3, html = "example3_3.html")

Examining the HTML output, sample 112 (yellow highlighting) matches two unique genotypes (unhighlighted). Sample 110 is very likely the same unique genotype as 112 because it differs mostly because of missing data, while sample 260 differs by both four mismatching alleles and missing data. In this case, sample 112 can be declared a match of sample 110 and a false match of sample 260. Using the same argumentation, sample 162 is likely a match of 161, and sample 218 a match of 213. This has resolved three of the four multipleMatch cases.

Sample 183 differs from two unique genotypes (samples 182 and 181) chiefly because of missing data. In this case, we should consider two unique genotype rows to be the same as the multipleMatch sample; we could think of sample 183 as an intermediate or “missing link” between two incorrectly designated unique genotypes. As a result, we must remove one of the unique genotypes. Because it has the least amount of missing data, let’s declare sample 181 the unique genotype, and sample 183 and 182 matches of this unique genotype.

Once again, we turn to the knownIndividual column for confirmation that we made the correct decisions for this simulated data. Indeed, the correct decision was to reduce the total number of unique genotypes by one (individual ACF appears in two unique genotypes).

Reviewing the unclassified samples added sample 251 to the unique genotype list, and reviewing the multipleMatch samples removed sample 183 from the unique genotype list. Finally, we produce a CSV file of the original unique analysis and use spreadsheet software to make the necessary changes in classification.

summary(uniqueExample3, csv = "example3_1.csv")

Example 4: Low quality data set

For this example, we have simulated6 a low quality data set where uncertainty created by genotyping error and missing data, combined with a lack of information in the form of allelic diversity across loci, will result in a low confidence assessment of the unique genotypes. At each locus, a random 6% of heterozygotes lost their second allele to simulate an allele dropout, and a random 20% of samples at each locus had alleles set to missing.

Load the data and create an amDataset object. Again, we retain the knownIndividual for instructional reasons. Then, we pre-screen the data using the amPreCheck() function.

data(amExample4)

example4 <- amDataset(
  amExample4, 
  indexColumn = "sampleId",
  metaDataColumn = "knownIndividual", 
  missingCode = "-99"
)

amPreCheck(example4)

The report for this dataset indicates low data quality: no single locus exceeds 19.9% missing data, the average individual missing data load is 19.9%, but the individual with the highest missing data load is at 70%. Here, at least one pair of samples is not called at any of the same loci. Given these metrics, the function provides a critical warning, indicating the dataset may not provide stable clusters, and suggests close inspection of the output, and consideration of pre-filtration of the data to remove individuals or loci with high missing data rates.

Find the optimal alleleMismatch parameter for this data set.

amUniqueProfile(example4, doPlot = TRUE)

The output profile (Figure 4) demonstrates a NoSecondMinimum profile, which is a sign that allelematch cannot make a confident assessment of the unique genotypes within the range of the alleleMismatch criterion examined (by default this is 0% to 40% of allele columns mismatching). We can also, therefore, disregard the optimal alleleMismatch criterion.