Question: Folding RNA on python
1
gravatar for gprezza
7 months ago by
gprezza10
gprezza10 wrote:

Hi all,

I'd like to write a script in python that, given a set of RNA sequences, finds the ones with the strongest predicted folding. A quick and dirty way like folding them and picking the ones with a - delta G below a set threshold would suffice in my case.

I was looking at two "classical" folding packages (mfold and RNAfold from the Vienna package), that I can run on command line from python. However, neither of them has an options in which just the delta G is output. They produce a series or files for each RNA, and it would be complicated to dive through each one just to extract the delta G value. RNA fold has an output that would fit better my needs, but is still not the ideal.

Does anyone have any suggestion on how to proceed?

EDIT: I'll add some examples, as some people suggested. I'm a bit of a python novice, so my apologies if I am missing any obvious ways to solve my issue.

the best I managed to get is through RNAfold. The output looks like this, repeated for each sequence:

>seq1
AUGGCUGUUCGCCAUUUAAAG
(((((.....)))))...... ( -5.20)

It's not too bad, and I could extract the delta G value (-5.20 here) by telling python to extract the value between the brackets. However, I'm not sure that's the most elegant way, and I was wondering if anyone knew of packages that output a tab delimited file, or of ways to have RNAfold (or mfold) do that.

Another "problem" is that both mfold and RNAfold also produce one file with an image of the folding for each sequence (which I don't need), so then I would have to delete them via python, and again I'm not sure this is ideal.

folding rna python • 537 views
ADD COMMENTlink modified 7 months ago by Joe14k • written 7 months ago by gprezza10
1

Are they just text files? Have you tried parsing the text output? Python is usually very good at processing text.

ADD REPLYlink written 7 months ago by curious140

Thanks for the comment, I updated my question.

ADD REPLYlink written 7 months ago by gprezza10
1

Have you looked at RNAsketch? This is a Python wrapper for the Vienna Package, which might give you just what you want.

ADD REPLYlink written 7 months ago by cschu1811.9k

Without examples of your output it’s impossible for us to say how to proceed

ADD REPLYlink written 7 months ago by Joe14k

Thanks for the comment, I updated my question.

ADD REPLYlink written 7 months ago by gprezza10

Your approach you described is a perfectly valid approach. Pull out a number between brackets via a regular expression (it wouldn’t be enough to look for brackets alone, since the dot bracket notation would mess it up.

See my answer below.

ADD REPLYlink modified 7 months ago • written 7 months ago by Joe14k

awk '{print $1"\t"$n}' fill in <n> as whatever the column the value is in, assuming it is a tab deliminated file. $1 is the sample name if not... substitute that for wherever the sample name is on that file.

like Brett said... can't really help you if we don't know what file type you are looking at.

ADD REPLYlink written 7 months ago by ggman80
0
gravatar for Joe
7 months ago by
Joe14k
United Kingdom
Joe14k wrote:

Here’s one approach:

import sys
import re

regex = re.compile(r'\( -?\d\.\d*\)|\( -?0\)')

with open(sys.argv[1], 'r') as fh:
    num = regex.findall("".join(fh.readlines()))[0].lstrip('( ').rstrip(')')
    print(num)

This is quite a cautious approach, as I pick up the ( ) as well as the leading space in your provided example (assuming thats not a typo), so consequently they need stripping off after to get the number.

I’ve also assumed that the minus symbol may be optional, and that 0 is a valid value (which may or may not have a minus sign).

This will break if, for some reason, your sequence name had a similarly formed decimal in it.

ADD COMMENTlink modified 7 months ago • written 7 months ago by Joe14k
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1711 users visited in the last hour