不同比例二进制数据的隐马尔可夫模型分割

机器算法验证 隐马尔可夫模型 变化点 时间序列分割
2022-04-11 17:48:17

我需要在相对较大的范围内按它们的比例分割一系列 0 和 1。例如,让我们定义 5 种不同的状态,代表 5 种不同的 1 和 0 比率。

Alphabet: 1 and 0

State      Definition                  emission prob.
state 0: 100% zeroes and 0% ones      0:0.999   1: 0.001
state 1: 75% zeroes and 25% ones      0:0.75     1: 0.25
state 2: 50% zeroes and 50% ones      0: 0.5     1: 0.5
state 3: 25% zeroes and 75% ones      0: 0.25    1: 0.25
state 4: 0% zeroes and 100% ones      0: 0.001   1: 0.999

尝试:使用到目前为止我尝试过的所有转换概率和每个状态的发射,我的模型状态序列的输出只是 state0或 state 4

例子:

数据(二进制):

00000000000001111111111110000000000101010101010101010000001000100010011001000010

无论我如何更改转换概率,我都会得到输出。(在各州):

00000000000004444444444440000000000404040404040404040000004000400040044004000040

我需要的输出(在各州):

00000000000004444444444440000000000222222222222222220000001111111111111111111111

我的印象是我缺少一些基本理论而不是实现问题。例如,我通过在任意定义的窗口中聚合获得1vs的比率来平滑数据0,这样我可以看到 state0和 state之间的中间状态4尽管如此,我不想平滑真实数据,因为我需要证明平滑窗口大小是合理的。

使用 HMM 是解决这个问题的好方法吗?

1个回答

我的回答分为两部分。首先,通过更改输入(初始)转换概率,您可以获得类似于您想要的东西。下面是一些 R 代码,为您的示例演示了这一点:

library(HMM)

States <- c("0","1","2","3","4")
Symbols <- c("0","1")
startProbs <- rep(0.2,5)
emissionProbs <- matrix(c(0.999,0.75,0.5,0.25,0.001,0.001,0.25,0.5,0.75,0.999),5,2)
transProbs <- matrix(0.025,5,5)
diag(transProbs) <- 0.9

hmm <- initHMM(States, Symbols, startProbs, transProbs, emissionProbs)
> print(hmm)
$States
[1] "0" "1" "2" "3" "4"

$Symbols
[1] "0" "1"

$startProbs
  0   1   2   3   4 
0.2 0.2 0.2 0.2 0.2 

$transProbs
    to
from     0     1     2     3     4
   0 0.900 0.025 0.025 0.025 0.025
   1 0.025 0.900 0.025 0.025 0.025
   2 0.025 0.025 0.900 0.025 0.025
   3 0.025 0.025 0.025 0.900 0.025
   4 0.025 0.025 0.025 0.025 0.900

$emissionProbs
      symbols
states     0     1
     0 0.999 0.001
     1 0.750 0.250
     2 0.500 0.500
     3 0.250 0.750
     4 0.001 0.999

使用这个初始转移矩阵,我们得到以下观测值 8、20、30 和 40 的概率,它们位于 0、1、0 和 0、1、0、1 序列的中间(大致)...分别:

obs <- as.character(c(0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,
         0,0,0,0,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,0,0,0,0,0,1,0,
         0,0,1,0,0,0,1,0,0,1,1,0,0,1,0,0,0,0,1,0))

post <- posterior(hmm, obs)
> post[,c(8,20,30,40)]
      index
states           8          20          30         40
     0 0.934764162 0.000001395 0.725475508 0.00004174
     1 0.059970011 0.000724501 0.244379742 0.31189082
     2 0.004632750 0.006383026 0.028836815 0.56445433
     3 0.000631774 0.082112681 0.001305885 0.11840354
     4 0.000001303 0.910778397 0.000002049 0.00520957

如您所见,最大值。如您所愿,概率状态分别为 0、4、0 和 2。

如果您不为状态 0 和 4 选择如此极端的概率,它也可能对您有所帮助,也许选择 0.95 / 0.05 而不是 0.999 / 0.001。这将更容易从给定状态中获得更高的转换概率,而不会一直处于状态 0 和 4 中。

如果您正在考虑 HMM 的替代方案,您可能会考虑连续状态空间模型,该模型可以表述为广义加法模型。使用mgcvR 中的包,可以如下设置:

library(mgcv)

obs <- c(0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,
                      0,0,0,0,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,0,0,0,0,0,1,0,
                      0,0,1,0,0,0,1,0,0,1,1,0,0,1,0,0,0,0,1,0)
time <- seq(1,length(obs))

foo <- gam(obs~s(time),family="binomial")

> predict(foo,type="response")[c(8,20,30,40)]
                   8                   20                   30 
0.000000000000000222 0.999999999999999778 0.000277113887323986 
                  40 
0.540166858432701846 

如您所见,概率与您想要的非常吻合。显然,可能需要对平滑项中的参数进行一些调整。