Friday, August 7, 2009

EM - algorithm: ABO


### ABO
Na<-25; Nab<-50; Nb<-25; No<-15; N<-Na+Nb+No+Nab

pa<-0.33; pb<-0.33; po<-0.34

for(i in 1:20){
Naa<-pa^2/(pa^2+2*pa*po)*Na
Nao<-2*pa*po/(pa^2+2*pa*po)*Na
Nbb<-pb^2/(pb^2+2*pb*po)*Nb
Nbo<-2*pb*po/(pb^2+2*pb*po)*Nb
#Nab<-2*pa*pb*Nab
#Noo<-po^2*No
pa<-(2*Naa+Nab+Nao)/(Naa+Nao+Nbb+Nbo+Nab+No)/2;
pb<-(2*Nbb+Nab+Nbo)/(Naa+Nao+Nbb+Nbo+Nab+No)/2;
po<-(2*No+Nao+Nbo)/(Naa+Nao+Nbb+Nbo+Nab+No)/2;
}

cbind(Naa,Nao,Nbb,Nbo,Nab,No)
cbind(pa,pb,po)

No comments: