#
#    USe try() and tryCatch() to loop going despite the error
#
#====================================================================

x <- rnorm(100, 2, 1)
x
which(x<0)

log(x)


#--- Define nonsensical function that gives error for neg x ---------
log.a <- function(x){   
  ifelse(x<0, a <- log("a"), a<- log(x))
  return(a)
}



#--- This loop stops where x<0, because log.a() gives error -----------
T <- 0
for (i in 1:100) {
  print(i)
  T[i] <- log.a(x[i])
}



#--- You can Use try() to keep the loop going --------------------------
T <- 0
for (i in 1:100) {
  print(i)
  T[i] <- try( log.a(x[i]), TRUE )    #returns TRUE when log.a() gives error
}
T             #loop went on despite of the error
mean(T)       #we still can't do this b/c of error error



#--- Use tryCatch() to get better output for error case ----------------
T <- 0
for (i in 1:100) {
  print(i)
  T[i] <- tryCatch( log.a(x[i]), error=function(ss){return(NA)}  )   # put NA if error is given
}
T
mean(T)                      #loop went on despite of the error
mean(  T[(is.na(T)==0)]  )   #now we can do this



#--- Or do this to see the error message -------------------------------
T <- 0
for (i in 1:100) {
  T[i] <- tryCatch( 
            log.a(x[i]),  
            error = function(error_condition) {
              #- See what error message was 
              cat(paste("\n  error msg when i=", i, ": \n  ", error_condition$message, "\n\n"))  
              return(NA)
	          }
          )
}
T
mean(T)                           #loop went on despite of the error
mean(  T[which(is.na(T)==0)]  )   #now we can do this




#--- While loop using tryCatch()  -------------------------------
n   = 100
itt = 50
MLE  = matrix(0, itt, 9)  
Vars = matrix(0, itt, 9)  

i=1
j=1
while (i <= itt) {
  X  <- arima.sim(n=n, 
                  list(ar = c(1.5051,  -0.7517,  -0.4294,  0.4672), 
                       ma = c(-2.6066,  2.8921,  -1.6503,  0.4296)  ) ) + 8
  Est <- tryCatch(
    Arima(X, order=c(4,0,4)),
    error = function(ss){return(NA)}  
  )  
  if ( length(is.na(Est))>1) {
    MLE[i,]  <- Est$coef
    Vars[i,] <- diag(Est$var.coef)
    i=i+1
  } else {
    j=j+1
  }
}

Result <- cbind( apply(MLE, 2, mean),
                 apply(MLE, 2, sd),
                 sqrt(apply(Vars, 2, mean)))
colnames(Result) <- c("Mean",  "SE by Sim",  "SE from Output")
Result