options = ['A', 'B'];
n = 10000;
a = 0;
b = 0;
countA = zeros(n, 1);
countB = zeros(n, 1);
for j = 1:n
newChoice = randsample(options, 1,true, [0.987, 0.013]);
if newChoice == 'A'
a = a + 1;
else
b = b + 1;
end
countB(j) = b;
end
Maybe this is much easier:
n = 10000;
countB = rand(n, 1) > 0.987;
countB = cumsum(countB);