Question CREM-A

[devetak-m]

Hello,

I am trying to implement CReM-A and I was checking your code. I think I am not understanding a small part in the sampling procedure. You implemented the sampling as:

q_ensemble = 1/(beta_i + beta_j) w_link = np.random.exponential(q_ensemble)

therefore you sample from a process with rate q_ensemble and mean 1/q_ensemble = beta_i + beta_j as far as I understand the documentation of numpy.random.exponential (https://numpy.org/doc/stable/reference/random/generated/numpy.random.exponential.html)

But in the paper, Equation III.14 you have that the weights are distributed according to an exponential distribution with rate beta_i + beta_j and mean 1/(beta_i + beta_j). Which makes sense since then the expected value of the row and column sumn match the observed value, equation III.18.

Is this correct?

And if so one could speed up the code by solving for 1/(beta_i + beta_j) directly by transforming it in a linear problem? I mean by solving a linear equation in terms of rates and not means.

[devetak-m]

I just realized that I was looking at the equation and not at the actual implementation. I apologise.