If the (say, $d$-dimensional) base manifold $M$ is *parallelizable* (i.e. $TM\to M$ is *trivial*), then the answer to both questions is *yes* even globally, provided we choose a (say, torsion-free) covariant derivative $\nabla^M$ on $TM$ and a covariant derivative $\nabla^E$ on $E\to M$. We may then combine these connections by means of the Leibniz rule to define a covariant derivative on $\otimes^k T^*M\otimes E\to M$ for all $k$. Let us denote this common extension of $\nabla^M$ and $\nabla^E$ simply by $\nabla$. One may then define iterated covariant derivatives $\nabla^k$ of any order $k\geq 0$ on $E\to M$ recursively as follows: $$\nabla^0\phi=\phi\ ,\,\nabla^1\phi=\nabla^E\phi\ ,\,\nabla^{k+1}\phi=\nabla(\nabla^k\phi)\ ,\quad\phi\in\Gamma(M,E)\ ,$$ so that $\nabla^k\phi\in\Gamma(M,\otimes^kT^*M\otimes E)$ for every $k\geq 0$. Let us denote the contraction of $\nabla^k\phi$ with $k$ smooth vector fields $Y_1,\ldots,Y_k\in\mathfrak{X}(M)$ by $\nabla^k_{Y_1,\ldots,Y_k}\phi\in\Gamma(M,E)$. The key fact we need (apparently due to R. Palais, but I am not sure), holding *even when $M$ is not parallelizable*, is:

Any linear partial differential operator of order $r$ $P:\Gamma(M,J^rE)\to\Gamma(M,F)$ may be written uniquely as $$\tag{1}\label{e1}P\phi=\sum^r_{k=0}A_k\nabla^k\phi\ ,\quad\phi\in\Gamma(M,E)\ ,$$ where $A_k\in\Gamma(M,S^k TM\otimes E^*\otimes F)$ for each $k=0,1,\ldots,r$, $S^kTM\to M$ is the symmetric, rank-$k$ contravariant tensor sub-bundle of $\otimes^kTM\to M$ and $E^*\to M$ is the dual bundle to $E\to M$. Particularly, one may replace $\nabla^k$ in \eqref{e1} by its symmetrized part.

The coefficients $A_k$ of the representation \eqref{e1} of $P$ may of course be seen as vector bundle morphisms from $\otimes^kT^*M\otimes E\to M$ to $F\to M$ covering $\mathrm{id}_M$. This is "almost" what you want, what gets in the way is precisely the potential non-triviality of $TM\to M$ - that is where the hypothesis of parallelizability of $M$ enters.

Recall now that the latter amounts to being able to choose $d$ smooth vector fields $X_1,\ldots,X_d\in\mathfrak{X}(M)$ such that $\{X_1(p),\ldots,X_d(p)\}$ is a basis of $T_pM$ for every $p\in M$ - that is, a (global) *linear frame* in $TM$. Let $\omega^1,\ldots,\omega^d\in\Omega^1(M)$ be the corresponding linear coframe (i.e. $\omega^i(X_j)=\delta^i_j$) - from \eqref{e1} one may write $$\tag{2}\label{e2}P\phi=\sum^r_{k=0}\sum^n_{j_1,\ldots,j_k=1}A_k(\omega^{j_1}\otimes\cdots\otimes\omega^{j_k})\nabla^k_{X_{j_1},\ldots,X_{j_k}}\phi\ ,\quad\phi\in\Gamma(M,E)\ ,$$ where $A_k(\omega^{j_1}\otimes\cdots\otimes\omega^{j_k})\in\Gamma(M,E^*\otimes F)$ is the contraction of $A_k$ with $\omega^{j_1}\otimes\cdots\otimes\omega^{j_k}\in\Gamma(M,\otimes^kT^*M)$. Formula \eqref{e2} is precisely the affirmative answer to your first question. From there to answering affirmatively your second question is just a(n essentially combinatorial) matter of rewriting $\nabla^k_{X_{j_1},\ldots,X_{j_k}}\phi$ in \eqref{e2} in terms of $\nabla^E_{X_{j_1}}\cdots\nabla^E_{X_{j_{m+1}}}\phi$ and $\nabla^M_{X_{l_1}}\cdots\nabla^M_{X_{l_m}}X_i$ for $i,j_1,\ldots,j_{m+1},l_1,\ldots,l_m=1,\ldots,d$, $m=1,\ldots,k-1$ using Leibniz's rule (the first-order operators you want are, of course, $\nabla^E_{X_j}$ with $j=1,\ldots,d$).

For most applications to global analysis, however, \eqref{e1} suffices, with the advantage that it holds in complete generality. There is, of course, the possibility of using a partition of unity subordinated to a finite atlas (an atlas with $d+1$ charts always exists whenever $M$ is second countable by topological dimension theory, see e.g. the Corollary to Theorem 1.2.I, pp. 20-21 of the book of W. Greub, S. Halperin and R. Vanstone, *Connections, Curvature, and Cohomology, Volume I* (Academic Press, 1972)) and write $P$ in each chart using coordinates, as suggested by you in your question and Dmitri Pavlov in his comment above, but another advantage of \eqref{e1} and \eqref{e2} is that they are both *coordinate-free* representations of $P$.

There is a middle course to circumvent the hypothesis of parallelizability which relies on the following remark: if $U\subset M$ is the (open) domain of a chart $x:U\to\mathbb{R}^d$, then $TU$ is *trivial* (just pick the coordinate vector fields $X_j(p)=(T_px)^{-1}\frac{\partial}{\partial x^j}$, $j=1,\ldots,d$). Let now $\{(U_\alpha,x_\alpha)\ |\ \alpha=1,\ldots,d+1\}$ be a finite atlas of $M$ as in the previous paragraph, $\{f_\alpha\ |\ \alpha=1,\ldots,d+1\}$ a partition of unity subordinate to the open cover $\{U_\alpha\ |\ \alpha=1,\ldots,d+1\}$ of $M$ and $X_{1,\alpha},\ldots,X_{d,\alpha}$ be the linear frame on $U_\alpha$ obtained from $x_\alpha$ as above for each $\alpha$, with corresponding linear coframe $\omega^{1,\alpha},\ldots,\omega^{d,\alpha}$. We can then write $$P\phi=\sum^{d+1}_{\alpha=1}\sum^r_{k=0}\sum^n_{j_1,\ldots,j_k=1}A_k(\omega^{j_1,\alpha}\otimes\cdots\otimes\omega^{j_k,\alpha})\nabla^k_{X_{j_1,\alpha},\ldots,X_{j_k,\alpha}}(f_\alpha\phi)\ ,$$ recalling that $$\nabla^k\phi=\sum^{d+1}_{\alpha=1}\nabla^k(f_\alpha\phi)$$ and, by Leibniz's rule, $$\nabla^k_{Y_1,\ldots,Y_k}(f_\alpha\phi)=\sum_{I\subset\{1,\ldots,k\}}\nabla^{M|I|}_{Y_I}f_\alpha\nabla^{n-|I|}_{Y_{I^c}}\phi\ ,\quad Y_I=(Y_j)_{j\in I}\ .$$ This leads us to the following, not-so-appetizing formula $$\tag{3}\label{e3}P\phi=\sum^{d+1}_{\alpha,\beta=1}\sum^r_{k=0}\sum^n_{j_1,\ldots,j_k=1}\sum_{I\subset\{1,\ldots,k\}}(\nabla^{M|I|}_{X_{j_I,\alpha}}f_\alpha)A_k(\omega^{j_1,\alpha}\otimes\cdots\otimes\omega^{j_k,\alpha})(f_\beta\nabla^{n-|I|}_{X_{j_{I^c},\alpha}}\phi)\ ,$$
where $Y_{j_I,\alpha}=(Y_{j_l,\alpha})_{l\in I}$. Since $\mathrm{supp}\nabla^kf_\alpha\subset\mathrm{supp} f_\alpha$ for all $k,\alpha$, formula \eqref{e3} is globally well defined and yields a positive answer to your first question even in the absence of parallelizability. To answer positively your second question, one proceeds as in the parallelizable case with each term $\nabla^{n-|I|}_{Y_{j_{I^c}}}\phi$ in \eqref{e3}. At this point, yet another advantage of working with \eqref{e1} instead should become clear - namely, economy and simplicity.

Differential analysis on complex manifoldslooks fine). Using a smooth partition of unity subordinate to such a finite cover to patch the local expressions indeed seems to solve the problem, but I will come back to this after I have checked the details (probably tomorrow; I am quite busy right now). $\endgroup$