Codebase list ghmm / HEAD ghmm / model.h
HEAD

Tree @HEAD (Download .tar.gz)

model.h @HEADraw · history · blame

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
/*******************************************************************************
*
*       This file is part of the General Hidden Markov Model Library,
*       GHMM version __VERSION__, see http://ghmm.org
*
*       Filename: ghmm/ghmm/model.h
*       Authors:  Benhard Knab, Bernd Wichern, Benjamin Georgi,
*                 Alexander Schliep, Janne Grunau
*
*       Copyright (C) 1998-2004 Alexander Schliep 
*       Copyright (C) 1998-2001 ZAIK/ZPR, Universitaet zu Koeln
*	Copyright (C) 2002-2004 Max-Planck-Institut fuer Molekulare Genetik, 
*                               Berlin
*                                   
*       Contact: schliep@ghmm.org             
*
*       This library is free software; you can redistribute it and/or
*       modify it under the terms of the GNU Library General Public
*       License as published by the Free Software Foundation; either
*       version 2 of the License, or (at your option) any later version.
*
*       This library is distributed in the hope that it will be useful,
*       but WITHOUT ANY WARRANTY; without even the implied warranty of
*       MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
*       Library General Public License for more details.
*
*       You should have received a copy of the GNU Library General Public
*       License along with this library; if not, write to the Free
*       Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
*
*
*       This file is version $Revision: 2277 $ 
*                       from $Date: 2009-04-28 08:44:31 -0400 (Tue, 28 Apr 2009) $
*             last change by $Author: grunau $.
*
*******************************************************************************/
#ifndef GHMM_MODEL_H
#define GHMM_MODEL_H

#ifdef __cplusplus
extern "C" {
#endif

#include <stdio.h>

#include "ghmm.h"

/**@name HMM-Modell */
/*@{ (Doc++-Group: model) */

/** State struct for discrete emissions HMMs.
 *
 * Keeps all parameters that belong to a state. 
 */
typedef struct {
  /** Initial probability */
  double pi;
  /** Output probability */
  double *b;

  /** IDs of the following states */
  int *out_id;
  /** IDs of the previous states */
  int *in_id;

  /** transition probabilities to successor states. */
  double *out_a;
  /** transition probabilities from predecessor states. */
  double *in_a;

  /** Number of successor states */
  int out_states;
  /** Number of precursor states */
  int in_states;
  /** if fix == 1 --> b stays fix during the training */
  int fix;
  /** contains a description of the state (null terminated utf-8)*/
  char *desc;
  /** x coordinate position for graph representation plotting **/
  int xPosition;
  /** y coordinate position for graph representation plotting **/
  int yPosition;
} ghmm_dstate;



/**
    Model struct for discrete emission HMMs.

    Contains all parameters, that define a HMM.
*/
typedef struct {
  /** Number of states */
  int N;
  /** Number of outputs */
  int M;
  /** Vector of the states */
  ghmm_dstate *s;
  /** The a priori probability for the model.
      A value of -1 indicates that no prior is defined. 
      Note: this is not to be confused with priors on emission
      distributions*/
  double prior;

  /* contains a arbitrary name for the model (null terminated utf-8) */
  char *name;

  /** Contains bit flags for varios model extensions such as
      kSilentStates, kTiedEmissions (see ghmm.h for a complete list)
  */
  int model_type;

  /** Flag variables for each state indicating whether it is emitting
      or not. 
      Note: silent != NULL iff (model_type & kSilentStates) == 1  */
  int *silent;
  /*AS*/
  /** Int variable for the maximum level of higher order emissions */
  int maxorder;
  /** saves the history of emissions as int, 
      the nth-last emission is (emission_history * |alphabet|^n+1) % |alphabet|
      see ...*/
  int emission_history;
  
  /** Flag variables for each state indicating whether the states emissions
      are tied to another state. Groups of tied states are represented
      by their tie group leader (the lowest numbered member of the group).
      
      tied_to[s] == kUntied  : s is not a tied state
      
      tied_to[s] == s        : s is a tie group leader
      
      tied_to[t] == s        : t is tied to state s (t>s)
      
      Note: tied_to != NULL iff (model_type & kTiedEmissions) != 0  */
  int *tied_to;

  /** Note: State store order information of the emissions.
      Classic HMMS have emission order 0, that is the emission probability
      is conditioned only on the state emitting the symbol.

      For higher order emissions, the emission are conditioned on the state s
      as well as the previous emission_order[s] observed symbols.

      The emissions are stored in the state's usual double* b. The order is
      set order.

      Note: order != NULL iff (model_type & kHigherOrderEmissions) != 0  */
  int * order;

  /** ghmm_dbackground is a pointer to a
      ghmm_dbackground structure, which holds (essentially) an
      array of background distributions (which are just vectors of floating
      point numbers like state.b).
      
      For each state the array background_id indicates which of the background
      distributions to use in parameter estimation. A value of kNoBackgroundDistribution
      indicates that none should be used.
      
      
      Note: background_id != NULL iff (model_type & kHasBackgroundDistributions) != 0  */
  int *background_id;
  ghmm_dbackground* bp;

  /** (WR) added these variables for topological ordering of silent states 
      Condition: topo_order != NULL iff (model_type & kSilentStates) != 0  */
  int *topo_order;
  int topo_order_length;
  
  /** pow_lookup is a array of precomputed powers
      It contains in the i-th entry M (alphabet size) to the power of i
      The last entry is maxorder+1 */
  int *pow_lookup;

  /** Store for each state a class label. Limits the possibly state sequence

      Note: label != NULL iff (model_type & kLabeledStates) != 0  */
  int* label;
  ghmm_alphabet* label_alphabet;

  ghmm_alphabet* alphabet;
} ghmm_dmodel;

#ifdef __cplusplus
}
#endif
/*
  Important: The inclusion of sequence.h ist not done before this point in
  order to avoid error by compiling.
*/
#include "sequence.h"
#include "scanner.h"
#ifdef __cplusplus
extern "C" {
#endif

/*----------------------------------------------------------------------------*/
/** 
    binary algorithm to compute powers of integers efficiently
    see Knuth, TAOCP, Vol 2, 4.6.3 
    uses if appropiate lookup table from struct ghmm_dmodel */
   int ghmm_ipow(const ghmm_dmodel * mo, int x, unsigned int n);


/*----------------------------------------------------------------------------*/
/** Allocates a ghmm_dmodel with following parameters:
    @return pointer to a ghmm_dmodel or NULL on error
    @param M           number of emissions
    @param N           number of states
    @param modeltype   model_type of the model, used to allocate optional arrays
    @param inDegVec    vector of in degrees of the states
    @param outDegVec   vector of out degrees of the states */
  ghmm_dmodel * ghmm_dmodel_calloc(int M, int N, int modeltype, int * inDegVec,
				 int * outDegVec);

/*----------------------------------------------------------------------------*/
/** Frees the memory of a model.
    @return 0 for succes; -1 for error
    @param mo:  address of pointer to a ghmm_dmodel */
  int ghmm_dmodel_free (ghmm_dmodel ** mo);

/** 
    Produces simple left-right models given sequences. The sequences
    are not read in from file, but exists already as a structur.
    @return vector of models
    @param sq:         sequence
    @param mo_number:  number of models to produce */
  ghmm_dmodel **ghmm_dmodel_from_sequence (ghmm_dseq * sq, long *mo_number);

/**
   Copies a given model. Allocates the necessary memory.
   @return copy of the model
   @param mo:  model to copy */
  ghmm_dmodel *ghmm_dmodel_copy (const ghmm_dmodel * mo);

/**
   Tests if all standardization requirements of model are fulfilled. 
   (That is, if the sum of the probabilities is 1).
   @return 0 for success; otherwise negated number of errors
   @param mo:  ghmm_dmodel to test */
  int ghmm_dmodel_check (const ghmm_dmodel * mo);

/**
   Tests if number of states and number of outputs in the models match.
   @return 0 for succes; -1 for error
   @param mo:           vector of models
   @param model_number: numbr of models */
  int ghmm_dmodel_check_compatibility (ghmm_dmodel ** mo, int model_number);

/**
   Test if to models are compatible. That means their states and outputs match.
   @return 0 for succes; -1 for error
   @param mo:     first ghmm_dmodel
   @param m2:     second ghmm_dmodel */
int ghmm_dmodel_check_compatibel_models (const ghmm_dmodel * mo, const ghmm_dmodel * m2);

/**
   Produces a model, which generates the given sequence with probability 1.
   The model is a strict left-right model with one state for each element 
   in the sequence and the output in state i is the i-th value in the sequence 
   with probability 1. The model also has a final state, a state with no output.
   @return         pointer to the produced ghmm_dmodel 
   @param seq:      sequence
   @param seq_len:  length of the sequence
   @param anz_symb: number of symbols in the sequence
*/
  ghmm_dmodel *ghmm_dmodel_generate_from_sequence (const int *seq, int seq_len,
                                       int anz_symb);

/** 
    Produces sequences to a given model. All memory that is needed for the 
    sequences is allocated inside the function. It is possible to define
    the length of the sequences global (global_len > 0) or it can be set 
    inside the function, when a final state in the model is reach (a state
    with no output). If the model has no final state, the sequences will
    have length MAX_SEQ_LEN.
    @return             pointer to an array of sequences
    @param mo:          model
    @param seed:        initial parameter for the random value generator
                        (an integer). If seed == 0, then the random value
			generator is not initialized.
    @param global_len:  length of sequences (=0: automatically via final states)
    @param seq_number:  number of sequences
    @param Tmax:        maximal sequence length, set to MAX_SEQ_LEN if -1
*/
  ghmm_dseq *ghmm_dmodel_generate_sequences (ghmm_dmodel * mo, int seed, int global_len,
                                        long seq_number, int Tmax);

/**
   Calculates the sum log( P( O | lambda ) ).
   Sequences, that can not be generated from the given ghmm_dmodel, are neglected.
   @return    log(P)
   @param mo model
   @param sq sequences       
*/
  double ghmm_dmodel_likelihood (ghmm_dmodel * mo, ghmm_dseq * sq);

/**
    Get transition probabality from state 'i' to state 'j' to value 'prob'.
    NOTE: No internal checks
    @return  transition probability from state i to state j
    @param mo model
    @param i  state index (source)
    @param j  state index (target)
*/
  double ghmm_dmodel_get_transition(ghmm_dmodel* mo, int i, int j);

/**
    Checks if a non zero transition exists between state 'i' to state 'j'.
    NOTE: No internal checks
    @return   1 if there is a transition, 0 otherwise
    @param mo model
    @param i  state index (source)
    @param j  state index (target)
*/
  int ghmm_dmodel_check_transition(ghmm_dmodel* mo, int i, int j);

/**
    Set transition from state 'i' to state 'j' to value 'prob'.
    NOTE: No internal checks - model might get broken if used without care.
    @param mo model
    @param i  state index (source)
    @param j  state index (target)
    @param prob probabilitys
    
*/
  void ghmm_dmodel_set_transition (ghmm_dmodel * mo, int i, int j, double prob);

/**
   Writes a model in matrix format.
   @param file: output file
   @param mo:   model
*/
  void ghmm_dmodel_print (FILE * file, ghmm_dmodel * mo);

/**
   Writes transition matrix of a model.
   @param file: output file
   @param mo:   model
   @param tab:  format: leading tabs
   @param separator: format: seperator for columns
   @param ending:    format: end of a row  
*/
  void ghmm_dmodel_A_print (FILE * file, ghmm_dmodel * mo, char *tab, char *separator,
                      char *ending);
/**
   Writes output matrix of a model.
   @param file: output file
   @param mo:   model
   @param tab:  format: leading tabs
   @param separator: format: seperator for columns
   @param ending:    format: end of a row  
*/
  void ghmm_dmodel_B_print (FILE * file, ghmm_dmodel * mo, char *tab, char *separator,
                      char *ending);
/**
   Writes initial allocation vector of a matrix.
   @param file: output file
   @param mo:   model
   @param tab:  format: leading Tabs
   @param separator: format: seperator for columns
   @param ending:    format: end of a row  
*/
  void ghmm_dmodel_Pi_print (FILE * file, ghmm_dmodel * mo, char *tab, char *separator,
                       char *ending);
/**
   Writes fix vector of a matrix.
   @param file: output file
   @param mo:   model
   @param tab:  format: leading Tabs
   @param separator: format: seperator for columns
   @param ending:    format: end of a row  
*/
  void ghmm_dmodel_fix_print (FILE * file, ghmm_dmodel * mo, char *tab, char *separator,
                        char *ending);
/**
   Writes transposed transition matrix of a model.
   @param file: output file
   @param mo:   model
   @param tab:  format: leading Tabs
   @param separator: format: seperator for columns
   @param ending:    format: end of a row  
*/
  void ghmm_dmodel_A_print_transp (FILE * file, ghmm_dmodel * mo, char *tab,
                             char *separator, char *ending);
/**
   Writes transposed output matrix of a model.
   @param file: output file
   @param mo:   model
   @param tab:  format: leading Tabs
   @param separator: format: seperator for columns
   @param ending:    format: end of a row  
*/
  void ghmm_dmodel_B_print_transp (FILE * file, ghmm_dmodel * mo, char *tab,
                             char *separator, char *ending);

/**
   Writes transposed initial allocation vector of a matrix.
   @param file: output file
   @param mo:   model
   @param tab:  format: leading Tabs
   @param ending:    format: end of a row  
*/
  void ghmm_dmodel_Pi_print_transp (FILE * file, ghmm_dmodel * mo, char *tab,
                              char *ending);

/** 
    Writes the parameters of a ghmm_dmodel sorted by states. 
    Is not very concise.   
    @param file: output file
    @param mo:   model
*/
  void ghmm_dmodel_states_print (FILE * file, ghmm_dmodel * mo);

/** Computes probabilistic distance of two models
    @return the distance
    @param m0  model used for generating random output
    @param m  model to compare with
    @param maxT  maximum output length (for HMMs with absorbing states multiple
                 sequences with a toal langth of at least maxT will be 
                 generated)
    @param symmetric  flag, whether to symmetrize distance (not implemented yet)
    @param verbose  flag, whether to monitor distance in 40 steps.
                    Prints to stdout (yuk!)
*/
  double ghmm_dmodel_prob_distance (ghmm_dmodel * m0, ghmm_dmodel * m, int maxT, int symmetric,
                              int verbose);

/** 
    Frees all memory from a ghmm_dstate, sets the pointers to NULL and 
    variables to zero.
    @author Peter Pipenbacher
    @param state  state to clean
*/
  void ghmm_dstate_clean (ghmm_dstate *state);


  ghmm_dseq *ghmm_dmodel_label_generate_sequences (ghmm_dmodel * mo, int seed,
                                              int global_len, long seq_number,
                                              int Tmax);

/** 
    Calculates the right index for emission array b of state j in model mo
    given an observation obs and taking the state order into account,
    returns -1 if state order exceeds number of so far emitted characters
    @param MO:  model
    @param  S:  state id 
    @param  O:  integer observation to be updated with
    @param  T:  position of obs in sequence (time)
*/
#define get_emission_index(MO, S, O, T)   (((MO)->model_type & GHMM_kHigherOrderEmissions) ?    \
                                            ((MO)->order[S] > (T)) ? -1 :                       \
                                              (((MO)->emission_history*(MO)->M) %               \
                                              ghmm_ipow((MO), (MO)->M, (MO)->order[S]+1)+(O)) \
                                            : (O))

/**
   Updates emission history of model mo, discarding the oldest and 'adding' the
   new observation by using modulo and multiplication
   left-shift the history, truncate to history length and
   add the last observation
   @param MO:  model to be updated
   @param  O:  integer observation to be updated with
*/
#define update_emission_history(MO, O)   if ((MO)->model_type & GHMM_kHigherOrderEmissions)               \
                                              (MO)->emission_history = ((MO)->emission_history*(MO)->M) % \
                                                ghmm_ipow((MO), (MO)->M, (MO)->maxorder) + (O)


/**
   Updates emission history of model mo for backward algorithm by 'adding'
   observation obs to the left,
   (example: obs = 3
   2 0 0 1 --> 3 2 0 0 )
   removes the most significant position (right-shift) and add the last seen
   observation (left-shifted with the length of history)
   @param MO:  model to be updated
   @param  O:  integer observation to be updated with
*/
#define update_emission_history_front(MO, O)   if ((MO)->model_type & GHMM_kHigherOrderEmissions)       \
                                                 (MO)->emission_history =                               \
                                                   ghmm_ipow((MO), (MO)->M, (MO)->maxorder-1) * (O) + \
                                                   (MO)->emission_history / (MO)->M


/**
    Uses ighmm_cvector_normalize in vector.h
    Normalizes the transition and output probs for each state
    in the given model
    @return 0 if normalization went through
    @param mo: model to be normalized
*/
  int ghmm_dmodel_normalize (ghmm_dmodel * mo);

/**
   Add a specific level of noise to the model parameters
   @return     :        -1 on error, 0 else
   @param mo   :        a pointer to the model
   @param level:        amount of noise to use,
                        a noise level of 0.0 doesn't change the model
   @param seed :        seed for ramdom number generator
*/
  int ghmm_dmodel_add_noise (ghmm_dmodel * mo, double level, int seed);


/** 
   Allocates a new ghmm_dbackground struct and assigs the arguments to
   the respective fields. Note: The arguments need allocation outside of this
   function.
   
   @return     :               0 on success, -1 on error
   @param mo   :               one model
   @param cur  :               a id of a state
   @param times:               number of times the state cur is at least evaluated
*/
  int ghmm_dmodel_duration_apply (ghmm_dmodel * mo, int cur, int times);


/**
   Apply the background distributions to the emission probabilities of states of
   the model which have one specified (background_id[state_id] != kNoBackgroundDistribution).

   @return    :                -1 on error, 0 else
   @param mo  :                a pointer to the ghmm_dmodel
   @param background_weight:   a parameter controlling the weight given to the
                               background. Note, should be between 0 and 1.
*/
  int ghmm_dmodel_background_apply (ghmm_dmodel * mo, double *background_weight);


/** 
   Allocates a new ghmm_dbackground struct and assigs the arguments to
   the respective fields. Note: The arguments need allocation outside of this
   function.
   
   @return                   new pointer to a ghmm_dbackground struct
   @param n                  number of distributions
   @param m                  number of symbols per distribution
   @param orders             orders of the distribtions
   @param B                  matrix of distribution parameters
*/
  ghmm_dbackground *ghmm_dbackground_alloc(int n, int m, int *orders,
                                           double **B);

  ghmm_dbackground *ghmm_dbackground_copy (ghmm_dbackground * bg);

  int ghmm_dbackground_free (ghmm_dbackground * bg);

/**
   Calculates the background distribution for a ghmm_dseq
   the calculated distribution is uniform, every ghmm_dstate with the same order
   has the same background distribution.
   Set more sophisticated background distribution from python.
   Caution: overwrites the pointer to previous defined background distributions!
   @return    : -1 on error, 0 else
   @param mo  :	a pointer to the ghmm_dmodel
   @param sq  : a pointer to a ghmm_dseq struct
   
*/
  int ghmm_dmodel_get_uniform_background (ghmm_dmodel * mo, ghmm_dseq * sq);

/**
   Calculates the squared distance between two compatible models.
   The distance is normalized by the number of parameters of the models.
   @return:      normalized squared distance
   @param mo:    first model
   @param m2:    second model
*/
  double ghmm_dmodel_distance(const ghmm_dmodel * mo, const ghmm_dmodel * m2);

/**
   Calculates a topological ordering of the silent states
   and saves it in the model. Detects cycles of silent states.
   @param mo:      model
*/
  void ghmm_dmodel_order_topological (ghmm_dmodel * mo);

/**
   Update the emissions according to the tie groups by computing the mean
   values within all groups.
   @param mo:      model
*/
  void ghmm_dmodel_update_tie_groups (ghmm_dmodel * mo);

/**
   Reads a xml file and returns an array of dmodel pointer
   @return           :array of dmodels, NULL on error
   @param filename   :filename of the xml file
   @param mo_number  :address of an int to store the length of the returned array
*/
  ghmm_dmodel** ghmm_dmodel_xml_read(const char *filename, int* mo_number);

/**
   Writes an array of dmodel pointer in a xml-file
   @return           :-1 on error, 0 else
   @param file       :filename of the xml file
   @param mo         :an array of pointers to ghmm_dmodel
   @param mo_number  :length of the array
*/
  int ghmm_dmodel_xml_write(ghmm_dmodel** mo, const char *file, int mo_number);

#ifdef __cplusplus
}
#endif
#endif
/*@} (Doc++-Group: model) */