# 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)
}


#exp 1 data
datafexp1<-read.table("datashortExp12AFC.txt")


#Experiment 1
newfit<- function(p,data,summary=0)
{
  if(any(p<0)) return(NA);

  a<-20
  sd<-10

  r1<-p[1] #rate for 1st letter of 8-letter string
  rA<-p[2] #rate for other letters of 8-letter string

  r1_7<-p[3]
  rA_7<-p[4]

  ee<-0.01 #motor error
  e<-0.01 #forgetting
  lambda<-1
  
  n<-length(data$V1)

  tars<-as.character(data$V5)
  fois<-as.character(data$V6)

  set8<-nchar(tars)==8
  set7<-nchar(tars)==7

  pk<-double(n)
  pkout<-double(sum(set8))
  pkout7<-double(sum(set7))

  dts<-rep(110,n)


  pars8<-c(a,sd,lambda,e,r1,rA,rA,rA,rA,rA,rA,rA)
  pars7<-c(a,sd,lambda,e,r1_7,rA_7,rA_7,rA_7,rA_7,rA_7,rA_7)

  
  pk[set8]<-.C("pk_vector",pk=pkout,pars8,tars[set8],fois[set8],dts,sum(set8))$pk #prob known
  pk[set7]<-.C("pk_vector",pk=pkout7,pars7,tars[set7],fois[set7],dts,sum(set7))$pk #prob known

  pg<-0.5 #prob guessed correctly

  pc<-(1-ee)*pk+pg*(1-(1-ee)*pk) #probability correct

  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
  }
}

optim(c(0.07,0.01,0.09,0.02),newfit,data=datafexp1,control=list(maxit=500,trace=T))->exp1


#exp1$par
#[1] 0.117941516 0.003657218 0.025701265 0.025088460


newfit(exp1$par,datafexp1,summary=1)$summary->predictionsexp1

with(subset(predictionsexp1,src=="model"),xtabs(V2~V3+V7)/xtabs(~V3+V7))
with(subset(predictionsexp1,src=="data"),xtabs(V2~V3+V7)/xtabs(~V3+V7))



##Experiment 1 with a repeated letter bias
newfitRL<- function(p,data,summary=0)
{
  if(any(p<0)) return(NA);

  a<-20
  sd<-10

  r1<-p[1] #rate for 1st letter of 8-letter string
  rA<-p[2] #rate for other letters of 8-letter string

  r1_7<-p[3]
  rA_7<-p[4]
  rl<-p[5]

  ee<-0.01 #motor error
  e<-0.01 #forgetting
  lambda<-1
  
  n<-length(data$V1)

  tars<-as.character(data$V5)
  fois<-as.character(data$V6)

  set8<-nchar(tars)==8
  set7<-nchar(tars)==7

  pk<-double(n)
  pkout<-double(sum(set8))
  pkout7<-double(sum(set7))

  dts<-rep(110,n)


  pars8<-c(a,sd,lambda,e,r1,rA,rA,rA,rA,rA,rA,rA)
  pars7<-c(a,sd,lambda,e,r1_7,rA_7,rA_7,rA_7,rA_7,rA_7,rA_7)

  
  pk[set8]<-.C("pk_vector",pk=pkout,pars8,tars[set8],fois[set8],dts,sum(set8))$pk #prob known
  pk[set7]<-.C("pk_vector",pk=pkout7,pars7,tars[set7],fois[set7],dts,sum(set7))$pk #prob known

  lapply(tars,function(x) table(strsplit(x,split="")))->tarIDs
  lapply(fois,function(x) table(strsplit(x,split="")))->foiIDs 
  mapply(function(x,y) names(x)[names(x)%in%names(y)],tarIDs,foiIDs,SIMPLIFY = FALSE)->commonlet
  mapply(function(x,y) x[y], tarIDs,commonlet,SIMPLIFY = FALSE)->tarcommon
  mapply(function(x,y) x[y], foiIDs,commonlet,SIMPLIFY = FALSE)->foicommon
  mapply(function(x,y) (x-y),tarcommon,foicommon,SIMPLIFY = FALSE)->diff
  ifelse(lapply(diff, function(x) sum(x)>0),"target",ifelse(lapply(diff, function(x) sum(x)<0),"foil","none"))->repin #repetition is in
  sapply(diff, function(x) abs(sum(x)))->nrl #number of exceeding repeated letters
  mapply(function(x,y) sum(as.numeric(abs(x-y))>0),tarcommon,foicommon)->nrID #number of repeated identities
  ifelse(repin!="none",ifelse(repin=="target",0.5+(nchar(tars)-nrl/nrID)/nchar(tars)*rl,
         0.5-(nchar(fois)-nrl/nrID)/nchar(fois)*rl),0.5)->pg #prob guessed correctly

  pc<-(1-ee)*pk+pg*(1-(1-ee)*pk) #probability correct

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

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

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


optim(c(0.07,0.01,0.09,0.02,0.2),newfitRL,data=datafexp1,control=list(maxit=500,trace=T))->exp1RL

#> exp1RL$par
#[1] 0.277747153 0.002676932 0.027027562 0.026540338 0.082483759


newfitRL(exp1RL$par,datafexp1,summary=1)$summary->predictionsexp1RL

with(subset(predictionsexp1RL,src=="modelRL"),xtabs(V2~V3+V7)/xtabs(~V3+V7))
with(subset(predictionsexp1RL,src=="data"),xtabs(V2~V3+V7)/xtabs(~V3+V7))


write.csv(rbind(predictionsexp1,predictionsexp1RL),file="../LTRSEXP1.txt",row.names=F,quote=F)



