Grand canonical partition function: factorization The grand canonical ensemble partition function is defined to be 
$$ \mathcal{Z} := \sum_{\forall |n\rangle} e^{\beta \mu N_{|n\rangle} - \beta E_{|n\rangle}} $$
being $|n\rangle$ a notation for each microstate (not necessarily quantum, just a microstate), $N_{|n\rangle}$ the number of particles of that microstate and $E_{|n\rangle}$ the energy corresponding to this microstate. This is the formal definition.
Now I have from my lecture notes that when we are dealing with non-interacting and indistinguishable particles (e.g. an ideal gas) the partition function can be written as 
$$ \mathcal{Z} = \prod_{\forall \epsilon_i} \sum_{\forall \text{ allowed } n}(z e^{-\beta \epsilon _i})^n $$
now being $\epsilon _i$ each "monopartiuclar state" and "$\forall \text{ allowed } n$" is $\{0,1\}$ for fermions and $\{0,1,2,\dots\}$ for bosons. The "monoparticular states" are the energy levels that each particle can be in.
The question: How do we go from the definition to the second formula?

My approach (wrong, or at least incomplete)
If particles do not interact between one another then the energy of the microstate $|n\rangle$ can be written as a summation 
$$ E_{|n\rangle} = \epsilon_1 + \epsilon_2 + \dots + \epsilon_{|n\rangle} = \sum_{i=1}^{N_{|n\rangle}} \epsilon_i $$
where $\epsilon_i$ is the energy of each particle. Thus the partition function writes 
$$ \mathcal{Z} = \sum_{\forall |n\rangle} e^{\beta \mu N_{|n\rangle}} \prod_{i=1}^{N_{|n\rangle}} e^{-\beta \epsilon _i}$$
where I have already expanded the exponential of a summation as a product  of exponentials.
Now, if particles are identical then the allowed values for $\epsilon _i$ are for all the particles the same, say $\epsilon _i \in \{\varepsilon_0, \varepsilon_1, \dots\}$. Using this notation, $\epsilon_8 = \varepsilon_3$ reads as "particle number 8 is in the third energy level". This allows to arrange the microstates of the system as follows (sorry for changing the language in the pic):

Now I don't know how to go on... Any help is appreciated.
 A: You are misunderstanding what is meant by a product over monopartiuclar states. This is not a product over the states of the $N$ particles in the system, it is a product over all possible single particles states. The sum is then over the occupation number of those states, i.e. the number of particles actually in that state (which may be $0$). The advantage of this occupation number approach to over keeping track of individual particles is that it automatically takes account of particle indistinguishably.
The total number of particles in a state $\gamma$ is clearly simply the sum of the number of particles in each single particle state
$$
N_\gamma = \sum_i n_i
$$
The total energy is the energy of each single particle state, times the number of particles in that state
$$
E_\gamma = \sum_i \epsilon_i n_i
$$
Now the formula follows from a simple manipulation 
\begin{align}
\mathcal{Z} &= \sum_{\gamma} e^{-\beta(E_\gamma - \mu N_\gamma)}\\
&= \sum_{\gamma} e^{-\beta\sum_i n_i(\epsilon_i-\mu)}\\
&= \sum_{n_0}\sum_{n_1}\ldots\; \left(\prod_i e^{\beta n_i(\mu - \epsilon_i)}\right)\\
&= \prod_i \sum_{n_i} \left(ze^{-\beta\epsilon_i}\right)^{n_i}
\end{align}
