R/03_RunIntLim.R

In Mathelab/IntLIM: Integration of Omics Data Using Linear Modeling

#' Run linear models and retrieve relevant statistics
#'
#' @include internalfunctions.R
#'
#' @param inputData MultiDataSet object (output of ReadData()) with gene expression,
#' metabolite abundances, and associated meta-data
#' @param stype column name that represents sample type (by default, it will be used 
#' in the interaction term). Only 2 categories are currently supported.
#' @param outcome 'metabolite' or 'gene' must be set as outcome/independent variable
#' (default is 'metabolite')
#' @param covar Additional variables from the phenotypic data that be integrated into linear model
#' @param class.covar Describing whether additional variables are 'numeric' or 'categorial'
#' @return IntLimModel object with model results
#'
#' @examples
#' dir <- system.file("extdata", package="IntLIM", mustWork=TRUE)
#' csvfile <- file.path(dir, "NCItestinput.csv")
#' mydata <- ReadData(csvfile,metabid='id',geneid='id')
#' \dontrun{
#' myres <- RunIntLim(mydata,stype="PBO_vs_Leukemia")
#' }
#' @export
RunIntLim <- function(inputData,stype=NULL,outcome="metabolite", covar=NULL, class.covar=NULL){

    if (class(inputData) != "MultiDataSet") {
        stop("input data is not a MultiDataSet class")
    }
    mytypes <- names(Biobase::assayData(inputData))
    if(!any(mytypes=="expression") || !any(mytypes=="metabolite")) {
        stop("input data must contain assayData of type 'metabolite' and 'expression.
        Try reading in the data with the ReadData function")
    }
    if (is.null(stype)) {stop("Please set the variable type (e.g. sample group)")}

#    incommon<-MultiDataSet::commonSamples(inputData)
#    mp <- as.character(Biobase::pData(incommon[["metabolite"]])[,stype])
#    gp <- as.character(Biobase::pData(incommon[["expression"]])[,stype])

    incommon <- getCommon(inputData,stype,covar,class.covar=class.covar)

    if(length(unique(stats::na.omit(incommon$p))) != 2) {
	stop(paste("IntLim currently requires only two categories.  Make sure the column",stype,"only has two unique values"))
    }

    print("Running the analysis on")
    print(table(incommon$p))

    ptm <- proc.time()
    myres <- RunLM(incommon,outcome=outcome,type=incommon$p,covar=covar)
    print(proc.time() - ptm)
    myres@stype=stype
    myres@outcome=outcome
    
    if(!is.null(covar)){
        covariate <- colnames(incommon$covar_matrix)
        class.var <- c()
        for(i in 1:length(covar)){
            class.var[i] <- class(incommon$covar_matrix[,i])
        }
        
        myres@covar <- data.frame(covariate,class.var)
    }
    return(myres)
}