ABSOLUTE肿瘤纯度分析测试

library(data.table)
library(foreach)
library(doMC)
library(ABSOLUTE)

setwd("/home/train/data/ABSOLUTE_test/Tumor_purity")

DoAbsolute <- function(scan,input) {
#  registerDoSEQ()  #多核心计算时候，需要注释掉
  library(ABSOLUTE)
  plate.name <- "test" #可以根据需要修改
  genome <- "hg18"
  platform <- c("SNP_6.0", "Illumina_WES", "SNP_250K_STY")
  sample.name <- scan
  primary.disease = "NA"
  sigma.p <- 0
  max.sigma.h <- 0.02
  min.ploidy <- 0.95
  max.ploidy <- 10
  max.as.seg.count <- 1500
  max.non.clonal <- 0
  max.neg.genome <- 0
  copy_num_type <- "total"
  seg.dat.fn <- input
  results.dir <- file.path(".", "output", scan, "absolute")
  print(paste("Starting scan", scan, "at", results.dir))
  log.dir <- file.path(".", "output", "abs_logs")
  if (!file.exists(log.dir)) {
     dir.create(log.dir, recursive=TRUE)
  }
  dz <- lapply(results.dir,function(x){
    if(!file.exists(x)){
      dir.create(x, recursive=TRUE)
    }
  })
  rm(dz)
  sink(file=file.path(log.dir, paste(scan, ".abs.out.txt", sep="")))
  RunAbsolute(seg.dat.fn, sigma.p, max.sigma.h, min.ploidy, max.ploidy, primary.disease, 
              platform, sample.name, results.dir, max.as.seg.count, max.non.clonal, 
              max.neg.genome, copy_num_type, verbose=TRUE)
  sink()
}


registerDoMC(20)    #对每组数据调用线程池

arrays.txt <- "SNP6_solid_tumor.seg.txt"

input_data <- "test.txt"
scans <- unique(as.data.frame(fread(arrays.txt))$Sample)
foreach (scan=scans, .combine=c) %dopar% {
  DoAbsolute(scan,input=input_data)
}

obj.name <- "test_summary"
results.dir <- file.path(".", "output", "abs_summary")
absolute.files <- file.path(".", "output",
                            scans, "absolute",
                            paste(scans, ".ABSOLUTE.RData", sep=""))


CreateReviewObject(obj.name, absolute.files, results.dir, "total", verbose=TRUE)

calls.path = file.path("output", "abs_summary", "test_summary.PP-calls_tab.txt")
modes.path = file.path("output", "abs_summary", "test_summary.PP-modes.data.RData")
output.path = file.path("output", "abs_extract")
ExtractReviewedResults(calls.path, "test", modes.path, output.path, "absolute", "total")


# doabsolute包更好使用，需要根据输入数据调整格式暂时分步处理，今后有时间写成pipline可以无缝对接doabsolute进行批量处理

mkdir T1_5
mkdir total
mv T1_5.cr.igv.seg T1_5
mv T1_5_absolute.maf T1_5
cd T1_5/
mv T1_5.cr.igv.seg tumor.seg
mv ../T1_3/step1_Xto23.sh .
bash step1_Xto23.sh 
cp ../T1_3/step2_MafFilebalance.pl .
perl step2_MafFilebalance.pl T1_5_absolute.maf
cp ../T1_3/step3_SegFilebalance.pl .
perl step3_SegFilebalance.pl new_tumor.seg
cp ../T1_3/step4_rebalance_seg_and_maf.R .
Rscript step4_rebalance_seg_and_maf.R
cp ../T1_3/step5_get_seperate_file.R .
Rscript step5_get_seperate_file.R
cp ../T1_3/step6_prepare_doabsolute.pl .
perl step6_prepare_doabsolute.pl ./seperate_file/
cp ../T1_3/step6.5_prepare_doabsolute.sh .
#修改step6.5中的文件名为T1_5
bash step6.5_prepare_doabsolute.sh

cd ..
mkdir total
cd total/
mkdir seperate_file
cp ../T1_3/seperate_file/*.maf seperate_file/
cp ../T1_3/seperate_file/new_T1_3.seg seperate_file/
cp ../T1_4/seperate_file/*.maf seperate_file/
cp ../T1_4/seperate_file/new_T1_4.seg seperate_file/
cp ../T1_5/seperate_file/*.maf seperate_file/
cp ../T1_5/seperate_file/new_T1_5.seg seperate_file/
cp ../T1_3/new_final.seg .
cat ../T1_4/new_final.seg | awk 'NR>1' >> new_final.seg
cat ../T1_5/new_final.seg | awk 'NR>1' >> new_final.seg
cp ../T1_3/step7_doabsolute.R .
#注意修改路径
Rscript step7_doabsolute.R