ContextSV 0.0.1
Integrative SV calling.
All Classes Files Functions Variables Typedefs Macros
khmm.h
Go to the documentation of this file.
1#ifndef _KHMM_H
2#define _KHMM_H
3
5#include <stdio.h>
6#include <stdlib.h>
7#include <math.h>
8#include <string>
9#include <vector>
10#include <map>
12
13typedef struct {
14 int N; /* number of states; Q={1,2,...,N} */
15 int M; /* number of observation symbols; V={1,2,...,M}*/
16 double **A; /* A[1..N][1..N]. a[i][j] is the transition prob
17 of going from state i at time t to state j
18 at time t+1 */
19 double **B; /* B[1..N][1..M]. b[j][k] is the probability of
20 of observing symbol k in state j */
21 double *pi; /* pi[1..N] pi[i] is the initial state distribution. */
22 double *B1_mean; /* B1_mean[1..N] mean of a continuous Gaussian distribution for state 1 through N*/
23 double *B1_sd; /*B1_sd standard deviation of B1 values, which is the same for all states*/
24 double B1_uf; /*B1_uniform_fraction: the contribution of uniform distribution to the finite mixture model */
25 double *B2_mean; /* B2_mean[1..4] is the average of B_allele_freq*/
26 double *B2_sd; /* B2_sd[1..4] is the standard deviation of four B_allele_freq, B2_sd[5] is specially for state1, where B is modelled as a wide normal distribution */
27 double B2_uf; /* B2_uniform_fraction: the fraction of uniform distribution in the finite mixture model */
28
29 int NP_flag; /*flag of 1 and 0 to indicate whether Non-Polymorhpic marker information is contained with HMM file*/
30 double *B3_mean; /* B3_mean[1..N] mean of non-polymorphic probe for state 1 through N*/
31 double *B3_sd; /* B3_sd[1..4] is the standard deviation of B3 values*/
32 double B3_uf; /* B3_uniform_fraction: */
33 int dist; /* new parameter to facilitate CNV calling from resequencing data (2009 April) */
34} CHMM;
35
36
37/************************************
38*subroutines for processing continuous HMM models
39************************************/
40
42CHMM ReadCHMM (const char *filename);
43
45void FreeCHMM(CHMM *phmm);
46
48std::vector<int> testVit_CHMM (CHMM hmm, int T, double *O1, double *O2, double *pfb, int *snpdist, double *plogproba);
49
51std::vector<int> ViterbiLogNP_CHMM(CHMM *phmm, int T, double *O1, double *O2, double *pfb, int *snpdist, double **delta, int **psi, double *pprob);
52
54double b1iot (int state, double *mean, double *sd, double uf, double o);
55
57double b2iot (int state, double *mean, double *sd, double uf, double pfb, double b);
58
59#endif // _HMM_H_
double mean(double *data, int n)
Definition kc.cpp:2789
double b1iot(int state, double *mean, double *sd, double uf, double o)
O1 emission probability.
Definition khmm.cpp:66
std::vector< int > ViterbiLogNP_CHMM(CHMM *phmm, int T, double *O1, double *O2, double *pfb, int *snpdist, double **delta, int **psi, double *pprob)
Viterbi algorithm.
Definition khmm.cpp:217
CHMM ReadCHMM(const char *filename)
Read an HMM from a file.
Definition khmm.cpp:371
double b2iot(int state, double *mean, double *sd, double uf, double pfb, double b)
O2 emission probability.
Definition khmm.cpp:86
std::vector< int > testVit_CHMM(CHMM hmm, int T, double *O1, double *O2, double *pfb, int *snpdist, double *plogproba)
Run the main HMM algorithm.
Definition khmm.cpp:23
void FreeCHMM(CHMM *phmm)
Free the memory allocated for an HMM.
Definition khmm.cpp:526
Definition khmm.h:13
double * B1_mean
Definition khmm.h:22
double * B2_sd
Definition khmm.h:26
int M
Definition khmm.h:15
double * B3_sd
Definition khmm.h:31
int dist
Definition khmm.h:33
double * pi
Definition khmm.h:21
double B2_uf
Definition khmm.h:27
int N
Definition khmm.h:14
double B1_uf
Definition khmm.h:24
double B3_uf
Definition khmm.h:32
double ** B
Definition khmm.h:19
int NP_flag
Definition khmm.h:29
double ** A
Definition khmm.h:16
double * B2_mean
Definition khmm.h:25
double * B3_mean
Definition khmm.h:30
double * B1_sd
Definition khmm.h:23