Identify working design model among numerous covariates

select.model

#Identify working design model among numerous covariates
# Returns design that can be used directly by DESeq2
#Usage:
#      select.model(covariates,main,region)
#      covariates: character list of covariates in meta data matrix (column names of meta data matrix)
#      main: main factor forumula (can include covariates)
#      region: brain region (temprorarily hard-coded for the immediate analysis)
#Example:
#       select.model(c("Antidepressant","Alcool","History.of.Abuse","Cause.of.death","PMI"),
#                      "RIN+Age.+Gender+Phenotype+Gender:Phenotype",
#                       "BA11")  
#



#TODO: Generalize "run.model" nested function. 


select.model<-function(covariates,main,region){

 run.model<-function(region){
    library(DESeq2)
    counts.1<-counts[,meta_data$region==region]
    meta<-meta_data[meta_data$region==region,]
    dds<-DESeqDataSetFromMatrix(countData= counts.1,colData= meta,design=as.formula(paste("~",design)))
    dds<-DESeq(dds)
 
  }

 next.design=main
 for(order in covariates){ 
   
     design<-paste(order,next.design,sep="+")
     
     checkModel<-tryCatch(run.model(region),error=function(e) e)
     
     if(!inherits(checkModel,"error")){
         next.design<-paste(order,next.design,sep="+")
    }else{
      print("error")
    }
  }
  return(as.formula(paste("~",next.design)))
}