Solving the Kolmogorov forward equation for transition probabilities Let $\lambda \mu > 0$ and let $X$ be a Markov chain on $\{1,2\}$ with generators
$$
Q = \begin{pmatrix} -\mu & \mu \\ \lambda & -\lambda \end{pmatrix}$$
Write down the forward equations and solve them fro the transition probabilities, $p_{i,j}(t), i,j = 1, 2$.
The answers say:
"We have that 
$$p_{11}' = -\mu p_{11} + \lambda p_{12}, \, \, \, p_{22}' = -\lambda p_{22} + \mu p_{21}$$
where $p_{12} = 1 - p_{11}$ and $p_{21} = 1 - p_{22}$. Solve these subject to $p_{ij}(t) = \delta_{ij}$, the Kronecker delta, to obtain that the matrix $P_t = \{p_{ij}(t)\}$."
The bits I don't understand are these
1) How do they get those equations in the first place? For example, why have they added $\lambda p_{12}$ when the coefficient for $p_{12} = \mu$?
2) How do you solve these equations? What do you do with that Kronecker delta to solve them?
 A: This is a continuous time Markov process; $\mathbb{Q}$ is an infinitesimal generator of the transition matrices $\mathbb{P}(t)$ giving the transition probabilities over a span of time $t \ge 0$, the primes denote differentiation with respect to $t$, and the "Kronecker delta" is the initial condition $\mathbb{P}(0) = \mathbb{Id}_2$ corresponding to making no transitions at all during zero time.
Think of the process (beginning in any state $\mathbb{x}_0$) as reaching a state $\mathbb{x}_t = (a_t, b_t)$ after time $t$ by means of many tiny jumps over short intervals $0=t_0 \lt t_1 \lt t_2 \lt \cdots \lt t_n = t$.  The infinitesimal generator is a first-order approximation to the changes.  To write this down, use the convenience notation $a(i) = a_{t_i}$ and $b(i) = b_{t_i}$:
$$(a(i+1), b(i+1)) = (a(i), b(i)) + (a(i), b(i)) \cdot \mathbb{Q} + O((t_{i+1}-t_{i})^2).$$
(The big-O notation refers to a term that is proportional to the square of the elapsed time during the interval from $t_i$ to $t_{i+1}$: for very small intervals, this square becomes negligible.)  In matrix notation this statement is
$$\mathbb{x}(i+1) \approx \mathbb{x}(i) \cdot (\mathbb{Id}_2 + \mathbb{Q}).$$
To focus on the change, subtract $\mathbb{x}(i)$ from both sides:
$$\mathbb{x}(i+1) - \mathbb{x}(i) \approx \mathbb{x}(i) \cdot \mathbb{Q}.$$
If we ignore the fact that the $\mathbb{x}$'s are vectors and $\mathbb{Q}$ is a matrix, and pretend they act like real-valued functions and a constant, respectively, and consider only "infinitesimal" changes in $\mathbb{x}$, we would write down a differential equation
$$\frac{d\mathbb{x}}{\mathbb{x}} = \mathbb{Q} dt$$
whose formal solution--paralleling the usual treatment in elementary Calculus--is
$$\mathbb{x}(t) = \mathbb{x}(0) \cdot \exp(\mathbb{Q} t),$$
indicating that
$$\mathbb{P}(t) = \exp(\mathbb{Q} t),$$
whatever that might mean!
It turns out that all this makes sense and can be justified rigorously.  Not only that, the exponential of a matrix can be found in various ways: you can use the series expansion for the exponential,
$$\exp(\mathbb{Q}t) = \sum_{i=0}^{\infty} \frac{(\mathbb{Q}t)^i}{i!},$$
or--this tends to be much easier in practice--when you can find a basis in which $\mathbb{Q}$ is diagonal, its exponential is obtained by exponentiating the diagonal entries (in the usual way: they are just numbers).  Specifically, I have computed that
$$\mathbb{Q}t = \pmatrix{-\mu t &\mu t \\ \lambda t & -\lambda t} = \pmatrix{1 & -\frac{\mu }{\lambda } \\ 1 & 1} \cdot \pmatrix{0 & 0 \\ 0 & -t(\lambda + \mu)} \cdot \pmatrix{1 & -\frac{\mu }{\lambda } \\ 1 & 1}^{-1}.$$
In this new basis, the infinitesimal generator is represented by the middle (diagonal) matrix and its exponential is
$$\exp{ \pmatrix{0 & 0 \\ 0 & -t(\lambda + \mu)}} =  \pmatrix{\exp(0) & 0 \\ 0 & \exp(-t(\lambda + \mu))} = \pmatrix{1 & 0 \\ 0 & \exp(-t(\lambda + \mu))}.$$
Changing back to the original basis gives
$$\mathbb{P}(t) = \exp(\mathbb{Q} t) = \frac{1}{\lambda + \mu}\pmatrix{\lambda + \mu \exp(-t(\lambda+\mu)) & \mu - \mu \exp(-t(\lambda+\mu)) \\ \lambda - \lambda\exp(-t(\lambda+\mu)) & \mu + \lambda\exp(-t(\lambda+\mu))}.$$
To check, it's easy to verify the rows sum to unity.  We had better assume $\mu$ and $\lambda$ have values that make all the entries non-negative, too: that's where the condition $\lambda \mu \gt 0$ comes in.  So at least we have a bona fide transition matrix.  To check that it has $\mathbb{Q}$ as its infinitesimal generator, let $t$ be such a small time interval that we can neglect all terms of order $t^2$ or greater in the calculations.  In particular, when computing the exponentials we can stop at the first term of the series: this is the approximation $\exp(x) \approx 1 + x$ for very small numbers $x$.  So, when $t$ is sufficiently small, we may estimate
$$\exp(-t(\lambda+\mu)) \approx 1 - t(\lambda+\mu)$$
and plugging this in to our formula for $\mathbb{P}(t)$ yields
$$\mathbb{P}(t) = \pmatrix{1 - t\mu & t \mu \\ t \lambda & 1 - t \lambda} + O(t^2) = \pmatrix{1 & 0 \\ 0 & 1} + t \mathbb{Q} + O(t^2).$$
It is also evident that $\mathbb{P}(0) = \mathbb{Id}_2$ (the identity matrix or "Kronecker delta"), exactly as intended.  Returning full circle to the initial setting, if we begin with the distribution $\mathbb{x}$, then after time $t$ the distribution will be
$$\mathbb{x} \cdot \mathbb{P}(t) = \mathbb{x} \cdot \left(\pmatrix{1 & 0 \\ 0 & 1} + t \mathbb{Q} + O(t^2)\right) = \mathbb{x} +  \mathbb{x Q} t + O(t^2).$$
That is, to first order in the time elapsed $t$, the change in $\mathbb{x}$ is proportional to $\mathbb{x Q}$ and to $t$: that's precisely what an infinitesimal generator tells us.  So even if you don't believe (or even understand) the manipulations with the matrix differential equation and the matrix exponential, this check fully justifies the answer.
