From 8f32f0a230fefef658dd21da9dd925600fcdbdd3 Mon Sep 17 00:00:00 2001 From: marcelcosta Date: Wed, 20 Jul 2022 16:19:33 +0200 Subject: [PATCH] =?UTF-8?q?Canviar=20threshold=20total=20per=20filtre=20lo?= =?UTF-8?q?cal=20i=20posar=20un=20m=C3=A0xim=20a=20la=202a=20agrupaci?= =?UTF-8?q?=C3=B3=20de=20pics.?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- chromatoR/app.R | 41 +++++++++++++++++++++++++++++++++-------- 1 file changed, 33 insertions(+), 8 deletions(-) diff --git a/chromatoR/app.R b/chromatoR/app.R index 29cb029..c06d3de 100644 --- a/chromatoR/app.R +++ b/chromatoR/app.R @@ -119,7 +119,7 @@ server <- function(input, output) { else return(c(max(region[,2], na.rm=TRUE), region[which.max(region[,2]),1])) } - redDiffs<-function(matrix, diff=5){ + redDiffs<-function(matrix, diff=5, max=0){ diffs<-diff(matrix[,1]) pos<-matrix(nrow=0, ncol=2) i_temp<-matrix(nrow = 0, ncol = 2) @@ -127,7 +127,7 @@ server <- function(input, output) { for(i in 1:nrow(matrix)){ i_temp<-rbind(i_temp, matrix[i,]) - if(diffs[i] >= diff | i == nrow(matrix)){ + if(diffs[i] >= diff | i == nrow(matrix) | nrow(i_temp) == max){ pos<-rbind(pos, i_temp[which(i_temp[,2] == max(i_temp[,2]))[1],]) i_temp<-matrix(nrow = 0, ncol = 2) } @@ -149,11 +149,30 @@ server <- function(input, output) { # Peaks joining and filtering by defined threshold peakthr<-input$thr - peakCusMatrix<-rbind(Gpeaks[Gpeaks[,2] > peakthr,], - Apeaks[Apeaks[,2] > peakthr,], - Tpeaks[Tpeaks[,2] > peakthr,], - Cpeaks[Cpeaks[,2] > peakthr,]) + filt_peaks<-function(Mpeaks, prob=0.75){ + Mpeaks_filt<-matrix(nrow = 0, ncol=2) + means<-matrix(nrow=0, ncol=2) + for (i in 1:nrow(Mpeaks)){ + pack<-Mpeaks[Mpeaks[,1] >= Mpeaks[i,1]-100 & Mpeaks[,1] <= Mpeaks[i,1]+100,] + if (Mpeaks[i,2] > quantile(pack[,2], na.rm=T, probs=prob)[[1]]){ + Mpeaks_filt<-rbind(Mpeaks_filt, Mpeaks[i,]) + } + means<-rbind(means, c(Mpeaks[i,1], quantile(pack[,2], na.rm=T, probs=prob)[[1]])) + } + return(list(Mpeaks_filt, means)) + } + + + # peakCusMatrix<-rbind(Gpeaks[Gpeaks[,2] > peakthr,], + # Apeaks[Apeaks[,2] > peakthr,], + # Tpeaks[Tpeaks[,2] > peakthr,], + # Cpeaks[Cpeaks[,2] > peakthr,]) + + peakCusMatrix<-filt_peaks(rbind(Gpeaks,Apeaks,Tpeaks,Cpeaks),0.75) + peakmeans<<-peakCusMatrix[[2]] + peakCusMatrix<-peakCusMatrix[[1]] + # In case there is more than one peak by position, we keep the higher peakCusMatrix<-peakCusMatrix %>% as.data.frame() %>% dplyr::rename(Int=V2) %>% group_by(indexes) %>% summarise(Int=max(Int)[1]) @@ -161,7 +180,7 @@ server <- function(input, output) { # Using the custom function, we aggregate peaks defining two distance thresholds pos1<-redDiffs(peakCusMatrix, input$dist1) - pos<-redDiffs(pos1, input$dist2)[,1] + pos<-redDiffs(pos1, input$dist2, max=2)[,1] progress$set(message = "Peaks Detected", value = 2) print("Peaks Detected") @@ -341,6 +360,9 @@ server <- function(input, output) { ranPeaks<-range(picks$rows) + # pmeans<-peakmeans %>% as.data.frame() %>% rename(index=1, value=2) + # pmeans<-pmeans[] + obj$plot<-obj_ab@traceMatrix %>% as.data.frame() %>% dplyr::rename(A=V1,C=V2,G=V3,T=V4) %>% add_column(rows=1:nrow(.)) %>% filter(rows >= ranPeaks[1] & rows <= ranPeaks[2]) %>% gather(Base, Value, -rows) %>% @@ -351,7 +373,10 @@ server <- function(input, output) { geom_text(data=picks, aes(label=Base, y=3800, color=factor(Base,levels=c("G","A","T","C"))))+ geom_text(data=picks, aes(label=num, y=4250), size=3, angle=90, hjust=0.5, vjust=0.5)+ scale_color_manual(values=c("G"="black", "A"="#66CC00", "T"="red","C"="blue"))+ - geom_hline(yintercept = input$thr)+ + # geom_hline(yintercept = input$thr)+ + geom_point(data=peakmeans %>% as.data.frame() %>% rename(index=1, value=2) %>% + filter(index >= ranPeaks[1] & index <= ranPeaks[2]), + aes(index, value))+ theme_classic() print("Final Imagen") }