@@ 44-82 (lines=39) @@ | ||
41 | def __init__(self): |
|
42 | pass |
|
43 | ||
44 | def cmscan(self, seq, verbose=False): |
|
45 | """Run cmscan on the seq. |
|
46 | ||
47 | Usage:: |
|
48 | ||
49 | >>> seq = RNASequence("GGCGCGGCACCGUCCGCGGAACAAACGG") |
|
50 | >>> rs = RfamSearch() |
|
51 | >>> hit = rs.cmscan(seq) |
|
52 | >>> print(hit) #doctest: +ELLIPSIS |
|
53 | # cmscan :: search sequence(s) against a CM database... |
|
54 | ||
55 | :param seq: string |
|
56 | :returns: result |
|
57 | :rtype: string """ |
|
58 | #if verbose: |
|
59 | print('RFAM_DB_PATH', RFAM_DB_PATH) |
|
60 | ||
61 | # make tmp file |
|
62 | tf = tempfile.NamedTemporaryFile(delete=False) |
|
63 | tf.name += '.fa' |
|
64 | with open(tf.name, 'w') as f: |
|
65 | f.write('>target\n') |
|
66 | f.write(seq.seq + '\n') |
|
67 | ||
68 | # make output file |
|
69 | of = tempfile.NamedTemporaryFile(delete=False) |
|
70 | ||
71 | # run cmscan |
|
72 | cmd = 'cmscan -E 1 ' + RFAM_DB_PATH + ' ' + tf.name + ' > ' + of.name |
|
73 | print(cmd) |
|
74 | o = subprocess.Popen( |
|
75 | cmd, shell=True, stdout=subprocess.PIPE, stderr=subprocess.PIPE) |
|
76 | # out = o.stdout.read().strip() |
|
77 | err = o.stderr.read().strip() |
|
78 | if err: |
|
79 | raise RfamSearchError(err) |
|
80 | self.output = open(of.name).read() |
|
81 | # os.chdir(old_pwd) |
|
82 | return self.output |
|
83 | ||
84 | ||
85 | # main |
@@ 1069-1104 (lines=36) @@ | ||
1066 | if outputfn: |
|
1067 | self.output = open(outputfn).read().strip().split('\n') |
|
1068 | ||
1069 | def run_cmalign(self, seq, cm, verbose=True): |
|
1070 | """Run cmalign and process the result. |
|
1071 | ||
1072 | :param seq: seq string |
|
1073 | :param cm: cm fn |
|
1074 | ||
1075 | Run:: |
|
1076 | ||
1077 | $ cmalign RF01831.cm 4lvv.seq |
|
1078 | # STOCKHOLM 1.0 |
|
1079 | #=GF AU Infernal 1.1.2 |
|
1080 | ||
1081 | 4lvv -GGAGAGUA-GAUGAUUCGCGUUAAGUGUGUGUGA-AUGGGAUGUCG-UCACACAACGAAGC---GAGA---GCGCGGUGAAUCAUU-GCAUCCGCUCCA |
|
1082 | #=GR 4lvv PP .********.******************9999998.***********.8999999******8...5555...8**************.************ |
|
1083 | #=GC SS_cons (((((----(((((((((((,,,,,<<-<<<<<<<<___________>>>>>>>>>>,,,<<<<______>>>>,,,)))))))))))-------))))) |
|
1084 | #=GC RF ggcaGAGUAGggugccgugcGUuAAGUGccggcgggAcGGGgaGUUGcccgccggACGAAgggcaaaauugcccGCGguacggcaccCGCAUcCgCugcc |
|
1085 | // |
|
1086 | ||
1087 | .. warning :: requires cmalign to be set in your shell |
|
1088 | """ |
|
1089 | tf = tempfile.NamedTemporaryFile(delete=False) |
|
1090 | tf.name += '.seq' |
|
1091 | ||
1092 | with open(tf.name, 'w') as f: |
|
1093 | f.write('>target\n') |
|
1094 | f.write(seq + '\n') |
|
1095 | ||
1096 | cmd = 'cmalign -g ' + cm + ' ' + tf.name # global |
|
1097 | if verbose: |
|
1098 | print('cmd' + cmd) |
|
1099 | o = subprocess.Popen(cmd, shell=True, stdout=subprocess.PIPE, stderr=subprocess.PIPE) |
|
1100 | stdout = o.stdout.read().strip() |
|
1101 | stderr = o.stderr.read().strip() |
|
1102 | if verbose: |
|
1103 | print(stdout) |
|
1104 | self.output = stdout.split('\n') |
|
1105 | ||
1106 | def get_gc_rf(self): |
|
1107 | """Get ``#=GC RF``. |