Created
June 17, 2015 12:56
-
-
Save gedankenstuecke/f237566c2a154e1795b5 to your computer and use it in GitHub Desktop.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
# this one works, plink is in cwd + conversion to BED files works flawlessly | |
for f in usable_files: | |
# gid is the genotype-ID | |
gid = f.split("/")[-1].split("_")[1].replace("file","") | |
# converts the genotyping file to plink format, using the gid as sample name | |
call = "./plink --23file "+ f + " F" + gid + "ID" + gid + "I 1" | |
call += " --out 23andme_plink/genotypeid_" + gid | |
print "convert gid " + gid | |
subprocess.call(call,shell=True) | |
# this is part of the same script, plink is still in cwd but it crashes, saying | |
# "Error: Failed to open /home/bastia.fam" | |
for i in list_bed: | |
# check that files have been processed and are working so far. | |
if os.path.isfile(i.replace(".bed",".fam")) and os.path.isfile(i.replace(".bed",".bim")): | |
listhandle.write(cwd + "/"+i.replace(".bed","")+"\n") | |
call_merge = "./plink --bfile " + cwd + "/" + start_bed + " --merge-list " + cwd + "/merge_list.txt --make-bed --out " + cwd + "/23andme_merged/merge-pass1" | |
print "merging 23andme data, first pass" | |
subprocess.call(["./plink","--bfile",cwd+"/"+start_bed,"--merge-list",cwd + "/merge_list.txt","--make-bed","--out",cwd+"/23andme_merged/merge-pass1"]) |
Nope, the problem remains. Actually first you have to change "--bfile" etc to " --bfile" to make sure it won't concatenate into a single whitespace-free string, but even then it won't work.
import subprocess
import os
from glob import glob
usable_files = glob('files/*')
# this one works, plink is in cwd + conversion to BED files works flawlessly
for f in usable_files:
gid = os.path.basename(f).split('.')[-1]
call = 'plink2 --23file %s F%sID%s I 1 --out 23andme_plugin/genotypeid_%s' \
% (f, gid, gid, gid)
subprocess.check_call(call, shell=True)
# this is part of the same script, plink is still in cwd but it crashes, saying
# "Error: Failed to open /home/bastia.fam"
cwd = os.path.dirname(os.path.abspath(__file__))
list_bed = glob('23andme_plugin/*bed')
with open('merge_list.txt', 'w') as listhandle:
for i in list_bed:
# check that bed and bim files exist,
# write bed file to merge_list
print i
if os.path.isfile(i.replace('.bed', '.fam')) \
and os.path.isfile(i.replace('.bed', '.bim')):
listhandle.write(cwd + '/' + i.replace('.bed', '') + '\n')
start_bed = list_bed[0].replace('.bed', '')
cmd = ' '.join([
'plink2',
'--bfile %s/%s' % (cwd, start_bed),
'--merge-list %s/%s' % (cwd, '/merge_list.txt'),
'--make-bed',
'--out %s/%s' % (cwd, '23andme_merged/merge-pass1')
])
# if subprocess.call doesn't work, copy paste this to debug
print cmd
subprocess.call(cmd, shell=True)
Put the input files in files
. Make sure you directories 23andme_merged
and 23andme_plugin
exist before running.
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
I refactored a bit so I could read it ;)
What happens when you change line 21 to the following snippet? If you copy/paste the command does it run?