helloworld_mfe.c
The following is an example showing the minimal requirements to compute the Minimum Free Energy (MFE) and corresponding secondary structure of an RNA sequence
#include <stdlib.h> #include <stdio.h> #include <string.h> #include <ViennaRNA/fold.h> #include <ViennaRNA/utils.h> int main() { /* The RNA sequence */ char *seq = "GAGUAGUGGAACCAGGCUAUGUUUGUGACUCGCAGACUAACA"; /* allocate memory for MFE structure (length + 1) */ char *structure = (char *)vrna_alloc(sizeof(char) * (strlen(seq) + 1)); /* predict Minmum Free Energy and corresponding secondary structure */ float mfe = vrna_fold(seq, structure); /* print sequence, structure and MFE */ printf("%s\n%s [ %6.2f ]\n", seq, structure, mfe); /* cleanup memory */ free(structure); return 0; }
helloworld_mfe_comparative.c
Instead of using a single sequence as done above, this example predicts a consensus structure for a multiple sequence alignment
#include <stdlib.h> #include <stdio.h> #include <string.h> #include <ViennaRNA/alifold.h> #include <ViennaRNA/aln_util.h> #include <ViennaRNA/utils.h> int main() { /* The RNA sequence alignment */ const char *sequences[] = { "CUGCCUCACAACGUUUGUGCCUCAGUUACCCGUAGAUGUAGUGAGGGU", "CUGCCUCACAACAUUUGUGCCUCAGUUACUCAUAGAUGUAGUGAGGGU", "---CUCGACACCACU---GCCUCGGUUACCCAUCGGUGCAGUGCGGGU", NULL /* indicates end of alignment */ }; /* compute the consensus sequence */ char *cons = consensus(sequences); /* allocate memory for MFE consensus structure (length + 1) */ char *structure = (char *)vrna_alloc(sizeof(char) * (strlen(sequences[0]) + 1)); /* predict Minmum Free Energy and corresponding secondary structure */ float mfe = vrna_alifold(sequences, structure); /* print consensus sequence, structure and MFE */ printf("%s\n%s [ %6.2f ]\n", cons, structure, mfe); /* cleanup memory */ free(cons); free(structure); return 0; }
helloworld_probabilities.c
This example shows how to compute the partition function and base pair probabilities with minimal implementation effort.
#include <stdlib.h> #include <stdio.h> #include <string.h> #include <ViennaRNA/fold.h> #include <ViennaRNA/part_func.h> #include <ViennaRNA/utils.h> int main() { /* The RNA sequence */ char *seq = "GAGUAGUGGAACCAGGCUAUGUUUGUGACUCGCAGACUAACA"; /* allocate memory for pairing propensity string (length + 1) */ char *propensity = (char *)vrna_alloc(sizeof(char) * (strlen(seq) + 1)); /* pointers for storing and navigating through base pair probabilities */ vrna_ep_t *ptr, *pair_probabilities = NULL; float en = vrna_pf_fold(seq, propensity, &pair_probabilities); /* print sequence, pairing propensity string and ensemble free energy */ printf("%s\n%s [ %6.2f ]\n", seq, propensity, en); /* print all base pairs with probability above 50% */ for (ptr = pair_probabilities; ptr->i != 0; ptr++) if (ptr->p > 0.5) printf("p(%d, %d) = %g\n", ptr->i, ptr->j, ptr->p); /* cleanup memory */ free(pair_probabilities); free(propensity); return 0; }
fold_compound_mfe.c
Instead of calling the simple MFE folding interface vrna_fold() , this example shows how to first create a vrna_fold_compound_t container with the RNA sequence to finally compute the MFE using this container. This is especially useful if non-default model settings are applied or the dynamic programming (DP) matrices of the MFE prediction are required for post-processing operations, or other tasks on the same sequence will be performed.
#include <stdlib.h> #include <stdio.h> #include <ViennaRNA/fold_compound.h> #include <ViennaRNA/mfe.h> #include <ViennaRNA/string_utils.h> #include <ViennaRNA/utils.h> int main() { /* initialize random number generator */ vrna_init_rand(); /* Generate a random sequence of 50 nucleotides */ char *seq = vrna_random_string(50, "ACGU"); /* Create a fold compound for the sequence */ vrna_fold_compound_t *fc = vrna_fold_compound(seq, NULL, VRNA_OPTION_DEFAULT); /* allocate memory for MFE structure (length + 1) */ char *structure = (char *)vrna_alloc(sizeof(char) * (strlen(seq) + 1)); /* predict Minmum Free Energy and corresponding secondary structure */ float mfe = vrna_mfe(fc, structure); /* print sequence, structure and MFE */ printf("%s\n%s [ %6.2f ]\n", seq, structure, mfe); /* cleanup memory */ free(seq); free(structure); vrna_fold_compound_free(fc); return 0; }
fold_compound_md.c
In the following, we change the model settings (model details) to a temperature of 25 Degree Celcius, and activate G-Quadruplex precition.
#include <stdlib.h> #include <stdio.h> #include <string.h> #include <ViennaRNA/mfe.h> #include <ViennaRNA/fold_compound.h> #include <ViennaRNA/model.h> #include <ViennaRNA/string_utils.h> #include <ViennaRNA/utils.h> int main() { /* initialize random number generator */ vrna_init_rand(); /* Generate a random sequence of 50 nucleotides */ char *seq = vrna_random_string(50, "ACGU"); /* allocate memory for MFE structure (length + 1) */ char *structure = (char *)vrna_alloc(sizeof(char) * (strlen(seq) + 1)); /* create a new model details structure to store the Model Settings */ vrna_md_t md; /* ALWAYS set default model settings first! */ vrna_md_set_default(&md); /* change temperature and activate G-Quadruplex prediction */ md.temperature = 25.0; /* 25 Deg Celcius */ md.gquad = 1; /* Turn-on G-Quadruples support */ /* create a fold compound */ vrna_fold_compound_t *fc = vrna_fold_compound(seq, &md, VRNA_OPTION_DEFAULT); /* predict Minmum Free Energy and corresponding secondary structure */ float mfe = vrna_mfe(fc, structure); /* print sequence, structure and MFE */ printf("%s\n%s [ %6.2f ]\n", seq, structure, mfe); /* cleanup memory */ free(structure); vrna_fold_compound_free(fc); return 0; }
callback_subopt.c
Here is a basic example how to use the callback mechanism in vrna_subopt_cb() . It simply defines a callback function (see interface definition for vrna_subopt_callback ) that prints the result and increases a counter variable.
#include <stdlib.h> #include <stdio.h> #include <ViennaRNA/fold_compound.h> #include <ViennaRNA/subopt.h> #include <ViennaRNA/string_utils.h> #include <ViennaRNA/utils.h> void subopt_callback(const char *structure, float energy, void *data) { /* simply print the result and increase the counter variable by 1 */ if (structure) printf("%d.\t%s\t%6.2f\n", (*((int *)data))++, structure, energy); } int main() { /* initialize random number generator */ vrna_init_rand(); /* Generate a random sequence of 50 nucleotides */ char *seq = vrna_random_string(50, "ACGU"); /* Create a fold compound for the sequence */ vrna_fold_compound_t *fc = vrna_fold_compound(seq, NULL, VRNA_OPTION_DEFAULT); int counter = 0; /* * call subopt to enumerate all secondary structures in an energy band of * 5 kcal/mol of the MFE and pass it the address of the callback and counter * variable */ vrna_subopt_cb(fc, 500, &subopt_callback, (void *)&counter); /* cleanup memory */ free(seq); vrna_fold_compound_free(fc); return 0; }
soft_constraints_up.c
In this example, a random RNA sequence is generated to predict its MFE under the constraint that a particular nucleotide receives an additional bonus energy if it remains unpaired.
#include <stdlib.h> #include <stdio.h> #include <ViennaRNA/fold_compound.h> #include <ViennaRNA/mfe.h> #include <ViennaRNA/string_utils.h> #include <ViennaRNA/constraints_soft.h> #include <ViennaRNA/utils.h> int main() { /* initialize random number generator */ vrna_init_rand(); /* Generate a random sequence of 50 nucleotides */ char *seq = vrna_random_string(50, "ACGU"); /* Create a fold compound for the sequence */ vrna_fold_compound_t *fc = vrna_fold_compound(seq, NULL, VRNA_OPTION_DEFAULT); /* Add soft constraint of -1.7 kcal/mol to nucleotide 5 whenever it appears in an unpaired context */ vrna_sc_add_up(fc, 5, -1.7, VRNA_OPTION_DEFAULT); /* allocate memory for MFE structure (length + 1) */ char *structure = (char *)vrna_alloc(sizeof(char) * 51); /* predict Minmum Free Energy and corresponding secondary structure */ float mfe = vrna_mfe(fc, structure); /* print seqeunce, structure and MFE */ printf("%s\n%s [ %6.2f ]\n", seq, structure, mfe); /* cleanup memory */ free(seq); free(structure); vrna_fold_compound_free(fc); return 0; }
example1.c
A more extensive example including MFE, Partition Function, and Centroid structure prediction.
#include <stdio.h> #include <stdlib.h> #include <string.h> #include <ViennaRNA/data_structures.h> #include <ViennaRNA/params.h> #include <ViennaRNA/utils.h> #include <ViennaRNA/eval.h> #include <ViennaRNA/fold.h> #include <ViennaRNA/part_func.h> int main(int argc, char *argv[]) { char *seq = "AGACGACAAGGUUGAAUCGCACCCACAGUCUAUGAGUCGGUGACAACAUUACGAAAGGCUGUAAAAUCAAUUAUUCACCACAGGGGGCCCCCGUGUCUAG"; char *mfe_structure = vrna_alloc(sizeof(char) * (strlen(seq) + 1)); char *prob_string = vrna_alloc(sizeof(char) * (strlen(seq) + 1)); /* get a vrna_fold_compound with default settings */ vrna_fold_compound_t *vc = vrna_fold_compound(seq, NULL, VRNA_OPTION_DEFAULT); /* call MFE function */ double mfe = (double)vrna_mfe(vc, mfe_structure); printf("%s\n%s (%6.2f)\n", seq, mfe_structure, mfe); /* rescale parameters for Boltzmann factors */ vrna_exp_params_rescale(vc, &mfe); /* call PF function */ FLT_OR_DBL en = vrna_pf(vc, prob_string); /* print probability string and free energy of ensemble */ printf("%s (%6.2f)\n", prob_string, en); /* compute centroid structure */ double dist; char *cent = vrna_centroid(vc, &dist); /* print centroid structure, its free energy and mean distance to the ensemble */ printf("%s (%6.2f d=%6.2f)\n", cent, vrna_eval_structure(vc, cent), dist); /* free centroid structure */ free(cent); /* free pseudo dot-bracket probability string */ free(prob_string); /* free mfe structure */ free(mfe_structure); /* free memory occupied by vrna_fold_compound */ vrna_fold_compound_free(vc); return EXIT_SUCCESS; }
#include <stdio.h> #include <stdlib.h> #include <math.h> #include <string.h> #include "utils.h" #include "fold_vars.h" #include "fold.h" #include "part_func.h" #include "inverse.h" #include "RNAstruct.h" #include "treedist.h" #include "stringdist.h" #include "profiledist.h" void main() { char *seq1 = "CGCAGGGAUACCCGCG", *seq2 = "GCGCCCAUAGGGACGC", *struct1, *struct2, *xstruc; float e1, e2, tree_dist, string_dist, profile_dist, kT; Tree *T1, *T2; swString *S1, *S2; float *pf1, *pf2; FLT_OR_DBL *bppm; /* fold at 30C instead of the default 37C */ temperature = 30.; /* must be set *before* initializing */ /* allocate memory for structure and fold */ struct1 = (char *)space(sizeof(char) * (strlen(seq1) + 1)); e1 = fold(seq1, struct1); struct2 = (char *)space(sizeof(char) * (strlen(seq2) + 1)); e2 = fold(seq2, struct2); free_arrays(); /* free arrays used in fold() */ /* produce tree and string representations for comparison */ xstruc = expand_Full(struct1); T1 = make_tree(xstruc); S1 = Make_swString(xstruc); free(xstruc); xstruc = expand_Full(struct2); T2 = make_tree(xstruc); S2 = Make_swString(xstruc); free(xstruc); /* calculate tree edit distance and aligned structures with gaps */ edit_backtrack = 1; tree_dist = tree_edit_distance(T1, T2); free_tree(T1); free_tree(T2); unexpand_aligned_F(aligned_line); printf("%s\n%s %3.2f\n", aligned_line[0], aligned_line[1], tree_dist); /* same thing using string edit (alignment) distance */ string_dist = string_edit_distance(S1, S2); free(S1); free(S2); printf("%s mfe=%5.2f\n%s mfe=%5.2f dist=%3.2f\n", aligned_line[0], e1, aligned_line[1], e2, string_dist); /* for longer sequences one should also set a scaling factor for * partition function folding, e.g: */ kT = (temperature + 273.15) * 1.98717 / 1000.; /* kT in kcal/mol */ pf_scale = exp(-e1 / kT / strlen(seq1)); /* calculate partition function and base pair probabilities */ e1 = pf_fold(seq1, struct1); /* get the base pair probability matrix for the previous run of pf_fold() */ bppm = export_bppm(); pf1 = Make_bp_profile_bppm(bppm, strlen(seq1)); e2 = pf_fold(seq2, struct2); /* get the base pair probability matrix for the previous run of pf_fold() */ bppm = export_bppm(); pf2 = Make_bp_profile_bppm(bppm, strlen(seq2)); free_pf_arrays(); /* free space allocated for pf_fold() */ profile_dist = profile_edit_distance(pf1, pf2); printf("%s free energy=%5.2f\n%s free energy=%5.2f dist=%3.2f\n", aligned_line[0], e1, aligned_line[1], e2, profile_dist); free_profile(pf1); free_profile(pf2); }
See also:
examples/helloworld_mfe.c
in the source code tarball
examples/helloworld_mfe_comparative.c
in the source code tarball
examples/helloworld_probabilities.c
in the source code tarball
examples/fold_compound_mfe.c
in the source code tarball
examples/fold_compound_md.c
in the source code tarball
examples/callback_subopt.c
in the source code tarball
examples/soft_constraints_up.c
in the source code tarball
examples/example1.c
in the source code tarball
examples/example_old.c
in the source code tarball