Yes, but you can get other quantities with this approach as well (though not in nearly as nice a form). To do this, consider the random matrix

t = Bernoulli(T) (i.e., T = E[t])

The distribution is

m_n = M_1 t_1 ... t_n,

which has expected value

M_n = E[m_n] = E[M_1 t_1 ... t_n] = M_1 T^n

Here, each of the t_i is an i.i.d. sample of t. You can similarly compute the variance:

Var[m_n] = Var[M_1 t_1 ... t_n]

I don't know a closed form for this off the top of my head, but at this point simulations are probably the best way to go.

edit: I completely thought I was answering a different topic. Please disregard the quoted, but I will leave it for completeness.

> I'm not quite sure what you mean by this question.

> What do you mean by "not sampling at each step"? The idea isn't to characterize the entirety of the sampling function along its support, but to draw a sequence of random variables from it.

Yes but that is not what you do. You calculate the expected proportion of each atom at each step, but:

1) this does not tell you what the distribution is. Decay is a random process... what is the probability that you have more than 90% of a left at step 1000?

2) this only works because you're taking a limit over an infinite number of atoms, which takes the variance to 0. If you have a discrete system, things become hairier.

Or take the following system for instance:

There are 100 cells, each of which has the same size. At every time step, a little bit of food falls on the cells and is eaten by a single cell. A cell gets the food with a probability proportional to its weight, it then grows by one unit.

If you use expectations, everyone grows at the same rate. If you simulate, you'll see that the results follow a Pareto distribution.

Whoops, I thought I was answering a question about my rejection sampling article. My mistake. I'll address this though.

1) No it doesn't. I don't make any distributional claims for these calculations. You're right. I haven't addressed questions like P(>90%|t=1000).

2) Actually no, it doesn't take the variance to 0. For any given t, the Markov approximation is equal to the expected value of the given atom at time t. There is a variance for every atom count at every time point. I didn't address calculating variance at all. We're speaking in expected values here, the interpretation of 54.8 atoms at time t is the same as the average american household having 2.5 children.

Something I addressed in a comment when I posted this on reddit was that the Markov chain approach is the theoretical mean. If you were to actually observe a decay process like this in the real world, it would essentially be one possible realization of the simulated example.

As far as the Pareto distribution is concerned, I didn't know that actually. That's really interesting.

Ah got it, you simulate individually per atom, I thought you were describing the entire population. Indeed, you can always use that vector as the parameter to a multinomial distribution and recover the real distribution.

Regarding the Pareto distribution, random growth model are probably the reason why power laws are everywhere.

The cell example is very cool. It's not a Markov chain though. The probability that a cell grows from size 2 to size 3 (in probabilistic notation, P(2->3|2) ), for instance, would be 2/(sum of all cell sizes), and the sum of all cell sizes is growing by 1 each period (assuming the unit of food given per period is 1, etc. etc.)

It's a markov chain if your state is the vector representing the size of each cell