I am trying to verify the following identity, where $M,N$ and $P$ are differentiable manifolds $\phi : M \rightarrow N$ and $\psi: N \rightarrow P$ are smooth maps between manifolds. Also, $f: P \rightarrow R$ is a $C^{\infty}$ function whose domain is a subset of the manifold P: $$(\psi \circ \phi)_{\star} = \psi_{\star} \circ \phi_{\star}. $$ However, I am getting the composition in reverse order.
Edit: The definition for the pushforward $\phi$ is given as:
$$ (\phi_{\star} X_{p})_{\phi(q)}(f) := X_{q}(\phi^{\star}f) $$
Let $q \in M$. By the definition of the pushforward:
EDIT: $$[(\psi \circ \phi)_{\star}X_{q}]_{\psi \circ \phi (q)} (f) : = X_{q} [(\psi \circ \phi)^{\star} (f)] = X_{q} [(\phi^{\star} \circ \psi^{\star}) (f)] = X_{q} [\phi^{\star} (\psi^{\star}(f))] = (\phi_{\star}X_{q})_{\phi(q)}[\psi^{\star} f] = [(\psi_{\star} \circ \phi_{\star})X_{q}]_{\psi \circ \phi (q)} (f)$$ Hence, $(\psi \circ \phi)_{\star} = \psi_{\star} \circ \phi_{\star}$.
EDIT: I understand that we can do it with the chain rule if M has a set of local coordinates $x^{\mu}$, $N$ has a set of coordinates $y^{\nu}$ and $P$ has a set of local coordinates $z^{\lambda}$. In that case,using Einstein summation notation
$$[(\psi\circ \phi)_{\star} X_{q}]_{\psi\circ\phi(q)} = X^{\mu} \frac{\partial y^{\nu}}{\partial x^{\mu}}\frac{\partial z^{\lambda}}{\partial y^{\nu}} \frac{\partial}{\partial z^{\lambda}} $$
The transformation $\phi_{\star}$ is represented locally by the Jacobian $\frac{\partial y^{\nu}}{\partial x^{\mu}}$ and $\psi_{\star}$ by the locally by the Jacobian $\frac{\partial z^{\lambda}}{\partial y^{\nu}}$.
Hence,
$$[(\psi\circ \phi)_{\star} X_{q}]_{\psi\circ\phi(q)} = (\psi_{\star})^{\lambda}_{\nu}(\phi_{\star})^{\nu}_{\mu}X^{\mu} \frac{\partial}{\partial z^{\lambda}}$$
But I don't see how this fixes the order of multiplication of the matrices.