Perhaps, the form given by [Lagrange inversion theorem](https://en.wikipedia.org/wiki/Lagrange_inversion_theorem) cannot be much simplified here. It expresses the $n$-th coefficient of a series reversion as the sum of $n-1$ values of exponential Bell polynomials. From the practical perspective, since $\theta_3-1$ contains nonzero coefficients only at square powers, computation of the reversion coefficients amounts to iterating over partitions of $n-1$ into squares ([OEIS A001156](https://oeis.org/A001156)).