• circos pipeline


    # /usr/bin/env python
    # coding=utf-8
    ###################################
    # Author : yunkeli
    # Version : 1.0(2015/6/20)
    # E-mail : 1316014512@qq.com
    ###################################
    import os
    import argparse
    import re
    import random
    def vcf_SNPdensity(snpvcffile,pathway):
    print "this step is vcf to SNPdensity "
    cmdvcftorate = "/home/liyunke/vcftools_0.1.12b/bin/vcftools --vcf "+snpvcffile+" --out " + pathway+"/SNPdensity100K --SNPdensity 1000000"
    result_analysis=os.popen(cmdvcftorate)
    print result_analysis.read()
    def density(SNPdensityfile,pathway):
    print "##############################"
    print "this step is vcf to densitysplit cat "
    fileopen = open(SNPdensityfile).readlines()[1:]
    savename = pathway+"/"+"SNPdensity50K.snpden.new.txt"
    filesave = open(savename,"w+")
    for i in fileopen:
    listi = i.split()
    filesave.write(listi[0].replace("chr","hs")+" "+listi[1]+" "+str(int(listi[1])+999999)+" "+str(float(listi[3])/10)+" ")
    filesave.close()
    def densitysplit(SNPdensityfile,pathway):
    print "##############################"
    print "this step is vcf to densitysplit "
    fileopen = open(SNPdensityfile).readlines()[1:]
    namelist = []
    for i in fileopen:
    if i.split()[0] not in namelist:
    namelist.append(i.split()[0])
    for j in namelist:
    savename = pathway+"/"+j.replace("chr","hs")+".snp.txt"
    filesave = open(savename,"w+")
    for x in fileopen:
    listx = x.split()
    if listx[0] == j:
    filesave.write(listx[0].replace("chr","hs")+" "+listx[1]+" "+str(int(listx[1])+499999)+" "+listx[3]+" ")
    filesave.close()
    print "densitysplit ok"
    def sv_split(svdensityfile,pathway):
    print "##############################"
    print "this step is vcf to sv_file split "
    fileopen = open(svdensityfile).readlines()[1:]
    namelist = []
    for i in fileopen:
    if i.split()[0] not in namelist:
    namelist.append(i.split()[0])
    for j in namelist:
    listrandom = []
    savename = pathway+"/"+j.replace("chr","hs")+".sv.txt"
    filesave = open(savename,"w+")
    for x in fileopen:
    listx = x.split()
    if listx[0] == j:
    if listx[0] != listx[5]:
    listrandom.append(x)
    if len(listrandom) > 10:
    slicelist = random.sample(listrandom, 10)
    for links in slicelist:
    listlinks = links.split()
    filesave.write("segdups"+str(listrandom.index(links))+" "+" ".join(listlinks[0:3]).replace("chr","hs")+" ")
    filesave.write("segdups"+str(listrandom.index(links))+" "+" ".join(listlinks[5:8]).replace("chr","hs")+" ")
    filesave.close()
    else:
    for links in listrandom:
    listlinks = links.split()
    filesave.write("segdups"+str(listrandom.index(links))+" "+" ".join(listlinks[0:3])+" ")
    filesave.write("segdups"+str(listrandom.index(links))+" "+" ".join(listlinks[5:8])+" ")
    filesave.close()
    def circos_config(npath,prefix):
    oldconfig = "/home/liyunke/circos/sof/circos-0.67-7/work/pipeline/etc/config"
    configopen = open(oldconfig).read()
    f1 = re.sub("pathway",npath,configopen)
    newconfig = "/home/liyunke/circos/sof/circos-0.67-7/work/pipeline/etc/"+prefix+".conf"
    newconfigsave = open(newconfig,"w+")
    newconfigsave.write(f1)
    newconfigsave.close()
    def main():
    p = argparse.ArgumentParser(usage='./circos.pipline.py [--vcf] [--sv] [--prefix] [--outdir] ', description='circos snp sv')
    p.add_argument('-v','--vcf', type=str, help='vcf file')
    p.add_argument('-s','--sv', type=str, help='sv file')
    p.add_argument('-p','--prefix', default="circostest",help='prefix or usrname')
    p.add_argument('-o','--outdir', default="./", help='document directory')
    args = p.parse_args()
    prefix = args.prefix
    vcffile = args.vcf
    outdir = args.outdir
    vcf_SNPdensity(vcffile,outdir)
    SNPdensityfile = outdir+"/SNPdensity100K.snpden"
    density(SNPdensityfile,outdir)
    densitysplit(SNPdensityfile,outdir)
    svdensityfile = args.sv
    sv_split(svdensityfile,outdir)
    circos_config(outdir,prefix)
    cmdstr = "/home/liyunke/circos/sof/circos-0.67-7/bin/circos -conf /home/liyunke/circos/sof/circos-0.67-7/work/pipeline/etc/"+prefix+".conf --outputdir "+ outdir+" -outputfile "+prefix
    result_analysis_circos =os.popen(cmdstr)
    print result_analysis_circos.read()
    rmcmd = "rm "+ outdir +"/hs*"
    result_analysis_rm =os.popen(rmcmd)
    print result_analysis_rm.read()
    if __name__ == '__main__':
    main()

  • 相关阅读:
    VS2008开发的MFC程序,静态连接的方法
    [delphi]参数带有默认值的函数
    __cplusplus的用处
    去掉输入法上的CH和EN
    Linux下Socket的简单使用及最简化封装
    VS2008 _CRT_SECURE_NO_WARNINGS 的问题
    VC:对话框中菜单的使用(WM_INITMENUPOPUP)
    VC:CFindReplaceDialog(非模态对话框、IsWindow()、m_fr、GetFindString())
    献给初学编程者
    VC:状态栏(AfxGetMainWnd()、GetDescendantWindow()、SetPaneInfo()、SetPaneText())
  • 原文地址:https://www.cnblogs.com/wipy/p/4600553.html
Copyright © 2020-2023  润新知