Matlab中的离散概率分布计算

Question

我给出了 P(x1...n) 个离散的独立概率值，它们代表例如发生 X 的可能性。

我想要一个通用代码来解决这个问题：X 发生的概率是多少次同时发生 0-n 次？

例如： 给定：每辆车（A，B，C）停放的 3 个概率 P（A），P（B），P（C）。问题是：无车、一车、两车和三车停放的概率有多大？

例如，两辆车同时停放的答案是：

P(A,B,~C) = P(A)*P(B)*(1-P(C))
P(A,~B,C) = P(A)*(1-P(B))*P(C)
P(~A,B,C) = (1-P(A))*P(B)*P(C)
P(2 of 3) = P(A,B,~C) + P(A,~B,C) + P(~A,B,C)

我已经为所有可能性编写了代码，但是我得到的值越多，由于可能的组合越多，它当然会越慢。

% probability: Vector with probabilities P1, P2, ... PN
% result: Vector with results as stated above.

% All possibilities:
result(1) = prod(probability);


shift_vector = zeros(anzahl_werte,1);
for i = 1:anzahl_werte
    % Shift Vector allocallization
    shift_vector(i) = 1;
    % Compute all unique permutations of the shift_vector
    mult_vectors = uperm(shift_vector);

    % Init Result Vector
    prob_vector = zeros(length(mult_vectors(:,1)), 1);

    % Calc Single Probabilities
    for k = 1:length(mult_vectors(:,1))
        prob_vector(k) = prod(abs(mult_vectors(k,:)'-probability));
    end

    % Sum of this Vector for one probability.
    result(i+1) = sum(prob_vector);
end



end


%%%%% Calculate Permutations
function p = uperm(a)
[u, ~, J] = unique(a);
p = u(up(J, length(a)));
end % uperm

function p = up(J, n)
ktab = histc(J,1:max(J));
l = n;
p = zeros(1, n);
s = 1;
for i=1:length(ktab)
    k = ktab(i);
    c = nchoosek(1:l, k);
    m = size(c,1);
    [t, ~] = find(~p.');
    t = reshape(t, [], s);
    c = t(c,:)';
    s = s*m;
    r = repmat((1:s)',[1 k]);
    q = accumarray([r(:) c(:)], i, [s n]);
    p = repmat(p, [m 1]) + q;
    l = l - k;
end
end
%%%%% Calculate Permutations End

有人知道加快这个功能的方法吗？或者也许 Matlab 有一个实现的功能？

---

我找到了计算的名称： 泊松二项分布

Answer 1

这个怎么样？

probability = [.3 .2 .4 .7];

n = numel(probability);
combs = dec2bin(0:2^n-1).'-'0'; %'// each column is a combination of n values,
%// where each value is either 0 or 1. A 1 value will represent an event
%// that happens; a 0 value will represent an event that doesn't happen.
result = NaN(1,n+1); %// preallocate
for k = 0:n; %// number of events that happen
    ind = sum(combs,1)==k; %// combinations with exactly k 1's
    result(k+1) = sum(prod(...
        bsxfun(@times, probability(:), combs(:,ind)) + ... %// events that happen
        bsxfun(@times, 1-probability(:), ~combs(:,ind)) )); %// don't happen
end