#> source("c:/Work07/ShortCourse07/genetic_em.txt",print.eval=TRUE)#
#This does EM for the genetic data#
nsim<-25;x<-c(125,18,20,34);t<-array(.3,dim=c(nsim,1));
std1<-array(.7,dim=c(nsim,1));std2<-array(.5,dim=c(nsim,1))
z2<-x[1]/2;
for (j in 2:nsim)  {
z2<-x[1]*t[j-1]/(2+t[j-1])
t[j]<- (z2+x[4])/(z2+x[2]+x[3]+x[4]);
std1[j]<-(x[1]*t[j]/((2+t[j])*t[j]^2))+(x[4]/(t[j]^2))+((x[2]+x[3])/((1-t[j])^2))
std1[j]<-1/sqrt(std1[j]);std1[j]<-std2[j]<-std1[j]
}
plot(t,type="l",lwd=2,xlim=c(1,nsim),ylim=c(.55,.7),xlab="iteration",main="",ylab="")
par(new=T)
plot(t+std1,type="l",lwd=2,lty=2,xlim=c(1,nsim),ylim=c(.55,.7),xlab="iteration",main="",ylab="")
par(new=T)
plot(t-std2,type="l",lwd=2,lty=2,xlim=c(1,nsim),ylim=c(.55,.7),xlab="iteration",main="",ylab="")