Ward identity derived from global symmetry and SDE, different from that derived from gauge symmetry? In QED, according to Schwinger-Dyson equation $^{[1]}$,
$$\left(\eta^{\mu\nu}(\partial ^2)-(1-\frac{1}{\xi})\partial^{\mu}\partial^{\nu}\right)\langle 0|\mathcal{T}A_{\nu}(x)...|0\rangle = e\,\langle 0|\mathcal{T}j^{\mu}(x)...|0\rangle + \text{contact terms}.$$
And the term $\left(\eta^{\mu\nu}(\partial ^2)-(1-\frac{1}{\xi})\partial^{\mu}\partial^{\nu}\right)$ is just the inverse bare photon propagator$^{[5]}$, so if we put the photon on shell, then the l.h.s will yield the complete $n$-point Green function with the complete photon propagator removed and also multiplied by a factor $Z_3$, the vector field renormalization constant.
But the r.h.s gives, according to the Ward Identity associated with global symmetry of the Lagrangian,
$$\partial_{\mu}\, \langle 0|\mathcal{T}j^{\mu}(x)...|0\rangle = \text{contact terms}^{[2]}$$
which is the sum of complete $(n-1)$-point complete Green functions each multiplied by a certain $\delta$ function.
So if we truncate all the $n-1$ external complete propagators, then we are left with the proper vertex Ward identity.  
The problem is, now the constant $Z_3$ appeared. 
For, example, if we apply this to $j^\mu=\bar\psi\gamma^\mu\psi$,
$$\partial_{\mu}\, \langle 0|\mathcal{T}j^{\mu}(x)\psi(x_1)\bar\psi(x_2)|0\rangle = -ie[\delta(x-x_1)\langle 0|\psi(x_1)\bar\psi(x_2)|0\rangle-\delta(x-x_2)\langle 0|\psi(x_1)\bar\psi(x_2)|0\rangle]^{[3]},$$
we will get an identity containing $Z_3$ RATHER THAN the well known proper vertex Ward identity (which is derived from the gauge or local symmetry of the Lagrangian), e.g. 
$$q_\mu\Gamma^\mu_P(p,q,p+q)=S^{-1}(p+q)-S^{-1}(p)^{[4]}$$
which doesn't contain $Z_3$. 
Where went wrong? Please help.
[1]. Peskin&Schroeder, An Introduction to Quantum Field Theory, page308, eq.(9.88)
[2]. Peskin&Schroeder, An Introduction to Quantum Field Theory, page310, eq.(9.97)
[3]. Peskin&Schroeder, An Introduction to Quantum Field Theory, page311, eq.(9.103)
[4]. Lewis H.Ryder, Quantum Field Theory, page 263-266, eq.(7.112)
[5]. Peskin&Schroeder, An Introduction to Quantum Field Theory, page297, eq.(9.58)
 A: We can write the Fourier transform of 
$\langle 0|\mathcal{T}A_{\nu}(x)\psi(x_1)\bar\psi(x_2)|0\rangle$
as $$S(p) D_{\nu\alpha}(q) \ e\,\Gamma^{\alpha}(p,q,p+q)S(p+q)$$ where $S(p)$ is the full fermion propagator, $D_{\nu\alpha}(q)$ is the full photon propagator, $\Gamma^{\alpha}(p,q,p+q)$ is the proper vertex function, and an overall momentum conservation delta function has been dropped. Similarly, we can write the Fourier transform of $\langle 0|\mathcal{T}j^{\mu}(x)\psi(x_1)\bar\psi(x_2)|0\rangle
$ as 
$$S(p)V^{\mu}(p,q,p+q)S(p+q)$$ where $V^{\mu}(p,q,p+q)$ is a vertex function that we want to relate to $\Gamma^{\mu}(p,q,p+q).$ The vertex function $V^{\mu}(p,q,p+q)$ enters into the derivation of the Ward-Takahashi identity in Peskin and Schroeder on page 311, but the Ward-Takahashi identity is normally stated in terms of $\Gamma^{\mu}(p,q,p+q)$. Your conundrum (as I understand it) is that according to your analysis of the Schwinger-Dyson equation, $V^{\mu}(p,q,p+q)$ and $\Gamma^{\mu}(p,q,p+q)$ ought to differ by a factor of $Z_3$, but this contradicts the usual statement of the Ward-Takahashi identity where no such factor of $Z_3$ appears. I will argue from the Schwinger-Dyson equation that the longitudinal parts (in $q^{\mu}$) of $V^{\mu}(p,q,p+q)$ and $\Gamma^{\mu}(p,q,p+q)$ are equal, but that the transverse parts differ by the factor of $Z_3$ that you have found. Since only the longitudinal part enters into the Ward-Takahashi identity, the factor of $Z_3$ does not inter into that identity. You may want to review page 246 of Peskin and Schroeder. There they show that only the transverse part of the photon propagator is modified by the self-energy, but that in calculating Feynman diagrams we can simplify the analysis by including the self-energy in the longitudinal part as well because the longitudinal part does not contribute to the Feynman diagrams due to the Ward identity. However, the Schwinger-Dyson equation involves an inverse propagator which does not arise in Feynman diagrams and we need to reevaluate where the self-energy does and does not enter.
Specializing the Schwinger-Dyson equation to the case of $\langle 0|\mathcal{T}A_{\nu}(x)\psi(x_1)\bar\psi(x_2)|0\rangle$ and Fourier transforming, we have
$$\tag{1} (D^{(0)\mu\nu}(q))^{-1} D_{\nu\alpha}(q) S(p) \ e\,\Gamma^{\alpha}(p,q,p+q)S(p+q) = \\ e\, S(p)V^{\mu}(p,q,p+q)S(p+q)$$
where $(D^{(0)\mu\nu}(q))^{-1}$ is the inverse of the non-interacting photon propagator. The Dyson equation for the photon propagator is 
$$\tag{2} D_{\nu\alpha}(q) = D^{(0)}_{\nu\alpha}(q) + D^{(0)}_{\nu\beta}(q) i\Pi^{\beta\gamma}(q)D_{\gamma\alpha}(q) ,$$
so 
$$\tag{3} (D^{(0)\mu\nu}(q))^{-1} D_{\nu\alpha}(q) = \delta^{\mu}_{\alpha} + i\Pi^{\mu\gamma}(q)D_{\gamma\alpha}(q).$$
Equation (1) then implies 
$$\tag{4}\Bigl(\delta^{\mu}_{\alpha} + i \Pi^{\mu\gamma}(q)D_{\gamma\alpha}(q)\Bigr) \Gamma^{\alpha}(p,q,p+q) = V^{\mu}(p,q,p+q).$$
The Ward identity forces the longitudinal part of $\Pi^{\mu\gamma}(q)$ to vanish; that is, $q_{\mu}\Pi^{\mu\gamma}(q) = 0.$
Contracting equation (4) with $q_{\mu}$,
we therefore have
$$\tag{5} q_{\alpha} \Gamma^{\alpha}(p,q,p+q) = q_{\mu}V^{\mu}(p,q,p+q)$$
so no factor of $Z_3$ appears between the longitudinal parts of $\Gamma^{\alpha}(p,q,p+q)$ and $V^{\mu}(p,q,p+q)$ and therefore no factor of $Z_3$ appears in the Ward-Takahashi identity.
The transverse component does not enter the Ward identity for the vertex function but it is useful to consider the transverse component to illustrate where the factor of $Z_3$ does arise. Define $\Pi(q^2)$ by the equation 
$\Pi^{\mu\nu}(q) = q^2(g^{\mu\nu} - q^{\mu}q^{\nu}/q^2)\Pi(q^2).$ The quantity $(g^{\mu\nu} - q^{\mu}q^{\nu}/q^2)$ can be described as a projection operator that projects out the transverse part of a vector. 
Contracting equation (4) with $(g_{\nu\mu} - q_{\nu}q_{\mu}/q^2)$ and using the fact that $\Pi^{\mu\gamma}(q)$ is already transverse, we have
$$\tag{7} \Bigl(g_{\nu\alpha} - q_{\nu}q_{\alpha}/q^2 + i q^2\Pi(q^2)D^T_{\nu\alpha}(q)\Bigr) \Gamma^{\alpha}(p,q,p+q) =\\ \bigl(g_{\nu\mu} - q_{\nu}q_{\mu}/q^2\bigr)V^{\mu}(p,q,p+q),$$
where $$D^T_{\nu\alpha}(q) = \frac{-i}{q^2 (1-\Pi(q^2))} \bigl(g_{\nu\alpha} - q_{\nu}q_{\alpha}/q^2\bigr)$$ is the transverse part of the photon propagator (see page 246, Peskin and Schroeder). Equation (7) can then be written
$$\tag{8} \bigl(g_{\nu\alpha} - q_{\nu}q_{\alpha}/q^2\bigr) \Bigl(1/\bigl(1-\Pi(q^2)\bigr)\Bigr) \Gamma^{\alpha}(p,q,p+q) =\\ \bigl(g_{\nu\mu} - q_{\nu}q_{\mu}/q^2\bigr)V^{\mu}(p,q,p+q).$$ 
Now consider $q^2$ small enough that $\Pi(q^2)\approx \Pi(0)$ and use the relation (Peskin and Schroeder, page 246)
$$Z_3 = \Bigl(1/\bigl(1-\Pi(0)\bigr)\Bigr).$$ We have
$$\tag{9} \bigl(g_{\nu\alpha} - q_{\nu}q_{\alpha}/q^2\bigr) Z_3 \Gamma^{\alpha}(p,q,p+q) = \bigl(g_{\nu\mu} - q_{\nu}q_{\mu}/q^2\bigr)V^{\mu}(p,q,p+q).$$ So we see that the transverse parts of $V^{\mu}(p,q,p+q)$ and $\Gamma^{\mu}(p,q,p+q)$ differ by a factor of $Z_3.$
A: The Schwinger-Dyson (SD) equation mentioned by OP reads
$$ \begin{align}
0~=~& \langle \Omega | T_{\rm cov} \frac{\delta S}{\delta A_{\mu}(x)}\psi(x_1) \bar{\psi}(x_2)|\Omega\rangle \cr
~=~&\left[\color{red}{Z_3}\delta^{\mu}_{\nu}\partial^2 +\left(\frac{1}{\xi}-\color{red}{Z_3}\right)\partial^{\mu}\partial_{\nu} \right] \langle \Omega | T  A^{\nu}(x)\psi(x_1) \bar{\psi}(x_2)|\Omega\rangle\cr
&\pm \color{red}{Z_1}\langle\Omega | T  j^{\mu}(x)\psi(x_1) \bar{\psi}(x_2)|\Omega\rangle
,\end{align}\tag{SD}$$
where the $\pm$ depends on sign conventions. If we contract eq. (SD) with $\partial_{\mu}$, the $Z_3$-factor is projected out, so we may write
$$ \begin{align}
&\partial_{\mu}\underbrace{\left[\delta^{\mu}_{\nu}\partial^2 +\left(\frac{1}{\xi}-1\right)\partial^{\mu}\partial_{\nu} \right]}_{\text{inverse free propagator}} \langle \Omega | T  A^{\nu}(x)\psi(x_1) \bar{\psi}(x_2)|\Omega\rangle\cr
~=~& \mp \color{red}{Z_1}\partial_{\mu}\langle \Omega | T  j^{\mu}(x)\psi(x_1) \bar{\psi}(x_2)|\Omega\rangle
.\end{align}\tag{SD'}$$
Eq. (SD') already essentially resolves OP's paradox, as we shall see.
At the classical level, the global gauge invariance leads via Noether's theorem to electric charge conservation, cf. e.g. this Phys.SE post. The Ward-Takahashi identity (WTI) can roughly speaking be thought of as the corresponding quantum version, cf. e.g. this Phys.SE post. In particular, we stress that the WTI is intimately tied to electric charge conservation.
OP's observation that the WTI can be derived in several ways is interesting, cf. e.g. this Phys.SE post.

*

*On one hand, the WTI for connected diagrams
$$\begin{align} & \color{red}{Z_2}\partial_{\mu}\langle \Omega | T j^{\mu}(x)\psi(x_1) \bar{\psi}(x_2)|\Omega\rangle  \cr 
~=~&e\left(\delta^4(x\!-\!x_2)-\delta^4(x\!-\!x_1)\right)
\langle \Omega |T\psi(x_1) \bar{\psi}(x_2)|\Omega\rangle\end{align}  
\tag{PS9.103}$$
is derived via the Schwinger-Dysons equations for a global symmetry, cf. e.g. Ref. [PS]. In detail, the change of integration variables in the path integral is a local gauge transformation of the matter fields without transforming the gauge field $A_{\mu}$. This transformation is just a global symmetry of the action. [The derivation simplifies because gauge-fixing terms do conveniently not transform. This makes this method a favourite of introductory QED textbooks, which often suppresses the role of gauge-fixing.] See also Ref. [W] for a derivation of the WTI
$$\begin{align} & \color{red}{Z_2} \frac{\partial}{\partial x^{\mu}} T\left\{ J^{\mu}(x)\Psi_n(y) \bar{\Psi}_m(z)\right\} \cr ~=~&e\left(\delta^4(x\!-\!z)-\delta^4(x\!-\!y)\right)
T\left\{ \Psi_n(y) \bar{\Psi}_m(z)\right\} \end{align}
\tag{W10.4.24}$$
using the operator formalism. Combined with eq. (SD'), we get
$$ \begin{align}
&\color{red}{\frac{Z_2}{Z_1}}\partial_{\mu}\underbrace{\left[\delta^{\mu}_{\nu}\partial^2 +\left(\frac{1}{\xi}-1\right)\partial^{\mu}\partial_{\nu} \right]}_{\text{inverse free propagator}} \langle \Omega | T  A^{\nu}(x)\psi(x_1) \bar{\psi}(x_2)|\Omega\rangle\cr
~=~& \mp e\left(\delta^4(x\!-\!x_2)-\delta^4(x\!-\!x_1)\right)
\langle \Omega |T\psi(x_1) \bar{\psi}(x_2)|\Omega\rangle
.\end{align}\tag{SD"}$$


*On the other hand, the WTI for proper diagrams
$$ \color{red}{\frac{Z_2}{Z_1}} q^{\mu} \Gamma_{\mu}(p,q,p+q)
~=~ e\left(S^{-1}(p+q)- S^{-1}(p) \right)\tag{R7.111}, $$
$$ \color{red}{\frac{Z_2}{Z_1}} (p^{\prime}-p)_{\mu} \Gamma^{\mu}(p^{\prime},p)
~=~ e\left(S^{-1}(p^{\prime})- S^{-1}(p)\right) \tag{B8.1.89}, $$
is derived via local gauge transformations in the path integral, cf. e.g. Refs. [R] and [B]. Equivalently, this WTI can be derived using BRST transformations, and an abelian Zinn-Justin equation.
Now let us compare the two methods.

*

*On one hand, the correlation functions in eq. (PS9.103) can be identified with non-amputated connected correlation functions. Diagrammatically, the non-amputated fermion/matter current $j^{\mu}= e:\bar{\psi}\gamma^{\mu}\psi:$ can be viewed as an amputated/striped free photon leg, cf. e.g. Ref. [K].


*On the other hand, eqs. (R7.111) and (B.1.89) use amputated proper 1PI 3-vertex $\Gamma_{\mu}$. Here $S$ denotes the dressed/renormalized fermion propagator.
The connected correlation functions and the proper 1PI correlation functions
are related via a Legendre transformation, cf. e.g. this Phys.SE post. In practice this means sewing a connected propagator on each external amputated leg of the proper 1PI diagram.
However, it turns out that it doesn't matter whether we amputate the connected WTI diagram with a connected photon propagator $G_c$ or free photon propagator $G_{(0)}$. This is because of the relation
$$ k_{\mu} G_c^{-1}(k)^{\mu}{}_{\nu}~G^{\nu\sigma}_{(0)}(k) ~=~k^{\sigma} \tag{B8.1.104},
$$
which, in turn, is a consequence of local gauge symmetry, see e.g. Ref. [B] and Jia Yiyang's Phys.SE answer here.
With the help of eq. (B8.1.104), the two WTI versions are identical up to trivial manipulations.
Eq. (B8.1.104) shows in particular, that gluing of the photon propagator to the WTI does not produce any renormalization $Z_3$ factor. See also Pointless' answer.
Philosophically, one may ponder the significance of the fact that method 1 uses global transformations, while method 2 uses local gauge transformations. The point is that the later method is able to trace the WTI at deeper more fundamental level of diagrams, namely the proper 1PI correlation functions, as opposed to the connected correlation functions.
Acknowledgements: We thank David Svoboda for discussions on the WTI.
References:

*

*[PS] M.E. Peskin & D.V. Schroeder, An Intro to QFT; Section 9.6. (For a diagrammatic treatment, see Section 7.4.) Note that in [PS] the gauge covariant derivative is defined with a plus: $D_{\mu}=\partial_{\mu} +ieA_{\mu}$ and $e=-|e|$, cf. e.g. eq. (PS4.5).


*[S] M. Srednicki, QFT; Chapters 67-68. Note that $e=-|e|$.


*[W] S. Weinberg, QFT, Vol. 1; Section 10.4.


*[R] L.H. Ryder, QFT; Section 7.4. Note that in [R] a factor of $e$ has been stripped from the 3-vertex function, cf. e.g. eq. (R7.109).


*[B] L.S. Brown, QFT; Section 8.1.3. Note that in [B] the gauge field $A_{\mu}$ contains a factor of $e$, cf. e.g. eq. (B8.1.90).


*[K] V. Kaplunovsky, WTI, lecture notes, 2012; p.17. The PDF file is available from the course homepage.
