R/mts_targetCoverageFromMutations.R

In MultiSEp: Predict Synthetic Lethality and Other Gene Dependency Relationships from Multiomics Data

mts_targetCoverageFromMutations <- function(mutData, edgeData, rowNum=nrow(edgeData), cores)
{
  
  if (missing(mutData)) stop("Please provide a mutation score matrix")
  if (missing(edgeData)) stop("Please provide edgeData")
  
  if (missing(cores)) {
    cores <- 1
    message("1 core selected")
  } else if (!is.numeric(cores) || cores < 1) {
    stop("cores should be >= 1")
  } else {
    message(cores, " cores selected")
  }
  
  # Parse mutation data
  mDet <- melt(as.matrix(mutData))
  colnames(mDet) <- c("gene_symbol", "sample_name", "mutation_status")
  
  # Parse Edge List
  if (is.data.frame(edgeData) && all(c("mRNA_gene","mutation_gene","tissue","chiSqPvalue","mutation_per_mode") %in% colnames(edgeData))) {
    mRes2 <- edgeData[,c(1:2,4)]
  } else {
    # If Edge List is parsed 
    mRes2 <- edgeData
  }
  
  # Subfunction Define
  genes <- function(edgeData){
    as_ids(V(graph_from_data_frame(edgeData)))
  } # gene IDs
  
  
  findNeighbour <- function(gene, edgeData, mutData){
    
    # Prepare empty vectors
    samples <- c()
    neighbor_gene <- c()
    graph  <- c()
    indexes <- c()
    num_mutation <- c()
    mut_frequency <- c()
    num_neighbor <- c()
    neighbor_list <- c()
    mut_samples <- c()
    
    #Create igraph object to quickly find the neighbors of the gene.
    graph <- graph_from_data_frame(edgeData, directed = F)
    neighbor_gene <- as_ids(neighbors(graph, gene))
    
    # Remove gene from the list of neighbors on the occasion it interacts with itself.
    neighbor_gene <- neighbor_gene[!neighbor_gene %in% gene]
    
    # Find rows in the sample data in which there is mutation in one or more of the neighbor genes.
    mutDataM <- mutData[which(mutData$mutation_status != "WT"),]
    
    indexes <- which(mutDataM$gene_symbol %in% neighbor_gene)
    
    # Find sample IDs of samples that have a mutation in one or more of the neighbor genes.
    if (length(indexes)) {   # BUT don't populate the samples object if there are no mutated samples
      for (i in 1:length(indexes)){
        samples <- append(samples,mutDataM[indexes[i],"sample_name"])
      }
    }
    # Number of samples with a mutation in one or more of the neighbor genes
    num_mutation <- length(unique(samples))
    
    # Frequency of samples with a mutation in one or more of the neighbor genes
    mut_frequency <- num_mutation / length(unique(mutData$sample_name)) # out of total samples
    
    # The mutated samples
    if (length(unique(samples)) > 1) {
      mut_samples <- paste(unique(samples), collapse=",")
      
    } else if (length(unique(samples)) == 1) {
      mut_samples <- unique(samples)
    } else {
      mut_samples <- "NA"
    }
    
    # Ensure mut_samples is treated as a character
    mut_samples <- as.character(mut_samples)
    
    # Count neighbors of the gene.
    num_neighbor <- length(neighbor_gene)
    
    if(length(neighbor_gene)) {
      # Name the neighbors of the gene.
      neighbor_list <- paste(neighbor_gene, collapse=",")
      
    } else {
      neighbor_gene = "NA"
      neighbor_list = "NA"
    }
    
    return(list(gene, num_neighbor, num_mutation, mut_frequency, neighbor_list, mut_samples))
  }
  
  # Resolve results
  makeNeighbourTable <- function(data){
    table <- as.data.frame(matrix(unlist(data, recursive = FALSE), ncol=6, byrow=T), stringsAsFactors = FALSE)
    colnames(table) <- c("Gene", "NeighbourCount", "NumSamplesWithMutation", "NeighbourMutationFrequency", "NeighbourGenes", "MutatedSamples")
    table[,2] <- as.numeric(table[,2]) 
    table[,3] <- as.numeric(table[,3])
    table[,4] <- as.numeric(table[,4])
    table[,1] <- as.character(table[,1]) 
    table[,5] <- as.character(table[,5]) 
    table[,6] <- as.character(table[,6])
    table <- table[order(-table$NeighbourMutationFrequency),]
    return(table)
  }
  # Implement functions
  neighbourAnalysis <- pbmclapply(genes(mRes2[1:rowNum,]), findNeighbour, edgeData = mRes2[1:rowNum,], mutData=mDet, mc.cores=cores)
  neighbourhoodTable <- makeNeighbourTable(neighbourAnalysis)
  
  return(neighbourhoodTable)
}