# Loads model.dll (32-bit Windows) or model.so (64-bit Linux)

if(!is.loaded("pk_vector")){
  if(.Platform$OS.type=="windows" && .Platform$r_arch!="i386) stop("On Windows, i386 version of R is needed for LTRS dll.")
  dry<-setwd("../LTRS")
  dyn.load(paste("model",sep="",.Platform$dynlib.ext))
  setwd(dry)
}



dataf<-read.table("datashortExp4SD.txt")#specify data file

  

samediffix<- function(p,data,summary=0)
{
  if(any(p<0)) return(NA);
  
  a<-38
  sd<-20
  r1<-p[1] #initial
  r2<-p[2] #all middle letters
  r3<-p[2] 
  r4<-p[2]
  r5<-p[2]
  r6<-p[2]
  r7<-p[2]
  r8<-p[3] #final

  ee<-0.01
  phi<-0.01
  lambda<-0.03
  ps<-p[4]

  n<-length(data$V1)
  pkout<-double(n)
  dts<-rep(120,n)
  tars<-as.character(data$V5)
  opts<-as.character(data$V6)

  pars<-c(a,sd,lambda,phi,r1,r2,r3,r4,r5,r6,r7,r8)
  
  pk<-.C("pk_vector",pk=pkout,pars,tars,opts,dts,n)$pk

  pc<-ifelse(tars==opts,    
      (1-ee)*pk+(ps)*(1-(1-ee)*pk),
      (1-ee)*pk+(1-ps)*(1-(1-ee)*pk))


  sse<-sum((data$V2-pc)^2)

  if(summary>0)
  {
    tempd<-data
    tempd$V2<-pc

    list(fit=sse,
      summary=rbind(data.frame(src="model",tempd),
        data.frame(src="data",data))
    )
  }
  else
  {
    sse
  }
}



pars<-c(0.02,0.001,0.2,0.78)

optim(pars,samediffix,data=dataf,control=list(maxit=1000,trace=T))->sdfix #0.2499375


sdfix

sdfix$par


samediffix(sdfix$par,dataf,summary=1)$summary->predictionsfix

with(subset(predictionsfix,src=="model"),xtabs(V2~V3+V4+V8)/xtabs(~V3+V4+V8))
with(subset(predictionsfix,src=="data"),xtabs(V2~V3+V4+V8)/xtabs(~V3+V4+V8))



write.csv(predictionsfix,file="../LTRSEXP4.txt",row.names=F,quote=F)



