Hi all, I appreciate the help this list has given me before. I have a
question which has been perplexing me. I have been working on doing a
Bayesian calculating inserting studies sequentially after using a
non-informative prior to get a meta-analysis type result. I created a
function using three iterations of this, my code is below. I insert prior
mean and precision (I add precision manually because I want it to be zero
but if I include zero as prior variance I get an error cant divide by zero.
Now my question is instead of having the three iterations below, is there
anyway to create a loop, the only problem is I want to have elements from
the previous posterior to be the new prior and now I cant figure out how to
do the code below in a loop. The data below is dummy data, I used a starting
mu of 1, and starting precision of 0.
bayes.analysis.treat<-function(mu0,p0){
n1 = 5
n2 = 10
n3 = 15
ybar1 = 12
ybar2 = 13
ybar3 = 14
sd1 = 2
sd2 = 3
sd3 = 4
#posterior 1
var1 = sd1^2 #sample variance
p1 = n1/var1 #sample precision
p1n = p0+p1
mu1 = ((p0)/(p1n)*mu0)+((p1)/(p1n))*ybar1
sigma1 = 1/sqrt(p1n)
#posterior 2
var2 = sd2^2 #sample variance
p2 = n2/var2 #sample precision
p2n = p1n+p2
mu2 = ((p1n)/(p2n)*mu1)+((p2)/(p2n))*ybar2
sigma2 = 1/sqrt(p2n)
#posterior 3
var3 = sd3^2 #sample variance
p3 = n3/var3 #sample precision
p3n = p2n+p3
mu3 = ((p2n)/(p3n)*mu2)+((p3)/(p3n))*ybar3
sigma3 = 1/sqrt(p3n)
posterior<-cbind(rbind(mu1,mu2,mu3),rbind(sigma1,sigma2,sigma3))
rownames(posterior)<-c("Posterior 1", "Posterior 2", "Posterior 3")
colnames(posterior)<-c("Mu", "SD")
return(posterior)}
-------------------------------------------
Joe King, M.A.
Ph.D. Student
University of Washington - Seattle
206-913-2912
<mailto:j...@joepking.com> j...@joepking.com
-------------------------------------------
"Never throughout history has a man who lived a life of ease left a name
worth remembering." --Theodore Roosevelt
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
