RNAlib-2.2.0-RC2
exterior_loops.h
Go to the documentation of this file.
1 #ifndef VIENNA_RNA_PACKAGE_EXTERIOR_LOOPS_H
2 #define VIENNA_RNA_PACKAGE_EXTERIOR_LOOPS_H
3 
4 #include <stdio.h>
5 #include <stdlib.h>
6 #include <math.h>
7 #include <ctype.h>
8 #include <string.h>
9 #include <ViennaRNA/fold_vars.h>
10 #include <ViennaRNA/energy_par.h>
11 #include <ViennaRNA/params.h>
12 #include <ViennaRNA/constraints.h>
13 #include <ViennaRNA/gquad.h>
14 
15 #ifdef __GNUC__
16 # define INLINE inline
17 #else
18 # define INLINE
19 #endif
20 
48 INLINE PRIVATE int E_ExtLoop(int type,
49  int si1,
50  int sj1,
51  vrna_param_t *P);
52 
58 INLINE PRIVATE double exp_E_ExtLoop( int type,
59  int si1,
60  int sj1,
61  vrna_exp_param_t *P);
62 
108 INLINE PRIVATE int E_Stem( int type,
109  int si1,
110  int sj1,
111  int extLoop,
112  vrna_param_t *P);
113 
122 INLINE PRIVATE double exp_E_Stem(int type,
123  int si1,
124  int sj1,
125  int extLoop,
126  vrna_exp_param_t *P);
127 
128 
129 /*
130 #################################
131 # BEGIN OF FUNCTION DEFINITIONS #
132 #################################
133 */
134 
135 INLINE PRIVATE int
136 E_ext_loop( int i,
137  int j,
138  vrna_fold_compound *vc){
139 
140  int ij, en, e, type;
141 
142  int cp = vc->cutpoint;
143  short *S = vc->sequence_encoding;
144  int *idx = vc->jindx;
145  char *ptype = vc->ptype;
146  vrna_param_t *P = vc->params;
147  vrna_md_t *md = &(P->model_details);
148  char *hard_constraints = vc->hc->matrix;
149 
150  e = INF;
151  ij = idx[j] + i;
152  type = ptype[ij];
153 
154  if((cp < 0) || (((i)>=cp)||((j)<cp))){ /* regular exterior loop */
155  switch(md->dangles){
156  case 0: if(hard_constraints[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP)
157  e = E_ExtLoop(type, -1, -1, P);
158  break;
159  case 2: if(hard_constraints[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP)
160  e = E_ExtLoop(type, S[i-1], S[j+1], P);
161  break;
162  default: if(hard_constraints[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP)
163  e = E_ExtLoop(type, -1, -1, P);
164  ij = idx[j - 1] + i;
165  if(hard_constraints[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
166  type = vc->ptype[ij];
167  en = E_ExtLoop(type, -1, S[j], P);
168  e = MIN2(e, en);
169  }
170  ij = idx[j] + i + 1;
171  if(hard_constraints[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
172  type = vc->ptype[ij];
173  en = E_ExtLoop(type, S[i], -1, P);
174  e = MIN2(e, en);
175  }
176  break;
177  }
178  }
179 
180  return e;
181 }
182 
183 
184 INLINE PRIVATE void
185 E_ext_loop_5( vrna_fold_compound *vc){
186 
187  int en, i, j, ij, type;
188  int length = (int)vc->length;
189  char *ptype = vc->ptype;
190  short *S = vc->sequence_encoding;
191  int *indx = vc->jindx;
192  char *hc = vc->hc->matrix;
193  int *hc_up = vc->hc->up_ext;
194  vrna_sc_t *sc = vc->sc;
195  int *f5 = vc->matrices->f5;
196  int *c = vc->matrices->c;
197  vrna_param_t *P = vc->params;
198  int dangle_model = P->model_details.dangles;
199  int *ggg = vc->matrices->ggg;
200  int with_gquad = P->model_details.gquad;
201  int turn = P->model_details.min_loop_size;
202 
203  f5[0] = 0;
204  for(i = 1; i <= turn + 1; i++){
205  if(hc_up[i]){
206  f5[i] = f5[i-1];
207  if(sc && (f5[i] != INF))
208  if(sc->free_energies)
209  f5[i] += sc->free_energies[i][1];
210  } else {
211  f5[i] = INF;
212  }
213  }
214 
215  /* duplicated code may be faster than conditions inside loop ;) */
216  switch(dangle_model){
217  /* dont use dangling end and mismatch contributions at all */
218  case 0: for(j=turn+2; j<=length; j++){
219  f5[j] = INF;
220  if(hc_up[j]){
221  f5[j] = f5[j-1];
222  if(sc && (f5[j] != INF))
223  if(sc->free_energies)
224  f5[j] += sc->free_energies[j][1];
225  }
226  for (i=j-turn-1; i>1; i--){
227  ij = indx[j]+i;
228  if(!(hc[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP))
229  continue;
230  if(f5[i-1] == INF)
231  continue;
232 
233  if(with_gquad){
234  f5[j] = MIN2(f5[j], f5[i-1] + ggg[indx[j]+i]);
235  }
236 
237  if(c[ij] != INF){
238  en = f5[i-1] + c[ij] + E_ExtLoop(ptype[ij], -1, -1, P);
239  f5[j] = MIN2(f5[j], en);
240  }
241  }
242 
243  ij = indx[j] + 1;
244  if(!(hc[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP))
245  continue;
246 
247  if(with_gquad){
248  f5[j] = MIN2(f5[j], ggg[indx[j]+1]);
249  }
250 
251  if(c[ij] != INF){
252  en = c[ij] + E_ExtLoop(ptype[ij], -1, -1, P);
253  f5[j] = MIN2(f5[j], en);
254  }
255  }
256  break;
257 
258  /* always use dangles on both sides */
259  case 2: for(j=turn+2; j<length; j++){
260  f5[j] = INF;
261  if(hc_up[j]){
262  f5[j] = f5[j-1];
263  if(sc && (f5[j] != INF))
264  if(sc->free_energies)
265  f5[j] += sc->free_energies[j][1];
266  }
267  for (i=j-turn-1; i>1; i--){
268  ij = indx[j] + i;
269  if(!(hc[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP))
270  continue;
271  if(f5[i-1] == INF)
272  continue;
273 
274  if(with_gquad){
275  f5[j] = MIN2(f5[j], f5[i-1] + ggg[indx[j]+i]);
276  }
277 
278  if(c[ij] != INF){
279  en = f5[i-1] + c[ij] + E_ExtLoop(ptype[ij], S[i-1], S[j+1], P);
280  f5[j] = MIN2(f5[j], en);
281  }
282  }
283  ij = indx[j] + 1;
284  if(!(hc[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP))
285  continue;
286 
287  if(with_gquad){
288  f5[j] = MIN2(f5[j], ggg[indx[j]+1]);
289  }
290 
291  if(c[ij] != INF){
292  en = c[ij] + E_ExtLoop(ptype[ij], -1, S[j+1], P);
293  f5[j] = MIN2(f5[j], en);
294  }
295  }
296  if(hc_up[length]){
297  f5[length] = f5[length-1];
298  if(sc && (f5[length] != INF))
299  if(sc->free_energies)
300  f5[length] += sc->free_energies[length][1];
301  }
302  for (i=length-turn-1; i>1; i--){
303  ij = indx[length] + i;
304  if(!(hc[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP))
305  continue;
306  if(f5[i-1] == INF)
307  continue;
308 
309  if(with_gquad){
310  f5[length] = MIN2(f5[length], f5[i-1] + ggg[indx[length]+i]);
311  }
312 
313  if(c[ij] != INF){
314  en = f5[i-1] + c[ij] + E_ExtLoop(ptype[ij], S[i-1], -1, P);
315  f5[length] = MIN2(f5[length], en);
316  }
317  }
318  ij = indx[length] + 1;
319  if(!(hc[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP))
320  break;
321 
322  if(with_gquad){
323  f5[length] = MIN2(f5[length], ggg[indx[length]+1]);
324  }
325 
326  if(c[ij] != INF){
327  en = c[ij] + E_ExtLoop(ptype[ij], -1, -1, P);
328  f5[length] = MIN2(f5[length], en);
329  }
330  break;
331 
332  /* normal dangles, aka dangle_model = 1 || 3 */
333  default: for(j=turn+2; j<=length; j++){
334  f5[j] = INF;
335  if(hc_up[j]){
336  f5[j] = f5[j-1];
337  if(sc && (f5[j] != INF))
338  if(sc->free_energies)
339  f5[j] += sc->free_energies[j][1];
340  }
341  for (i=j-turn-1; i>1; i--){
342  ij = indx[j] + i;
343  if(hc[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
344 
345  if(f5[i-1] != INF){
346  if(with_gquad){
347  f5[j] = MIN2(f5[j], f5[i-1] + ggg[indx[j]+i]);
348  }
349 
350  if(c[ij] != INF){
351  type = ptype[ij];
352  en = f5[i-1] + c[ij] + E_ExtLoop(type, -1, -1, P);
353  f5[j] = MIN2(f5[j], en);
354  }
355  }
356  if(hc_up[i-1]){
357  if((f5[i-2] != INF) && c[ij] != INF){
358  en = f5[i-2] + c[ij] + E_ExtLoop(type, S[i-1], -1, P);
359 
360  if(sc)
361  if(sc->free_energies)
362  en += sc->free_energies[i-1][1];
363 
364  f5[j] = MIN2(f5[j], en);
365  }
366  }
367  }
368  ij = indx[j-1] + i;
369  if(hc[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
370  if(hc_up[j]){
371  if(c[ij] != INF){
372  type = ptype[ij];
373  if(f5[i - 1] != INF){
374  en = f5[i-1] + c[ij] + E_ExtLoop(type, -1, S[j], P);
375 
376  if(sc)
377  if(sc->free_energies)
378  en += sc->free_energies[j][1];
379 
380  f5[j] = MIN2(f5[j], en);
381  }
382 
383  if(hc_up[i-1]){
384  if(f5[i - 2] != INF){
385  en = f5[i-2] + c[ij] + E_ExtLoop(type, S[i-1], S[j], P);
386 
387  if(sc)
388  if(sc->free_energies)
389  en += sc->free_energies[i-1][1] + sc->free_energies[j][1];
390 
391  f5[j] = MIN2(f5[j], en);
392  }
393  }
394  }
395  }
396  }
397  }
398  ij = indx[j] + 1;
399  if(hc[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
400 
401  if(with_gquad){
402  f5[j] = MIN2(f5[j], ggg[indx[j]+1]);
403  }
404 
405  if(c[ij] != INF){
406  type = ptype[ij];
407  en = c[ij] + E_ExtLoop(type, -1, -1, P);
408  f5[j] = MIN2(f5[j], en);
409  }
410  }
411  ij = indx[j-1] + 1;
412  if(hc[ij] & VRNA_CONSTRAINT_CONTEXT_EXT_LOOP){
413  if(hc_up[j]){
414  if(c[ij] != INF){
415  type = ptype[ij];
416  en = c[ij] + E_ExtLoop(type, -1, S[j], P);
417 
418  if(sc)
419  if(sc->free_energies)
420  en += sc->free_energies[j][1];
421 
422  f5[j] = MIN2(f5[j], en);
423  }
424  }
425  }
426  }
427  }
428 }
429 
430 INLINE PRIVATE int
431 E_Stem( int type,
432  int si1,
433  int sj1,
434  int extLoop,
435  vrna_param_t *P){
436 
437  int energy = 0;
438  int d5 = (si1 >= 0) ? P->dangle5[type][si1] : 0;
439  int d3 = (sj1 >= 0) ? P->dangle3[type][sj1] : 0;
440 
441  if(type > 2)
442  energy += P->TerminalAU;
443 
444  if(si1 >= 0 && sj1 >= 0)
445  energy += (extLoop) ? P->mismatchExt[type][si1][sj1] : P->mismatchM[type][si1][sj1];
446  else
447  energy += d5 + d3;
448 
449  if(!extLoop) energy += P->MLintern[type];
450  return energy;
451 }
452 
453 INLINE PRIVATE int
454 E_ExtLoop(int type,
455  int si1,
456  int sj1,
457  vrna_param_t *P){
458 
459  int energy = 0;
460  if(si1 >= 0 && sj1 >= 0){
461  energy += P->mismatchExt[type][si1][sj1];
462  }
463  else if (si1 >= 0){
464  energy += P->dangle5[type][si1];
465  }
466  else if (sj1 >= 0){
467  energy += P->dangle3[type][sj1];
468  }
469 
470  if(type > 2)
471  energy += P->TerminalAU;
472 
473  return energy;
474 }
475 
476 INLINE PRIVATE double
477 exp_E_Stem( int type,
478  int si1,
479  int sj1,
480  int extLoop,
481  vrna_exp_param_t *P){
482 
483  double energy = 1.0;
484  double d5 = (si1 >= 0) ? P->expdangle5[type][si1] : 1.;
485  double d3 = (sj1 >= 0) ? P->expdangle3[type][sj1] : 1.;
486 
487  if(si1 >= 0 && sj1 >= 0)
488  energy = (extLoop) ? P->expmismatchExt[type][si1][sj1] : P->expmismatchM[type][si1][sj1];
489  else
490  energy = d5 * d3;
491 
492  if(type > 2)
493  energy *= P->expTermAU;
494 
495  if(!extLoop) energy *= P->expMLintern[type];
496  return energy;
497 }
498 
499 INLINE PRIVATE double
500 exp_E_ExtLoop(int type,
501  int si1,
502  int sj1,
503  vrna_exp_param_t *P){
504 
505  double energy = 1.0;
506  if(si1 >= 0 && sj1 >= 0){
507  energy = P->expmismatchExt[type][si1][sj1];
508  }
509  else if(si1 >= 0){
510  energy = P->expdangle5[type][si1];
511  }
512  else if(sj1 >= 0){
513  energy = P->expdangle3[type][sj1];
514  }
515 
516  if(type > 2)
517  energy *= P->expTermAU;
518 
519  return energy;
520 }
521 
527 #endif
PRIVATE double exp_E_Stem(int type, int si1, int sj1, int extLoop, vrna_exp_param_t *P)
Definition: exterior_loops.h:477
int * ggg
Energies of g-quadruplexes.
Definition: data_structures.h:429
int * c
Energy array, given that i-j pair.
Definition: data_structures.h:422
PRIVATE double exp_E_ExtLoop(int type, int si1, int sj1, vrna_exp_param_t *P)
Definition: exterior_loops.h:500
struct vrna_hc_t * hc
The hard constraints data structure used for structure prediction.
Definition: data_structures.h:707
int cutpoint
The position of the (cofold) cutpoint within the provided sequence. If there is no cutpoint...
Definition: data_structures.h:703
char * matrix
Upper triangular matrix encoding where a base pair or unpaired nucleotide is allowed.
Definition: constraints.h:379
struct vrna_param_t * params
The precomputed free energy contributions for each type of loop.
Definition: data_structures.h:712
int * jindx
DP matrix accessor.
Definition: data_structures.h:716
int gquad
Include G-quadruplexes in structure prediction.
Definition: model.h:184
PRIVATE int E_Stem(int type, int si1, int sj1, int extLoop, vrna_param_t *P)
Definition: exterior_loops.h:431
int * up_ext
A linear array that holds the number of allowed unpaired nucleotides in an exterior loop...
Definition: constraints.h:382
#define MIN2(A, B)
Get the minimum of two comparable values.
Definition: utils.h:118
vrna_md_t model_details
Model details to be used in the recursions.
Definition: params.h:79
The most basic data structure required by many functions throughout the RNAlib.
Definition: data_structures.h:689
PRIVATE int E_ExtLoop(int type, int si1, int sj1, vrna_param_t *P)
Definition: exterior_loops.h:454
#define VRNA_CONSTRAINT_CONTEXT_EXT_LOOP
Hard constraints flag, base pair in the exterior loop.
Definition: constraints.h:189
#define INF
Definition: energy_const.h:16
vrna_mx_mfe_t * matrices
The MFE DP matrices.
Definition: data_structures.h:709
int min_loop_size
Minimum size of hairpin loops.
Definition: model.h:194
int dangles
Specifies the dangle model used in any energy evaluation (0,1,2 or 3)
Definition: model.h:172
unsigned int length
The length of the sequence (or sequence alignment)
Definition: data_structures.h:702
Various functions related to G-quadruplex computations.
char * ptype
Pair type array.
Definition: data_structures.h:740
The datastructure that contains temperature scaled energy parameters.
Definition: params.h:41
The soft constraints data structure.
Definition: constraints.h:401
The data structure that contains the complete model details used throughout the calculations.
Definition: model.h:169
Here all all declarations of the global variables used throughout RNAlib.
short * sequence_encoding
Numerical encoding of the input sequence.
Definition: data_structures.h:735
struct vrna_sc_t * sc
The soft constraints for usage in structure prediction and evaluation.
Definition: data_structures.h:754
The datastructure that contains temperature scaled Boltzmann weights of the energy parameters...
Definition: params.h:86
int * f5
Energy of 5' end.
Definition: data_structures.h:423