matlab - 估计马尔可夫转移矩阵的置信区间

Question

我有一系列包含字母 ACGTE 的不同长度的n=400序列。例如，在A之后出现C的概率为：

并且可以从一组经验序列中计算出来，因此

假设：

然后我得到一个转换矩阵：

但我对计算 Phat 的置信区间很感兴趣，有什么想法可以解决这个问题吗？

Answer 1

您可以使用自举来估计置信区间。MATLABbootci在统计工具箱中提供了函数。这是一个例子：

%# generate a random cell array of 400 sequences of varying length
%# each containing indices from 1 to 5 corresponding to ACGTE
sequences = arrayfun(@(~) randi([1 5], [1 randi([500 1000])]), 1:400, ...
    'UniformOutput',false)';

%# compute transition matrix from all sequences
trans = countFcn(sequences);

%# number of bootstrap samples to draw
Nboot = 1000;

%# estimate 95% confidence interval using bootstrapping
ci = bootci(Nboot, {@countFcn, sequences}, 'alpha',0.05);
ci = permute(ci, [2 3 1]);

我们得到：

>> trans         %# 5x5 transition matrix: P_hat
trans =
      0.19747       0.2019      0.19849       0.2049      0.19724
      0.20068      0.19959      0.19811      0.20233      0.19928
      0.19841      0.19798       0.2021       0.2012      0.20031
      0.20077      0.19926      0.20084      0.19988      0.19926
      0.19895      0.19915      0.19963      0.20139      0.20088

以及另外两个包含置信区间下限和上限的类似矩阵：

>> ci(:,:,1)     %# CI lower bound
>> ci(:,:,2)     %# CI upper bound

我正在使用以下函数从一组序列中计算转换矩阵：

function trans = countFcn(seqs)
    %# accumulate transition matrix from all sequences
    trans = zeros(5,5);
    for i=1:numel(seqs)
        trans = trans + sparse(seqs{i}(1:end-1), seqs{i}(2:end), 1, 5,5);
    end

    %# normalize into proper probabilities
    trans = bsxfun(@rdivide, trans, sum(trans,2));
end

作为奖励，我们可以使用bootstrp函数来获取从每个引导样本计算的统计量，我们用它来显示转换矩阵中每个条目的直方图：

%# compute multiple transition matrices using bootstrapping
stat = bootstrp(Nboot, @countFcn, sequences);

%# display histogram for each entry in the transition matrix
sub = reshape(1:5*5,5,5);
figure
for i=1:size(stat,2)
    subplot(5,5,sub(i))
    hist(stat(:,i))
end

Answer 2

不确定它在统计上是否合理，但一种获得指示性上限和下限的简单方法：

将您的样本切成 n 等份（例如 1:40,41:80,...,361:400）并计算每个子样本的概率矩阵。

通过查看子样本之间的概率分布，您应该对方差是什么有一个很好的了解。

这种方法的缺点是可能无法实际计算出具有所需给定概率的区间。优点是它应该让您对系列的行为有很好的感觉，并且它可能会捕获由于其他方法（例如自举）所基于的假设而在其他方法中可能丢失的一些信息。