20170210

## for the ask of https://www.ptt.cc/bbs/R_Language/M.1486659291.A.337.html

# t
          [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
 [1,] 2.947500    2    3    1   NA  2.5    1    0    0
 [2,] 4.084530    4    5    1   NA  4.5    1    0    0
 [3,] 8.808208    8    9    1   NA  8.5    1    0    0
 [4,] 6.712504    6    7    1   NA  6.5    1    0    0
 [5,] 8.503359    8    9    1   NA  8.5    1    0    0
 [6,] 6.304814    6    7    1   NA  6.5    1    0    0
 [7,] 8.348371    8    9    1   NA  8.5    1    0    0
 [8,] 6.647538    6    7    1   NA  6.5    1    0    0
 [9,] 4.227048    4    5    1   NA  4.5    1    0    0
[10,] 4.270762    4    5    1   NA  4.5    1    0    0
[11,] 6.573022    6    7    1   NA  6.5    1    0    0
[12,] 6.164170    6    7    1   NA  6.5    1    0    0
[13,] 4.439435    4    5    1   NA  4.5    1    0    0
[14,] 3.662537    3    4    1   NA  3.5    1    0    0
[15,] 1.979395    1    2    1   NA  1.5    1    0    0
[16,] 6.143496    6    7    1   NA  6.5    1    0    0
[17,] 7.382536    7    8    1   NA  7.5    1    0    0
[18,] 7.260993    7    8    1   NA  7.5    1    0    0
[19,] 4.436384    4    5    1   NA  4.5    1    0    0
[20,] 6.629455    6    7    1   NA  6.5    1    0    0

surlikelihood<-function(x,t){
  p<-x[1]
  alpha<-x[2]
  beta<-x[3]
  logL<-numeric(5)
  for(i in 0:4){
    if(i==0){
      sright<-(p+(1-p)*(1+(t[t[,4]==i,2]/alpha)^beta)^-1)
      logL[i+1]<-sum(log(sright))
    }
    if(i==1){
      sinterval<-((1-p)*((1+(t[t[,4]==i,2]/alpha)^beta)^-1-(1+(t[t[,4]==i,3]/alpha)^beta)^-1))
      logL[i+1]<-sum(log(sinterval))
    }
    if(i==2){
      sleft<-(p+(1-p)*(1-(1+(t[t[,4]==i,3]/alpha)^beta)^-1))
      logL[i+1]<-sum(log(sleft))
    }
    if(i==3){
      suncensored<-(1-p)*(((beta/alpha)*(t[t[,4]==i,1]/alpha)^(beta-1))/(1+(t[t[,4]==i,1]/alpha)^beta)^2)
      logL[i+1]<-sum(log(suncensored))
    }
    if(i==4){
      scure<-(p+(1-p)*(1+(t[t[,4]==i,2]/alpha)^beta)^-1)
      logL[i+1]<-sum(log(scure))
    }
  }
  return(-1*sum(logL))
}