DFE parametrizations

DFE parametrizations#

The DFE needs to be parametrized in some way in order to be amenable to maximum likelihood estimation. GammaExpParametrization is used by default. Other parametrizations are also implemented (cf. parametrization), and custom parametrizations can be created by subclassing Parametrization.

To see how the parametrization affects the shape of the DFE, we use our example data for B. pendula.

import fastdfe as fd

# only show very important log messages
fd.logger.setLevel('FATAL')

parametrizations = [
    fd.GammaExpParametrization(),
    fd.DiscreteFractionalParametrization(),
    fd.GammaDiscreteParametrization(),
    fd.DisplacedGammaParametrization()
]

inferences = []
for model in parametrizations:
    inf = fd.BaseInference(
        sfs_neut=fd.Spectrum([177130, 997, 441, 228, 156, 117, 114, 83, 105, 109, 0]),
        sfs_sel=fd.Spectrum([797939, 1329, 499, 265, 162, 104, 117, 90, 94, 119, 0]),
        fixed_params=dict(all=dict(h=0.5, eps=0)),
        model=model
    )

    inf.run()

    inferences.append(inf)
Discretization>Precomputing: 100%|██████████| 9/9 [--:--<--:--, --it/s]
BaseInference>Performing inference: 100%|██████████| 10/10 [--:--<--:--, --it/s]
BaseInference>Bootstrapping (2 runs each): 100%|██████████| 100/100 [--:--<--:--, --it/s]
Discretization>Precomputing: 100%|██████████| 9/9 [--:--<--:--, --it/s]
BaseInference>Performing inference: 100%|██████████| 10/10 [--:--<--:--, --it/s]
BaseInference>Bootstrapping (2 runs each): 100%|██████████| 100/100 [--:--<--:--, --it/s]
Discretization>Precomputing: 100%|██████████| 9/9 [--:--<--:--, --it/s]
BaseInference>Performing inference: 100%|██████████| 10/10 [--:--<--:--, --it/s]
BaseInference>Bootstrapping (2 runs each): 100%|██████████| 100/100 [--:--<--:--, --it/s]
Discretization>Precomputing: 100%|██████████| 9/9 [--:--<--:--, --it/s]
BaseInference>Performing inference: 100%|██████████| 10/10 [--:--<--:--, --it/s]
BaseInference>Bootstrapping (2 runs each): 100%|██████████| 100/100 [--:--<--:--, --it/s]
library(fastdfe)
fd <- load_fastdfe()

# only show very important log messages
fd$logger$setLevel('FATAL')

parametrizations <- list(
  fd$GammaExpParametrization(),
  fd$DiscreteFractionalParametrization(),
  fd$GammaDiscreteParametrization(),
  fd$DisplacedGammaParametrization()
)

inferences <- list()
for (model in parametrizations) {
  inf <- fd$BaseInference(
    sfs_neut = fd$Spectrum(c(177130, 997, 441, 228, 156, 117, 114, 83, 105, 109, 0)),
    sfs_sel = fd$Spectrum(c(797939, 1329, 499, 265, 162, 104, 117, 90, 94, 119, 0)),
    fixed_params = list(all = list(h = 0.5, eps = 0)),
    model = model
  )

  sfs_modelled <- inf$run()

  inferences <- c(inferences, inf)
}
Discretization>Precomputing: 100%|██████████| 9/9 [--:--<--:--, --it/s]
BaseInference>Performing inference: 100%|██████████| 10/10 [--:--<--:--, --it/s]
BaseInference>Bootstrapping (2 runs each): 100%|██████████| 100/100 [--:--<--:--, --it/s]
Discretization>Precomputing: 100%|██████████| 9/9 [--:--<--:--, --it/s]
BaseInference>Performing inference: 100%|██████████| 10/10 [--:--<--:--, --it/s]
BaseInference>Bootstrapping (2 runs each): 100%|██████████| 100/100 [--:--<--:--, --it/s]
Discretization>Precomputing: 100%|██████████| 9/9 [--:--<--:--, --it/s]
BaseInference>Performing inference: 100%|██████████| 10/10 [--:--<--:--, --it/s]
BaseInference>Bootstrapping (2 runs each): 100%|██████████| 100/100 [--:--<--:--, --it/s]
Discretization>Precomputing: 100%|██████████| 9/9 [--:--<--:--, --it/s]
BaseInference>Performing inference: 100%|██████████| 10/10 [--:--<--:--, --it/s]
BaseInference>Bootstrapping (2 runs each): 100%|██████████| 100/100 [--:--<--:--, --it/s]

The inferred DFEs are plotted in discretized form.

import numpy as np

fd.Inference.plot_discretized(
    inferences=inferences,
    labels=['GammaExp', 'DiscreteFractional', 'GammaDiscrete', 'DisplacedGamma'],
    intervals=[-np.inf, -100, -10, -1, 1, np.inf]
);
../_images/2687055207bf04c476680cba531849c5f8e3cf1bebe7bcdf8cff909ed9556751.png
p <- fd$Inference$plot_discretized(
  inferences = inferences,
  labels = c('GammaExp', 'DiscreteFractional', 'GammaDiscrete', 'DisplacedGamma'),
  intervals = c(-Inf, -100, -10, -1, 1, Inf)
)
../_images/294349f7631ec37f8d4e3b5f6d579efcc851e6e55fdd95798952bc1642bd1ac3.png

The overall shape is similar, but DiscreteFractionalParametrization shows noticeably wider confidence intervals. In general, estimating the full DFE with a sample size of 10 and the limited SNP count used here leads to substantial uncertainty.