R/gen_result.R

In ChenyueLiu/ivqr: Perform inverse quantile regression by solving mixed integer linear programming

#' Generate a table containing the estimation of \eqn{\beta}_D using various \eqn{\tau} and plot \eqn{\beta} against \eqn{\tau};
#' 
#' @param data a list containing four variables Y, D, Z, X
#' @param lpsolver "gurobi","cplexapi","lpsolveapi"
#' @return a named vector containing the estimation of \eqn{\beta}_D 
#' @examples 
#' angrist_data=angrist(n=100,pD=3)
#' data=list()
#' data$Y=angrist_data$Y
#' data$Z=angrist_data$Z
#' data$D=angrist_data$D
#' data$X=angrist_data$X[,1:3,drop=FALSE]
#' reg_out=gen_result(data,"gurobi")
 
gen_result<-function(data,lpsolver=NULL){
  #### error checking for dataset TBA ####
if(is.null(data)){
  stop("please specify a dataset",call.=FALSE)
}
  
Y = data$Y; D = data$D; Z = data$Z; X = data$X
results<-c()
for(tau in seq(0.3,0.7,by=0.025)){
  fit_full = iqr_milp(Y, D, X, Z, tau=tau, O_pos= NULL, O_neg= NULL, TimeLimit = 1440,lpsolver=lpsolver)
  if(is.null(fit_full$B_D)){
    stop(paste("no solution for tau = ",tau),call.=FALSE)
  }
  results<-cbind(results,fit_full$B_D)
}
if(is.null(results)){
  stop("no solution for this dataset",call.=FALSE)
}
row.names(results)<-paste0("B_D",1:pD)
colnames(results)<-paste0("tau",seq(0.3,0.7,by=0.025))

par(mfrow=c(1,pD))
for(i in 1:pD){
  plot(seq(0.3,0.7,by=0.025),results[i,],type="l",xlab=expression(tau),ylab=bquote('B_D('~tau~')_'~.(i)))
  #expression(paste("B_D(",tau,")_"))
}
return(results)
}