Benefit of using Matsubara Green function Physicists often calculate Matsubara Green function and then perform an analytic continuation $i\omega_n \rightarrow \omega +i\eta$ to obtain the retarded Green function.
Why is doing so better than directly computing the retarded Green function?
 A: That's a good question. 
There is no perturbation series expansion for the retarded Green function because to apply Wicks theorem you need time-ordered expressions. You could of course just write down the time-ordered Green functions on the real axis and try to obtain the retarded Green function from the corresponding perturbation expansion. However, this path gets significantly cumbersome for non-zero temperatures. Because then you need an extra perturbation expansion for the Boltzmann weights in your Green function. Essentially its possible though.
Here is where the power of Matsubara comes in. 
(i) Matsubara merges time and temperature into a single variable. The perturbation series is then carried out for this imaginary time parameter. 
(ii) While this is already a great simplification, it turns out, that the Matsubara Green function and the retarded Green function are related by just the simple analytic continuation $i\omega_n \rightarrow \omega + i\delta$.
Hence the answer to your question is, Matsubara makes calculations easier.
A: Matsubara formalism is a clever trick to simplify the calculations. Indeed, when we are dealing with perturbations that are homogeneous in time and space (like particle-particle interactions), the calculations are significantly simplified by performing a Fourier transform in time and space. This works very well for zero temperature formalism, but fails at finite temperatures. Matsubara trick saves the advantages of the Fourier transform, by considering imaginary time interval $[0, -i\beta]$ (or, equivalently, imaginary frequencies). Note that in either case one deals with the time ordered Green's function, for which the Wick's theorem can be applied (it is true that some problems can be solved using only the retarded Green's, but these are usually solvable by even simpler methods).
Another approach is using the Keldysh formalism (also associated with the names of Kadanoff and Baym), where the time contour is taken to run to $t_0+\infty$, then back to $t_0$ and then down to $t_0-i\beta$. Then one can essentially neglect the $[0, -i\beta]$ part, and avoid using the Matsubara frequencies. There is however a price to pay: one now has to keep track of whether the time variables are on the forward or backward branches of the time contour, which necessitates using three different Green's functions  - e.g., retarded, advanced and the Keldysh ($G^R, G^A, G^K$) or retarded, advance and smaller ($G^R, G^A, G^<$). You can find the references here (Review by Rammer&Smith could serve as a crash course on Keldysh).
