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
|
// $Id: betaOmegaDistribution.cpp 962 2006-11-07 15:13:34Z privmane $
#include "betaOmegaDistribution.h"
#include "gammaUtilities.h"
#include "betaUtilities.h"
#include "errorMsg.h"
#include "logFile.h"
#include <cmath>
betaOmegaDistribution::betaOmegaDistribution()
{
_omega=1;
_betaProb = 0.5;
}
// note that the order of initalization makes a diffrence.
betaOmegaDistribution::betaOmegaDistribution(const betaOmegaDistribution& other) :
_betaDistr(other._betaDistr),
_omega(other._omega),
_betaProb(other._betaProb){
}
betaOmegaDistribution::betaOmegaDistribution(MDOUBLE alpha,MDOUBLE beta,int in_number_of_categories,MDOUBLE betaProb,MDOUBLE omega) :distribution(){
_omega = omega;
_betaProb = betaProb;
_betaDistr.setGlobalRate(1.0);
_betaDistr.setBetaParameters(in_number_of_categories,alpha,beta);
}
betaOmegaDistribution::~betaOmegaDistribution() {}
void betaOmegaDistribution::setBetaOmegaParameters(int in_number_of_categories,MDOUBLE alpha, MDOUBLE beta,MDOUBLE betaProb,MDOUBLE omega){
_omega = omega;
_betaProb = betaProb;
_betaDistr.setBetaParameters(in_number_of_categories, alpha, beta);
}
const MDOUBLE betaOmegaDistribution::ratesProb(const int i) const {
if (i < _betaDistr.categories())
return _betaDistr.ratesProb(i)*_betaProb;
else return (1-_betaProb); //omega prob
}
const MDOUBLE betaOmegaDistribution::rates(const int i) const {
if (i < _betaDistr.categories())
return _betaDistr.rates(i);
else return _omega; //omega
}
const MDOUBLE betaOmegaDistribution::getCumulativeProb(const MDOUBLE x) const
{ return _betaDistr.getCumulativeProb(x);
}
|
