Сlosed formula for $(g\partial)^n$ The objective is to obtain a closed formula for:
$$
\boxed{A(n)=\big(g(z)\,\partial_z\big)^n,\qquad n=1,2,\dots}
$$
where $g(z)$ is smooth in $z$ and $\partial_z$ is a derivative with respect to $z$. 
I think the first few terms are,
\begin{equation}
\begin{aligned}
A(1) &= g\,\partial\\
A(2)&= g\,(\partial g)\,\partial+g^2\,\partial^2\\
A(3)&= \big[(\partial^2g)g^2+(\partial g)^2g\big]\partial+3(\partial g)g^2\,\partial^2+g^3\partial^3\\
A(4) &= \big[(\partial^3g)g^3+4(\partial^2g)(\partial g)g^2+(\partial g)^3g\big]\partial\\
&\quad +\big[4(\partial^2g)g^3+7(\partial g)^2g^2\big]\partial^2+6(\partial g)g^3\partial^3+g^4\partial^4\\
&\,\,\vdots
\end{aligned}
\end{equation}
and perhaps there is a simple pattern that I'm failing to see. 
The partitionings of the $\partial$ and $g$ are reminiscent of Bell polynomials but the coefficients are more complicated. Perhaps it is useful to make explicit that the general expansion is of the form:
$$
(g\,\partial)^n=g^n\sum_{p=0}^{n-1}a_{n,p}(g)\,\partial^{\,n-p}
$$
with,
$$
a_{n,p}(g)=\sum_{m_1+2m_2+\dots+pm_{p}=p} C_{n,p}(m_1,\dots,m_{p})\Big(\frac{\partial g}{g}\Big)^{m_1}\Big(\frac{\partial^2 g}{g}\Big)^{m_2}\dots \Big(\frac{\partial^{p} g}{g}\Big)^{m_{p}}\qquad (*)
$$
and the latter sum is over all non-negative integers, $\{m_1,\dots,m_{p}\}$, subject to:
$$
m_1+2m_2+\dots+pm_{p}=p
$$ 
From this viewpoint the objective is to determine the coefficients $C_{n,p}(m_1,\dots,m_{p})$, which in turn depend on all integers, $n$, $p$ and $\{m_1,\dots,m_p\}$. 
Any ideas?
Many thanks in advance.
 A: The Ihara reference "Derivations and automorphisms on non-commutative power series" (open archive now) in OEIS A139605 contains an explicit formula for the coefficients you are looking for, obtained from the Comtet ref. "Une formule explicite pour les puissances successives de l'opérateur de dérivation de Lie."
See A139605 (also related OEIS A145271) for simple matrix computations for these partition polynomials and numerous other references.  
The formula section of A139605 contains the matrix formula. Multiply the $n$-th diagonal (with $n=0$  the main diagonal) of the lower triangular Pascal matrix A007318 by $g_n = D_x^n g(x)$  to obtain the matrix $VP$  with $VP_{n,k} = \binom{n}{k}g_{n-k} $. Then 
$$(g(x)D_x)^n = (1, 0, 0,..) [VP \dot \; S]^n (1, D, D^2, ..)^T,$$
where S is the shift matrix A129185, representing differentiation in the divided powers basis $x^n/n!$. 
Example:
$$(g(x)D_x)^3$$
$$= (1, 0, 0, 0) [VP \dot \; S]^3 (1, D, D^2, D^3)^T$$
$$= \begin{pmatrix}
1 &  0 &  0 & 0
\end{pmatrix} \begin{pmatrix}
0 & g_0 & 0 & 0 \\ 
0 & g_1 & g_0 & 0\\ 
0 & g_2 & 2g_1 & g_0 \\ 
0 & g_3 & 3g_2 & 3g_1
\end{pmatrix}^3 \begin{pmatrix}
1 \\ 
D \\ 
D^2 \\ 
D^3 \end{pmatrix} $$
$$ = [g_0g_1^2  + g_0^2 g_2] D + 3 g_0^2g_1 D^2 + g_0^3D^3 $$
And, the pdf Mathemagical Forests gives a diagrammatic method for creating forests of trees through "natural growth" that represent the partition polynomials.
