@@ -23,6 +23,7 @@

 ^src/data_structures/Vector.o

 ^src/file_parser/FileParser.o

 ^src/file_parser/CsvParser.o

+^src/file_parser/GctParser.o

 ^src/file_parser/TsvParser.o

 ^src/file_parser/MtxParser.o

 ^src/math/Algorithms.o

@@ -1,5 +1,5 @@

 Package: CoGAPS

-Version: 3.3.5

+Version: 3.3.6

 Date: 2018-04-24

 Title: Coordinated Gene Activity in Pattern Sets

 Author: Thomas Sherman, Wai-shing Lee, Conor Kelton, Ondrej Maxian, Jacob Carey,

@@ -49,8 +49,10 @@ Collate:

     'class-CogapsParams.R'

     'CoGAPS.R'

     'DistributedCogaps.R'

+    'HelperFunctions.R'

     'Package.R'

     'RcppExports.R'

+    'SubsetData.R'

     'class-CogapsResult.R'

     'methods-CogapsParams.R'

     'methods-CogapsResult.R'

@@ -9,15 +9,20 @@ export(calcCoGAPSStat)

 export(calcGeneGSStat)

 export(calcZ)

 export(computeGeneGSProb)

+export(getClusteredPatterns)

+export(getCorrelationToMeanPattern)

 export(getMeanChiSq)

 export(getOriginalParameters)

 export(getParam)

+export(getSubsets)

+export(getUnmatchedPatterns)

 export(getVersion)

 export(patternMarkers)

 export(plotPatternMarkers)

 export(plotResiduals)

 export(reconstructGene)

 export(scCoGAPS)

+export(setAnnotationWeights)

 export(setDistributedParams)

 export(setParam)

 exportClasses(CogapsParams)

@@ -45,6 +50,7 @@ importFrom(methods,"slot<-")

 importFrom(methods,callNextMethod)

 importFrom(methods,is)

 importFrom(methods,new)

+importFrom(methods,show)

 importFrom(methods,slot)

 importFrom(methods,slotNames)

 importFrom(methods,validObject)

@@ -1,51 +1,17 @@

 #' @include class-CogapsParams.R

 NULL

-#' Checks if file is supported

-#' @param file path to file

-#' @return TRUE if file is supported, FALSE if not

-#' @importFrom tools file_ext

-supported <- function(file)

-{

-    if (!is(file, "character"))

-        return(FALSE)

-    return(tools::file_ext(file) %in% c("tsv", "csv", "mtx"))

-}

-

-#' get number of rows from supported file name or matrix

-#' @param data either a file name or a matrix

-#' @return number of rows

-#' @importFrom data.table fread

-#' @importFrom tools file_ext

-nrow_helper <- function(data)

-{

-    if (is(data, "character"))

-    {

-        return(switch(tools::file_ext(data),

-            "csv" = nrow(data.table::fread(data, select=1)),

-            "tsv" = nrow(data.table::fread(data, select=1)),

-            "mtx" = as.numeric(data.table::fread(data, nrows=1, fill=TRUE)[1,1])

-        ))

-    }

-    return(nrow(data))

-}

-

-#' get number of columns from supported file name or matrix

-#' @param data either a file name or a matrix

-#' @return number of columns

-#' @importFrom data.table fread

-#' @importFrom tools file_ext

-ncol_helper <- function(data)

+#' Information About Package Compilation

+#' @export

+#'

+#' @details returns information about how the package was compiled, i.e. which

+#'  compiler/version was used, which compile time options were enabled, etc...

+#' @return string containing build report

+#' @examples

+#' CoGAPS::buildReport()

+buildReport <- function()

 {

-    if (is(data, "character"))

-    {

-        return(switch(tools::file_ext(data),

-            "csv" = ncol(data.table::fread(data, nrows=1)) - 1,

-            "tsv" = ncol(data.table::fread(data, nrows=1)) - 1,

-            "mtx" = as.numeric(data.table::fread(data, nrows=1, fill=TRUE)[1,2])

-        ))

-    }

-    return(ncol(data))

+    getBuildReport_cpp()

 }

 #' CoGAPS Matrix Factorization Algorithm

@@ -73,6 +39,11 @@ ncol_helper <- function(data)

 #' for data that is stored as samples x genes since CoGAPS requires data to be

 #' genes x samples

 #' @param BPPARAM BiocParallel backend 

+#' @param geneNames vector of names of genes in data

+#' @param sampleNames vector of names of samples in data

+#' @param matchedPatterns manually matched patterns for distributed CoGAPS

+#' @param outputToFile name of a file to save the output to, will create 4 files

+#' of the form "filename_nPatterns_[Amean, Asd, Pmean, Psd].extension"

 #' @param ... allows for overwriting parameters in params

 #' @return CogapsResult object

 #' @examples

@@ -94,9 +65,12 @@ ncol_helper <- function(data)

 CoGAPS <- function(data, params=new("CogapsParams"), nThreads=1,

 messages=TRUE, outputFrequency=500, uncertainty=NULL,

 checkpointOutFile="gaps_checkpoint.out", checkpointInterval=1000,

-checkpointInFile=NULL, transposeData=FALSE, BPPARAM=NULL, ...)

+checkpointInFile=NULL, transposeData=FALSE, BPPARAM=NULL,

+geneNames=NULL, sampleNames=NULL, matchedPatterns=NULL,

+outputToFile=NULL, ...)

 {

     # store all parameters in a list and parse parameters from ...

+    validObject(params)

     allParams <- list("gaps"=params,

         "nThreads"=nThreads,

         "messages"=messages,

@@ -106,12 +80,14 @@ checkpointInFile=NULL, transposeData=FALSE, BPPARAM=NULL, ...)

         "checkpointInFile"=checkpointInFile,

         "transposeData"=transposeData,

         "bpBackend"=BPPARAM,

+        "matchedPatterns"=matchedPatterns,

+        "outputToFile"=outputToFile,

         "whichMatrixFixed"=NULL # internal parameter

     )

     allParams <- parseExtraParams(allParams, list(...))

     # display start up message for the user

-    startupMessage(data, allParams$transposeData, allParams$gaps@distributed)

+    startupMessage(data, allParams)

     # check file extension

     if (is(data, "character") & !supported(data))

@@ -125,7 +101,7 @@ checkpointInFile=NULL, transposeData=FALSE, BPPARAM=NULL, ...)

     if (!is(data, "character") & !is.null(uncertainty) & !is(uncertainty, "matrix"))

         stop("uncertainty must be a matrix unless data is a file path")

     if (!is(data, "character"))

-        checkDataMatrix(data, uncertainty, params)

+        checkDataMatrix(data, uncertainty, allParams$gaps)

     # check single cell parameter

     if (!is.null(allParams$gaps@distributed))

@@ -136,16 +112,15 @@ checkpointInFile=NULL, transposeData=FALSE, BPPARAM=NULL, ...)

         stop("can't run multi-threaded and distributed CoGAPS at the same time")

     # convert data to matrix

+    if (is(data, "matrix"))

+        data <- data

     if (is(data, "data.frame"))

         data <- data.matrix(data)

     else if (is(data, "SummarizedExperiment"))

         data <- SummarizedExperiment::assay(data, "counts")

     else if (is(data, "SingleCellExperiment"))

         data <- SummarizedExperiment::assay(data, "counts")

-

-    # label matrix

-    

-

+   

     # determine which function to call cogaps algorithm

     if (!is.null(allParams$gaps@distributed))

         dispatchFunc <- distributedCogaps # genome-wide or single-cell cogaps

@@ -160,9 +135,24 @@ checkpointInFile=NULL, transposeData=FALSE, BPPARAM=NULL, ...)

         if (!is.null(allParams$gaps@distributed))

             stop("checkpoints not supported for distributed cogaps")

         else

-            message("Running CoGAPS from a checkpoint")

+            cat("Running CoGAPS from a checkpoint")

     }

+    # get gene/sample names

+    if (is.null(geneNames)) geneNames <- getGeneNames(data, allParams$transposeData)

+    if (is.null(sampleNames)) sampleNames <- getSampleNames(data, allParams$transposeData)

+

+    nGenes <- ifelse(allParams$transposeData, ncolHelper(data), nrowHelper(data))

+    nSamples <- ifelse(allParams$transposeData, nrowHelper(data), ncolHelper(data))

+

+    if (length(geneNames) != nGenes)

+        stop("incorrect number of gene names given")

+    if (length(sampleNames) != nSamples)

+        stop("incorrect number of sample names given")

+

+    allParams$geneNames <- geneNames

+    allParams$sampleNames <- sampleNames

+

     # run cogaps

     gapsReturnList <- dispatchFunc(data, allParams, uncertainty)

@@ -172,8 +162,9 @@ checkpointInFile=NULL, transposeData=FALSE, BPPARAM=NULL, ...)

         Asd         = gapsReturnList$Asd,

         Pmean       = gapsReturnList$Pmean,

         Psd         = gapsReturnList$Psd,

-        seed        = gapsReturnList$seed,

         meanChiSq   = gapsReturnList$meanChiSq,

+        geneNames   = gapsReturnList$geneNames,

+        sampleNames = gapsReturnList$sampleNames,

         diagnostics = append(gapsReturnList$diagnostics,

             list("params"=allParams$gaps, "version"=utils::packageVersion("CoGAPS")))

     ))

@@ -189,13 +180,14 @@ checkpointInFile=NULL, transposeData=FALSE, BPPARAM=NULL, ...)

 scCoGAPS <- function(data, params=new("CogapsParams"), nThreads=1,

 messages=TRUE, outputFrequency=500, uncertainty=NULL,

 checkpointOutFile="gaps_checkpoint.out", checkpointInterval=1000,

-checkpointInFile=NULL, transposeData=FALSE, BPPARAM=NULL, ...)

+checkpointInFile=NULL, transposeData=FALSE, BPPARAM=NULL,

+geneNames=NULL, sampleNames=NULL, matchedPatterns=NULL, ...)

 {

     params@distributed <- "single-cell"

     params@singleCell <- TRUE

     CoGAPS(data, params, nThreads, messages, outputFrequency, uncertainty,

         checkpointOutFile, checkpointInterval, checkpointInFile, transposeData,

-        BPPARAM, ...)

+        BPPARAM, geneNames, sampleNames, matchedPatterns, ...)

 }

 #' Genome Wide CoGAPS

@@ -208,83 +200,11 @@ checkpointInFile=NULL, transposeData=FALSE, BPPARAM=NULL, ...)

 GWCoGAPS <- function(data, params=new("CogapsParams"), nThreads=1,

 messages=TRUE, outputFrequency=500, uncertainty=NULL,

 checkpointOutFile="gaps_checkpoint.out", checkpointInterval=1000,

-checkpointInFile=NULL, transposeData=FALSE, BPPARAM=NULL, ...)

+checkpointInFile=NULL, transposeData=FALSE, BPPARAM=NULL,

+geneNames=NULL, sampleNames=NULL, matchedPatterns=NULL, ...)

 {

     params@distributed <- "genome-wide"

     CoGAPS(data, params, nThreads, messages, outputFrequency, uncertainty,

         checkpointOutFile, checkpointInterval, checkpointInFile, transposeData,

-        BPPARAM, ...)

-}   

-

-#' write start up message

-#'

-#' @param data data set

-#' @param transpose if we are transposing the data set

-#' @param distributed if we are running distributed CoGAPS

-#' @return message displayed to screen

-startupMessage <- function(data, transpose, distributed)

-{

-    nGenes <- ifelse(transpose, ncol_helper(data), nrow_helper(data))

-    nSamples <- ifelse(transpose, nrow_helper(data), ncol_helper(data))

-

-    dist_message <- "Standard"

-    if (!is.null(distributed))

-        dist_message <- distributed

-    message(paste("Running", dist_message, "CoGAPS on", nGenes, "genes and",

-        nSamples, "samples"))

-}

-

-#' parse parameters passed through the ... variable

-#'

-#' @param allParams list of all parameters

-#' @param extraParams list of parameters in ...

-#' @return allParams with any valid parameters in extraParams added

-#' @note will halt with an error if any parameters in extraParams are invalid

-#' @importFrom methods slotNames

-parseExtraParams <- function(allParams, extraParams)

-{

-    # parse direct params

-    for (s in slotNames(allParams$gaps))

-    {

-        if (!is.null(extraParams[[s]]))

-        {

-            allParams$gaps <- setParam(allParams$gaps, s, extraParams[[s]])

-            extraParams[[s]] <- NULL

-        }

-    }

-

-    # check for unrecognized options

-    if (length(extraParams) > 0)

-        stop(paste("unrecognized argument:", names(extraParams)[1]))

-

-    return(allParams)

-}

-

-#' check that provided data is valid

-#'

-#' @param data data matrix

-#' @param uncertainty uncertainty matrix, can be null

-#' @param params CogapsParams object

-#' @return throws an error if data has problems

-checkDataMatrix <- function(data, uncertainty, params)

-{

-    if (sum(data < 0) > 0 | sum(uncertainty < 0) > 0)

-        stop("negative values in data and/or uncertainty matrix")

-    if (nrow(data) <= params@nPatterns | ncol(data) <= params@nPatterns)

-        stop("nPatterns must be less than dimensions of data")

-    if (sum(uncertainty < 1e-5) > 0)

-        warning("small values in uncertainty matrix detected")

-}

-

-#' Information About Package Compilation

-#' @export

-#'

-#' @details returns information about how the package was compiled, i.e. which

-#'  compiler/version was used, which compile time options were enabled, etc...

-#' @return string containing build report

-#' @examples

-#' CoGAPS::buildReport()

-buildReport <- function()

-{

-    getBuildReport_cpp()

-}

+        BPPARAM, geneNames, sampleNames, matchedPatterns, ...)

+}

\ No newline at end of file

@@ -1,4 +1,5 @@

 #' CoGAPS Distributed Matrix Factorization Algorithm

+#' @keywords internal

 #'

 #' @description runs CoGAPS over subsets of the data and stitches the results

 #' back together

@@ -10,13 +11,20 @@

 #' @importFrom BiocParallel bplapply MulticoreParam

 distributedCogaps <- function(data, allParams, uncertainty)

 {

-    FUN <- function(index, sets, data, allParams, uncertainty, fixedMatrix=NULL)

+    FUN <- function(index, sets, data, allParams, uncertainty, geneNames,

+    sampleNames, fixedMatrix=NULL)

     {

+        if (allParams$gaps@distributed == "genome-wide")

+            geneNames <- geneNames[sets[[index]]]

+        else

+            sampleNames <- sampleNames[sets[[index]]]

+

         internal <- ifelse(is(data, "character"), cogaps_cpp_from_file, cogaps_cpp)

         raw <- internal(data, allParams, uncertainty, sets[[index]],

             fixedMatrix, index == 1)

         new("CogapsResult", Amean=raw$Amean, Asd=raw$Asd, Pmean=raw$Pmean,

-            Psd=raw$Psd, seed=raw$seed, meanChiSq=raw$meanChiSq)

+            Psd=raw$Psd, meanChiSq=raw$meanChiSq, geneNames=geneNames,

+            sampleNames=sampleNames)

     }

     # randomly sample either rows or columns into subsets to break the data up

@@ -25,71 +33,78 @@ distributedCogaps <- function(data, allParams, uncertainty)

     if (is.null(allParams$bpBackend))

         allParams$bpBackend <- BiocParallel::MulticoreParam(workers=length(sets))

-    # run Cogaps normally on each subset of the data

-    if (allParams$messages)

-        message("Running Across Subsets...")

-    initialResult <- bplapply(1:length(sets), FUN, BPPARAM=allParams$bpBackend,

-        sets=sets, data=data, allParams=allParams, uncertainty=uncertainty)

-

-    # get all unmatched patterns

-    if (allParams$gaps@distributed == "genome-wide")

-        unmatchedPatterns <- lapply(initialResult, function(x) x@sampleFactors)

+    if (is.null(allParams$matchedPatterns))

+    {

+        # run Cogaps normally on each subset of the data

+        if (allParams$messages)

+            cat("Running Across Subsets...\n\n")

+        initialResult <- bplapply(1:length(sets), FUN, BPPARAM=allParams$bpBackend,

+            sets=sets, data=data, allParams=allParams, uncertainty=uncertainty,

+            geneNames=allParams$geneNames, sampleNames=allParams$sampleNames)

+

+        # get all unmatched patterns

+        if (allParams$gaps@distributed == "genome-wide")

+            unmatchedPatterns <- lapply(initialResult, function(x) x@sampleFactors)

+        else

+            unmatchedPatterns <- lapply(initialResult, function(x) x@featureLoadings)

+

+        # match patterns in either A or P matrix

+        if (allParams$messages)

+            cat("\nMatching Patterns Across Subsets...\n")

+        matchedPatterns <- findConsensusMatrix(unmatchedPatterns, allParams)

+        allParams$gaps@nPatterns <- ncol(matchedPatterns$consensus)

+

+        # set fixed matrix

+        allParams$whichMatrixFixed <- ifelse(allParams$gaps@distributed

+            == "genome-wide", "P", "A")

+    }

     else

-        unmatchedPatterns <- lapply(initialResult, function(x) x@featureLoadings)

-

-    # match patterns in either A or P matrix

-    if (allParams$messages)

-        message("Matching Patterns Across Subsets...")

-    matchedPatterns <- findConsensusMatrix(unmatchedPatterns, allParams)

-    allParams$gaps@nPatterns <- ncol(matchedPatterns$consensus)

-

-    # set fixed matrix

-    allParams$whichMatrixFixed <- ifelse(allParams$gaps@distributed

-        == "genome-wide", "P", "A")

-

+    {

+        matchedPatterns <- list(consensus=allParams$matchedPatterns)

+        allParams$gaps@nPatterns <- ncol(matchedPatterns$consensus)

+        allParams$whichMatrixFixed <- ifelse(allParams$gaps@distributed

+            == "genome-wide", "P", "A")

+    }

+        

     # run final phase with fixed matrix

     if (allParams$messages)

-        message("Running Final Stage...")

+        cat("Running Final Stage...\n\n")

     finalResult <- bplapply(1:length(sets), FUN, BPPARAM=allParams$bpBackend,

         sets=sets, data=data, allParams=allParams, uncertainty=uncertainty,

+        geneNames=allParams$geneNames, sampleNames=allParams$sampleNames,    

         fixedMatrix=matchedPatterns$consensus)

     # concatenate final result

     fullResult <- stitchTogether(finalResult, allParams)

     # add diagnostic information before returning

-    fullResult$diagnostics$unmatchedPatterns <- unmatchedPatterns

-    fullResult$diagnostics$clusteredPatterns <- matchedPatterns$clusteredPatterns

-    fullResult$diagnostics$RtoMeanPattern <- lapply(matchedPatterns$clusteredPatterns, correlationToMeanPattern)

-    fullResult$diagnostics$subsets <- sets

-    return(fullResult)

-}

-

-#' partition genes/samples into subsets

-#' @description either genes or samples or partitioned depending on the type

-#' of distributed CoGAPS (i.e. genome-wide or single-cell)

-#' @param data either file name or matrix

-#' @param allParams list of all CoGAPS parameters

-#' @return list of sorted subsets of either genes or samples

-createSets <- function(data, allParams)

-{

-    total <- ifelse(xor(allParams$transposeData, allParams$gaps@distributed == "genome-wide"),

-        nrow_helper(data), ncol_helper(data))

-    setSize <- floor(total / allParams$gaps@nSets)

-

-    sets <- list()

-    remaining <- 1:total

-    for (n in 1:(allParams$gaps@nSets - 1))

+    if (is.null(allParams$matchedPatterns))

     {

-        selected <- sample(remaining, setSize, replace=FALSE)

-        sets[[n]] <- sort(selected)

-        remaining <- setdiff(remaining, selected)

+        fullResult$diagnostics$unmatchedPatterns <- unmatchedPatterns

+        fullResult$diagnostics$clusteredPatterns <- matchedPatterns$clusteredPatterns

+        fullResult$diagnostics$CorrToMeanPattern <- lapply(matchedPatterns$clusteredPatterns, corrToMeanPattern)

     }

-    sets[[allParams$gaps@nSets]] <- sort(remaining)

-    return(sets)

+

+    if (allParams$gaps@distributed == "genome-wide")

+        allNames <- allParams$geneNames

+    else

+        allNames <- allParams$sampleNames

+    fullResult$diagnostics$subsets <- lapply(sets, function(s) allNames[s])

+

+    # rename genes/samples if dimension was subsetted incompletely

+    allUsedIndices <- sort(unlist(sets))

+    if (allParams$gaps@distributed == "genome-wide")

+        fullResult$geneNames <- allParams$geneNames[allUsedIndices]

+    else

+        fullResult$sampleNames <- allParams$sampleNames[allUsedIndices]

+

+    return(fullResult)

 }

+

 #' find the consensus pattern matrix across all subsets

+#' @keywords internal

+#'

 #' @param unmatchedPatterns list of all unmatched pattern matrices from initial

 #' run of CoGAPS

 #' @param allParams list of all CoGAPS parameters

@@ -103,73 +118,84 @@ findConsensusMatrix <- function(unmatchedPatterns, allParams)

 }

 #' Match Patterns Across Multiple Runs

+#' @keywords internal

+#'

 #' @param allPatterns matrix of patterns stored in the columns

 #' @param allParams list of all CoGAPS parameters

 #' @return a matrix of consensus patterns

 #' @importFrom stats weighted.mean

 patternMatch <- function(allPatterns, allParams)

 {

-    PatsByClust <- corcut(allPatterns, allParams)

+    # cluster patterns

+    clusters <- corcut(allPatterns, allParams$gaps@cut, allParams$gaps@minNS)

-    # split by maxNS

-    indx <- which(sapply(PatsByClust, function(x) ncol(x) > allParams$gaps@maxNS))

-    while (length(indx) > 0)

-    { 

-        allParams$gaps@cut <- 2

-        internalPatsByClust <- corcut(PatsByClust[[indx[1]]], allParams)

-

-        PatsByClust[[indx[1]]] <- internalPatsByClust[[1]]

-        if (length(internalPatsByClust) > 1)

-        {

-            PatsByClust <- append(PatsByClust, internalPatsByClust[2])

-        }

-        indx <- which(sapply(PatsByClust, function(x) ncol(x) > allParams$gaps@maxNS))

+    # function to split a cluster in two (might fail to do so)

+    splitCluster <- function(list, index, minNS)

+    {

+        split <- corcut(list[[index]], 2, minNS)

+        list[[index]] <- split[[1]]

+        if (length(split) > 1)

+            list <- append(list, split[2])

+        return(list)

     }

-    # create matrix of mean patterns - weighted by coefficient of determination

-    PatsByCDSWavg <- sapply(PatsByClust, function(clust)

-        apply(clust, 1, function(row) weighted.mean(row, correlationToMeanPattern(clust)^3)))

-    colnames(PatsByCDSWavg) <- paste("Pattern", 1:length(PatsByClust))

+    # split large clusters into two

+    tooLarge <- function(x) ncol(x) > allParams$gaps@maxNS

+    indx <- which(sapply(clusters, tooLarge))

+    while (length(indx) > 0)

+    {

+        clusters <- splitCluster(clusters, indx[1], allParams$gaps@minNS)

+        indx <- which(sapply(clusters, tooLarge))

+    }

+    names(clusters) <- as.character(1:length(clusters))

-    # scale

-    return(list("clusteredPatterns"=PatsByClust,

-        "consensus"=apply(PatsByCDSWavg, 2, function(col) col / max(col))))

+    # create matrix of mean patterns - weighted by correlation to mean pattern

+    meanPatterns <- sapply(clusters, function(clust) apply(clust, 1,

+        function(row)  weighted.mean(row, corrToMeanPattern(clust)^3)))

+    colnames(meanPatterns) <- paste("Pattern", 1:length(clusters))

+    # returned patterns after scaling max to 1

+    return(list("clusteredPatterns"=clusters,

+        "consensus"=apply(meanPatterns, 2, function(col) col / max(col))))

 }

-correlationToMeanPattern <- function(cluster)

+#' calculate correlation of each pattern in a cluster to the cluster mean

+#' @keywords internal

+corrToMeanPattern <- function(cluster)

 {

     meanPat <- rowMeans(cluster)

     sapply(1:ncol(cluster), function(j) round(cor(x=cluster[,j], y=meanPat), 3))

 }

 #' cluster patterns together

+#' @keywords internal

+#'

 #' @param allPatterns matrix of all patterns across subsets

-#' @param allParams list of all CoGAPS parameters

+#' @param cut number of branches at which to cut dendrogram

+#' @param minNS minimum of individual set contributions a cluster must contain

 #' @return patterns listed by which cluster they belong to

 #' @importFrom cluster agnes

 #' @importFrom stats cutree as.hclust cor

-corcut <- function(allPatterns, allParams)

+corcut <- function(allPatterns, cut, minNS)

 {

     corr.dist <- cor(allPatterns)

     corr.dist <- 1 - corr.dist

-    clust <- cluster::agnes(x=corr.dist, diss=TRUE, "complete")

-    patternIds <- stats::cutree(stats::as.hclust(clust), k=allParams$gaps@cut)

+    clusterSummary <- cluster::agnes(x=corr.dist, diss=TRUE, "complete")

+    clusterIds <- stats::cutree(stats::as.hclust(clusterSummary), k=cut)

-    PatsByClust <- list()

-    for (cluster in unique(patternIds))

+    clusters <- list()

+    for (id in unique(clusterIds))

     {

-        if (sum(patternIds==cluster) >= allParams$gaps@minNS)

-        {

-            clusterPats <- allPatterns[,patternIds==cluster]

-            PatsByClust[[as.character(cluster)]] <- clusterPats

-        }

+        if (sum(clusterIds==id) >= minNS)

+            clusters[[as.character(id)]] <- allPatterns[,clusterIds==id,drop=FALSE]

     }

-    return(PatsByClust)

+    return(clusters)

 }

 #' concatenate final results across subsets

+#' @keywords internal

+#'

 #' @param result list of CogapsResult object from all runs across subsets

 #' @param allParams list of all CoGAPS parameters

 #' @return list with all CoGAPS output

new file mode 100644

@@ -0,0 +1,152 @@

+#' wrapper around cat

+#' @keywords

+#'

+#' @description cleans up message printing

+#' @param allParams all cogaps parameters

+#' @param ... arguments forwarded to cat

+#' @return displays text

+gapsCat <- function(allParams, ...)

+{

+    if (allParams$messages)

+        cat(...)

+}

+

+#' checks if file is supported

+#' @keywords internal

+#'

+#' @param file path to file

+#' @return TRUE if file is supported, FALSE if not

+#' @importFrom tools file_ext

+supported <- function(file)

+{

+    if (!is(file, "character"))

+        return(FALSE)

+    return(tools::file_ext(file) %in% c("tsv", "csv", "mtx"))

+}

+

+#' get number of rows from supported file name or matrix

+#' @keywords internal

+#'

+#' @param data either a file name or a matrix

+#' @return number of rows

+#' @importFrom data.table fread

+#' @importFrom tools file_ext

+nrowHelper <- function(data)

+{

+    if (is(data, "character"))

+    {

+        return(switch(tools::file_ext(data),

+            "csv" = nrow(data.table::fread(data, select=1)),

+            "tsv" = nrow(data.table::fread(data, select=1)),

+            "mtx" = as.numeric(data.table::fread(data, nrows=1, fill=TRUE)[1,1])

+        ))

+    }

+    return(nrow(data))

+}

+

+#' get number of columns from supported file name or matrix

+#' @keywords internal

+#'

+#' @param data either a file name or a matrix

+#' @return number of columns

+#' @importFrom data.table fread

+#' @importFrom tools file_ext

+ncolHelper <- function(data)

+{

+    if (is(data, "character"))

+    {

+        return(switch(tools::file_ext(data),

+            "csv" = ncol(data.table::fread(data, nrows=1)) - 1,

+            "tsv" = ncol(data.table::fread(data, nrows=1)) - 1,

+            "mtx" = as.numeric(data.table::fread(data, nrows=1, fill=TRUE)[1,2])

+        ))

+    }

+    return(ncol(data))

+}

+

+#' extract gene names from data

+#' @keywords internal

+getGeneNames <- function(data, transpose)

+{

+    nGenes <- ifelse(transpose, ncolHelper(data), nrowHelper(data))

+    return(paste("Gene", 1:nGenes, sep="_"))

+}

+

+#' extract sample names from data

+#' @keywords internal

+getSampleNames <- function(data, transpose)

+{

+    nSamples <- ifelse(transpose, nrowHelper(data), ncolHelper(data))

+    return(paste("Sample", 1:nSamples, sep="_"))

+}

+

+#' write start up message

+#' @keywords internal

+#'

+#' @param data data set

+#' @param allParams list of all parameters

+#' @return message displayed to screen

+#' @importFrom methods show

+startupMessage <- function(data, allParams)

+{

+    nGenes <- ifelse(allParams$transposeData, ncolHelper(data), nrowHelper(data))

+    nSamples <- ifelse(allParams$transposeData, nrowHelper(data), ncolHelper(data))

+

+    dist_message <- "Standard"

+    if (!is.null(allParams$gaps@distributed))

+        dist_message <- allParams$gaps@distributed

+

+    cat("Running", dist_message, "CoGAPS on", nGenes, "genes and",

+        nSamples, "samples")

+

+    if (allParams$messages)

+    {

+        cat(" with parameters:\n\n")

+        methods::show(allParams$gaps)

+    }

+    cat("\n")

+}

+

+#' parse parameters passed through the ... variable

+#' @keywords internal

+#'

+#' @param allParams list of all parameters

+#' @param extraParams list of parameters in ...

+#' @return allParams with any valid parameters in extraParams added

+#' @note will halt with an error if any parameters in extraParams are invalid

+#' @importFrom methods slotNames

+parseExtraParams <- function(allParams, extraParams)

+{

+    # parse direct params

+    for (s in slotNames(allParams$gaps))

+    {

+        if (!is.null(extraParams[[s]]))

+        {

+            allParams$gaps <- setParam(allParams$gaps, s, extraParams[[s]])

+            extraParams[[s]] <- NULL

+        }

+    }

+

+    # check for unrecognized options

+    if (length(extraParams) > 0)

+        stop(paste("unrecognized argument:", names(extraParams)[1]))

+

+    return(allParams)

+}

+

+#' check that provided data is valid

+#' @keywords internal

+#'

+#' @param data data matrix

+#' @param uncertainty uncertainty matrix, can be null

+#' @param params CogapsParams object

+#' @return throws an error if data has problems

+checkDataMatrix <- function(data, uncertainty, params)

+{

+    if (sum(data < 0) > 0 | sum(uncertainty < 0) > 0)

+        stop("negative values in data and/or uncertainty matrix")

+    if (nrow(data) <= params@nPatterns | ncol(data) <= params@nPatterns)

+        stop("nPatterns must be less than dimensions of data")

+    if (sum(uncertainty < 1e-5) > 0)

+        warning("small values in uncertainty matrix detected")

+}

\ No newline at end of file

new file mode 100644

@@ -0,0 +1,95 @@

+#' use user provided subsets

+#' @keywords internal

+#'

+#' @param allParams list of all CoGAPS parameters

+#' @param total total number of rows (cols) that are being paritioned

+#' @return list of subsets

+sampleWithExplictSets <- function(allParams, total)

+{

+    if (all(sapply(allParams$gaps@explicitSets, function(s) is(s, "numeric"))))

+    {

+        gapsCat(allParams, "using provided indexed subsets\n")

+        return(allParams$gaps@explicitSets)

+    }

+    else if (all(sapply(allParams$gaps@explicitSets, function(s) is(s, "character"))))

+    {

+        gapsCat(allParams, "using provided named subsets\n")

+        if (allParams$gaps@distributed == "genome-wide")

+            allNames <- allParams$geneNames

+        else

+            allNames <- allParams$sampleNames

+        return(lapply(allParams$gaps@explicitSets, function(set) which(allNames %in% set)))

+    }

+}

+

+#' subset rows (cols) proportional to the user provided weights

+#' @keywords internal

+#'

+#' @param allParams list of all CoGAPS parameters

+#' @param setSize the size of each subset of the total

+#' @return list of subsets

+sampleWithAnnotationWeights <- function(allParams, setSize)

+{

+    groups <- unique(allParams$gaps@samplingAnnotation)

+    return(lapply(1:allParams$gaps@nSets, function(i)

+    {

+        groupCount <- sample(groups, size=setSize, replace=TRUE,

+            prob=allParams$gaps@samplingWeight)

+        sort(unlist(sapply(groups, function(g)

+        {

+            groupNdx <- which(allParams$gaps@samplingAnnotation == g)

+            sample(groupNdx, size=sum(groupCount == g), replace=TRUE)

+        })))

+    }))

+}

+

+#' subset data by uniformly partioning rows (cols)

+#' @keywords internal

+#'

+#' @param allParams list of all CoGAPS parameters

+#' @param total total number of rows (cols) that are being paritioned

+#' @param setSize the size of each subset of the total

+#' @return list of subsets

+sampleUniformly <- function(allParams, total, setSize)

+{

+    sets <- list()

+    remaining <- 1:total

+    for (n in 1:(allParams$gaps@nSets - 1))

+    {

+        selected <- sample(remaining, setSize, replace=FALSE)

+        sets[[n]] <- sort(selected)

+        remaining <- setdiff(remaining, selected)

+    }

+    sets[[allParams$gaps@nSets]] <- sort(remaining)

+    return(sets)

+}

+

+#' partition genes/samples into subsets

+#' @keywords internal

+#'

+#' @description either genes or samples or partitioned depending on the type

+#' of distributed CoGAPS (i.e. genome-wide or single-cell)

+#' @param data either file name or matrix

+#' @param allParams list of all CoGAPS parameters

+#' @return list of sorted subsets of either genes or samples

+createSets <- function(data, allParams)

+{

+    subsetRows <- xor(allParams$transposeData, allParams$gaps@distributed == "genome-wide")

+    total <- ifelse(subsetRows, nrowHelper(data), ncolHelper(data))

+    setSize <- floor(total / allParams$gaps@nSets)

+

+    gapsCat(allParams, "Creating subsets...")

+

+    if (!is.null(allParams$gaps@explicitSets))

+    {

+        return(sampleWithExplictSets(allParams, total))

+    }

+

+    if (!is.null(allParams$gaps@samplingAnnotation))

+    {

+        gapsCat(allParams, "sampling with annotation weights\n")

+        return(sampleWithAnnotationWeights(allParams, setSize))

+    }

+    gapsCat(allParams, "\n")

+    return(sampleUniformly(allParams, total, setSize))

+}