varEst {DNAmixturesLite}

Estimated asymptotic variance matrix for MLE

Description

Provided that the user specifies the MLE as well as the constraints used in the maximisation, this function computes an estimate of the variance of the MLE based on the observed information. The observed information is obtained by numerically deriving the Hessian of the log-likelihood function.

IMPORTANT: This is the DNAmixturesLite package, which is intended as a service to enable users to try DNAmixtures without purchasing a commercial licence for Hugin. When at all possible, we strongly recommend the use of DNAmixtures rather than this lite-version. See https://dnamixtures.r-forge.r-project.org/ for details on both packages.

While the lite-version seeks to provide the full functionality of DNAmixtures, note that computations are much less efficient and that there are some differences in available functionality. Be aware that the present documentation is copied from DNAmixtures and thus may not accurately describe the implementation of this lite-version.

Usage

varEst(mixture, mle, npars, method = "Richardson", ...)

## S3 method for class 'mixVarEst'
summary(object, transform = FALSE, ...)

## S3 method for class 'summary.mixVarEst'
print(
  x,
  digits = max(3, getOption("digits") - 3),
  scientific = FALSE,
  print.gap = 3L,
  ...
)

Arguments

mixture

A DNAmixture.

mle

A mixpar, typically obtained by mixML.

npars

A list of integers specifying the number of each of the four parameters rho, eta, xi and phi. Allowed values are

0

The parameter is fixed, but might differ across mixtures.

1

The parameter is equal across mixtures.

N

There is one parameter for each of the N mixtures in the model.

method

Method for numeric differentiation used in hessian

...

Arguments to be passed on to other methods.

object

An object of class mixVarEst, typically obtained by a call to varEst.

transform

Should the parameterisation (\mu, \sigma) be used? Defaults to FALSE.

x

An object of class "summary.mixVarEst".

digits

Number of significant digits to print

scientific

Should scientific notation be used?

print.gap

Distance between columns in the printing of the summary.

Details

As the user can apply highly customized constraints to the model parameters when maximising with mixML, it is a complicated matter to write a generic function for computing the asymptotic variance matrix. We have thus restricted attention to the case where each of the (multi-dimensional) parameters rho, eta, xi and phi can be either

• fixed at known values

• unknown, but common across traces

• unconstrained

Value

The mle and the estimated covariance matrix in different parametrisations

cov and mle

\rho, \eta, \xi, \phi

cov.trans and mle.trans

\mu, \sigma, \xi, \phi

cov.res

A non-singular covariance matrix for a reparametrisation of \rho, \eta, \xi, \phi, collapsing the parameters according to the specified constrains and removing one contributor from \phi.

An integer suffix is used to indicate which mixture the parameter is associated with. In the restricted covariance matrix, all fixed parameters are left out. If parameters are equal accross mixtures, the suffix for this parameter will be .1. If parameters are unconstrained and there are N mixtures in the model, suffixes are .1, ..., .N

Examples

data(MC18, USCaucasian)
mixHp <- DNAmixture(list(MC18), k = 3, K = c("K1", "K2", "K3"), C = list(50),
                    database = USCaucasian)
p <- mixpar(rho = list(30), eta = list(34), xi = list(0.08),
            phi = list(c(K1 = 0.71, K3 = 0.1, K2 = 0.19)))
mlHp <- mixML(mixHp, pars = p)
## Find the estimated covariance matrix of the MLE
V.Hp <- varEst(mixHp, mlHp$mle, npars = list(rho=1,eta=1,xi=1,phi=1))
V.Hp$cov ## using (rho, eta)
V.Hp$cov.trans ## using (mu, sigma)
## The summary is a table containing the MLE and their standard errors
summary(V.Hp)


data(MC18, USCaucasian)
mixmult <- DNAmixture(list(MC18), C = list(50), k = 3, K = c("K1", "K2"), database = USCaucasian)
startpar <- mixpar(rho = list(30), eta = list(28), xi = list(0.08),
                   phi = list(c(U1 = 0.2, K1 = 0.7, K2 = 0.1)))
ml.mult <- mixML(mixmult, startpar)
Vmult <- varEst(mixmult, ml.mult$mle, list(rho=1,eta=1,xi=1,phi=1))
summary(Vmult)



## Be aware that the following two advanced examples are computationally demanding and
## typically have a runtime of several minutes with the lite-version of DNAmixtures.

data(MC15, MC18, USCaucasian)
mix <- DNAmixture(list(MC15, MC18), C = list(50, 38), k = 3, K = c("K1", "K2"),
database = USCaucasian)
startpar <- mixpar(rho = list(30, 30), eta = list(28, 35), xi = list(0.08, 0.1),
                   phi = list(c(U1 = 0.05, K1 = 0.7, K2 = 0.25),
                                c(K1 = 0.7, K2 = 0.1, U1 = 0.2)))
eqxis <- function(x){ diff(unlist(x[,"xi"])) }
## Here we set stutter equal for all traces
ml.diff <- mixML(mix, startpar, eqxis, val = 0, phi.eq = FALSE)
V.diff <- varEst(mix, ml.diff$mle, list(rho=2,eta=2,xi=1,phi=2))
summary(V.diff)

## Fixing stutter to 0.07
xival <- function(x){unlist(x[,"xi"])}
ml.eq <- mixML(mix, startpar, xival, val = c(0.07, 0.07), phi.eq = FALSE)
V.eq <- varEst(mix, ml.eq$mle, list(rho=2,eta=2,xi=0,phi=2))
summary(V.eq)