1 #ifndef VIENNA_RNA_PACKAGE_LOOPS_HAIRPIN_H
2 #define VIENNA_RNA_PACKAGE_LOOPS_HAIRPIN_H
3
4 #include <math.h>
5 #include <string.h>
6 #include <ViennaRNA/utils/basic.h>
7 #include <ViennaRNA/datastructures/basic.h>
8 #include <ViennaRNA/fold_compound.h>
9 #include <ViennaRNA/params/basic.h>
10
11 #ifdef VRNA_WARN_DEPRECATED
12 # if defined(DEPRECATED)
13 # undef DEPRECATED
14 # endif
15 # if defined(__clang__)
16 # define DEPRECATED(func, msg) func __attribute__ ((deprecated("", msg)))
17 # elif defined(__GNUC__)
18 # define DEPRECATED(func, msg) func __attribute__ ((deprecated(msg)))
19 # else
20 # define DEPRECATED(func, msg) func
21 # endif
22 #else
23 # define DEPRECATED(func, msg) func
24 #endif
25
26 #ifdef __GNUC__
27 # define INLINE inline
28 #else
29 # define INLINE
30 #endif
31
32 /**
33 *
34 * @file ViennaRNA/loops/hairpin.h
35 * @ingroup eval, eval_loops, eval_loops_hp
36 * @brief Energy evaluation of hairpin loops for MFE and partition function calculations
37 */
38
39 /**
40 *
41 * @addtogroup eval_loops_hp
42 * @{
43 * @brief Functions to evaluate the free energy contributions for hairpin loops
44 */
45
46
47 /**
48 * @name Basic free energy interface
49 * @{
50 */
51
52
53 /**
54 * @brief Evaluate the free energy of a hairpin loop
55 * and consider hard constraints if they apply
56 *
57 * This function evaluates the free energy of a hairpin loop
58 *
59 * In case the base pair is not allowed due to a constraint
60 * conflict, this function returns #INF.
61 *
62 * @note This function is polymorphic! The provided #vrna_fold_compound_t may be of type
63 * #VRNA_FC_TYPE_SINGLE or #VRNA_FC_TYPE_COMPARATIVE
64 *
65 * @param fc The #vrna_fold_compound_t that stores all relevant model settings
66 * @param i The 5' nucleotide of the base pair (3' to evaluate the pair as exterior hairpin loop)
67 * @param j The 3' nucleotide of the base pair (5' to evaluate the pair as exterior hairpin loop)
68 * @returns The free energy of the hairpin loop in 10cal/mol
69 */
70 int
71 vrna_E_hp_loop(vrna_fold_compound_t *fc,
72 int i,
73 int j);
74
75
76 /**
77 * @brief Evaluate the free energy of an exterior hairpin loop
78 * and consider possible hard constraints
79 *
80 * @note This function is polymorphic! The provided #vrna_fold_compound_t may be of type
81 * #VRNA_FC_TYPE_SINGLE or #VRNA_FC_TYPE_COMPARATIVE
82 *
83 */
84 int
85 vrna_E_ext_hp_loop(vrna_fold_compound_t *fc,
86 int i,
87 int j);
88
89
90 /**
91 * @brief Evaluate free energy of an exterior hairpin loop
92 */
93 int
94 vrna_eval_ext_hp_loop(vrna_fold_compound_t *fc,
95 int i,
96 int j);
97
98
99 /**
100 * @brief Evaluate free energy of a hairpin loop
101 *
102 * @note This function is polymorphic! The provided #vrna_fold_compound_t may be of type
103 * #VRNA_FC_TYPE_SINGLE or #VRNA_FC_TYPE_COMPARATIVE
104 *
105 * @param fc The #vrna_fold_compound_t for the particular energy evaluation
106 * @param i 5'-position of the base pair
107 * @param j 3'-position of the base pair
108 * @returns Free energy of the hairpin loop closed by @f$ (i,j) @f$ in deka-kal/mol
109 */
110 int
111 vrna_eval_hp_loop(vrna_fold_compound_t *fc,
112 int i,
113 int j);
114
115
116 /**
117 * @brief Compute the Energy of a hairpin-loop
118 *
119 * To evaluate the free energy of a hairpin-loop, several parameters have to be known.
120 * A general hairpin-loop has this structure:<BR>
121 * <PRE>
122 * a3 a4
123 * a2 a5
124 * a1 a6
125 * X - Y
126 * | |
127 * 5' 3'
128 * </PRE>
129 * where X-Y marks the closing pair [e.g. a <B>(G,C)</B> pair]. The length of this loop is 6 as there are
130 * six unpaired nucleotides (a1-a6) enclosed by (X,Y). The 5' mismatching nucleotide is
131 * a1 while the 3' mismatch is a6. The nucleotide sequence of this loop is "a1.a2.a3.a4.a5.a6" <BR>
132 * @note The parameter sequence should contain the sequence of the loop in capital letters of the nucleic acid
133 * alphabet if the loop size is below 7. This is useful for unusually stable tri-, tetra- and hexa-loops
134 * which are treated differently (based on experimental data) if they are tabulated.
135 * @see scale_parameters()
136 * @see vrna_param_t
137 * @warning Not (really) thread safe! A threadsafe implementation will replace this function in a future release!\n
138 * Energy evaluation may change due to updates in global variable "tetra_loop"
139 *
140 * @param size The size of the loop (number of unpaired nucleotides)
141 * @param type The pair type of the base pair closing the hairpin
142 * @param si1 The 5'-mismatching nucleotide
143 * @param sj1 The 3'-mismatching nucleotide
144 * @param string The sequence of the loop (May be @p NULL, otherwise mst be at least @f$size + 2@f$ long)
145 * @param P The datastructure containing scaled energy parameters
146 * @return The Free energy of the Hairpin-loop in dcal/mol
147 */
148 PRIVATE INLINE int
E_Hairpin(int size,int type,int si1,int sj1,const char * string,vrna_param_t * P)149 E_Hairpin(int size,
150 int type,
151 int si1,
152 int sj1,
153 const char *string,
154 vrna_param_t *P)
155 {
156 int energy;
157
158 if (size <= 30)
159 energy = P->hairpin[size];
160 else
161 energy = P->hairpin[30] + (int)(P->lxc * log((size) / 30.));
162
163 if (size < 3)
164 return energy; /* should only be the case when folding alignments */
165
166 if ((string) && (P->model_details.special_hp)) {
167 if (size == 4) {
168 /* check for tetraloop bonus */
169 char tl[7] = {
170 0
171 }, *ts;
172 memcpy(tl, string, sizeof(char) * 6);
173 tl[6] = '\0';
174 if ((ts = strstr(P->Tetraloops, tl)))
175 return P->Tetraloop_E[(ts - P->Tetraloops) / 7];
176 } else if (size == 6) {
177 char tl[9] = {
178 0
179 }, *ts;
180 memcpy(tl, string, sizeof(char) * 8);
181 tl[8] = '\0';
182 if ((ts = strstr(P->Hexaloops, tl)))
183 return P->Hexaloop_E[(ts - P->Hexaloops) / 9];
184 } else if (size == 3) {
185 char tl[6] = {
186 0
187 }, *ts;
188 memcpy(tl, string, sizeof(char) * 5);
189 tl[5] = '\0';
190 if ((ts = strstr(P->Triloops, tl)))
191 return P->Triloop_E[(ts - P->Triloops) / 6];
192
193 return energy + (type > 2 ? P->TerminalAU : 0);
194 }
195 }
196
197 energy += P->mismatchH[type][si1][sj1];
198
199 return energy;
200 }
201
202
203 /* End basic interface */
204 /**@}*/
205
206
207 /**
208 * @name Boltzmann weight (partition function) interface
209 * @{
210 */
211
212
213 /**
214 * @brief Compute Boltzmann weight @f$e^{-\Delta G/kT} @f$ of a hairpin loop
215 *
216 * multiply by scale[u+2]
217 * @see get_scaled_pf_parameters()
218 * @see vrna_exp_param_t
219 * @see E_Hairpin()
220 * @warning Not (really) thread safe! A threadsafe implementation will replace this function in a future release!\n
221 * Energy evaluation may change due to updates in global variable "tetra_loop"
222 *
223 * @param u The size of the loop (number of unpaired nucleotides)
224 * @param type The pair type of the base pair closing the hairpin
225 * @param si1 The 5'-mismatching nucleotide
226 * @param sj1 The 3'-mismatching nucleotide
227 * @param string The sequence of the loop (May be @p NULL, otherwise mst be at least @f$size + 2@f$ long)
228 * @param P The datastructure containing scaled Boltzmann weights of the energy parameters
229 * @return The Boltzmann weight of the Hairpin-loop
230 */
231 PRIVATE INLINE FLT_OR_DBL
exp_E_Hairpin(int u,int type,short si1,short sj1,const char * string,vrna_exp_param_t * P)232 exp_E_Hairpin(int u,
233 int type,
234 short si1,
235 short sj1,
236 const char *string,
237 vrna_exp_param_t *P)
238 {
239 double q, kT;
240
241 kT = P->kT; /* kT in cal/mol */
242
243 if (u <= 30)
244 q = P->exphairpin[u];
245 else
246 q = P->exphairpin[30] * exp(-(P->lxc * log(u / 30.)) * 10. / kT);
247
248 if (u < 3)
249 return (FLT_OR_DBL)q; /* should only be the case when folding alignments */
250
251 if ((string) && (P->model_details.special_hp)) {
252 if (u == 4) {
253 char tl[7] = {
254 0
255 }, *ts;
256 memcpy(tl, string, sizeof(char) * 6);
257 tl[6] = '\0';
258 if ((ts = strstr(P->Tetraloops, tl))) {
259 if (type != 7)
260 return (FLT_OR_DBL)(P->exptetra[(ts - P->Tetraloops) / 7]);
261 else
262 q *= P->exptetra[(ts - P->Tetraloops) / 7];
263 }
264 } else if (u == 6) {
265 char tl[9] = {
266 0
267 }, *ts;
268 memcpy(tl, string, sizeof(char) * 8);
269 tl[8] = '\0';
270 if ((ts = strstr(P->Hexaloops, tl)))
271 return (FLT_OR_DBL)(P->exphex[(ts - P->Hexaloops) / 9]);
272 } else if (u == 3) {
273 char tl[6] = {
274 0
275 }, *ts;
276 memcpy(tl, string, sizeof(char) * 5);
277 tl[5] = '\0';
278 if ((ts = strstr(P->Triloops, tl)))
279 return (FLT_OR_DBL)(P->exptri[(ts - P->Triloops) / 6]);
280
281 if (type > 2)
282 return (FLT_OR_DBL)(q * P->expTermAU);
283 else
284 return (FLT_OR_DBL)q;
285 }
286 }
287
288 q *= P->expmismatchH[type][si1][sj1];
289
290 return (FLT_OR_DBL)q;
291 }
292
293
294 /**
295 * @brief High-Level function for hairpin loop energy evaluation (partition function variant)
296 *
297 * @see vrna_E_hp_loop() for it's free energy counterpart
298 *
299 * @note This function is polymorphic! The provided #vrna_fold_compound_t may be of type
300 * #VRNA_FC_TYPE_SINGLE or #VRNA_FC_TYPE_COMPARATIVE
301 *
302 */
303 FLT_OR_DBL
304 vrna_exp_E_hp_loop(vrna_fold_compound_t *fc,
305 int i,
306 int j);
307
308
309 /* End partition function interface */
310 /**@}*/
311
312
313 /**
314 * @}
315 */
316
317
318 /**
319 * @addtogroup mfe_backtracking
320 * @{
321 */
322
323 /**
324 * @brief Backtrack a hairpin loop closed by @f$ (i,j) @f$
325 *
326 * @note This function is polymorphic! The provided #vrna_fold_compound_t may be of type
327 * #VRNA_FC_TYPE_SINGLE or #VRNA_FC_TYPE_COMPARATIVE
328 *
329 */
330 int
331 vrna_BT_hp_loop(vrna_fold_compound_t *fc,
332 int i,
333 int j,
334 int en,
335 vrna_bp_stack_t *bp_stack,
336 int *stack_count);
337
338
339 /**
340 * @}
341 */
342
343 #ifndef VRNA_DISABLE_BACKWARD_COMPATIBILITY
344
345 /**
346 * @addtogroup eval_deprecated
347 * @{
348 */
349
350 /**
351 * @}
352 */
353
354 #endif
355
356 #endif
357