  We give a full account of the Numerical Stochastic Perturbation Theory method
for Lattice Gauge Theories. Particular relevance is given to the inclusion of
dynamical fermions, which turns out to be surprisingly cheap in this context.
We analyse the underlying stochastic process and discuss the convergence
properties. We perform some benchmark calculations and - as a byproduct - we
present original results for Wilson loops and the 3-loop critical mass for
Wilson fermions.
