#' Determine relationship between pairwise sample
#'
#' This function identifies relationship of pairwise sample
#' [00: same patient | 10: X contaminate Y | 01 Y contaminate X | 11: both way contamination]
#'
#' @param sample_pairs data frame, pairwise samples information
#' @param center_cutoff numeric, cutoff to be classified as same subject samples (default: 0.5)
#' @param source_cutoff numeric, cutoff to be classififed as source sample (default: 0.4)
#' @param target_cutoff numeric, cutoff to be classified as target sample (default: 0.1)
#' @param localPcomm_cutoff numeric, cutoff to manage high or low SNPs sharing between a pair of samples
#' @param region_cutoff numeric, cutoff for case with low SNP sharing (default: 0.75)
#' @param num_round_digit integer, rounding numeric up to this decimal point (default: 3)
#' @param output_path character, output directory
#'
#' @usage pairRelation(sample_pairs,center_cutoff,source_cutoff,target_cutoff,
#' localPcomm_cutoff,region_cutoff,num_round_digit,output_path)
#'
#' @return matrix, 2 columns containing pcommon and relation
#' @references
#' {TBA}
#'
#' @export
#'
pairRelation <- function(sample_pairs,center_cutoff=0.5,source_cutoff=0.4,target_cutoff=0.1,localPcomm_cutoff=0.1,region_cutoff=0.75,num_round_digit=3,output_path){
##
logfile <- paste0(output_path,"/tmp/log.txt")
out <- t(c('pid1','pid2',"ncommon","pcommon","pn1","pn2","pn3","qn1","qn3","rn1","rn3"))
write.table(out,quote=F,sep='\t',logfile,append=F,row.names = F,col.names = F)
relation <- apply(sample_pairs, 1, function(i) {
pn1=pn2=pn3=qn1=qn3=rn1=rn3=0
rel <-"-1"
ncommon <- as.numeric(i["x2y"])+as.numeric(i["y2x"])+as.numeric(i["center"])
pcommon <- ncommon /(ncommon + as.numeric(i["xonly"])+as.numeric(i["yonly"]) )
pcommon <- round(pcommon, digits = num_round_digit)
pn1 <- as.numeric(i["x2y"])/ncommon
pn1 <- round(pn1,digits =num_round_digit)
pn2 <- as.numeric(i["center"])/ncommon
pn2 <- round(pn2,digits =num_round_digit)
pn3 <- as.numeric(i["y2x"])/ncommon
pn3 <- round(pn3,digits =num_round_digit)
qn1 <- as.numeric(i["x2y"])/(as.numeric(i["x2y"])+as.numeric(i["xonly"]))
qn1 <- round(qn1,digits =num_round_digit)
qn3 <- as.numeric(i["y2x"])/(as.numeric(i["y2x"])+as.numeric(i["yonly"]))
qn3 <- round(qn3,digits =num_round_digit)
rn1 <- as.numeric(i["x2y"])/(as.numeric(i["x2y"])+as.numeric(i["center"]))
rn1 <- round(rn1,digits =num_round_digit)
rn3 <- as.numeric(i["y2x"])/(as.numeric(i["y2x"])+as.numeric(i["center"]))
rn3 <- round(rn3,digits =num_round_digit)
if(is.na(qn1)){
qn1=0
}
if(is.na(qn3)){
qn3=0
}
if(is.na(rn1)){
rn1=0
}
if(is.na(rn3)){
rn3=0
}
if(is.na(pn1)){
pn1=0
}
if(is.na(pn2)){
pn2=0
}
if(is.na(pn3)){
pn3=0
}
if(pcommon > localPcomm_cutoff){
if(pn2>center_cutoff){
rel <-"00"
}else if(pn1>source_cutoff & pn3<target_cutoff){ #x2y
rel<-"10"
}else if(pn3>source_cutoff & pn1<target_cutoff){ #y2x
rel<-"01"
}else{ #bothway
rel<-"11"
}
}else{## Mar 24 2017
cutoff <- 20
if((qn1>target_cutoff&rn1>region_cutoff) & (qn3<=target_cutoff&rn3<=region_cutoff)){
rel <- "10" #x2y
}else if((qn3>target_cutoff&rn3>region_cutoff) & (qn1<=target_cutoff&rn1<=region_cutoff)){
rel <- "01" #y2x
}else if((qn1>target_cutoff&rn1>region_cutoff) & (qn3>target_cutoff&rn3>region_cutoff)){
rel <- "11" #x<->y
}else if(rn1>region_cutoff|rn3>region_cutoff){ # little portion of center sharing
if(qn1>qn3&qn3>0){
ratio <- qn1/qn3
if(ratio>=cutoff){
rel <- "10"
}
}else if(qn3>qn1&qn1>0){
ratio <- qn3/qn1
if(ratio>=cutoff){
rel <- "01"
}
}
}else{ # center sharing larger portion
if(qn1>qn3&qn3>0){
ratio <- qn1/qn3
if(ratio>=cutoff){
rel <- "10"
}
}else if(qn3>qn1&qn1>0){
ratio <- qn3/qn1
if(ratio>=cutoff){
rel <- "01"
}
}
}
}
out <- t(c(i['pid1'],i['pid2'],ncommon,pcommon,pn1,pn2,pn3,qn1,qn3,rn1,rn3))
write.table(out,quote=F,sep='\t',logfile,append=T,row.names = F,col.names = F)
c(pcommon,rel)
})
relation = t(relation)
colnames(relation) <- c("pcommon","rel")
relation
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.