It's based on Burnside's lemma.
The rotation group $G$ of the necklace is a cyclic group of order $n$. Let $\alpha$ be a generator of $G$, i.e. a rotation of order $n$, such as the rotation by one bead in the positive direction. Thus $G=\{\alpha^1,\alpha^2,\dots,\alpha^n\}$ where of course $\alpha^n$ is the identity permutation.
For $k\in\{1,2,\dots,n\}$, the rotation $\alpha^k$ is a permutation of order $\frac n{(k,\,n)}$ where $(k,\,n)$ is the greatest common divisor of $k$ and $n$; it partitions the set of $n$ beads into $(k,\,n)$ orbits, each of size $\frac n{(k,\,n)}$.
A coloring is invariant under $\alpha^k$ if and only if it's constant on each orbit; thus, with $2$ colors, the number of invariant colorings for $\alpha^k$ is $2^{(k,\,n)}$. According to Burnside's lemma, the number of indistinguishable colorings is obtained by averaging the number of invariant colorings over all elements of the group; thus
$$Z_n=\frac1n\sum_{k=1}^n2^{(k,\,n)}.$$
All that remains is to verify that
$$\sum_{k=1}^n2^{(k,\,n)}=\sum_{d|n}\phi(d)2^{n/d}.$$
This is true because
$$\phi(d)=|\{k\in\{1,\dots,n\}:(k,\,n)=\frac nd\}|$$
or, equivalently,
$$\phi(\frac nd)=|\{k\in\{1,\dots,n\}:(k,\,n)=d\}|.$$
For an arbitrary number $c$ of colors, just replace $2$ with $c$ in all the formulas; the number of indistinguishable colorings is
$$\frac1n\sum_{d|n}\phi(d)c^{n/d}.$$
When $c=1$ this reduces to the familiar identity
$$\sum_{d|n}\phi(d)=n.$$
