unafold

Function

Description

unafold can find structures containing GC, AU (or AT) and GU (or GT) base pairs. It relies on Gibbs free enthalpy of formation calculations and has tables with energy values for base stacking and the effect of the formation of stems, hairpin loops, bulge loops and multibranch loops. It does however not take into account tertiary interactions, modified bases and the possible formation of "pseudo-knots".

By default unafold considers the molecule to be linear RNA at 37°C. You can however precise that it must consider the molecule to be circular, (single-stranded) DNA and/or at another temperature. The current version (3.6) of UNAFOLD has three different built-in sets of energy tables :

• for RNA at 37°C (from Turner, introduced in version 3.0)
• for RNA at other temperatures (also from Turner, older data introduced in version 2.3)
• for DNA (from SantaLucia, introduced in version 3.1). The model for DNA takes the salt concentration into account.

unafold can operate in 3 different modes :

1. energy minimization : in this mode unafold searches the folding with the lowest Gibbs free enthalpy.
2. energy minimization with suboptimal foldings : this corresponds to the mode of action of the earlier version of the suite called MFOLD. unafold reports besides the optimal folding also a series of suboptimal foldings with free enthalpies of formation within a certain percentage above the optimal. It also creates an "energy dotplot" with dots indicating base pairs in these foldings.
3. partition function calculation : instead of mimimizing energy it computes the complete partition function and from it base pair probabilities. It samples stochastically a number of foldings and draws a "basepair probability plot".

The chance of finding the correct folding can be greatly improved by adding a number of "constaints", that is provividing directives (based on experimental evidence) about which bases cannot or must pair (see input).

Algorithm

• Pseudoknots are prohibited. That is, if r_i-r_j and r_i'-r_j' are part of the same folding and assuming i < i', either i < i' < j' < j or i < j < i' < j' but never i < i' < j < j'. The reason is that when pseudoknots are included, the loop decomposition of a secondary structure breaks down and the energy rules break down. Although we can assign reasonable free energies to the helices in a pseudoknot, and even to possible coaxial stacking between them, it is not possible to estimate the effects of the new kinds of loops that are created.
• In the same folding base pairs r_i-r_j and r_i'-r_j' are only allowed if and . Thus base triples are deliberately excluded from the definition of secondary structure. The reason is that base triples pose an even greater challenge, because the exact nature of the triple cannot be predicted in advance, and even if it could, we have no data for assigning free energies.
• Sharp U-turns are prohibited. A U-turn, called a hairpin loop, must contain at least 3 bases.
• By default only stems of size at least 2 are considered, with other words lone base pairs (base pairs that do not stack on other base pairs) are not considered. You can allow lone base pairs by setting -minstem=1.
• By default only loops up to 30 bases long are considered. You can increase this with the parameter maxlp=n.
• You can decrease the computing time by setting a maximum distance between pairing bases with -maxbp=n

If UNAFOLD is working in partition function calculation mode it stochastically samples a number of likely foldings, as set by the -maxfold=n parameter ; note that you can get several times the same folding. If UNAFOLD is working in MFOLD mode it will compute a number of foldings as set by the -suboptimality=P, -window=W and -maxfold=MAX parameters. W may be thought of as a distance parameter. The distance between 2 base pairs i.j and i'.j' may be defined as max{|i-i'|,|j-j'|}. Then if k-1 foldings have already been predicted, the k^th will have at leas W base pairs that are at least a distance W from any of the base pairs in the first k-1 foldings. UNAFOLD continues sampling foldings until k = MAX or until there are no more foldings with a Gibbs free enthalpy of formation above %P of the minimum.

If W is not specified, UNAFOLD will choose its value from this table based on sequence length. The user is encouraged to experiment with this parameter.

Sequence length	Default window size
1-29	0
30-49	1
50-119	2
120-199	3
200-299	5
300-399	7
400-499	8
500-599	10
600-699	11
700-799	12
800-1199	15
1200-1999	20
2000	25

Usage

> unafold -storeall Prediction of optimal and suboptimal RNA or DNA secondary structure Input nucleotide sequence: embl:k00152 1 : energy minimization 2 : energy minimization with suboptimal foldings 3 : partition function calculation Mode [2]: Molecule is single-stranded DNA rather than RNA [N]: Molecule is circular [N]: Base name for output files [K00152]: png : PNG jpg : JPEG Graphic output format [png]:

Go to the input files for this example
Go to the output files for this example

Command line arguments

   Standard (Mandatory) qualifiers (* if not always prompted):
  [-sequence]          sequence   Nucleotide sequence filename and optional
                                  format, or reference (input USA)
   -mode               menu       [2] Mode. The default mode 2 corresponds to
                                  the old MFOLD software. (Values: 1 (energy
                                  minimization); 2 (energy minimization with
                                  suboptimal foldings); 3 (partition function
                                  calculation))
   -dna                toggle     Molecule is single-stranded DNA rather than
                                  RNA
*  -circular           boolean    Molecule is circular
   -basename           string     [$(sequence.name)] Base name for output
                                  files (Any string is accepted)
   -graph              menu       [png] Graphic output format (Values: png
                                  (PNG); jpg (JPEG))

   Additional (Optional) qualifiers (* if not always prompted):
   -temperature        integer    [37] Temperature (Integer from -5 to 105)
*  -na                 float      [1.0] Na+ molar concentration (for DNA only)
                                  (Number 0.000 or more)
*  -mg                 float      [0.0] Mg++ molar concentration (for DNA
                                  only) (Number 0.000 or more)
   -minstem            integer    [2] Minimum stem size. The default implies
                                  that a helix of size 1, that is a base pair
                                  that does not stack on another base pair on
                                  either side, is not considered (Integer 1 or
                                  more)
   -maxloop            integer    [30] Maximum interior or bulge loop size
                                  (Integer 0 or more)
   -maxbp              integer    [No limit] Maximum distance between pairing
                                  bases (Integer 0 or more)
*  -suboptimality      integer    [5] Maximum percent free energy
                                  suboptimality for computing suboptimal
                                  structures when UNAFOLD is working in MFOLD
                                  mode (Integer from 0 to 100)
*  -window             integer    [depends on sequence length] Base pair
                                  distance window for suboptimal foldings when
                                  UNAFOLD is working in MFOLD mode. A
                                  negative value sets sequence length
                                  dependent value (see on-line manual). (Any
                                  integer value)
*  -maxfold            integer    [maximum 100 suboptimal foldings or 10
                                  stochastically sampled foldings] Number of
                                  predicted foldings. If UNAFOLD is working in
                                  MFOLD mode this corresponds to the maximum
                                  number of foldings within the limits set by
                                  the -suboptimality and -window parameters,
                                  if UNAFOLD is working in partition function
                                  calculation mode this corresponds to the
                                  number of stochastically sampled foldings.
                                  (Integer 1 or more)
*  -annotation         menu       [none] Structure annotation mode. If UNAFOLD
                                  is working in MFOLD mode 'ann' and
                                  'ss-count' stand for colouring the bases
                                  according to the number of alternative base
                                  pairings respectively the number of
                                  single-stranded foldings. If UNAFOLD is
                                  working in partition function calculation
                                  mode 'ann' and 'ss-count' stand for
                                  colouring the bases according to base
                                  pairing probability respectively
                                  single-stranded probability. (Values: none
                                  (No annotation); ann (Base pairing
                                  annotation); ss-count (Single stranded
                                  annotation))

   Advanced (Unprompted) qualifiers:
   -constraints        infile     File with constraints that force or prohibit
                                  base pairings. See on-line manual for
                                  syntax
   -cutoff             float      [1e-6] Cutoff for display in basepair
                                  probability plot (Number from 0.000 to
                                  1.000)
   -display            menu       [auto] Structure display mode (Values: auto
                                  (write base symbols if sequence is shorter
                                  than 800); bases (always write base
                                  symbols); lines (do not write base symbols))
   -numbering          integer    [depends on sequence length] Base numbering
                                  frequency. A value of 0 triggers the default
                                  (10 for seq. shorter than 50, 50 for seq.
                                  longer than 300, otherwise 20). (Integer 0
                                  or more)
   -storeall           boolean    Store all UNAFOLD output files

   Associated qualifiers:

   "-sequence" associated qualifiers
   -sbegin1            integer    Start of the sequence to be used
   -send1              integer    End of the sequence to be used
   -sreverse1          boolean    Reverse (if DNA)
   -sask1              boolean    Ask for begin/end/reverse
   -snucleotide1       boolean    Sequence is nucleotide
   -sprotein1          boolean    Sequence is protein
   -slower1            boolean    Make lower case
   -supper1            boolean    Make upper case
   -sformat1           string     Input sequence format
   -sdbname1           string     Database name
   -sid1               string     Entryname
   -ufo1               string     UFO features
   -fformat1           string     Features format
   -fopenfile1         string     Features file name

   General qualifiers:
   -auto               boolean    Turn off prompts
   -stdout             boolean    Write standard output
   -filter             boolean    Read standard input, write standard output
   -options            boolean    Prompt for standard and additional values
   -debug              boolean    Write debug output to program.dbg
   -verbose            boolean    Report some/full command line options
   -help               boolean    Report command line options. More
                                  information on associated and general
                                  qualifiers can be found with -help -verbose
   -warning            boolean    Report warnings
   -error              boolean    Report errors
   -fatal              boolean    Report fatal errors
   -die                boolean    Report dying program messages

Standard (Mandatory) qualifiers		Allowed values	Default
[-sequence] (Parameter 1)	Nucleotide sequence filename and optional format, or reference (input USA)	Readable sequence	Required
-mode	Mode. The default mode 2 corresponds to the old MFOLD software.	1 (energy minimization) 2 (energy minimization with suboptimal foldings) 3 (partition function calculation)	2
-dna	Molecule is single-stranded DNA rather than RNA	Toggle value Yes/No	No
-circular	Molecule is circular	Boolean value Yes/No	No
-basename	Base name for output files	Any string is accepted	<sequence name>
-graph	Graphic output format	png (PNG) jpg (JPEG)	png
Additional (Optional) qualifiers		Allowed values	Default
-temperature	Temperature	Integer from -5 to 105	37
-na	Na+ molar concentration (for DNA only)	Number 0.000 or more	1.0
-mg	Mg++ molar concentration (for DNA only)	Number 0.000 or more	0.0
-minstem	Minimum stem size. The default implies that a helix of size 1, that is a base pair that does not stack on another base pair on either side, is not considered	Integer 1 or more	2
-maxloop	Maximum interior or bulge loop size	Integer 0 or more	30
-maxbp	Maximum distance between pairing bases	Integer 0 or more	No limit
-suboptimality	Maximum percent free energy suboptimality for computing suboptimal structures when UNAFOLD is working in MFOLD mode	Integer from 0 to 100	5
-window	Base pair distance window for suboptimal foldings when UNAFOLD is working in MFOLD mode. A negative value sets sequence length dependent value (see on-line manual).	Any integer value	depends on sequence length
-maxfold	Number of predicted foldings. If UNAFOLD is working in MFOLD mode this corresponds to the maximum number of foldings within the limits set by the -suboptimality and -window parameters, if UNAFOLD is working in partition function calculation mode this corresponds to the number of stochastically sampled foldings.	Integer 1 or more	maximum 100 suboptimal foldings or 10 stochastically sampled foldings
-annotation	Structure annotation mode. If UNAFOLD is working in MFOLD mode 'ann' and 'ss-count' stand for colouring the bases according to the number of alternative base pairings respectively the number of single-stranded foldings. If UNAFOLD is working in partition function calculation mode 'ann' and 'ss-count' stand for colouring the bases according to base pairing probability respectively single-stranded probability.	none (No annotation) ann (Base pairing annotation) ss-count (Single stranded annotation)	none
Advanced (Unprompted) qualifiers		Allowed values	Default
-constraints	File with constraints that force or prohibit base pairings. See on-line manual for syntax	Input file	Required
-cutoff	Cutoff for display in basepair probability plot	Number from 0.000 to 1.000	1e-6
-display	Structure display mode	auto (write base symbols if sequence is shorter than 800) bases (always write base symbols) lines (do not write base symbols)	auto
-numbering	Base numbering frequency. A value of 0 triggers the default (10 for seq. shorter than 50, 50 for seq. longer than 300, otherwise 20).	Integer 0 or more	depends on sequence length
-storeall	Store all UNAFOLD output files	Boolean value Yes/No	No

Input file format

You can optionally (-constraints=infile) give to mfold a "constraints" file with information (based on laboratory evidence) about which bases cannot or must pair. This can greatly improve the chance of finding the correct structure. An example of a file, that works for the usage example (it makes mfold find the correct structure intstead of what you get by default), is :

constraint file for usage example

P 16 0 1
P 35 0 3
P 48 0 1
P 73 0 1
P 75 0 2

The complete set of syntax rules for the "constraints" file is as follows :

1. Forcing a string of consecutive bases to pair :
   Syntax : F   i   0   k
   Ribonucleotides r_i,r_i+1,...,r_i+k-1 are forced to be double stranded. The partners for these bases are chosen by the program. As an example, the command :
   F   23   0   5
   would cause bases 23, 24, 25, 26 and 27 to pair.

2. Forcing a string of consecutive base pairs :
   Syntax : F   i   j   k
   Base pairs r_i-r_j,r_i+1-r_j-1,r_i+2-r_j-2,...,r_i+k-1-r_j-k+1 are forced to occur. This is the same thing as forcing a helix to form. The helix is designated by its (external) closing base pair, i.j. As an example, the command :
   F   2   110   3
   would force base pairs 2.110, 3.109 and 4.108. Note that these ase pairs must be able to form ! Be aware also that unafold by default filters out isolated base pairs.

3. Prohibiting a string of consecutive bases from pairing :
   Syntax : P   i   0   k
   Ribonucleotides r_i,r_i+1,...,r_i+k-1 are prevented from pairing.

4. Prohibiting a string of consecutive base pairs :
   Syntax : P   i   j   k
   r_i-r_j,r_i+1-r_j-1,r_i+2-r_j-2,...,r_i+k-1-r_j-k+1 are not allowed to form. This is equivalent to prohibiting a helix.

5. Prohibiting 1 segment of a sequence from pairing with another :
   Syntax : P   i-j   k-l
   where and . In this case, no base pairs are allowed between r_i,r_i+1,...,r_j and r_k,r_k+1,...,r_l. Note that the 2 segments need not be distinct. For example, the command :
   P   i-j   i-j
   will not allow r_i,r_i+1,...,r_j to pair with itself.

Output file format

Output files for usage example

File: K00152.png

File: K00152.plot

level	length	i	j	energy
1	5	49	65	-297
1	6	35	57	-287
1	3	38	47	-292
1	1	29	49	-293
1	6	30	48	-296
1	7	1	72	-297
1	4	10	44	-297
1	1	13	34	-294
1	5	14	33	-297
1	1	19	28	-291
1	3	19	27	-297
1	2	10	35	-284
1	4	10	25	-296

This file contains instructions for drawing the energy dotplot. The first record is a header, and each subsequent record describes a single helix, with the position of the first base, the position of the last base and the length. The energy is the smallest free energy change (dG), in kcal/mol, of the folding that still contains the helix.

File: K00152_1.png

File: K00152_1.ct

76	dG = -29.7	K00152 K00152.1 E.coli Arg-tRNA-1.
1	G	0	2	72	1	0	2
2	C	1	3	71	2	1	3
3	A	2	4	70	3	2	4
4	T	3	5	69	4	3	5
5	C	4	6	68	5	4	6
6	C	5	7	67	6	5	7
7	G	6	8	66	7	6	8

  [Part of this file has been deleted for brevity]

70	T	69	71	3	70	69	71
71	G	70	72	2	71	70	72
72	C	71	73	1	72	71	73
73	A	72	74	0	73	72	0
74	C	73	75	0	74	0	0
75	C	74	76	0	75	0	0
76	A	75	0	0	76	0	0

This file contains instructions for drawing the secondary structure above (see below for exact syntax).

[The files K00152.ann, K00152.dG, and the structure 2, 3 and 4 have been omitted for brevity]

Prediction reliability annotation

1. none : secondary structures are drawn without any special annotation. Letters or outline are in black, while base pairs are coloured lines (red for GC, blue for AU and green for GU).
2. ann : Base characters are replaced or covered by coloured dots that show how well-determined the folding is. If unafold in running in MFOLD mode a "p-num" file with the number of alternative pairings is used, bases that are never paired or are not paired with more than one partner are in red, the other bases are in a scale from red to black according to the number of alternative partners. If unafold in running in partition function calculation mode a "prob" file with the probabilities is used and the probabilities are expressed in a scale from black to red. Also see below for the colour scale and the .ann file syntax.
3. ss-count : Base characters are replaced or covered by coloured dots that show how likely they are to be single-stranded. An "ss-count" file is used. If unafold in running in MFOLD mode bases that are never paired are in red, bases that are paired in all alternative foldings are in black and the others in an intermediate colour. If unafold in running in partition function calculation mode a scale from black to red indicates the prabability of the base being single-stranded.

The full list of the possible text output files

• <basename>.ann : an "annotation" file used for adding annotation to the graphical representation of the foldings. If unafold is run in MFOLD mode there is a single file with for every base the number of alternative partners with which it can form a base-pair in the different alternative foldings (or a 0 if it is always single-stranded). If unafold is run in partition function calculation mode there is a file for each of the stochastically sampled foldings, with for each paired base the probability of formation of that particular base pair and for each unpaired base the probability that that base is single-stranded.
• <basename>.ct : a file that defines a nucleic acid sequence together with its structure. There is one for each folding. Besides the header it contains lines of type :
  

  19      G       18      20      27      19      18      20

  There are as many lines as bases in the sequence and the 8 colums have the following meaning :
  1. the position number of the base (from 1 to N)
  2. the base symbol (A, C, G, T or U)
  3. the position of the 5'-upstream base. For the first base this number is 0 if the sequence is linear and N if the sequence is circular.
  4. the position of the 3'-downstream base. For the last base this number is 0 if the sequence is linear and 1 if the sequence is circular.
  5. If the base form a base pair this is the position number of its partner ; if the base is single-stranded it is 0.
  6. (redundant with column 1, but could be used to indicate the "historical numbering" if a program "remembers" the numbering of a larger sequence of which the handled sequence is a fragment)
  7. if the base is not involved in stacking this is 0, otherwise it is the position number of the base that stacks 5'
  8. if the base is not involved in stacking this is 0, otherwise it is the position number of the base that stacks 3'
• <basename>.dG : if unafold is run in energy minimization mode or MFOLD mode this file contains the dG (free enthalpy change of formation) for the most stable folding, expressed in kcal/mol, as well as a Z that could be interpreted as the equilibrium constant for folding. If unafold is run in partition function calculation mode this file contains the partition function Z and the corresponding Gibbs free enthalpy (dG=-RTlnZ).
• <basename>.ext : if unafold is run in partition function calculation mode this file contains for each base the probability that the base is single-stranded and the joint probability that the base as well as the 3'-downstream base are both single-stranded (in energy minimization mode a quite useless file with 0's and 1's is created).
• <basename>.plot : if unafold is run in partition function calculation mode this file contains for each possible base pair the probability of formation of that base pair. if unafold is run in energy minimization mode this file contains just the list of base pairs in the most stable folding (with 1's). If unafold is run in MFOLD mode this file contains for each possible helix a line of type :
  

  1       6       35      57      -287

  with :
  1. (deprecated, just the number 1)
  2. the length of the helix
  3. the position of the first base in the helix
  4. the position of the last base in the helix (which forms a base-pair with the first base, so the pairs go from 35-57 to 40-52)
  5. the Gibss free energy of formation, not of the helix itself (!) but of the least stable folding that still contains it
• <basename>.ss-count : only created when 'ss-count' annotation is requested, contains for each base the number of alternative foldings in which the base is single-stranded

Data files

Notes

References

1. The "classic" article on the suboptimal algorithm :
   M. Zuker
   On Finding All Suboptimal Foldings of an RNA Molecule. Science 244, 48-52 (1989)
2. The reference for the MFOLD software :
   M. Zuker, D.H. Mathews & D.H. Turner
   Algorithms and Thermodynamics for RNA Secondary Structure Prediction: A Practical Guide in RNA Biochemistry and Biotechnology, J. Barciszewski & B.F.C. Clark, eds., NATO ASI Series, Kluwer Academic Publishers, (1999)
3. The reference for the UNAFOLD software :
   N.R. Markham & M. Zuker
   UNAFold: Software for Nucleic Acid Folding and Hybridization. Methods Mol Biol. 453:3-31, 2008.
4. The version 2.3 free energy rules for folding RNA :
   A.E. Walter, D.H. Turner, J. Kim, M.H. Lyttle, P. Müller, D.H. Mathews & M. Zuker
   Coaxial stacking of helixes enhances binding of oligoribonucleotides and improves predictions of RNA folding. Proc. Natl. Acad. Sci. USA 91, 9218-9222 (1994)
5. The version 3.1 free energy rules for folding RNA :
   D.H. Mathews, J. Sabina, M. Zuker & D.H. Turner
   Expanded Sequence Dependence of Thermodynamic Parameters Provides Robust Prediction of RNA Secondary Structure J. Mol. Biol. 288:911-940 (1999)
6. The free energy rules for folding DNA :
   J.Jr SantaLucia
   "A unified view of polymer, dumbbell, and oligonucleotide DNA nearest-neighbor thermodynamics." Proc. Natl. Acad. Sci. USA 95, 1460-1465 (1998)
7. Annotation of base pair prediction reliability for MFOLD :
   M. Zker & A.B. Jacobson
   Using reliability information to annotate RNA secondary structures. RNA, 4:669-679, 1998.
8. Annotation of base pair probabilities for partition function calculation mode :
   D.H. Matthews
   Using an RNA secondary structure partition function to determine confidence in base pairs predicted by free energy minimizatin. RNA, 10:1174-1177, 2004.

Warnings

Diagnostic Error Messages

      No foldings were found  

If you add -annotation=ann or -annotation=ss-count to the command line, while the program is running in energy minimization mode, it issues the following warning :

      annotation reset to 'none' ! If only one folding is computed,
            counting alternative base pairings is pointless.

Exit status

Known bugs

See also

Author(s)

The UNAFOLD suite itself was written by Michael Zuker (zukerm@rpi.edu) and Nicholas Markham (markhn@rpi.edu) at Rensselaer Polytechnic Institute (Troy, New York, USA).

History

Completed 23 October 2008.

Target users