Hi to everyone.
I should compute the spectral radius of a Next Generation Matrix, which is non-symmetric.
Do you know the Rstan function able to compute the eigenvalues of a non-symmetric matrix?
Thank you !
No. Depends on your model a possible alternative is to write out the close form of the reproduction number.
1 Like
Oh. Do you know why Rstan doesn’t have a function, able to compute the eigenvalues of a non-symmetric matrix?
There’s no such a function out of box. See Stan’s matrix ops here. In general finding non-sym matrix eigenvalues is more difficult. What’s the size of your NGM?
The size is 5x5. However we can split the NGM in the product of a matrix C that doesn’t depend from the parameters, and a diagonal matrix D. Since NGM=CD=DC we have that (spectral radius = sr ) sr(NGM)=sr(CD)<=sr(C)*sr(D)=sr(C)*max(D_ii). Using the same argument we have also sr(C)*min(D_ii)<=sr(NGM).
At the end we have an estimate: lower=sr(C)*min(D_ii) <= sr(NGM) <= sr(C)*max(D_ii)=upper.
Since C doesn’t depend from the parameters, we can compute the spectral radius in R and we can give to Rstan this value as a data. Rstan is able to compute the maximum and minimum of the components of a vector. So we can compute the lower and the upper values.
Since upper-lower is very little we can take for instance R0=(upper-lower)/2.
Does it make sense? What is your opinion?
It sounds fine except I’m missing the “upper-lower is very little” part. On the other hand, in general I’d probably take one of the three alternatives, assuming not much special nature of the matrix can be exploited:
- Anoldi method
- brute-force Gelfand formula
- write out the eigenvalue explicitly as it’s not too hard for a 5x5 matrix.
We have that R0 belongs to the interval [lower, upper], where lower = sr(C)*min(D_ii) and upper=sr(C)max(D_ii). If upper-lower= sr(C)(max(D_ii)-min(D_ii)) is little we have a distribution very concentrated, it means that if wi e randomly extract a value from the interval [lower, upper], we have a very good approximation of R0.
Thank you for the suggestions, I didn’t know Anoldi method.
The eigenvalues for a 5x5 matrix are difficult to obtain explicitaly, we have to solve a 5th degree polynomial, I know explicit formula at most 4th degree.
If you know the formula for 5th degree polynomial, please let me know.
You are right that I shouldn’t propose close-form eigenvalue as a general alternative, instead it should be used along with exploiting the structure of the matrix. IIRC NGM tend to be sparse or blocked, unless each compartment is directly tied to the rest, which would be uncommon. If that’s the case we can probably construct the eigenvalues from eigenvalues from the blocks.
One can actually use lower-upper bound to put R_0 into a truncated distribution, so that we turn a numerical issue into an uncertainty issue. This applies regardless if the gap is small or not.
Sure, this Theorem is the base of Galois’s Theory!
1 Like
Yeah, however the original NGM is 20x20, using permutation matrices and Binet’t Theorem I obtained that R0, is the spectral radius of a 5x5 matrix with every strictly positive components.
It a age-structured SIR-type model, with 5 age groups, exposed, presymptomatic, symptomatic, asymptomatic.
Thank you for all suggestions, I very appreciated!
1 Like