生成確率の計算

• このときのつづきで
• ディプロタイプからどんなハプロタイプが作られるか、その生成確率を求める
• rがその生成確率
  • 最後の確率の計算を行列で一気に計算するように直せるだろうか
• ソースが正確かどうかまだ確かめていない
• 間違っているので修正が必要
  • 修正したので改良する

i1<-sample(1:2,1);i2<-sample(1:2,1);i3<-sample(1:2,1);j1<-sample(1:2,1);j2<-sample(1:2,1);j3<-sample(1:2,1);k1<-sample(1:2,1);k2<-sample(1:2,1);k3<-sample(1:2,1)

print(paste(i1,i2,i3,j1,j2,j3,k1,k2,k3))

hap1<-c(i1,i2,i3)
hap2<-c(j1,j2,j3)
hap<-c(k1,k2,k3)

#s12<-0.1
#s23<-0.2

v1<-hap1==hap
v2<-hap2==hap

D<-matrix(c(v1,v2),byrow="TRUE",nrow=2)
#D<-matrix(c(v1,v2),byrow="FALSE",ncol=2)

s<-rep(0.1,length=(length(v1)-1))     #各ローカス間の組換え確率
if(min(v1 + v2) > 0){
	r1<-v1[1]/2
	r2<-v2[1]/2
	rs<-c(r1,r2)
for ( i in 1:(length(v1)-1)){
	
	M<-matrix(c(v1[i]*(1-s[i])*v1[i+1],v2[i]*s[i]*v1[i+1],v1[i]*s[i]*v2[i+1],v2[i]*(1-s[i])*v2[i+1]),byrow="TRUE",nrow=2)
	rs<-M%*%rs
	}
	
	r<-sum(rs)
}else{
r<-0
}
print(v1)
print(v2)
print(r)