/viennarna-structure-prediction

name: "viennarna-structure-prediction"
description: "Predict RNA secondary structure, MFE folding, base-pair probabilities, RNA-RNA interactions via ViennaRNA Python bindings. Pipeline: sequence → MFE → partition function and pair-probability matrix → dot-bracket → duplex. Use for siRNA/sgRNA targeting, ribozyme design, RNA accessibility. Use RNAfold CLI for batch use without Python."
license: "MIT"

ViennaRNA Structure Prediction

Overview

ViennaRNA is the gold-standard toolkit for RNA secondary structure prediction based on thermodynamic nearest-neighbor parameters. It predicts the minimum free energy (MFE) structure and dot-bracket notation for a given RNA sequence, computes the full partition function to obtain base pair probabilities, and models RNA-RNA interactions via co-folding and duplex prediction. The Python bindings (`import RNA`) expose the full ViennaRNA C library with sequence-level and fold-compound APIs. Command-line programs (`RNAfold`, `RNAalifold`, `RNAduplex`) are also available and demonstrated here.

When to Use

• Predicting the minimum free energy secondary structure of an RNA sequence (mRNA, lncRNA, miRNA precursor, aptamer)
• Computing base pair probability matrices to assess structural uncertainty and identify well-defined stem-loops
• Designing or evaluating siRNA accessibility by folding the target mRNA region and checking for double-stranded structure
• Assessing sgRNA targeting efficiency by predicting guide RNA secondary structure that may reduce on-target activity
• Modeling RNA-RNA interactions (co-folding or duplex prediction) for miRNA-target binding or antisense oligonucleotide design
• Calculating folding free energies for a set of sequences to compare thermodynamic stability
• Use `mfold` (web server) or `RNAstructure` instead when you need Mfold algorithm predictions specifically or need the Efold partition function; ViennaRNA uses the Turner 2004 nearest-neighbor parameters and is the standard for research-grade thermodynamic prediction

Prerequisites

• **Python packages**: `ViennaRNA` (Python bindings), `matplotlib`, `numpy`
• **Data requirements**: RNA sequences as strings (ACGU alphabet; T is auto-converted to U by ViennaRNA)
• **Environment**: Python 3.8+; conda installation strongly recommended (handles C library dependencies)

# Install via conda (recommended)
conda install -c conda-forge -c bioconda viennarna

# Verify installation
python -c "import RNA; print(RNA.__version__)"
# 2.6.4

# Install additional Python dependencies
pip install matplotlib numpy pandas

# Optional: verify CLI tools are available
RNAfold --version
# RNAfold 2.6.4

Quick Start

import RNA

# Predict MFE structure for an RNA sequence
sequence = "GCGGAUUUAGCUCAGUUGGGAGAGCGCCAGACUGAAGAUCUGGAGGUCCUGUGUUCGAUCCACAGAAUUCGCACCA"
structure, mfe = RNA.fold(sequence)

print(f"Sequence:  {sequence}")
print(f"Structure: {structure}")
print(f"MFE:       {mfe:.2f} kcal/mol")
# Sequence:  GCGGAUUUAGCUCAGUUGGGAGAGCGCCAGACUGAAGAUCUGGAGGUCCUGUGUUCGAUCCACAGAAUUCGCACCA
# Structure: (((((((..((((........)))).(((((.......))))).....(((((.......))))))))))))....
# MFE:       -31.30 kcal/mol

Workflow

Step 1: Sequence Preparation and MFE Folding

Load an RNA sequence and compute its minimum free energy secondary structure using `RNA.fold()`. Validate the input and inspect the dot-bracket output.

import RNA

def prepare_sequence(seq: str) -> str:
    """Normalize sequence: uppercase, replace T→U, validate alphabet."""
    seq = seq.upper().replace("T", "U").strip()
    invalid = set(seq) - set("ACGUNX")
    if invalid:
        raise ValueError(f"Invalid characters in sequence: {invalid}")
    return seq

# E. coli tRNA-Phe (GenBank: M10217)
raw_seq = "GCGGAUUUAGCUCAGUUGGGAGAGCGCCAGACUGAAGAUCUGGAGGUCCUGUGUUCGAUCCACAGAAUUCGCACCA"
sequence = prepare_sequence(raw_seq)

structure, mfe = RNA.fold(sequence)

print(f"Sequence length: {len(sequence)} nt")
print(f"Structure:       {structure}")
print(f"MFE:             {mfe:.2f} kcal/mol")

# Validate: structure length must equal sequence length
assert len(structure) == len(sequence), "Structure and sequence length mismatch"

# Count stems (paired bases)
n_paired   = structure.count("(") + structure.count(")")
n_unpaired = structure.count(".")
print(f"Paired bases: {n_paired}  |  Unpaired bases: {n_unpaired}")
print(f"Stem fraction: {n_paired/len(sequence):.2f}")

Step 2: Create a Fold Compound for Advanced Analysis

The `RNA.fold_compound` object is the central API for partition function, base pair probabilities, and constrained folding.

import RNA

sequence = "GCGGAUUUAGCUCAGUUGGGAGAGCGCCAGACUGAAGAUCUGGAGGUCCUGUGUUCGAUCCACAGAAUUCGCACCA"

# Create fold compound (wraps the sequence with model parameters)
fc = RNA.fold_compound(sequence)

# Compute MFE structure via the fold compound API
structure, mfe = fc.mfe()
print(f"MFE structure: {structure}")
print(f"MFE:           {mfe:.2f} kcal/mol")

# Evaluate free energy of an alternative structure
alt_structure = "." * len(sequence)   # fully unfolded
energy = fc.eval_structure(alt_structure)
print(f"Fully unfolded energy: {energy:.2f} kcal/mol")
print(f"Folding stabilization:  {energy - mfe:.2f} kcal/mol")

Step 3: Partition Function and Base Pair Probabilities

Compute the thermodynamic partition function to obtain ensemble-level base pair probabilities. High-probability pairs indicate well-defined structural elements.

import RNA
import numpy as np

sequence = "GCGGAUUUAGCUCAGUUGGGAGAGCGCCAGACUGAAGAUCUGGAGGUCCUGUGUUCGAUCCACAGAAUUCGCACCA"
n = len(sequence)

fc = RNA.fold_compound(sequence)

# Step 1: MFE folding (required before pf for proper initialization)
structure_mfe, mfe = fc.mfe()

# Step 2: Rescale Boltzmann factors for numerical stability (optional but recommended)
fc.exp_params_rescale(mfe)

# Step 3: Compute partition function
structure_pf, gibbs_free_energy = fc.pf()
print(f"Gibbs free ene