 Solution

For transition-intensity vector $q$ and generator $Q(q)$, the expected infectious occupancy over three months from $S$ is the Markov reward model quantityInsert the fitted intensities to obtain $m(\hat{q})$, evaluating the matrix exponential and integral numerically if necessary.

If $\eta =(\log q_{SE},\log q_{EI},\log q_{IR})$ has estimated covariance $V$, calculate the numerical gradient $g={\nabla }_{\eta }m(e^{\eta })$ at $\hat{\eta }$. The delta method gives estimated variance $g^{T}Vg$ and the approximate confidence intervalSimulation from $N(\hat{\eta },V)$ followed by transformation through $m$ gives a useful asymmetric alternative.