the inverse of $f(x)=x-x^p$

As written, this works only for integer $p>1$. You can use (binomial expansion) $$(1-w^{p-1})^{-k}=\sum_{n\geqslant 0}\binom{n+k-1}{n}w^{n(p-1)},$$ see that $g_k\neq 0$ only when $k-1=n(p-1)$ for an integer $\color{red}{n\geqslant 0}$, and get $$g_{n(p-1)+1}=\binom{np}{n}\big(n(p-1)\big)!$$ as expected. But there is a simplification that also handles non-integer $p>1$. The solution $w(z)$ of $w-w^p=z$ we look for is of the form $w(z)=zy(z)$, where $y-z^{p-1}y^p=1$, so that $y$, as a function of $x=z^{p-1}$, satisfies $$x=\frac{y-1}{y^p},$$ and we apply the theorem to this equation (at $y=1$). This gives $y=1+\sum\limits_{k\geqslant 1}g_k x^k$, where $$g_k=\frac{1}{k!}\lim_{y\to 1}\left(\frac{d}{dy}\right)^{k-1}\color{blue}{y^{pk}}=\frac{1}{k}\binom{pk}{k-1}=\frac{1}{(p-1)k+1}\binom{pk}{k}.$$

Suppose we have

$$z = q(z) - q(z)^p$$

with $p\ge 2$ an integer and we seek

$$q(z) = \sum_{n\ge 0} Q_n z^n.$$

Start with some basic observations namely that

$$[z^0] (q(z)-q(z)^p) = [z^0] z = 0 = Q_0 - Q_0^p.$$

We will choose the branch that has $Q_0 = 0.$ Furthermore we have

$$[z^1] (q(z)-q(z)^p) = [z^1] z = 1 = [z^1] (Q_1 z + \cdots - Q_1^p z^p - \cdots) = Q_1$$

and hence $Q_1 = 1.$ Using the Cauchy Coefficient Formula we write

$$n Q_n = [z^{n-1}] q'(z) = \frac{1}{2\pi i} \int_{|z|=\epsilon} \frac{1}{z^n} q'(z) \; dz.$$

We put $q(z) = w$ to that $q'(z) \; dz = dw.$ With the chosen branch we obtain

$$\frac{1}{2\pi i} \int_{|w|=\gamma} \frac{1}{(w-w^p)^n} \; dw = \frac{1}{2\pi i} \int_{|w|=\gamma} \frac{1}{w^n} \frac{1}{(1-w^{p-1})^n} \; dw.$$

This yields with the factor in front

$$\frac{1}{n} [w^{n-1}] \frac{1}{(1-w^{p-1})^n}$$

so that we must have $n=(p-1)k+1$ where $k\ge 0.$ We find

$$\frac{1}{(p-1)k+1} [w^{(p-1)k}] \frac{1}{(1-w^{p-1})^{(p-1)k+1}} \\= \frac{1}{(p-1)k+1} [w^{k}] \frac{1}{(1-w)^{(p-1)k+1}} = \frac{1}{(p-1)k+1} {k+(p-1)k\choose k}.$$

This finally yields

$$\bbox[5px,border:2px solid #00A000]{ q(z) = \sum_{k\ge 0} \frac{z^{(p-1)k+1}}{(p-1)k+1} {pk\choose k}.}$$