admix (version 2.1-3)

IBM_hessian_contrast: Hessian matrix of the contrast function in the Inversion - Best Matching (IBM) method

Description

Compute the hessian matrix of the contrast as defined in the IBM approach, at point (p1,p2). Here, based on two samples following admixture models, where we recall that admixture models have probability distribution function (pdf) given by l where l = p*f + (1-p)*g, where g represents the only known quantity and l is the pdf of the observed sample. See 'Details' below for further information about the definition of the contrast.

Usage

IBM_hessian_contrast(
  par,
  fixed.p1 = NULL,
  known.p = NULL,
  sample1,
  sample2,
  G,
  comp.dist = NULL,
  comp.param = NULL
)

Value

the hessian matrix of the contrast.

Arguments

par

Numeric vector with two elements (corresponding to the two unknown component weights) at which the hessian is computed.

fixed.p1

(optional, NULL by default) Arbitrary value chosen by the user for the component weight related to the unknown component of the first admixture model. Only useful for optimization when the known components of the two models are identical (G1=G2, leading to unidimensional optimization).

known.p

(optional, NULL by default) Numeric vector with two elements, the known (true) mixture weights.

sample1

Observations of the first sample under study.

sample2

Observations of the second sample under study.

G

Distribution on which to integrate when calculating the contrast.

comp.dist

A list with four elements corresponding to the component distributions (specified with R native names for these distributions) involved in the two admixture models. The two first elements refer to the unknown and known components of the 1st admixture model, and the last two ones to those of the second admixture model. If there are unknown elements, they must be specified as 'NULL' objects. For instance, 'comp.dist' could be specified as follows: list(f1=NULL, g1='norm', f2=NULL, g2='norm').

comp.param

A list with four elements corresponding to the parameters of the component distributions, each element being a list itself. The names used in this list must correspond to the native R argument names for these distributions. The two first elements refer to the parameters of unknown and known components of the 1st admixture model, and the last two ones to those of the second admixture model. If there are unknown elements, they must be specified as 'NULL' objects. For instance, 'comp.param' could be specified as follows: : list(f1=NULL, g1=list(mean=0,sd=1), f2=NULL, g2=list(mean=3,sd=1.1)).

Author

Xavier Milhaud

Details

See the paper presenting the IBM approach at the following HAL weblink: https://hal.science/hal-03201760

Examples

## Simulate data:
list.comp <- list(f1 = 'norm', g1 = 'norm',
                  f2 = 'norm', g2 = 'norm')
list.param <- list(f1 = list(mean = 3, sd = 0.5), g1 = list(mean = 0, sd = 1),
                   f2 = list(mean = 1, sd = 0.1), g2 = list(mean = 5, sd = 2))
sample1 <- rsimmix(n=1500, unknownComp_weight=0.5, comp.dist = list(list.comp$f1,list.comp$g1),
                                                   comp.param=list(list.param$f1,list.param$g1))
sample2 <- rsimmix(n=2000, unknownComp_weight=0.7, comp.dist = list(list.comp$f2,list.comp$g2),
                                                   comp.param=list(list.param$f2,list.param$g2))
## Define the distribution over which to integrate:
fit.all <- stats::density(x = c(sample1[['mixt.data']],sample2[['mixt.data']]))
G <- stats::rnorm(n = 1000, mean = sample(c(sample1[['mixt.data']], sample2[['mixt.data']]),
                                          size = 1000, replace = TRUE), sd = fit.all$bw)
## Evaluate the hessian matrix at point (p1,p2) = (0.3,0.6):
list.comp <- list(f1 = NULL, g1 = 'norm',
                  f2 = NULL, g2 = 'norm')
list.param <- list(f1 = NULL, g1 = list(mean = 0, sd = 1),
                   f2 = NULL, g2 = list(mean = 5, sd = 2))
IBM_hessian_contrast(par = c(0.3,0.6), fixed.p1 = NULL, known.p = NULL,
                     sample1 = sample1[['mixt.data']],  sample2 = sample2[['mixt.data']], G = G,
                     comp.dist = list.comp, comp.param = list.param)