Takes as input (i) the parent genotypes, (ii) the offspring genotype frequency, (iii) sequencing error rate, (iv) read depth, (v) bias, and (vi) overdispersion. It returns genotype likelihoods.
Usage
po_gl(
genovec,
ploidy,
p1_geno = NULL,
p2_geno = NULL,
rd = 10,
seq = 0.01,
bias = 1,
od = 0.01
)Arguments
- genovec
Offspring genotypes.
genovec[i]is the dosage for individual i.- ploidy
Ploidy
- p1_geno
Parent 1 genotype
- p2_geno
Parent 2 genotype
- rd
Read depth. Lower is more uncertain.
- seq
Sequencing error rate. Higher means more uncertain.
- bias
Bias. 1 means no bias.
- od
Overdispersion. Typical value is like 0.01. Higher means more uncertain.
Examples
set.seed(1)
po_gl(genovec = c(1, 2, 1, 1, 3), p1_geno = 2, p2_geno = 2, ploidy = 4)
#> $bias
#> [1] 1
#>
#> $seq
#> [1] 0.01
#>
#> $od
#> [1] 0.01
#>
#> $num_iter
#> [1] 8
#>
#> $llike
#> [1] -3.773022
#>
#> $postmat
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 1.858060e-03 0.67185769 0.3242752 0.002009102 3.191680e-11
#> [2,] 5.463499e-06 0.21609824 0.7526822 0.031214040 1.554090e-08
#> [3,] 1.858060e-03 0.67185769 0.3242752 0.002009102 3.191680e-11
#> [4,] 1.858060e-03 0.67185769 0.3242752 0.002009102 3.191680e-11
#> [5,] 1.554090e-08 0.03121404 0.7526822 0.216098243 5.463499e-06
#>
#> $gene_dist
#> [1] 0.02138753 0.20642301 0.54437892 0.20642301 0.02138753
#>
#> $par
#> $par$ell1
#> [1] 2
#>
#> $par$ell2
#> [1] 2
#>
#> $par$tau1
#> [1] 0.03424972
#>
#> $par$tau2
#> [1] 0.03424972
#>
#> $par$gamma1
#> [1] 0.4677409
#>
#> $par$gamma2
#> [1] 0.4677409
#>
#> $par$alpha
#> [1] 0.001
#>
#>
#> $genologlike
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] -4.947896 -1.315725 -3.006407 -7.128083 -22.82756
#> [2,] -10.226929 -1.899873 -1.614191 -3.834738 -16.08930
#> [3,] -4.947896 -1.315725 -3.006407 -7.128083 -22.82756
#> [4,] -4.947896 -1.315725 -3.006407 -7.128083 -22.82756
#> [5,] -16.089299 -3.834738 -1.614191 -1.899873 -10.22693
#>
#> $geno
#> [1] 1 2 1 1 2
#>
#> $maxpostprob
#> [1] 0.6718577 0.7526822 0.6718577 0.6718577 0.7526822
#>
#> $postmean
#> [1] 1.326435 1.815105 1.326435 1.326435 2.184895
#>
#> $input
#> $input$refvec
#> [1] 2 4 2 2 6
#>
#> $input$sizevec
#> [1] 10 10 10 10 10
#>
#> $input$ploidy
#> [1] 4
#>
#> $input$model
#> [1] "f1pp"
#>
#> $input$p1ref
#> [1] 7
#>
#> $input$p1size
#> [1] 10
#>
#> $input$p2ref
#> [1] 6
#>
#> $input$p2size
#> [1] 10
#>
#>
#> $prop_mis
#> [1] 0.2958125
#>
#> attr(,"class")
#> [1] "flexdog"