David Huen 编写

我们将讲解Biojava源代码中的例子，目录为demos/dp,文件为Dice.java(Samiul Hasan贡献)

介绍

这个程序实现了R.Durbin,S.Eddy,A.Krogh,G.Mitchison编著的“生物序列分析”一书中的“作弊的骰子”的例子。
这个游戏使用两个骰子，一个是正常的，另一个是做过弊的骰子。正常的骰子六个面朝上的概率一样，做过弊的骰子
点数六朝上的概率是一半，其他面具有相同的概率。这些概率分别代表正常骰子和做过弊的骰子的生成概率。
玩家来回切换这两个骰子投掷。当使用正常骰子时，下一次投掷仍使用正常骰子的概率是0.95。同样，当使用做过弊
的骰子时，下一次仍使用作弊骰子的概率是0.90。这些概率生成了状态的转移分布。（叫做转移概率－－译者注）。
隐马模型由上述再加上魔术状态(MagicalState)建模。魔术状态代表模型的开始和结束的状态。魔术状态转移到正常
骰子的概率是0.8。也就是说魔术状态转移到作弊骰子的概率是0.2。这里还有个终止条件，从正常或作弊的骰子转移
到魔术状态的概率是0.01。
隐马模型展示如下：

源代码：

程序的关键是createCasino()方法。这个方法创建了实现模型的马尔可夫模型实例。

public static MarkovModel createCasino() {
Symbol[] rolls = new Symbol[6];

// 设置骰子字母表
SimpleAlphabet diceAlphabet = new SimpleAlphabet();
diceAlphavet.setName("DiceAlphabet");

for (int i = 1; i<7; i++){
try{
rolls[i-1] = AlphabetManager.createSymbol((char)('0'+i),""+i,Annotation.EMPTY_ANNOTATION);
diceAlphabet.addSymbol(rolls[i-1]);
} catch(Exception e){
throw new NestedError(
e,"Can't create symbols to represent dice rolls");
}
}
将骰子的点数作为标志，字母表就是1－6点，由字母表管理器管理。
下一步，代表正常骰子生成概率的分布和作弊骰子生成概率的分布被创建(分别叫fairD,loadedD)。骰子状态分别被创建SimpleEmissionStates,fairS,loadedS。

你将看到一个以单个值1开始的整数数组.在单向隐马模型中(就像我们创建的),只有一条生成序列.我们沿着这个唯一的序列上的每前进一个位置就使用一次模型的转移概率.在多向隐马模型中,拥有多条生成序列,不同的序列延长可能不同.

例如蛋白质序列每增加一个残基,相对应的DNA序列就增加三个碱基.

int [] advance = { 1 };
Distribution fairD;
Distribution loadedD;
try {
fairD = DistributionFactory.DEFAULT.createDistribution(diceAlphabet);
loadedD = DistributionFactory.DEFAULT.createDistribution(diceAlphabet);
} catch (Exception e) {
throw new NestedError(e, "Can't create distributions");
}
EmissionState fairS = new SimpleEmissionState("fair", Annotation.EMPTY_ANNOTATION, advance, fairD);
EmissionState loadedS = new SimpleEmissionState("loaded", Annotation.EMPTY_ANNOTATION, advance, loade

// 使用这些状态,隐马模型被创建
SimpleMarkovModel casino = new SimpleMarkovModel(1, diceAlphabet, "Casino");
try {
casino.addState(fairS);
casino.addState(loadedS);
} catch (Exception e) {
throw new NestedError(e, "Can't add states to model");
}

下一步我们需要对状态间的转移建模:

try {
casino.createTransition(casino.magicalState(),fairS);
casino.createTransition(casino.magicalState(),loadedS);
casino.createTransition(fairS,casino.magicalState());
casino.createTransition(loadedS,casino.magicalState());
casino.createTransition(fairS,loadedS);
casino.createTransition(loadedS,fairS);
casino.createTransition(fairS,fairS);
casino.createTransition(loadedS,loadedS);
} catch (Exception e) {
throw new NestedError(e, "Can't create transitions");
}

注意casino.magicalState()方法返回MagicalState对象.这个由SimpleMarkovModel类负责,用户不用自己创建.
接下来我们以前创建的fairD,loadedD的生成分布需要初始化.
try {
for(int i=0;i<rolls.length;i++) {
fairD.setWeight(rolls[i],1.0/6.0);
loadedD.setWeight(rolls[i], 0.1);
}
loadedD.setWeight(rolls[5],0.5);
} catch (Exception e) {
throw new NestedError(e, "Can't set emission probabilities");
}

我们同时需要初始化转移分布.注意这是如何实现的: 每个状态的转移分布可以对模型调用getWeights()方法获得.然后将需要的值通过调用分布的getWeight()方法更新.

从这以后就不必调用setWeights()方法将状态分布传递回模型.这看起来有些奇怪.但采用这种方法是因为模型对象可能只使用唯一的分布类,不能被通用分布类替换,从而提高内部效率.每个状态必须正确初始化自己的转移概率.



//设置转移分值.
try {
Distribution dist;

dist = casino.getWeights(casino.magicalState());
dist.setWeight(fairS, 0.8);
dist.setWeight(loadedS, 0.2);

dist = casino.getWeights(fairS);
dist.setWeight(loadedS, 0.04);
dist.setWeight(fairS, 0.95);
dist.setWeight(casino.magicalState(), 0.01);

dist = casino.getWeights(loadedS);
dist.setWeight(fairS, 0.09);
dist.setWeight(loadedS, 0.90);
dist.setWeight(casino.magicalState(), 0.01);
} catch (Exception e) {
throw new NestedError(e, "Can't set transition probabilities");
}

完成创建隐马模型后,所需要做的就是将这个模型返回给调用者.

return casino;

使用隐马模型

创建了隐马模型后,我们创建相应的动态规划对象.

DP dp = DPFactory.DEFAULT.createDP(casino);

现在, 至少我们有些东西可用了.用这个模型生成一条投掷骰子的序列:

StatePath obs_rolls = dp.generate(300);

SymbolList roll_sequence = obs_rolls.symbolListForLabel(StatePath.SEQUENCE);

generate()方法通过这个模型生成了300个标志的一个路径,在第二行我们将这个路径转化成
标志链.

我们先简单介绍一下StatePath对象.这个对象包含一条DP对象生成的比对(Alignment).在一个多向
对象中,将生成多条比对序列.在我们这个例子中,只生成一条让状态路径SEQUENCE读取的序列.这条
序列是用掷骰子游戏投掷出来的.

下面,我们用DP对象检测一种DP算法.我们想处理我们刚生成的投掷序列标志链,使用Viterbi方法预测
每次投掷使用的是哪一个骰子.

为了做到这点,我们创建包含投掷序列(roll_sequence)的多条标志链数组.我们使用单向隐马模型.
应用Viterbi算法,使用模型概率(ScoreType.PROBABILITY)来计算(你也可以使用空模型或者基于对数
几率的概率).这将生成模型支持的状态路径,这个路径就是通过这个模型预测的每次投掷使用哪个骰子
的结果.

SymbolList[] res_array = {roll_sequence};
StatePath v = dp.viterbi(res_array, ScoreType.PROBABILITY);

接下来所做的就是打印出生成序列,指出序列实际的状态路径(使用的是哪个骰子).还有隐马模型对状态
路径的估计:v (隐马模型预测的每次投掷使用的骰子).

//打印 obs_sequence, output, state symbols.
for(int i = 1; i <= obs_rolls.length()/60; i++) {
for(int j=i*60; j

头两个将分别用StatePath.SEQUENCE和StatePath.STATES来访问状态路径obs_rolls产生的生成序列和实际状态路径打印出来.第三个访问v状态序列打印出预测的状态路径.

输出的结果象这样:

544552213525245666363632432522253566166546666666533666543261
fffffffffffflllllllllllfffffffffffflllllllllllllllllllffffff
ffffffffffffffffffffffffffffffffffllllllllllllllllllllffffff

363546253252546524422555242223224344432423341365415551632161
ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff
ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff

144212242323456563652263346116214136666156616666566421456123
fffffflllfffffffffffffffffffffffflfllllllllllllllllfffffffff
fffffffffffffffffffffffffffffffffffllllllllllllllllfffffffff

346313546514332164351242356166641344615135266642261112465663
ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff
ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff

最上面的一行是隐马模型生成的300个投掷结果(数字代表点数).第二行是投掷时使用的状态(f代表正常的骰子,l代表作弊的骰子,我们用它们来创建SimpleEmissionState对象,这个对象代表骰子).

最后一行和上一行很相似,但是是从Viterbi算法的结果--状态路径v生成的.这个结果很符合我们使用的场景,但是也有可能变化很大.

--BACK TO TOP

Maintainted by Wu Xin, CBI, Peking University, China, 2003