From f94a879fad391a01df1d09798cf7837f83b7b82e Mon Sep 17 00:00:00 2001 From: ClarissaFeuersteinAkgoz <48914703+ClarissaFeuersteinAkgoz@users.noreply.github.com> Date: Tue, 2 Jun 2020 20:17:20 +0200 Subject: [PATCH 1/4] trial :( --- .gitignore | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.gitignore b/.gitignore index 46b6fe1..f0929d0 100644 --- a/.gitignore +++ b/.gitignore @@ -1,5 +1,5 @@ Meta -.DS_Store +.DS_Store .Rproj.user .Rproj .Rhistory From 903deba06cb3c60f35b3589260a4f737d3ac6117 Mon Sep 17 00:00:00 2001 From: ClarissaFeuersteinAkgoz <48914703+ClarissaFeuersteinAkgoz@users.noreply.github.com> Date: Thu, 3 Jun 2021 22:07:15 +0200 Subject: [PATCH 2/4] Update methrix_operations.R Function to calculate hydroxymethylation proportions from BS and oxBS data returns a methrix object with proportions for mC (~oxBS) hmC C mC+hmC (~BS) dependencies to add?: MLML2R, rlist, manuall is still missing --- R/methrix_operations.R | 80 ++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 80 insertions(+) diff --git a/R/methrix_operations.R b/R/methrix_operations.R index ea84464..0e4e4b6 100644 --- a/R/methrix_operations.R +++ b/R/methrix_operations.R @@ -871,3 +871,83 @@ convert_methrix <- function(m = NULL) { desc = metadata(m)$descriptive_stats) return(m) } + +#-------------------------------------------------------------------------------------------------------------------------- +#' Calculates hydroxymethylation from BS and oxBS data +#' @details Takes a \code{\link{methrix}} object and returns with the same object with delayed array assay slots +#' @param m An object of class \code{\link{methrix}} +#' @return An object of class \code{\link{methrix}} +calculate_proportions<- function (m){ + + # check if input data set is ok + # for each "Sample_ID" there must be a oxBS and BS data set, defined in "Sequencing_type" + anno <- data.frame(colData(m)) + anno <- as.data.frame(lapply(anno, factor)) + for(i in levels(anno$Sample_ID)){ + anno_sub <- anno[anno$Sample_ID==i,] + print(i) + if (nrow(anno_sub)!=2 | sum(is.element(c("oxBS", "BS"), anno_sub$Sequencing_type))!=2 ) + { stop("oxBS or BS sample is missing or in excess for a sample type.", + call. = T) + } + } + + # extract data for MLML2R + methylation <- get_matrix(m, type = 'M') + covarage <- get_matrix(m, type = 'C') + + + id_BS <- which(colData(m)$Sequencing_type=="BS") + id_oxBS <- which(colData(m)$Sequencing_type=="oxBS") + + #match samples of BS and oxBS to same order + id_oxBS <- id_oxBS[match(colData(m)$Sample_ID[id_BS], colData(m)$Sample_ID[id_oxBS])] + + colData(m)$Sample_ID[id_BS]==colData(m)$Sample_ID[id_oxBS] + + + #total count + TotalBS <- covarage[,id_BS] + TotalOxBS <- covarage[,id_oxBS] + + #methylation count + C_BS <- data.frame(methylation[,id_BS])*data.frame(covarage[,id_BS]) + C_OxBS<- data.frame(methylation[,id_oxBS])*data.frame(covarage[,id_oxBS]) + + #unmethylation count + T_BS<- data.frame(TotalBS) - data.frame(C_BS) + T_OxBS<- data.frame(TotalOxBS) - data.frame(C_OxBS) + + Tm = as.matrix(C_BS) + Um = as.matrix(T_BS) + Lm = as.matrix(T_OxBS) + Mm = as.matrix(C_OxBS) + + results_exact <- MLML2R::MLML(T.matrix = Tm, + U.matrix = Um, + L.matrix = Lm, + M.matrix = Mm) + + + #create new methrix object + + # get the lowest covarage among both sequencing types + + min_cov <- rlist::list.cbind(lapply(1:length(id_BS), function(x) rowMins(covarage[,c(id_BS[x],id_oxBS[x])]))) + + #prepare methrix dataset for mC, hmC and C data + m_MLML2R <- m[,c(id_BS,id_BS,id_BS,id_BS)] + m_MLML2R@colData$Sequencing_type <- NULL + m_MLML2R@colData$Mark <-c(rep("mC",length(id_BS)),rep("hmC",length(id_BS)),rep("C",length(id_BS)),rep("mChmC",length(id_BS))) + mChmC<- as.matrix(results_exact$mC)+as.matrix(results_exact$hmC) + m_MLML2R@ assays@ data@ listData$ beta <- cbind(results_exact$mC,results_exact$hmC,results_exact$C,mChmC) + m_MLML2R@ assays@ data@ listData$ cov <- cbind(min_cov,min_cov,min_cov,min_cov) + + #change names in methrix object + names <- gsub("bs", "",m_MLML2R@colData@ listData$ Sample_Name) + m_MLML2R <-replace_names(m_MLML2R,paste0(names, m_MLML2R@colData$Mark)) + + return(m_MLML2R) + +} + From 2342b52efa62f133f2746154d8ccc0d55f847fac Mon Sep 17 00:00:00 2001 From: ClarissaFeuersteinAkgoz <48914703+ClarissaFeuersteinAkgoz@users.noreply.github.com> Date: Thu, 3 Jun 2021 23:35:01 +0200 Subject: [PATCH 3/4] Update accessory_funcs.R function needed when calculating hmC from oxBS and BS --- R/accessory_funcs.R | 11 +++++++++++ 1 file changed, 11 insertions(+) diff --git a/R/accessory_funcs.R b/R/accessory_funcs.R index e7eced1..a93d579 100644 --- a/R/accessory_funcs.R +++ b/R/accessory_funcs.R @@ -612,3 +612,14 @@ get_y_lims <- function(vec) { list(y_lims = y_lims, y_at = y_at) } + + +#-------------------------------------------------------------------------------------------------------------------------- +# function to replace names when working with oxBS and BS data +replace_names = function(m, value) { + colData(m)$rownames = value + colnames(assays(m)$beta) = value + colnames(assays(m)$cov) = value + colData(m)$Sample_Name = value + return(m) +} From 0db7e75f692e04c61bcd7a660a7687cdcc1443ef Mon Sep 17 00:00:00 2001 From: ClarissaFeuersteinAkgoz <48914703+ClarissaFeuersteinAkgoz@users.noreply.github.com> Date: Thu, 3 Jun 2021 23:40:52 +0200 Subject: [PATCH 4/4] Update methrix_operations.R use of accessor functions --- R/methrix_operations.R | 16 +++++++++------- 1 file changed, 9 insertions(+), 7 deletions(-) diff --git a/R/methrix_operations.R b/R/methrix_operations.R index 0e4e4b6..8954e68 100644 --- a/R/methrix_operations.R +++ b/R/methrix_operations.R @@ -937,15 +937,17 @@ calculate_proportions<- function (m){ #prepare methrix dataset for mC, hmC and C data m_MLML2R <- m[,c(id_BS,id_BS,id_BS,id_BS)] - m_MLML2R@colData$Sequencing_type <- NULL - m_MLML2R@colData$Mark <-c(rep("mC",length(id_BS)),rep("hmC",length(id_BS)),rep("C",length(id_BS)),rep("mChmC",length(id_BS))) + colData(m_MLML2R)$Sequencing_type <- NULL + colData(m_MLML2R)$Mark <-c(rep("mC",length(id_BS)),rep("hmC",length(id_BS)),rep("C",length(id_BS)),rep("mChmC",length(id_BS))) mChmC<- as.matrix(results_exact$mC)+as.matrix(results_exact$hmC) - m_MLML2R@ assays@ data@ listData$ beta <- cbind(results_exact$mC,results_exact$hmC,results_exact$C,mChmC) - m_MLML2R@ assays@ data@ listData$ cov <- cbind(min_cov,min_cov,min_cov,min_cov) - + assays(m_MLML2R)$beta <- cbind(results_exact$mC,results_exact$hmC,results_exact$C,mChmC) + min_cov_mat <- cbind(min_cov,min_cov,min_cov,min_cov) + colnames(min_cov_mat) <- colnames(assays(m_MLML2R)$cov) + assays(m_MLML2R)$cov <- min_cov_mat + #change names in methrix object - names <- gsub("bs", "",m_MLML2R@colData@ listData$ Sample_Name) - m_MLML2R <-replace_names(m_MLML2R,paste0(names, m_MLML2R@colData$Mark)) + names <- gsub("bs", "", colData(m_MLML2R)$Sample_Name) + m_MLML2R <-replace_names(m_MLML2R,paste0(names, colData(m_MLML2R)$Mark)) return(m_MLML2R)