# # Model for Hardy-Weinberg data in which we have genotypes corresponding to k alleles - # we have k parameters, k-1 distinct allele probs q and an inbreeding coefficient f. # We have a dirichlet prior on q. f is constrained to be > fmin = -qmin/(1-qmin) where qmin # is the minimum of q[1],...,q[k]. # # model{ # genotype frequencies, q's are allele freqs, f is the common inbreeding coeff for (i in 1:k){ for(j in 1:i){ p[i*(i-1)/2+j] <- equals(i,j)*(q[i]*q[i] + f*q[i]*(1-q[i])) + (1-equals(i,j))*2*q[i]*q[j]*(1-f) } } # likelihood N <- sum(n[]) # ncell <- k*(k+1)/2 n[1:ncell] ~ dmulti(p[],N) # priors for (i in 1:k){ q[i] <- delta[i]/sum(delta[]) delta[i] ~ dgamma(1,1) } qmin <- ranked(q[],1) fmin <- -qmin/(1-qmin) lambda ~ dnorm(lambdamu,lambdasd) f <- (exp(lambda)+fmin)/(exp(lambda)+1) for (i in 1:k-1){ theta[i] <- log(q[i]/q[k]) } theta[k] <- log((f-fmin)/(1-f)) } # Data is in the order n11, n21, n22, n31, n32, ..., nkk list(n=c(5,6,1,7,11,2,8,19,26,15),k=4,ncell=10,lambdamu=-1,lambdasd=1) list(n=c(49,45,6),k=2,ncell=3,lambdamu=-1,lambdasd=1) list(delta=c(.25,.25),lambda=0) list(n=c(0,3,1,5,18,1,3,7,5,2),k=4,ncell=10,lambdamu=-2.95,lambdasd=1.07) list(delta=c(.25,.25,.25,.25),lambda=0.5) list(n=c(1236,120,3,18,0,0,982,55,7,249,32,1,0,12,0,2582,132,20,1162,29,1312,6,0,0,4,0,4,0,2,0,0,0,0,0,0,0,115,5,2,53,1,149,0,0,4),k=9,ncell=45,lambdamu=-4.65,lambdasd=0.21) # orig for 0,0.27 were lambdamu=-4.95, lambdasd=2.09 list(delta=c(1,1,1,1,1,1,1,1,1),lambda=0.0)