/lib/galaxy/datatypes/converters/pbed_ldreduced_converter.py
Python | 123 lines | 116 code | 4 blank | 3 comment | 0 complexity | 1edc8ea713b989f572d2d04a7ce2249e MD5 | raw file
1# converter for ldreduced rgenetics datatype 2# used for grr and eigenstrat - shellfish if we get around to it 3# 4 5import os,sys,tempfile,subprocess,time 6 7from galaxy import eggs 8 9prog="pbed_ldreduced_converter.py" 10 11galhtmlprefix = """<?xml version="1.0" encoding="utf-8" ?> 12<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN" "http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd"> 13<html xmlns="http://www.w3.org/1999/xhtml" xml:lang="en" lang="en"> 14<head> 15<meta http-equiv="Content-Type" content="text/html; charset=utf-8" /> 16<meta name="generator" content="Galaxy %s tool output - see http://getgalaxy.org" /> 17<title></title> 18<link rel="stylesheet" href="/static/style/base.css" type="text/css" /> 19</head> 20<body> 21<div class="document"> 22""" 23 24plinke = 'plink' 25 26 27def timenow(): 28 """return current time as a string 29 """ 30 return time.strftime('%d/%m/%Y %H:%M:%S', time.localtime(time.time())) 31 32 33def pruneLD(plinktasks=[],cd='./',vclbase = []): 34 """ 35 """ 36 fplog,plog = tempfile.mkstemp() 37 alog = [] 38 alog.append('## Rgenetics: http://rgenetics.org Galaxy Tools rgQC.py Plink pruneLD runner\n') 39 for task in plinktasks: # each is a list 40 vcl = vclbase + task 41 sto = file(plog,'w') 42 x = subprocess.Popen(' '.join(vcl),shell=True,stdout=sto,stderr=sto,cwd=cd) 43 retval = x.wait() 44 sto.close() 45 try: 46 lplog = file(plog,'r').readlines() 47 lplog = [x for x in lplog if x.find('Pruning SNP') == -1] 48 alog += lplog 49 alog.append('\n') 50 os.unlink(plog) # no longer needed 51 except: 52 alog.append('### %s Strange - no std out from plink when running command line\n%s\n' % (timenow(),' '.join(vcl))) 53 return alog 54 55 56def makeLDreduced(basename,infpath=None,outfpath=None,plinke='plink',forcerebuild=False,returnFname=False, 57 winsize="60", winmove="40", r2thresh="0.1" ): 58 """ not there so make and leave in output dir for post job hook to copy back into input extra files path for next time 59 """ 60 ldr = basename # we store ld reduced and thinned data 61 ldreduced = os.path.join(outfpath,ldr) # note where this is going 62 outbase = os.path.join(outfpath,basename) 63 inbase = os.path.join(infpath) 64 loglines = [] 65 ldbedname = '%s.bed' % ldreduced 66 bedname = '%s.bed' % basename 67 ldbedfn = os.path.join(infpath,ldbedname) 68 bedfn = os.path.join(infpath,bedname) 69 bmap = os.path.join(infpath,'%s.bim' % basename) 70 plinktasks = [] 71 vclbase = [plinke,'--noweb'] 72 plinktasks += [['--bfile',inbase,'--indep-pairwise %s %s %s' % (winsize,winmove,r2thresh),'--out %s' % outbase], 73 ['--bfile',inbase,'--extract %s.prune.in --make-bed --out %s' % (outbase, outbase)]] 74 vclbase = [plinke,'--noweb'] 75 loglines = pruneLD(plinktasks=plinktasks,cd=outfpath,vclbase = vclbase) 76 77def main(): 78 """ 79 need to work with rgenetics composite datatypes 80 so in and out are html files with data in extrafiles path 81 82 .. raw:: xml 83 84 <command interpreter="python"> 85 pbed_ldreduced_converter.py '$input1.extra_files_path/$input1.metadata.base_name' '$winsize' '$winmove' '$r2thresh' 86 '$output1' '$output1.files_path' 'plink' 87 </command> 88 89 """ 90 nparm = 7 91 if len(sys.argv) < nparm: 92 sys.stderr.write('## %s called with %s - needs %d parameters \n' % (prog,sys.argv,nparm)) 93 sys.exit(1) 94 inpedfilepath = sys.argv[1] 95 base_name = os.path.split(inpedfilepath)[-1] 96 winsize = sys.argv[2] 97 winmove = sys.argv[3] 98 r2thresh = sys.argv[4] 99 outhtmlname = sys.argv[5] 100 outfilepath = sys.argv[6] 101 try: 102 os.makedirs(outfilepath) 103 except: 104 pass 105 plink = sys.argv[7] 106 makeLDreduced(base_name,infpath=inpedfilepath,outfpath=outfilepath,plinke=plink,forcerebuild=False,returnFname=False, 107 winsize=winsize,winmove=winmove,r2thresh=r2thresh) 108 f = file(outhtmlname,'w') 109 f.write(galhtmlprefix % prog) 110 flist = os.listdir(outfilepath) 111 s1 = '## Rgenetics: http://rgenetics.org Galaxy Tools %s %s' % (prog,timenow()) # becomes info 112 s2 = 'Input %s, winsize=%s, winmove=%s, r2thresh=%s' % (base_name,winsize,winmove,r2thresh) 113 print '%s %s' % (s1,s2) 114 f.write('<div>%s\n%s\n<ol>' % (s1,s2)) 115 for i, data in enumerate( flist ): 116 f.write('<li><a href="%s">%s</a></li>\n' % (os.path.split(data)[-1],os.path.split(data)[-1])) 117 f.write("</div></body></html>") 118 f.close() 119 120 121if __name__ == "__main__": 122 main() 123