Why can I replace a gauge field $A^{\mu}$ by the current $j^{\mu}$ it couples to in the calculation of a Greens function? I am reading about anomalies in QFT at the moment and have a related question.
Often people calculate the time ordered expectation value of some fields (in QED for example) by replacing the field $A_\mu$ with the current $j_\mu = \bar{\psi}\gamma_\mu\psi$.
For example in Peskin and Schröder (p. 311) they derive the Ward identities in QED and state that the object $$\langle T j^\mu(x)\psi(y)\overline{\psi}(z)\rangle$$ corresponds to the electron-photon vertex. But should the vertex not correspond to $$\langle T A^\mu(x)\psi(y)\overline{\psi}(z)\rangle$$
where $A^\mu$ is the EM potential?
The same thing happens in the calculation of the triangle anomaly. Often you see people calculate something like $\langle j_1^\mu(x)j_2^\nu(y)j_3^\rho(z)\rangle$ where the $j$s are (axial or vector) currents of fermions coupled to some vector fields. This gives than the diagram with three external vector bosons and a circulating fermion.
It is clear to me that $\langle T j^\mu \ldots \rangle$ and $\langle T A^\mu \ldots \rangle$ give me the same diagrams when I start to calculate. I am just wondering if I can prove in any way (from a generating functional perhaps) that these Green's functions are really the same. (I use Green's function in the sense of time ordered vacuum expectation value of fields.)
Hopefully somebody can clarify this a little.
 A: There's a minor mistake in your question, $\langle T j^\mu(x)\psi(y)\overline{\psi}(z)\rangle$ is not the same diagram as $\langle T A^\mu(x)\psi(y)\overline{\psi}(z)\rangle$, the latter contains one extra bare photon propagator, this can be seen by considering possible contractions of the fields. So your question should be modified as: how to prove the equation
$\langle T A^\mu(x)\psi(y)\overline{\psi}(z)\rangle=\int d^4s G^\mu_{\ \ \nu}(x,s)\langle T j^\nu(s)\psi(y)\overline{\psi}(z)\rangle\cdots\cdots(1),$
where $G(x,s)$ is the bare photon propagator?
Like you said, this can be justified just by inspecting the corresponding Feynman diagrams, thus you already have a perturbative proof! So I guess what you really want is a non-perturbative proof, here it goes:
Recall the bare photon propagator is really just the inverse operator(alias Green's function) of the differential operator $d^{\mu}_{\ \ \nu} $in the field equation acting on $A^\nu$, for example, with a Lorentz gauged Lagrangian $\mathcal{L}=-\frac{1}{4}F_{\mu\nu}F^{\mu\nu}+\frac{1}{2\xi}(\partial_\mu A^\mu)^2-j^\mu A_\mu+\mathcal{L}_{matter}$, the field equation is 
$d^{\mu}_{\ \ \nu}A^\nu=[-\delta^\mu_{\ \ \nu}\partial^2+(1+\frac{1}{\xi})\partial^\mu\partial_\nu]A^\nu=j^\mu\cdots\cdots(2),$
and the propagator is just $G^\mu_{\ \ \nu}=(d^{-1})^\mu_{\ \ \nu}$, that is,
$[-\delta^\mu_{\ \ \nu}\partial^2+(1+\frac{1}{\xi})\partial^\mu\partial_\nu]G^\nu_{\ \ \rho}(x,s)=\delta^\mu_{\ \ \rho}\delta^4(x-s)$. So if we "de-convolute"(i.e. act $d^{\mu}_{\ \  \nu}$ on) both sides of equation (1), we would get an equivalent expression as
$d^{\mu}_{x\ \nu}\langle T A^\nu(x)\psi(y)\overline{\psi}(z)\rangle=\langle T j^\mu(x)\psi(y)\overline{\psi}(z)\rangle\cdots\cdots(3),$
where the subscript $x$ in $d^{\mu}_{x\ \nu}$ indicates the space-time variable it's acting on. So we might as well just prove (3). Now (3) looks like (2) a lot, but one has to be careful because $d^{\mu}_{x\ \nu}$ involves time derivatives hence does not commute with the time-ordering symbol. One can expand time-ordering in terms of step functions and differentiate term by term carefully, then summing up together the resulting terms will show (3) is indeed true, as if $d^{\mu}_{x\ \nu}$ commutes with time ordering. This is a lot of work, but one can use a much more efficient method called
Dyson-Schwinger(D-S) equation(although there's still  something subtle about it that I don't understand, but it has never disappointed me in practice). Applying D-S equation to $\langle T\psi(y)\overline{\psi}(z)\rangle$ by considering variations in gauge fields will immediately lead to equation (3).
Last but not least(maybe OP is aware of the following, but I can't tell from the question, so I'll write it anyway), I must emphasize in general a simple relation like equation (1) doesn't hold if you have more than one gauge field in the time-ordered product. For example, you can't get $\langle T A^\mu(x) A^\rho(y) j^\sigma(z) \rangle$ just by attaching two bare photon propagators to $\langle T j^\mu(x) j^\rho(y) j^\sigma(z) \rangle$. Again you can convince yourself in two different ways: 
In terms of Feynman diagrams, now you can have a Wick contraction between external points $A^\mu(x)$ and $A^\rho(y)$(unlike when you only have one gauge field operator), this is just a bare photon propagator. The differences can be visually captured by the following figures(note that a blob is not necessarily a connected part):    

The second way is again in terms of D-S equation, you start from $\langle j^\sigma(z)\rangle$ and apply D-S equation iteratively twice, you'll get exactly the right terms corresponding to the third diagram above, and as it turns out the contraction between $A^\mu(x)$ and $A^\rho(y)$ corresponds to what's called a "contact term" in D-S equation. This time you can't pretend $d^{\mu}_{\ \ \nu}$ commutes with time-ordering symbol!
