批量运行弦图,封装函数

一共30个样本,要出弦图,先在jupyter里批量运行了PlotFancyVJUsage:

def circos():
    import os
    for i in (1,7,35,189,217):
        for j in (1003,1501,1503,2001,2501,2503):
            x = "day" + str(i) + "-" + str(j)
            print(x)
            cmd_string = "java -jar ./vdjtools-1.2.1.jar PlotFancyVJUsage "+x+".txt "+ x
            print('x:{}'.format(cmd_string))
            print(os.popen(cmd_string).read())
                
circos()

之前是在产生的错误文件(vj_pairing_plot.r)之里修改内容后一个个样本文件进行运行,想到这样太麻烦,就想写成一个函数,变化的是输入文件的名字
输入的文件名称的规律:

代码里输入文件的地方在:

args <- c("day1-2501.fancyvj.wt.txt", "day1-2501.fancyvj.wt.pdf")
temp <- read.table("day1-2501.fancyvj.wt.txt", sep="\t", comment="")

本来这个循环我是写成这样的:

但是要么报这个错

要么报这个错

因为手动输入文件名是可以运行的,所以觉得200G这个不存在……一直在找文件名哪里输的不对
然后改成只运行一个循环,就成功了
可能真的是需要那么多运行内存吧...

circle_plot <- function(){
  for(i in c(1,7,35,189,217)){
      x = paste('day',i,'-2503',sep = "")
      y = '.fancyvj.wt.txt'
      z = '.fancyvj.wt.pdf'
      a = paste(x,y,sep = "")
      b = paste(x,z,sep = "")
      args <- c(a,b)
      #args<-commandArgs(TRUE)
    
      file_in  <- args[1]
      file_out <- args[2]
  
      require(circlize); require(RColorBrewer)

      # load data and preproc to fit formats
    
      temp <- read.table(a, sep="\t", comment="")
      n <- nrow(temp)
      m <- ncol(temp)
      rn = as.character(temp[2:n,1])
      cn = apply(temp[1,2:m], 2 , as.character)
      mat <- matrix(apply(temp[2:n, 2:m], 1:2, as.numeric), n - 1, m-1) * 100
      
      n <- nrow(temp)
      m <- ncol(temp)
      
      # Here columns and rows correspond to V and J segments respectively
      # Also replace possible duplicates (undef, '.', ...)
      
      duplicates <- intersect(rn, cn)
      
      rownames(mat) <- replace(rn, rn==duplicates, paste("V", duplicates, sep=""))
      colnames(mat) <- replace(cn, cn==duplicates, paste("J", duplicates, sep=""))
      
      # sort
      
      col_sum = apply(mat, 2, sum)
      row_sum = apply(mat, 1, sum)
      
      mat <- mat[order(row_sum), order(col_sum)]
      
      # equal number of characters for visualizaiton
      
      rn <- rownames(mat)
      cn <- colnames(mat)
      
      maxrn <- max(nchar(rn))
      maxcn <- max(nchar(cn))
      
      for(i in seq_len(length(rn))) {
            rn[i] <- paste(rn[i], paste(rep(" ", maxrn - nchar(rn[i])), collapse = ''))
      }
      
      for(i in seq_len(length(cn))) {
            cn[i] <- paste(cn[i], paste(rep(" ", maxcn - nchar(cn[i])), collapse = ''))
      }
      
      rownames(mat) <- rn
      colnames(mat) <- cn
      
      # viz using circlize
      
      if (grepl("\\.pdf$",file_out)){
         pdf(file_out)
      } else if (grepl("\\.png$",file_out)) {
         png(file_out, width     = 3.25,
                       height    = 3.25,
                       units     = "in",
                       res       = 1200,
                       pointsize = 4)
      } else {
         stop('Unknown plotting format')
      }
      
      circos.par(gap.degree = c(rep(1, nrow(mat)-1), 10, rep(1, ncol(mat)-1), 15), start.degree = 5)
      
      rcols <- rep(brewer.pal(12, "Paired"), nrow(mat)/12 + 1)[1:nrow(mat)]
      ccols <- rep(brewer.pal(12, "Paired"), ncol(mat)/12 + 1)[1:ncol(mat)]
      
      names(rcols) <- sort(rownames(mat))
      names(ccols) <- sort(colnames(mat))
      
      chordDiagram(mat, annotationTrack = "grid",
                   reduce = 0,
                   grid.col = c(rcols, ccols),
                   preAllocateTracks = list(track.height = 0.2), transparency = 0.5)
      
      circos.trackPlotRegion(track.index = 1, bg.border = NA,
             panel.fun = function(x, y) {
                         sector.name = get.cell.meta.data("sector.index")
                         xlim = get.cell.meta.data("xlim")
                         ylim = get.cell.meta.data("ylim")
                         circos.text(mean(xlim), ylim[1], cex = 0.5, sector.name, facing = "clockwise", adj = c(0, 0.5))
                         }
             )
      
      circos.clear()
      
      dev.off()
    }
}
  
circle_plot()

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 203,937评论 6 478
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 85,503评论 2 381
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 150,712评论 0 337
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 54,668评论 1 276
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 63,677评论 5 366
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,601评论 1 281
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 37,975评论 3 396
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,637评论 0 258
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 40,881评论 1 298
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,621评论 2 321
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,710评论 1 329
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,387评论 4 319
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 38,971评论 3 307
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,947评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,189评论 1 260
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 44,805评论 2 349
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,449评论 2 342

推荐阅读更多精彩内容