Description Usage Arguments Author(s) Examples
The (incomplete data) loglikelihood function for the negative binomial mixture model
1 | loglikelihood(theta, x.syn, x.non, cvec)
|
theta |
c(alpha, beta, p1, p2, ..., pk-1) |
x.syn |
vector of synonymous mutation count |
x.non |
vector of non-synonymous mutation count |
cvec |
vector of k positive values c(c1, ..., ck) |
Johanna Bertl
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 | c1 = 0.5; c2 = 10
p1 = 0.2; p2 = 0.8
alpha = 10
beta = 5
mutations = rnegbinmix_syn_nonsyn(n=3000, alpha=alpha, beta=beta, c=c(c1, c2), p = c(p1, p2))
# log likelihood at the true value
loglikelihood(theta=c(alpha,beta,p1), x.syn = mutations$Syn, x.non = mutations$Non, cvec=c(c1, c2))
# likelihood surface
alphavec = 1:20
betavec = seq(0.5, 10, by=0.5)
p1vec = seq(0.11, 0.3, by=0.01)
llalpha = numeric(20)
llbeta = numeric(20)
llp1 = numeric(20)
for(i in 1:20){
llalpha[i] = loglikelihood(theta=c(alphavec[i],beta,p1), x.syn = mutations$Syn, x.non = mutations$Non, cvec=c(c1, c2))
llbeta[i] = loglikelihood(theta=c(alpha,betavec[i],p1), x.syn = mutations$Syn, x.non = mutations$Non, cvec=c(c1, c2))
llp1[i] = loglikelihood(theta=c(alpha,beta,p1vec[i]), x.syn = mutations$Syn, x.non = mutations$Non, cvec=c(c1, c2))
}
plot(alphavec, llalpha, t="b")
abline(v=alpha, col="red")
plot(betavec, llbeta, t="b")
abline(v=beta, col="red")
plot(p1vec, llp1, t="b")
abline(v=p1, col="red")
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.