Guessing of the vertex factor of derivative interactions I want to ask if my "empirical" method of guessing the vertex factor is correct. As far as I understand in canonical QFT the Feynman rules are only a way to rearrange the important parts of the perturbative expansion. So, if I have an interaction lagrangian (possibly also involving derivative terms) like scalar QED: $$\mathcal{L} = -i e A_{\mu}(\psi^{\dagger} \partial^{\mu} \psi - \psi \partial^{\mu} \psi^{\dagger})$$ one has to guess what is the vertex term in order to avoid doing the explicit calculation. 

In a vertex like the one above I would consider an initial state like $b^{\dagger}(p) |0 \rangle $ and a final state $a^{\dagger}(k) b^{\dagger}(q) | 0 \rangle $ (one photon and one scalar). The vertex should bring me from the initial state to the final state. Only the first part of the interaction is $-i e A_{\mu} \psi^{\dagger} \partial^{\mu} \psi$. The $\psi$ field annihilate the scalar in the initial state through $b(p) e^{-i p x}$ so the derivative brings down a factor $-ip$. The second part has $ b^{\dagger} (p) e^{i p x}$ acting on the initial state that brings down a factor $ip$. So due to the minus sign the two momentum get summed and up to an overall phase the vertex is $e (p+q)^{\mu}$. As $\psi$ and $\psi^{\dagger}$ are different fields I don't have to worry about combinatorial factors.
I would apply a similar reasoning to every vertex (even if the legs are propagator) for every interacting field theory. Is it sensible?
If not how can I guess every vertex without having to do long calculations?
I want to avoid the use of the functional formulation of QFT. Anyone?
 A: Your derivation is indeed sometimes used as a shortcut to argue that a derivative interaction corresponds to a factor of momentum for each derivative. Unfortunately, it is far from convincing for several reasons.
The most important one is that the momentum in $b^\dagger(p)|0\rangle$ is on-shell, $p^2=m^2$, while the momentum in the Feynman rules is off-shell, $p^2\neq m^2$. This distinction is very important because Feynman intergrals are over all $p$, not over the mass-shell $\delta(p^2-m^2)$, so the $p$ in $b^\dagger(p)|0\rangle$ is very different from the $p$ in $-iep^\mu$. Even if you use the same symbol, they are fundamentally different objects. The first one is the eigenvalue of $\hat P^\mu$ and the second one is a Fourier variable.
That being said, your derivation is the best you can do in the canonical formalism. For a better justification you need to use the functional formulation of QFT, something you don't want to do. For completeness, let us mention that in the functional formalism a Feynman vertex is just a functional derivative,
$$
\frac{\delta}{\delta\psi}\frac{\delta}{\delta\bar\psi}\frac{\delta}{\delta A}\mathcal L(\psi,\bar\psi,A)
$$
from which the factor $-iep^\mu$ readily follows. It is an immediate result, and much more easily derived than your heuristic calculation. The functional formalism is more rigorous and more convenient -- although I can see why you want to avoid it at first.
How can we offer a slightly better derivation of the Feynman rules of derivative interactions within the canonical formalism? Well, you can introduce the covariant time-ordering, as discussed in this PSE post. It is kind of tricky, so it may not be along the lines of what you are looking for. But it is there, and it is comforting to know that it can be done.
