The most interesting part, I think, is that you can write the equation quickly by inspection.
For the impedance matrix, the terms on the diagonal are just the sum of the impedances (resistances in this case) in the given loop. The off-diagonal terms are just the negative of the impedances that are shared by the loops indexed by the matrix entry. So in this case, for example, R2 is shared by loops 1 and 2, so z12 and z21 are both -R2. Note that the impedance matrix is symmetrical about the diagonal, so it's very quick to fill in.
The voltage vector is just the sum of the voltage source rises in the given loop (where a rise is taken to mean that a given voltage source goes from - to + in the direction of the loop current).