Tools for analysing forensic genetic DNA databases

Description

Computational efficient tools for comparing all pairs of profiles in a DNA database. The expectation and covariance of the summary statistic is implemented for fast computing. Routines for estimating proportions of close related individuals are available. The use of wildcards (also called F-designation) is implemented. Dedicated functions ease plotting the results.

Details

dbCompare: Compares make all n(n-1)/2 pairwise comparisons between profiles of a database with n DNA profiles. dbExpect: Computes the expected number of matching and partial matching loci for a given number of profiles in a database. dbVariance: Calculates the associated covariance matrix.

Author(s)

Torben Tvedebrink <tvede@math.aau.dk>, James Curran <j.curran@auckland.ac.nz> and Mikkel Meyer Andersen <mikl@math.aau.dk>.

References

Tvedebrink T, JM Curran, PS Eriksen, HS Mogensen and N Morling (2012). Analysis of matches and partial-matches in a Danish STR data set. Forensic Science International: Genetics, 6(3): 387-392.

Read the vignette: vigette('DNAtools')

Examples

  ## Not run: 
  data(dbExample)
  dbCompare(dbExample,hit=5,trace=TRUE)
  
## End(Not run)

The exact distribution of the number of alleles in a m-person DNA mixture

Description

Computes the exact distribution of the number of alleles in a m-person DNA mixture typed with STR loci. For a m-person DNA mixture it is possible to observe 1,\ldots,2\times m \times L alleles, where L is the total number of typed STR loci. The method allows incorporation of the subpopulation correction, the so-called \theta-correction, to adjust for shared ancestry. If needed, the locus-specific probabilities can be obtained using the locuswise argument.

Usage

Pnm_all(m, theta, probs, locuswise = FALSE)
Pnm_locus(m, theta, alleleProbs)

Arguments

Details

Computes the exact distribution of the number of alleles for a m-person DNA mixture.

Value

Returns a vector of probabilities, or a matrix of locuswise probability vectors.

Author(s)

Torben Tvedebrink, James Curran, Mikkel Andersen

References

T. Tvedebrink (2014). 'On the exact distribution of the number of alleles in DNA mixtures', International Journal of Legal Medicine; 128(3):427–37. <https://doi.org/10.1007/s00414-013-0951-3>

Examples

  ## Simulate some allele frequencies:
  freqs <-  structure(replicate(10, { g = rgamma(n = 10, scale = 4, shape = 3); 
                                      g/sum(g)
                                    },
              simplify = FALSE), .Names = paste('locus', 1:10, sep = '.'))

  ## Compute \eqn{\Pr(N(m = 3) = n)}, \eqn{n = 1,\ldots,2 * L *m}, where \eqn{L = 10}
  ## here
  Pnm_all(m = 2, theta = 0, freqs)
  ## Same, but locuswise results
  Pnm_all(m = 2, theta = 0, freqs, locuswise = TRUE)
  

Collapse m/p output to vector

Description

Collapse a m/p-matrix from dbCompare/dbExpect to a vector.

Usage

dbCollapse(x)

Arguments

Details

Collapse a m/p-matrix from dbCompare/dbExpect to a vector with entry i being the sum of all entries from m/p-matrix satisfying 2*m+p=i.

Value

A vector of length 2*max(m)+1 with entries begin the sum of entries i in m/p-matrix satisfying i=2*m+p.

Author(s)

Torben Tvedebrink

Examples

  ## Not run: 
  data(dbExample)
  res <- dbCompare(dbExample, hit=5, trace=TRUE)
  dbCollapse(res) ## same as dbCompare(dbExample, hit=5, trace=TRUE, collapse=TRUE)
  
## End(Not run)

Compare DNA profiles

Description

Compare DNA profiles

Usage

dbCompare(
  x,
  profiles = NULL,
  hit = 7,
  trace = TRUE,
  vector = FALSE,
  collapse = FALSE,
  wildcard = FALSE,
  wildcard.effect = FALSE,
  wildcard.impose = FALSE,
  Rallele = FALSE,
  threads = 2
)

Arguments

Details

Computes the distance between DNA profiles in terms of matching and partially-matching STR loci.

Value

Returns a matrix with the number of pairs mathcing/partially-matching at (i,j)-loci.

Author(s)

James Curran and Torben Tvedebrink. The multicore/CPU implementation was provided by Mikkel Meyer Andersen.

Examples

  ## Not run: 
  data(dbExample)
  dbCompare(dbExample,hit=5,trace=TRUE)
  
## End(Not run)

Simulated database with 1,000 individuals

Description

Database containing 1,000 simulated DNA profiles typed on ten autosomal markers.

Format

A data frame with each row being a DNA profile and each column a part of a genetic marker. Note that homozygote profiles has the same allelic value in the two columns associated to the same marker.

Expected value of cell counts in DNA database comparison

Description

Computes the expected number of cell counts when comparing DNA profiles in a DNA database. For every pair of DNA profiles in a database the number of matching and partial matching loci is recorded. A match is declared if the two DNA profiles coincide for both alleles in a locus and a partial-match is recorded if only one allele is shared between the profiles. With a total of L loci the number of matching loci is 0,...,L and partial number of matches is 0,...,L-m, where m is the number of matching loci.

Usage

dbExpect(
  probs,
  theta = 0,
  k = c(0, 0, 1),
  n = 1,
  r = 0,
  R = 0,
  round = FALSE,
  na = TRUE,
  vector = FALSE,
  collapse = FALSE,
  wildcard = FALSE,
  no.wildcard = NULL,
  rare.allele = FALSE,
  no.rare.allele = NULL
)

Arguments

Details

Computes the expected cell counts using a recursion formula. See Tvedebrink et al (2011) for details.

Value

Returns a matrix (or vector, see above) of expected cell counts.

Author(s)

James Curran and Torben Tvedebrink

References

T Tvedebrink, PS Eriksen, J Curran, HS Mogensen, N Morling. 'Analysis of matches and partial-matches in Danish DNA reference profile database'. Forensic Science International: Genetics, 2011.

Examples

  ## Not run: 
  ## Simulate some allele frequencies:
  freqs <-  replicate(10, { g = rgamma(n=10,scale=4,shape=3); g/sum(g)},
              simplify=FALSE)
  ## Compute the expected number for a DB with 10000 profiles:
  dbExpect(freqs,theta=0,n=10000)
  
## End(Not run) 

Simulate a DNA database

Description

Simulates a DNA database given a set of allele probabilities and theta value. It is possible to have close relatives in the database simulated in pairs, such that within each pair the profiles are higher correlated due to close familial relationship, but between pairs of profiles the correlation is only modelled by theta.

Usage

dbSimulate(probs, theta = 0, n = 1000, relatives = NULL)

Arguments

Details

Simulates a DNA database with a given number of DNA profiles (and possibly relatives) with a correlation between profiles governed by theta.

Value

A data frame where each row represents a DNA profile. The first column is a profile identifier (id) and the next 2*L columns contains the simulated genotype for each of the L loci. L is determined by the length of the list 'probs' with allele probabilities

Author(s)

James Curran and Torben Tvedebrink

Examples

  ## Not run: 
  ## Simulate some allele frequencies:
                                                                                
  freq <-  replicate(10, { g = rgamma(n=10,scale=4,shape=3); g/sum(g)},
             simplify=FALSE)
  ## Simulate a single database with 5000 DNA profiles:
  simdb <- dbSimulate(freq,theta=0,n=5000)
  ## Simulate a number of databases, say N=50. For each database compute
  ## the summary statistic using dbCompare:
  N <- 50
  Msummary <- matrix(0,N,(length(freq)+1)*(length(freq)+2)/2)
  for(i in 1:N)
    Msummary[i,] <- dbCompare(dbSimulate(freq,theta=0,n=1000),
                      vector=TRUE,trace=FALSE)$m
  ## Give the columns  representative names:
  dimnames(Msummary)[[2]] <- DNAtools:::dbCats(length(freq),vector=TRUE)
  ## Plot the simulations using a boxplot
  boxplot(log10(Msummary))
  ## There might come some warnings due to taking log10 to zero-values (no counts)
  ## Add the expected number to the plot:
  points(1:ncol(Msummary),log10(dbExpect(freq,theta=0,n=1000,vector=TRUE)),
         col=2,pch=16)
  
## End(Not run)

Covariance matrix of cell counts in DNA database comparison

Description

Computes the covariance matrix for the cell counts when comparing DNA profiles in a DNA database. For every pair of DNA profiles in a database the number of matching and partial matching loci is recorded. A match is declared if the two DNA profiles coincide for both alleles in a locus and a partial-match is recorded if only one allele is shared between the profiles. With a total of L loci the number of matching loci is 0,...,L and partial number of matches is 0,...,L-m, where m is the number of matching loci. The expression is given by:

latex

Usage

dbVariance(probs, theta = 0, n = 1, collapse = FALSE)

Arguments

Details

Computes the covariance matrix of the cell counts using a recursion formula. See Tvedebrink et al (2011) for details.

Value

Returns a covariance matrix for the cell counts.

Author(s)

James Curran and Torben Tvedebrink

References

T Tvedebrink, PS Eriksen, J Curran, HS Mogensen, N Morling. 'Analysis of matches and partial-matches in Danish DNA reference profile database'. Forensic Science International: Genetics, 2011.

Examples

  ## Not run: 
  ## Simulate some allele frequencies:
  freqs <-  replicate(10, { g = rgamma(n=10,scale=4,shape=3); g/sum(g)}, simplify=FALSE)
  ## List of elements needed to compute the covariance matrix.
  ## Useful option when the covariance needs to be computed for varying
  ## database sizes but for identical theta-value.
  comps <- dbVariance(freqs,theta=0,n=1)
  ## Covariance for a DB with 1000 DNA profiles
  cov1000 <- dbVariance(freqs,theta=0,n=1000)
  ## The result is the same as:
  comps1000 <- choose(1000,2)*comps$V1 + 6*choose(1000,3)*comps$V2 + 6*choose(1000,4)*comps$V3
  
## End(Not run)

Estimate the drop-out probability based on number of alleles

Description

An inferior may to estimate the drop-out probability compared to using the peak heights from the electropherogram. However, to compare the performance with Gill et al. (2007) this implements a theoretical approach based on their line of arguments.

Usage

estimatePD(n0, m, pnoa = NULL, probs = NULL, theta = 0, locuswise = FALSE)

Arguments

Details

Computes the \Pr(D) that maximises equation (10) in Tvedebrink (2014).

Value

Returns the MLE of \Pr(D) based on equation (10) in Tvedebrink (2014)

Author(s)

Torben Tvedebrink

References

Gill, P., A. Kirkham, and J. Curran (2007). LoComatioN: A software tool for the analysis of low copy number DNA profiles. Forensic Science International 166(2-3): 128 - 138.

T. Tvedebrink (2014). 'On the exact distribution of the number of alleles in DNA mixtures', International Journal of Legal Medicine; 128(3):427–37. <https://doi.org/10.1007/s00414-013-0951-3>

Examples

  ## Simulate some allele frequencies:
  freqs <-  simAlleleFreqs()
  ## Assume 15 alleles are observed in a 2-person DNA mixture with 10 loci:
  estimatePD(n0 = 15, m = 2, probs = freqs)

Simple allele frequency estimation

Description

Estimates allele frequencies from a database with DNA profiles

Usage

freqEst(x)

Arguments

Details

Computes the allele frequencies for a given database.

Value

Returns a list of probability vectors - one vector for each locus.

Author(s)

James Curran and Torben Tvedebrink

Examples

  data(dbExample)
  freqEst(dbExample)

Generates DNA profiles of n individuals.

Description

These are formed as n/2 pairs for relatives with a IDB-vector given by k. I.e. the profiles are mutually unrelated between pairs.

Usage

genRypeRec(x, t, k, n, print = FALSE)

Arguments

Generates DNA profiles of n unrelated individuals for a locus

Description

Generates DNA profiles of n unrelated individuals for a locus

Usage

genTypeRec(x, t, n, z = rep(0, lx <- length(x)))

Arguments

Estimate theta and the fraction of comparisons between close relatives

Description

Estimates the fraction of comparisons between pairs of close relatives while fitting the theta parameter minimising the object function. The function makes use of the R-package 'Rsolnp' which is an implementation of an solver for non-linear minimisation problems with parameter constraints.

Usage

optim.relatedness(
  obs,
  theta0 = 0,
  theta1 = 0.03,
  theta.tol = 10^(-7),
  theta.step = NULL,
  max.bisect = 15,
  probs,
  var.list = NULL,
  init.alpha = 10^c(-4, -6, -8, -10),
  init.keep = FALSE,
  objFunction = c("T2", "T1", "C3", "C2", "C1"),
  collapse = FALSE,
  trace = FALSE,
  solnp.ctrl = list(tol = 10^(-9), rho = 10, delta = min(init.alpha) * 0.01, trace =
    FALSE)
)

Arguments

Details

Computes the proportion of comparisons between close relatives in a database matching exercise for each theta value under investigation.

Value

Returns a list of three components: value, solution and var.list. The first element, value, is a dataframe with the value of the objection function for each of the theta values investigated. Solution is the estimated alpha-vector where the objection function was minimised. Finally, var.list is a names list of components for computing variances. May be reused in later computations for increased speed in some iterations.

Author(s)

James Curran and Torben Tvedebrink

References

T Tvedebrink, PS Eriksen, J Curran, HS Mogensen, N Morling. 'Analysis of matches and partial-matches in Danish DNA reference profile database'. Forensic Science International: Genetics, 2011.

Examples

  ## Not run: 
  ## Simulate some allele frequencies:
  freqs <-  replicate(10, { g = rgamma(n=10,scale=4,shape=3); g/sum(g)},
              simplify=FALSE)
  ## Load the sample database:
  data(dbExample)
  obs <- dbCompare(dbExample,trace=FALSE)$m
  C3 <- optim.relatedness(obs,theta0=0.0,theta1=0.03,probs=freqs,
          objFunction='C3',max.bisect=30,trace=TRUE)
  
## End(Not run)

Compute the posterior probabilities for P(m|n0) for a given prior P(m) and observed vector n0 of locus counts

Description

where m ranges from 1 to m_{\max} and n_0 is the observed locus counts.

Usage

pContrib(n0, probs = NULL, m.prior = rep(1/m.max, m.max), m.max = 8, theta = 0)

Arguments

Details

Computes a vector P(m|n0) evaluated over the plausible range 1,...,m.max.

Value

Returns a vector P(m|n0) for m=1,...,m.max

Author(s)

Torben Tvedebrink, James Curran

References

T. Tvedebrink (2014). 'On the exact distribution of the number of alleles in DNA mixtures', International Journal of Legal Medicine; 128(3):427–37. <https://doi.org/10.1007/s00414-013-0951-3>

Examples

  ## Simulate some allele frequencies:
  freqs <-  simAlleleFreqs()
  m <- 2
  n0 <- sapply(freqs, function(px){
                            peaks = unique(sample(length(px),
                                             size = 2 * m,
                                             replace = TRUE,
                                             prob = px))
                            return(length(peaks))
                       })
  ## Compute P(m|n0) for m=1,...,4 and the sampled n0
  pContrib(n0=n0,probs=freqs,m.max=4)

Compute the posterior probabilities for \Pr(m|n_0) for a given prior \Pr(m).

Description

Compute a matrix of posterior probabilties \Pr(m|n_0) where m ranges from 1 to m_{\max}, and n_0 is 0,\ldots,2m_{\max}. This is done by evaluating \Pr(m|n_0)=Pr(n_0|m)Pr(m)/Pr(n), where \Pr(n_0|m) is evaluated by pNoA.

Usage

pContrib_locus(
  prob = NULL,
  m.prior = NULL,
  m.max = 8,
  pnoa.locus = NULL,
  theta = 0
)

Arguments

Details

Computes a matrix of \Pr(m|n_0) values for a specific locus.

Value

Returns a matrix [\Pr(m|n_0)] for m = 1,\ldots,m.max and n_0 = 1,\ldots,2m.max.

Author(s)

Torben Tvedebrink, James Curran

References

T. Tvedebrink (2014). 'On the exact distribution of the number of alleles in DNA mixtures', International Journal of Legal Medicine; 128(3):427–37. <https://doi.org/10.1007/s00414-013-0951-3>

Examples

  ## Simulate some allele frequencies:
  freqs <-  simAlleleFreqs()
  
  ## Compute Pr(m|n0) for m = 1, ..., 5 and n0 = 1, ..., 10 for the first locus:
  pContrib_locus(prob = freqs[[1]], m.max = 5)

Plots the fitted object function for estimated familial relationships in the database and theta.

Description

Plots the minimised object function for included values of theta

Usage

## S3 method for class 'dbOptim'
plot(x, type = "l", ...)

Arguments

Details

Plots the object function

Value

A plot of the object function

Author(s)

James Curran and Torben Tvedebrink

See Also

optim.relatedness

Examples

  ## Not run: 
  ## Simulate some allele frequencies:
  freqs <-  replicate(10, { g = rgamma(n=10,scale=4,shape=3); g/sum(g)},
              simplify=FALSE)
  ## Load the sample database:
  data(dbExample)
  obs <- dbCompare(dbExample,trace=FALSE)$m
  C3 <- optim.relatedness(obs,theta0=0.0,theta1=0.03,probs=freqs,
          objFunction='C3',max.bisect=30,trace=TRUE)
  plot(C3)
  
## End(Not run)

Plots the summary matrix

Description

Plots the summary matrix with counts on y-axis and classification on x-axis.

Usage

## S3 method for class 'dbcompare'
plot(x, log = "y", las = 3, xlab = "Match/Partial", ylab = "Counts", ...)

Arguments

Value

A plot of the summary matrix. The counts are on log10 scale and the x-axis is labeled by appropriate matching/partially-matching levels.

Author(s)

James Curran and Torben Tvedebrink

See Also

dbCompare,print.dbcompare

Examples

  ## Not run: 
  data(dbExample)
  M = dbCompare(dbExample,hit=5)
  plot(M)
  
## End(Not run)

Prints the results from optim.relatedness()

Description

Prints the evaluated functions for the object function, best estimate of alpha and possibly list of variances.

Usage

## S3 method for class 'dbOptim'
print(x, var.list = FALSE, ...)

Arguments

Details

Prints the summary details of the fit

Value

A dataframe with [theta,value] and a vector of fitted alpha parameters

Author(s)

James Curran and Torben Tvedebrink

See Also

optim.relatedness

Examples

  ## Not run: 
  ## Simulate some allele frequencies:
  freqs <-  replicate(10, { g = rgamma(n=10,scale=4,shape=3); g/sum(g)},
              simplify=FALSE)
  ## Load the sample database:
  data(dbExample)
  obs <- dbCompare(dbExample,trace=FALSE)$m
  C3 <- optim.relatedness(obs,theta0=0.0,theta1=0.03,probs=freqs,
          objFunction='C3',max.bisect=30,trace=TRUE)
  print(C3)
  
## End(Not run)

Prints the summary matrix

Description

Prints the summary matrix and possible 'big hits'.

Usage

## S3 method for class 'dbcompare'
print(x, ...)

Arguments

Details

Prints the summary matrix

Value

Prints the summary matrix and data frame with 'big hits'

Author(s)

James Curran and Torben Tvedebrink

See Also

dbCompare,plot.dbcompare

Examples

  ## Not run: 
  data(dbExample)
  M = dbCompare(dbExample,hit=5)
  M
  
## End(Not run)

Simulate Allele Frequencies

Description

Simulate some allele frequencies using Dirichlet Random variables

Usage

simAlleleFreqs(
  nLoci = 10,
  allelesPerLocus = rep(10, nLoci),
  shape = rep(3, nLoci)
)

Arguments

Value

a list with elements locus.l where l=1,\ldots,L, each of which are vectors of length allelesPerLocus[l], consisting of allele frequencies for that locus

Examples

set.seed(123)
simAlleleFreqs()