RNAlib-2.2.0-RC2
hairpin_loops.h
Go to the documentation of this file.
1 #ifndef VIENNA_RNA_PACKAGE_HAIRPIN_LOOPS_H
2 #define VIENNA_RNA_PACKAGE_HAIRPIN_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 
63 INLINE PRIVATE int E_Hairpin(int size,
64  int type,
65  int si1,
66  int sj1,
67  const char *string,
68  vrna_param_t *P);
69 
88 INLINE PRIVATE double exp_E_Hairpin( int u,
89  int type,
90  short si1,
91  short sj1,
92  const char *string,
93  vrna_exp_param_t *P);
94 
95 
96 INLINE PRIVATE int
98  int i,
99  int j);
100 
101 INLINE PRIVATE int
103  int i,
104  int j);
105 
106 /*
107 #################################
108 # BEGIN OF FUNCTION DEFINITIONS #
109 #################################
110 */
111 
112 INLINE PRIVATE int
113 E_Hairpin(int size,
114  int type,
115  int si1,
116  int sj1,
117  const char *string,
118  vrna_param_t *P){
119 
120  int energy;
121 
122  if(size <= 30)
123  energy = P->hairpin[size];
124  else
125  energy = P->hairpin[30] + (int)(P->lxc*log((size)/30.));
126 
127  if(size < 3) return energy; /* should only be the case when folding alignments */
128 
129  if(P->model_details.special_hp){
130  if(size == 4){ /* check for tetraloop bonus */
131  char tl[7]={0}, *ts;
132  strncpy(tl, string, 6);
133  if ((ts=strstr(P->Tetraloops, tl)))
134  return (P->Tetraloop_E[(ts - P->Tetraloops)/7]);
135  }
136  else if(size == 6){
137  char tl[9]={0}, *ts;
138  strncpy(tl, string, 8);
139  if ((ts=strstr(P->Hexaloops, tl)))
140  return (energy = P->Hexaloop_E[(ts - P->Hexaloops)/9]);
141  }
142  else if(size == 3){
143  char tl[6]={0,0,0,0,0,0}, *ts;
144  strncpy(tl, string, 5);
145  if ((ts=strstr(P->Triloops, tl))) {
146  return (P->Triloop_E[(ts - P->Triloops)/6]);
147  }
148  return (energy + (type>2 ? P->TerminalAU : 0));
149  }
150  }
151  energy += P->mismatchH[type][si1][sj1];
152 
153  return energy;
154 }
155 
160 INLINE PRIVATE int
162  int i,
163  int j){
164 
165  int u, *hc_up;
166  char hc;
167 
168  u = j - i - 1;
169  hc_up = vc->hc->up_hp;
170  hc = vc->hc->matrix[vc->jindx[j] + i];
171 
172  /* is this base pair allowed to close a hairpin (like) loop ? */
174  /* are all nucleotides in the loop allowed to be unpaired ? */
175  if(hc_up[i+1] >= u){
176  return vrna_eval_hp_loop(vc, i, j);
177  }
178  }
179  return INF;
180 }
181 
186 INLINE PRIVATE int
188  int i,
189  int j){
190 
191  int u, *hc_up;
192  char hc;
193 
194  u = vc->length - j + i - 1;
195  hc_up = vc->hc->up_hp;
196  hc = vc->hc->matrix[vc->jindx[j] + i];
197 
198  /* is this base pair allowed to close a hairpin (like) loop ? */
200  /* are all nucleotides in the loop allowed to be unpaired ? */
201  if(hc_up[j+1] >= u){
202  return vrna_eval_ext_hp_loop(vc, i, j);
203  }
204  }
205  return INF;
206 }
207 
214 INLINE PRIVATE int
216  int i,
217  int j){
218 
219  int u, e, type;
220  char loopseq[10];
221  short *S;
222  vrna_param_t *P;
223  vrna_sc_t *sc;
224  vrna_md_t *md;
225 
226  S = vc->sequence_encoding;
227  P = vc->params;
228  sc = vc->sc;
229  md = &(P->model_details);
230  u = vc->length - j + i - 1;
231  type = md->pair[S[j]][S[i]];
232 
233  if(type == 0)
234  type = 7;
235 
236  if (u<7) {
237  strcpy(loopseq , vc->sequence + j - 1);
238  strncat(loopseq, vc->sequence, i);
239  }
240 
241  e = E_Hairpin(u, type, S[j + 1], S[i - 1], loopseq, P);
242 
243  if(sc){
244  if(sc->free_energies)
245  e += sc->free_energies[j + 1][vc->length - j]
246  + sc->free_energies[1][i - 1];
247 
248  if(sc->f)
249  e += sc->f(j, i, i, j, VRNA_DECOMP_PAIR_HP, sc->data);
250  }
251 
252  return e;
253 }
254 
261 INLINE PRIVATE int
263  int i,
264  int j){
265 
266  int u, e, ij, type;
267 
268  int cp = vc->cutpoint;
269  short *S = vc->sequence_encoding;
270  int *idx = vc->jindx;
271  vrna_param_t *P = vc->params;
272  vrna_sc_t *sc = vc->sc;
273  vrna_md_t *md = &(P->model_details);
274 
275 
276  u = j - i - 1;
277  ij = idx[j] + i;
278  type = md->pair[S[i]][S[j]];
279 
280  if(type == 0)
281  type = 7;
282 
283  if((cp < 0) || ((i >= cp) || (j < cp))){ /* regular hairpin loop */
284  e = E_Hairpin(u, type, S[i+1], S[j-1], vc->sequence+i-1, P);
285  } else { /* hairpin-like exterior loop (for cofolding) */
286  short si, sj;
287  si = ((i >= cp) || ((i + 1) < cp)) ? S[i+1] : -1;
288  sj = (((j - 1) >= cp) || (j < cp)) ? S[j-1] : -1;
289  if (md->dangles)
290  e = E_ExtLoop(md->rtype[type], sj, si, P);
291  else
292  e = E_ExtLoop(md->rtype[type], -1, -1, P);
293  }
294 
295  /* add soft constraints */
296  if(sc){
297  if(sc->free_energies)
298  e += sc->free_energies[i+1][u];
299 
300  if(sc->en_basepair){
301  e += sc->en_basepair[ij];
302  }
303  if(sc->f)
304  e += sc->f(i, j, i, j, VRNA_DECOMP_PAIR_HP, sc->data);
305  }
306 
307  return e;
308 }
309 
310 INLINE PRIVATE int
311 E_hp_loop_ali(int i,
312  int j,
313  vrna_fold_compound *vc){
314 
315  int u, e, s, *type;
316  int n_seq = vc->n_seq;
317  int *idx = vc->jindx;
318  vrna_param_t *P = vc->params;
319  vrna_md_t *md = &(P->model_details);
320  short **S = vc->S;
321  short **S5 = vc->S5; /*S5[s][i] holds next base 5' of i in sequence s*/
322  short **S3 = vc->S3; /*Sl[s][i] holds next base 3' of i in sequence s*/
323  char **Ss = vc->Ss;
324  unsigned short **a2s = vc->a2s;
325  char hc = vc->hc->matrix[idx[j]+i];
326  int *hc_up = vc->hc->up_hp;
327  vrna_sc_t **sc = vc->scs;
328 
329  type = (int *)vrna_alloc(sizeof(int) * n_seq);
330 
331  for (s=0; s<n_seq; s++) {
332  type[s] = md->pair[S[s][i]][S[s][j]];
333  if (type[s]==0) type[s]=7;
334  }
335 
336  /* is this base pair allowed to close a hairpin loop ? */
338  if(hc_up[i+1] >= j - i - 1){
339  for (e=s=0; s<n_seq; s++) {
340  u = a2s[s][j-1] - a2s[s][i];
341  if (u < 3) e += 600;
342  else e += E_Hairpin(u, type[s], S3[s][i], S5[s][j], Ss[s]+(a2s[s][i-1]), P);
343  }
344 
345  if(sc)
346  for(s = 0; s < n_seq; s++){
347  if(sc[s]){
348  u = a2s[s][j-1]-a2s[s][i];
349 
350  if(sc[s]->free_energies)
351  e += sc[s]->free_energies[a2s[s][i]+1][u];
352 
353  if(sc[s]->en_basepair)
354  e += sc[s]->en_basepair[idx[j] + i];
355 
356  if(sc[s]->f)
357  e += sc[s]->f(i, j, i, j, VRNA_DECOMP_PAIR_HP, sc[s]->data);
358  }
359  }
360 
361  free(type);
362  return e;
363  }
364  }
365 
366  free(type);
367  return INF;
368 }
369 
370 /*
371 *************************************
372 * Partition function variants below *
373 *************************************
374 */
375 
376 
377 INLINE PRIVATE double
379  int type,
380  short si1,
381  short sj1,
382  const char *string,
383  vrna_exp_param_t *P){
384 
385  double q, kT;
386  kT = P->kT; /* kT in cal/mol */
387 
388  if(u <= 30)
389  q = P->exphairpin[u];
390  else
391  q = P->exphairpin[30] * exp( -(P->lxc*log( u/30.))*10./kT);
392 
393  if(u < 3) return q; /* should only be the case when folding alignments */
394 
395  if(P->model_details.special_hp){
396  if(u==4){
397  char tl[7]={0,0,0,0,0,0,0}, *ts;
398  strncpy(tl, string, 6);
399  if ((ts=strstr(P->Tetraloops, tl))){
400  if(type != 7)
401  return (P->exptetra[(ts-P->Tetraloops)/7]);
402  else
403  q *= P->exptetra[(ts-P->Tetraloops)/7];
404  }
405  }
406  else if(u==6){
407  char tl[9]={0,0,0,0,0,0,0,0,0}, *ts;
408  strncpy(tl, string, 8);
409  if ((ts=strstr(P->Hexaloops, tl)))
410  return (P->exphex[(ts-P->Hexaloops)/9]);
411  }
412  else if(u==3){
413  char tl[6]={0,0,0,0,0,0}, *ts;
414  strncpy(tl, string, 5);
415  if ((ts=strstr(P->Triloops, tl)))
416  return (P->exptri[(ts-P->Triloops)/6]);
417  if (type>2)
418  return q * P->expTermAU;
419  else
420  return q;
421  }
422  }
423  q *= P->expmismatchH[type][si1][sj1];
424 
425  return q;
426 }
427 
433 INLINE PRIVATE double
435  int i,
436  int j){
437 
438  int u, ij, type;
439  char hc;
440  double q;
441 
442  int cp = vc->cutpoint;
443  short *S = vc->sequence_encoding;
444  int *idx = vc->jindx;
445  vrna_exp_param_t *P = vc->exp_params;
446  int *hc_up = vc->hc->up_hp;
447  vrna_sc_t *sc = vc->sc;
448  FLT_OR_DBL *scale = vc->exp_matrices->scale;
449 
450  q = 0.;
451  u = j - i - 1;
452  ij = idx[j] + i;
453  type = vc->ptype[ij];
454  hc = vc->hc->matrix[ij];
455 
456  /* is this base pair allowed to close a hairpin (like) loop ? */
457  if(hc & VRNA_CONSTRAINT_CONTEXT_HP_LOOP){
458  /* are all nucleotides in the loop allowed to be unpaired ? */
459  if(hc_up[i+1] >= u){
460 
461  if((cp < 0) || ((i >= cp) || (j < cp))){ /* regular hairpin loop */
462  q = exp_E_Hairpin(u, type, S[i+1], S[j-1], vc->sequence+i-1, P)
463  * scale[u+2];
464  } else { /* hairpin-like exterior loop (for cofolding) */
465 
466  }
467 
468  /* add soft constraints */
469  if(sc){
470  if(sc->boltzmann_factors)
471  q *= sc->boltzmann_factors[i+1][u];
472 
473  if(sc->exp_en_basepair)
474  q *= sc->exp_en_basepair[ij];
475 
476  if(sc->exp_f)
477  q *= sc->exp_f(i, j, i, j, VRNA_DECOMP_PAIR_HP, sc->data);
478  }
479  }
480  }
481  return q;
482 }
483 
489 #endif
void * vrna_alloc(unsigned size)
Allocate space safely.
int special_hp
Include special hairpin contributions for tri, tetra and hexaloops.
Definition: model.h:178
FLT_OR_DBL ** boltzmann_factors
Boltzmann Factors of the energy contributions for unpaired sequence stretches.
Definition: constraints.h:404
int * up_hp
A linear array that holds the number of allowed unpaired nucleotides in a hairpin loop...
Definition: constraints.h:385
PRIVATE int E_Hairpin(int size, int type, int si1, int sj1, const char *string, vrna_param_t *P)
Compute the Energy of a hairpin-loop.
Definition: hairpin_loops.h:113
struct vrna_hc_t * hc
The hard constraints data structure used for structure prediction.
Definition: data_structures.h:707
#define VRNA_DECOMP_PAIR_HP
Generalized constraint folding flag indicating hairpin loop decomposition step.
Definition: constraints.h:254
int cutpoint
The position of the (cofold) cutpoint within the provided sequence. If there is no cutpoint...
Definition: data_structures.h:703
short ** S5
S5[s][i] holds next base 5' of i in sequence s.
Definition: data_structures.h:788
PRIVATE int vrna_eval_ext_hp_loop(vrna_fold_compound *vc, int i, int j)
Evaluate free energy of an exterior hairpin loop.
Definition: hairpin_loops.h:215
PRIVATE int vrna_eval_hp_loop(vrna_fold_compound *vc, int i, int j)
Evaluate free energy of a hairpin loop.
Definition: hairpin_loops.h:262
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 ** free_energies
Energy contribution for unpaired sequence stretches.
Definition: constraints.h:402
vrna_md_t model_details
Model details to be used in the recursions.
Definition: params.h:79
short ** S3
Sl[s][i] holds next base 3' of i in sequence s.
Definition: data_structures.h:791
unsigned int n_seq
The number of sequences in the alignment.
Definition: data_structures.h:776
vrna_mx_pf_t * exp_matrices
The PF DP matrices.
Definition: data_structures.h:710
FLT_OR_DBL * exp_en_basepair
Boltzmann Factors of the energy contribution for base pairs.
Definition: constraints.h:405
PRIVATE int vrna_E_hp_loop(vrna_fold_compound *vc, int i, int j)
Evaluate the free energy of a hairpin loop and consider possible hard constraints.
Definition: hairpin_loops.h:161
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 INF
Definition: energy_const.h:16
vrna_md_t model_details
Model details to be used in the recursions.
Definition: params.h:136
short ** S
Numerical encoding of the sequences in the alignment.
Definition: data_structures.h:785
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
PRIVATE double vrna_exp_E_hp_loop(vrna_fold_compound *vc, int i, int j)
High-Level function for hairpin loop energy evaluation (partition function variant) ...
Definition: hairpin_loops.h:434
Various functions related to G-quadruplex computations.
char * ptype
Pair type array.
Definition: data_structures.h:740
struct vrna_exp_param_t * exp_params
The precomputed free energy contributions as Boltzmann factors.
Definition: data_structures.h:713
#define VRNA_CONSTRAINT_CONTEXT_HP_LOOP
Hard constraints flag, base pair encloses hairpin loop.
Definition: constraints.h:197
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
void * data
A pointer to the data object provided for for pseudo energy contribution functions of the generalized...
Definition: constraints.h:447
Here all all declarations of the global variables used throughout RNAlib.
PRIVATE double exp_E_Hairpin(int u, int type, short si1, short sj1, const char *string, vrna_exp_param_t *P)
Compute Boltzmann weight of a hairpin loop.
Definition: hairpin_loops.h:378
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
PRIVATE int vrna_E_ext_hp_loop(vrna_fold_compound *vc, int i, int j)
Evaluate the free energy of an exterior hairpin loop and consider possible hard constraints.
Definition: hairpin_loops.h:187
int * en_basepair
Energy contribution for base pairs.
Definition: constraints.h:403
The datastructure that contains temperature scaled Boltzmann weights of the energy parameters...
Definition: params.h:86
char * sequence
The input sequence string.
Definition: data_structures.h:732
FLT_OR_DBL(* exp_f)(int, int, int, int, char, void *)
A function pointer used for pseudo energy contribution boltzmann factors in PF calculations.
Definition: constraints.h:421
int(* f)(int, int, int, int, char, void *)
A function pointer used for pseudo energy contribution in MFE calculations.
Definition: constraints.h:411
struct vrna_sc_t ** scs
A set of soft constraints (for each sequence in the alignment)
Definition: data_structures.h:799