R/Pt.R

In DNAseqtest: Generating and Testing DNA Sequences

Pt <-
function(S,Pi,t){
ax<-matrix(0,nrow=4,ncol=4)
fx<-matrix(0,nrow=4,ncol=4)
Rx<-S%*%Pi
a1<-eigen((sqrt(Pi))%*%Rx%*%(solve(sqrt(Pi))),symmetric=T)$values
a2<-eigen((sqrt(Pi))%*%Rx%*%(solve(sqrt(Pi))),symmetric=T)$vectors
 for(j in 1:4){
 ax<-exp(a1[j]*(t))*a2[,j]%*%t(a2[,j])
 fx<-fx+ax
 }
ptx<-(solve(sqrt(Pi)))%*%fx%*%(sqrt(Pi))
ptx
}